Use NetCDF Real, Parameter :: Miss=-99.0 Integer, Parameter :: IMiss=-99 ! , Nmem=11 !Character (Len=5), Parameter :: memdir(Nmem)=(/"mem1","mem2","mem3","mem4","mem5","mem6","mem7","mem8","mem9","mem10","cf"/) Integer, Allocatable :: time(:) Integer :: Nst, Ntime, i, LenSt, sys Character (Len=255) :: Fname, t_unit, DirIn, Str, time_unit, Fname1, Fname2 Character (Len=255), Allocatable :: Files(:) Character (Len=12) :: date0, timeres, idate0, idate1 Character (Len=1) :: Ch Integer :: yy1,mm1,dd1,hh1,yy2,mm2,dd2,hh2 Real*8 :: factor, offset, Missing Real, Allocatable :: Prec(:,:), lon(:), lat(:), Reg(:) Character (Len=20), Allocatable :: Stname(:) Integer, Allocatable, Dimension (:) :: dates ! (Nfcs) Integer :: Nlon, Nlat, Nfcs, date_fcs, Nlon0, Nlat0 ! Nfcs = So lan du bao, date_fcs = yyyymmdd cua ngay DB Integer :: imem, ifcs, iy, im, id DirIn="/home/tanpv/Progs/S2S/Data_Stn/" Fname=trim(DirIn)//"R_Daily.nc" Call Read_Prec (Fname) !Fname=trim(DirIn)//"R_Stn_Daily.nc" Fname="/work/projects/G09/S2S/OBS_Stn/R_Stn_Daily.nc" Call Write_obs_to_file (Fname) CONTAINS !---------------------------------------------------------------------------------------------- Subroutine Write_obs_to_file (Fname1) ! ,X4d,X1d,Nlon,Nlat,Nt,Nfcs) Character (Len=*) :: Fname1 Integer, Allocatable :: StnID(:) Integer :: it Character (Len=50) :: t_unit, var_ln, var_unit, lon_nm, lat_nm, z_nm, t_nm Allocate (StnID(Nst)) print*,Nst,Ntime Do it=1,Nst; StnID(it)=it; Enddo var_ln="Daily precipitation Obs"; var_unit="mm/day"; lon_nm="stnID"; lat_nm="date" Call wr_ncdf_2d_new (Fname1,StnID,dates,Nst,Ntime,Prec,"pr", var_ln, var_unit, lon_nm, lat_nm, Miss) Call wr_ncdf_1d_apd_char_new (Fname1,Stname,"stname", Nst, "Name of Stations", "-", lon_nm, Miss) Deallocate (StnID) End Subroutine Write_obs_to_file !---------------------------------------------------------------------------------------------- Subroutine Read_Prec (Fname) Character (Len=*) :: Fname Real, Allocatable :: tmp(:,:,:,:), stnID(:), xyear(:), xmon(:), xday(:), xlon(:), xlat(:), xreg(:) Real, Allocatable :: tmp1(:,:) ! Prec (ist, ymd) since 2000 Integer, Allocatable :: ymd(:) Character (Len=20), Allocatable :: Stnm(:) Integer :: Nt, it, i, j, ilon,ilat, Nyr, Nmon, Nday, ist, iyr, imon, iday Integer :: Yr1, Yr2, Nyear, date, DayInMon !Call rd_ncdf_3d (Fname,lon,lat,xtime,nlon,nlat,Nt, Prec,"tp") Call rd_ncdf_4d (Fname,stnID,xyear,xmon,xday, Nst, Nyr, Nmon, Nday, tmp, "R") Deallocate (stnID) Call rd_ncdf_1d (Fname, stnID, Nst, xlon, "Lon") Deallocate (stnID) Call rd_ncdf_1d (Fname, stnID, Nst, xlat, "Lat") Deallocate (stnID) Call rd_ncdf_1d (Fname, stnID, Nst, xreg, "reg") Deallocate (stnID) Call rd_ncdf_1d_Char (Fname, stnID, Nst, Stnm, "Stname") Yr1=2000; Yr2=2020; Nyear=Yr2-Yr1+1 Nt=Nyear*Nmon*Nday Allocate (tmp1(Nst,Nt),ymd(Nt)) tmp1=Miss; ymd=IMiss Do ist=1,Nst it=0 Do i=1,Nyr iyr=Int(xyear(i)) if (iyr>=Yr1.AND.iyr<=Yr2) then Do imon=1,Nmon DayInMon=DayOfMon(imon,iyr) Do iday=1,DayInMon it=it+1 tmp1(ist,it)=tmp(ist,i,imon,iday) ymd(it)=iyr*10000+imon*100+iday !print*,it,ymd(it),tmp1(ist,it) Enddo Enddo endif Enddo Enddo print*,Nt,it, Nst Nt=it Allocate (Prec(Nst,Nt),dates(Nt),Lon(Nst),Lat(Nst),Reg(Nst),Stname(Nst)) Prec(:,1:Nt)=tmp1(:,1:Nt); dates(1:Nt)=ymd(:) Lon(:)=xlon(:); Lat(:)=xlat(:); Reg(:)=xreg(:); Stname(:)=Stnm(:) Ntime=Nt Deallocate (tmp,tmp1,ymd,stnID, xyear, xmon, xday, xlon, xlat, xreg) End Subroutine Read_Prec !---------------------------------------------------------------------------------- SUBROUTINE wr_ncdf_1d_apd_char_new (filenm,var,var_nm, nlon, var_long_name, var_unit, lon_nm, flag) character (len = *), intent (in) :: filenm, var_nm, var_long_name, var_unit character (len = *) :: lon_nm integer, parameter :: NDIMS = 1 ! 2 integer :: ncid, var_id, charID, dims_id(NDIMS) integer, intent (in) :: nlon character (len = *), intent (in) :: var(nlon) real, intent (in) :: flag integer :: istat ! nx print*," Writing append to netCDF file: ",filenm,"..." call hulk ( nf90_open (filenm, nf90_write, ncid) ) call hulk ( nf90_inq_dimid (ncid, lon_nm, dims_id(1)) ) istat = ( nf90_inq_dimid (ncid, "charlen", charID )) if (istat .ne. nf90_noerr) then print *,"new dimension" call hulk ( nf90_redef (ncid) ) call hulk (nf90_def_dim (ncID, "charlen", 20, charID) ) else call hulk ( nf90_redef (ncid) ) endif call hulk ( nf90_def_var (ncid, var_nm, nf90_char, (/charID,dims_id(1:1)/), var_id) ) call hulk ( nf90_put_att (ncid, var_id, "long_name",trim(var_long_name))) call hulk ( nf90_put_att (ncid, var_id, "units",trim(var_unit))) call hulk ( nf90_put_att (ncid, var_id, "missing_value", flag) ) call hulk ( nf90_enddef (ncid) ) call hulk ( nf90_put_var (ncid, var_id, var) ) call hulk ( nf90_close (ncid) ) print*," Completed writing append to file ",filenm return END SUBROUTINE wr_ncdf_1d_apd_char_new !---------------------------------------------------------------------------------- SUBROUTINE wr_ncdf_2d_new (filenm,lon,lat,nlon,nlat,var,var_nm, var_long_name, var_unit, lon_nm, lat_nm, flag) character (len = *), intent (in) :: filenm, var_nm, var_long_name, var_unit ! character (len = *), parameter :: lon_nm = "lon", lat_nm = "lat" character (len = *) :: lon_nm, lat_nm ! Names of Dim_1, Dim_2 integer, parameter :: NDIMS = 2 integer :: ncid, var_id, dims_id(NDIMS) integer :: lon_id, lat_id integer, intent (in) :: nlon, nlat real, dimension (nlon,nlat), intent (in) :: var Integer, intent (in) :: lon(nlon), lat(nlat) Real :: flag print*," Writing netCDF file: ",filenm,"..." call hulk ( nf90_create (filenm, nf90_write, ncid) ) call hulk ( nf90_def_dim (ncid, lon_nm, nlon, dims_id(1)) ) call hulk ( nf90_def_dim (ncid, lat_nm, nlat, dims_id(2)) ) call hulk ( nf90_def_var (ncid, lon_nm, nf90_Int, dims_id(1), lon_id) ) call hulk ( nf90_def_var (ncid, lat_nm, nf90_Int, dims_id(2), lat_id) ) call hulk ( nf90_def_var (ncid, var_nm, nf90_float, dims_id(1:2), var_id) ) call hulk ( nf90_put_att (ncid, lon_id, "Station_ID", "ID") ) call hulk ( nf90_put_att (ncid, lon_id, "units", "1") ) call hulk ( nf90_put_att (ncid, lat_id, "Date", "Date") ) call hulk ( nf90_put_att (ncid, lat_id, "units", "yyyymmdd") ) call hulk ( nf90_put_att (ncid, var_id, "long_name", trim(var_long_name))) !"What the hell (!@.@) !!!") ) call hulk ( nf90_put_att (ncid, var_id, "units", trim(var_unit))) ! "OMG (@.@!)") ) call hulk ( nf90_put_att (ncid, var_id, "missing_value", flag) ) call hulk ( nf90_enddef(ncid) ) call hulk ( nf90_put_var (ncid, lon_id, lon) ) call hulk ( nf90_put_var (ncid, lat_id, lat) ) call hulk ( nf90_put_var (ncid, var_id, var) ) call hulk ( nf90_close (ncid) ) print*," Completed writing file ",filenm return END SUBROUTINE wr_ncdf_2d_new !---------------------------------------------------------------------------------- Subroutine wr_ncdf_4d_new (filenm,lon,lat,z,t,nlon,nlat,nz,nt,var,var_nm, & t_unit, var_long_name, var_unit, lon_nm, lat_nm, z_nm, t_nm, flag) character (len = *), intent (in) :: filenm, var_nm, t_unit, var_long_name, var_unit character (len = *) :: lon_nm, lat_nm, z_nm, t_nm integer, parameter :: NDIMS = 4 integer :: ncid, var_id, dims_id(NDIMS) integer :: lon_id, lat_id, z_id, t_id integer, intent (in) :: nlon, nlat, nz, nt real, dimension (nlon,nlat,nz,nt), intent (in) :: var real, intent (in) :: lon(nlon), lat(nlat), flag !, z(nz), t(nt), flag integer, intent (in) :: z(nz), t(nt) print*," Writing netCDF file: ",filenm,"..." call hulk ( nf90_create (filenm, nf90_write, ncid) ) call hulk ( nf90_def_dim (ncid, lon_nm, nlon, dims_id(1)) ) call hulk ( nf90_def_dim (ncid, lat_nm, nlat, dims_id(2)) ) call hulk ( nf90_def_dim (ncid, z_nm, nz, dims_id(3)) ) call hulk ( nf90_def_dim (ncid, t_nm, nt, dims_id(4)) ) call hulk ( nf90_def_var (ncid,lon_nm,nf90_float,dims_id(1),lon_id) ) call hulk ( nf90_def_var (ncid,lat_nm,nf90_float,dims_id(2),lat_id) ) call hulk ( nf90_def_var (ncid,z_nm,NF90_INT,dims_id(3),z_id) ) call hulk ( nf90_def_var (ncid,t_nm,NF90_INT,dims_id(4),t_id) ) call hulk ( nf90_def_var (ncid,var_nm,nf90_float,dims_id(1:4),var_id) ) call hulk ( nf90_put_att (ncid,lon_id,"long_name","Longitude") ) call hulk ( nf90_put_att (ncid,lon_id,"units","degree_east") ) call hulk ( nf90_put_att (ncid,lon_id,"axis","X") ) call hulk ( nf90_put_att (ncid,lat_id,"long_name","Latitude") ) call hulk ( nf90_put_att (ncid,lat_id,"units","degree_north") ) call hulk ( nf90_put_att (ncid,lat_id,"axis","Y") ) call hulk ( nf90_put_att (ncid,z_id, "long_name","z") ) call hulk ( nf90_put_att (ncid,z_id, "units","Lev") ) call hulk ( nf90_put_att (ncid,z_id, "axis","Z") ) call hulk ( nf90_put_att (ncid,t_id,"long_name","Time") ) call hulk ( nf90_put_att (ncid,t_id,"units",trim(t_unit)//" 00:00:00") ) call hulk ( nf90_put_att (ncid,t_id,"calendar","proleptic_gregorian") ) call hulk ( nf90_put_att (ncid,t_id,"axis","T") ) call hulk ( nf90_put_att (ncid,var_id,"long_name",trim(var_long_name))) call hulk ( nf90_put_att (ncid,var_id,"units", trim(var_unit))) call hulk ( nf90_put_att (ncid,var_id,"missing_value",flag) ) call hulk ( nf90_enddef (ncid) ) call hulk ( nf90_put_var (ncid,lon_id,lon) ) call hulk ( nf90_put_var (ncid,lat_id,lat) ) call hulk ( nf90_put_var (ncid,z_id, z) ) call hulk ( nf90_put_var (ncid,t_id, t) ) call hulk ( nf90_put_var (ncid,var_id,var) ) call hulk ( nf90_close (ncid) ) print*," Completed writing file ",filenm End Subroutine wr_ncdf_4d_new !---------------------------------------------------------------------------------------------- ! <-~~~=| Nf90 checking errors ... SUBROUTINE hulk (stat) use netcdf integer :: stat if ( stat .ne. nf90_noerr) then print*, trim(nf90_strerror(stat)) stop endif END SUBROUTINE hulk !---------------------------------------------------------------------------------- Include "/home/tanpv/Mylibs/DayOfMon.f90" Include "/home/tanpv/Mylibs/sub_rd_ncdf_Tan.F90" ! Include "/home/tanpv/Mylibs/sub_wr_ncdf_Tan.F90" END