Use NetCDF
      Implicit None
      Real, Parameter :: Miss=-99., Miss_Clu=-9999.0

      Character (Len=255) :: Fname
      Real, Allocatable :: Reg2D(:,:), X2d(:,:), lon(:), lat(:)
      Integer :: Nlon, Nlat, ilon, ilat
      Logical :: status

      Fname="/home/tanpv/Progs/2021/NC/kmeans_4.nc"

      Call rd_ncdf_2d (Fname,lon,lat,Nlon,Nlat, Reg2D, "cluster")
      Where (Reg2D==Miss_Clu) Reg2D=Miss

      Allocate (X2d(Nlon,Nlat))
      X2d = Reg2D

      Call Replace (96.,115.,7.,20., 3.)  ! VN
      Call Replace (120.,128.,5.,20., 3.)   ! Phi
      Call Replace (116.,120.,4.,7., 1.)    ! Ind
      Call Replace (94.,100.,0.,6., 1.)     ! Sum
      Call Replace (100.,106.,0.,5., 1.)     ! Mas
      Call Replace (94.,146.5,-11.,4., 1.)   ! Ind
      Call Replace (120.,123.,22.,26., 2.)  ! TW
      Call Replace (95.,100.,8.,10., 3.)  ! Tha

      Where (Reg2D==1) Reg2D=11
      Where (Reg2D==3) Reg2D=12
      Where (Reg2D==4) Reg2D=13
      Where (Reg2D==2) Reg2D=14

      Where (Reg2D==11) Reg2D=1
      Where (Reg2D==12) Reg2D=2
      Where (Reg2D==13) Reg2D=3
      Where (Reg2D==14) Reg2D=4

      Fname="/home/tanpv/Progs/2021/NC/kmeans_4_modified.nc"
      Call wr_ncdf_2d (Fname,lon,lat,nlon,nlat,Reg2D,"cluster", &
           "Sub_reg","1","lon","lat", Miss) 

      CONTAINS

      Subroutine Replace (lon1,lon2,lat1,lat2, new)
         Integer :: ilon, ilat
         Real :: lon1,lon2,lat1,lat2, new

         Do ilon=1,Nlon
            Do ilat=1,Nlat
               if (lon(ilon)>=lon1.AND.lon(ilon)<=lon2.AND. &
                   lat(ilat)>=lat1.AND.lat(ilat)<=lat2) then
                   if (Reg2D(ilon,ilat)/=Miss) Reg2D(ilon,ilat)=New
               endif
            Enddo
         Enddo
      End Subroutine Replace

      Include "/home/tanpv/Mylibs/sub_rd_ncdf_Tan.F90"
      Include "/home/tanpv/Mylibs/sub_wr_ncdf_Tan.F90"

      END