Use NetCDF
Implicit None
Integer :: Nlon, Nlat, Nlon1, Nlat1, ilon, ilat, ncid
   integer, parameter :: NDIMS = 2
   integer :: var_id, dims_id(NDIMS)
   integer :: lon_id, lat_id, varid
   character (len = 15) :: lon_nm, lat_nm, var_nm

Real*8, Allocatable :: Lon(:), Lat(:)
Real, Allocatable :: Hgt(:,:)
Real*8, Allocatable :: Lon1(:), Lat1(:) 
Real, Allocatable :: Hgt1(:,:)
Real*8, Allocatable :: Lon0(:), Lat0(:)
Real, Allocatable :: Hgt0(:,:)
Character (Len=255) :: Fname
Logical :: Status
Real :: Miss
Miss=32767.

Fname="VN_Topo.nc"

status = nf90_open (trim(Fname), nf90_Nowrite,ncid)

!status =  nf90_inq_varid (ncid,var_nm,varid)
!status =  nf90_inquire_variable (ncid,varid, dimids = dims_id(1:2))
!status =  nf90_inquire_dimension (ncid, dims_id(1),name = lon_nm, len = nlon)
!status =  nf90_inquire_dimension (ncid, dims_id(2),name = lat_nm, len = nlat)

status = nf90_inq_dimid(ncid,"lat", lat_id)
status = nf90_inquire_dimension(ncid, lat_id, len=Nlat)
status = nf90_inq_dimid(ncid,"lon", lon_id)
status = nf90_inquire_dimension(ncid, lon_id, len=Nlon)

print*, "nlon, nlat, lat_nm, lon_nm, lat_id, lon_id="
print*, nlon, nlat, lat_nm, lon_nm, lat_id, lon_id
pause
ALLOCATE (lon(nlon),lat(nlat), Hgt(Nlon,Nlat))

print*,"OK"

status = nf90_inq_varid(ncid,"lat", varid)
status = nf90_get_var(ncid, varid, Lat)

print*,"Reading Lat ... OK"

status = nf90_inq_varid(ncid,"lon", varid)
status = nf90_get_var(ncid, varid, Lon)

print*,"Reading Lon ... OK"

status = nf90_inq_varid(ncid,"Band1", varid)
status = nf90_get_var(ncid, varid, Hgt)

print*, "Reading Hgt ... OK"

pause

!Call rd_ncdf_2d (trim(Fname), Lon, Lat, Nlon, Nlat, Hgt, "Band1")

Allocate (Hgt0(Nlon,Nlat),Lon0(Nlon),Lat0(Nlat))
print*, "Allocate Hgt0 ... OK"

Nlon1=0;
Do ilon=1, Nlon, 5
   Nlon1=Nlon1+1; Lon0(Nlon1)=Lon(ilon); Nlat1=0
   Do ilat=1, Nlat, 5
      Nlat1=Nlat1+1; Lat0(Nlat1)=Lat(ilat)
      Hgt0(Nlon1,Nlat1)=Hgt(ilon,ilat)
   Enddo
Enddo
Allocate(Hgt1(Nlon1,Nlat1),Lon1(Nlon1),Lat1(Nlat1))
Hgt1(1:Nlon1,1:Nlat1)=Hgt0(1:Nlon1,1:Nlat1)
Lon1(1:Nlon1)=Lon0(1:Nlon1); Lat1(1:Nlat1)=Lat0(1:Nlat1)

Fname="Vietnam_1.nc"
Call wr_ncdf_2d_dbl (trim(Fname), lon1 ,lat1 ,nlon1, nlat1, Hgt1, "Hgt", "Elevation", "m", "Lon", "Lat", Miss)

CONTAINS

Include "sub_rd_ncdf.F90"
Include "sub_wr_ncdf_Tan.F90"

END