HRLDAS_forcing for GLDAS NetCDF

HRLDAS/HRLDAS_forcing/run/examples/GLDAS/Makefile

@@ -0,0 +1,14 @@
+
+include ../../../../user_build_options
+
+all:	
+	$(COMPILERF90) $(NETCDFMOD) -o combine_precips_netcdf.exe $(F90FLAGS) combine_precips_netcdf.f90 $(NETCDFLIB)
+	$(COMPILERF90) $(NETCDFMOD) -o create_UV_netcdf.exe $(F90FLAGS) create_UV_netcdf.f90 $(NETCDFLIB)
+clean:
+	$(RM) *.o *~ *.exe 
+#
+
+
+
+
+

HRLDAS/HRLDAS_forcing/run/examples/GLDAS/combine_precips_netcdf.f90

@@ -0,0 +1,107 @@
+use netcdf
+implicit none
+integer            :: ncid,ierr
+integer, parameter :: NDIMS = 3, NRECS = 1
+integer, parameter :: NLATS = 600, NLONS = 1440
+character (len = *), parameter :: LAT_NAME = "latitude"
+character (len = *), parameter :: LON_NAME = "longitude"
+character (len = *), parameter :: REC_NAME = "time"
+integer :: lon_dimid, lat_dimid, rec_dimid
+
+integer :: start(NDIMS), count(NDIMS)
+real :: lats(NLATS), lons(NLONS)
+integer :: lon_varid, lat_varid
+
+character (len = *), parameter :: RAIN_NAME = "Rainf_tavg" 
+character (len = *), parameter :: SNOW_NAME = "Snowf_tavg"
+character (len = *), parameter :: PREC_NAME = "Precf_tavg"
+integer :: rain_varid,snow_varid,prec_varid
+integer :: dimids(NDIMS)
+
+! We recommend that each variable carry a "units" attribute.
+character (len = *), parameter :: UNITS = "units"
+character (len = *), parameter :: RAIN_UNITS = "kg m-2 s-1"
+character (len = *), parameter :: SNOW_UNITS = "kg m-2 s-1"
+character (len = *), parameter :: PREC_UNITS = "kg m-2 s-1"
+character (len = *), parameter :: LAT_UNITS = "degrees_north"
+character (len = *), parameter :: LON_UNITS = "degrees_east"
+real    :: rain_in(NLONS, NLATS), snow_in(NLONS, NLATS), prec_out(NLONS, NLATS)
+
+! Loop indices
+integer :: lat, lon, rec, i
+character*100   :: infile1, infile2, outfile  
+
+call getarg(1,infile1)
+call getarg(2,infile2)
+call getarg(3,outfile)
+
+! READ IN RAIN VARIABLE
+print *,"Reading in infile: ",infile1
+ierr = nf90_open(infile1,nf90_nowrite, ncid)
+ierr = nf90_inq_varid(ncid,RAIN_NAME,rain_varid)
+ierr = nf90_get_var(ncid,rain_varid,rain_in)
+ierr = nf90_inq_varid(ncid,"lat",lat_varid)
+ierr = nf90_inq_varid(ncid,"lon",lon_varid)
+ierr = nf90_get_var(ncid,lat_varid,lats)
+ierr = nf90_get_var(ncid,lon_varid,lons)
+ierr = nf90_close(ncid)
+! READ IN SNOW VARIABLE
+print *,"Reading in infile: ",infile2
+ierr = nf90_open(infile2,nf90_nowrite, ncid)
+ierr = nf90_inq_varid(ncid,SNOW_NAME,snow_varid)
+ierr = nf90_get_var(ncid,snow_varid,snow_in)
+ierr = nf90_get_var(ncid,lat_varid,lats)
+ierr = nf90_get_var(ncid,lon_varid,lons)
+ierr = nf90_close(ncid)
+! WRITE WIND SPEED to U component
+prec_out = rain_in + snow_in
+where (prec_out<0) prec_out = -9999 
+! Create the file. 
+ierr = nf90_create(outfile, nf90_clobber, ncid) 
+! Define the dimensions. The record dimension is defined to have
+! unlimited length - it can grow as needed. In this example it is
+! the time dimension.
+ierr = nf90_def_dim(ncid, LAT_NAME, NLATS, lat_dimid) 
+ierr = nf90_def_dim(ncid, LON_NAME, NLONS, lon_dimid)
+ierr = nf90_def_dim(ncid, REC_NAME, NF90_UNLIMITED, rec_dimid)
+! Define the coordinate variables. We will only define coordinate
+! variables for lat and lon.  Ordinarily we would need to provide
+! an array of dimension IDs for each variable's dimensions, but
+! since coordinate variables only have one dimension, we can
+! simply provide the address of that dimension ID (lat_dimid) and
+! similarly for (lon_dimid).
+ierr = nf90_def_var(ncid, LAT_NAME, NF90_REAL, lat_dimid, lat_varid) 
+ierr = nf90_def_var(ncid, LON_NAME, NF90_REAL, lon_dimid, lon_varid)
+! Assign units attributes to coordinate variables.
+ierr = nf90_put_att(ncid,lat_varid,UNITS,LAT_UNITS)
+ierr = nf90_put_att(ncid,lon_varid,UNITS,LON_UNITS)
+! The dimids array is used to pass the dimids of the dimensions of
+! the netCDF variables. Both of the netCDF variables we are creating
+! share the same four dimensions. In Fortran, the unlimited
+! dimension must come last on the list of dimids.
+dimids = (/ lon_dimid, lat_dimid, rec_dimid /)
+! Define the netCDF variables for the u variable
+ierr = nf90_def_var(ncid, PREC_NAME, NF90_REAL, dimids, prec_varid)
+! Assign units attributes to the netCDF variables.
+ierr = nf90_put_att(ncid, prec_varid, UNITS, PREC_UNITS) 
+ierr = nf90_put_att(ncid, prec_varid, "missing_value",-9999)
+! End define mode
+ierr = nf90_enddef(ncid)
+! Write the coordinate variable data. This will put the latitudes
+! and longitudes of our data grid into the netCDF file.
+ierr = nf90_put_var(ncid, lat_varid, lats) 
+ierr = nf90_put_var(ncid, lon_varid, lons) 
+! These settings tell netcdf to write one timestep of data. (The
+! setting of start(4) inside the loop below tells netCDF which
+! timestep to write.)
+count = (/ NLONS, NLATS, 1 /)
+start = (/ 1, 1, 1 /)
+do rec = 1, NRECS
+   start(3) = rec
+   ierr =  nf90_put_var(ncid, prec_varid, prec_out, start = start, &
+                            count = count) 
+end do
+ierr = nf90_close(ncid)
+
+print *, "*** SUCCESS combining precipitation file ", outfile
+end

HRLDAS/HRLDAS_forcing/run/examples/GLDAS/combine_precips_netcdf.perl

@@ -0,0 +1,57 @@
+#!/usr/bin/perl
+
+$filename = "combine_precips_netcdf.exe";
+system ("make clean");
+system ("make");
+@nums = ("00","01","02","03","04","05","06","07","08","09", 
+         "10","11","12","13","14","15","16","17","18","19", 
+	 "20","21","22","23","24","25","26","27","28","29", 
+	 "30","31");
+
+@yrs = ("00");
+
+$day_start = 275;
+$day_end = 275;
+
+@hrs = ("00","03","06","09","12","15","18","21");
+
+@noleap_days = (0,31,59,90,120,151,181,212,243,273,304,334,365);
+@leap_days   = (0,31,60,91,121,152,182,213,244,274,305,335,366);
+
+$data_dir = "/glade/work/zhezhang/GLDAS/extracted/";
+$results_dir = "/glade/work/zhezhang/GLDAS/extracted/";
+
+for $yy (@yrs)
+ {
+# This will be the jday time loop
+
+for($julday=$day_start;$julday<=$day_end;$julday++)
+
+ { 
+
+ # This little section finds the text month and day
+
+ @modays = @noleap_days;
+ if($yy == "00" || $yy == "04" || $yy == "08" || $yy == "12") {@modays = @leap_days}
+
+ for($mo=1;$mo<=12;$mo++)
+  {
+    if($julday>$modays[$mo-1] && $julday<=$modays[$mo]) 
+     {
+       $mon = $mo;
+       $day = $julday - $modays[$mo-1];
+     }
+  }
+
+for $hr (@hrs)
+ {
+  $file1_in  =    "$data_dir/Rainf/GLDAS_Rainf_tavg_20$yy$nums[$mon]$nums[$day]$hr";
+  $file2_in  =    "$data_dir/Snowf/GLDAS_Snowf_tavg_20$yy$nums[$mon]$nums[$day]$hr";
+  $file_out  =  "$results_dir/Precip/GLDAS_Precip_20$yy$nums[$mon]$nums[$day]$hr";
+  print ("$file_out \n");
+  system ("./$filename $file1_in $file2_in $file_out");
+ }
+ }
+ }
+
+system("rm $filename");

HRLDAS/HRLDAS_forcing/run/examples/GLDAS/create_UV_netcdf.f90

@@ -0,0 +1,126 @@
+use netcdf
+implicit none
+integer            :: ncid,ierr
+integer, parameter :: NDIMS = 3, NRECS = 1
+integer, parameter :: NLATS = 600, NLONS = 1440
+character (len = *), parameter :: LAT_NAME = "latitude"
+character (len = *), parameter :: LON_NAME = "longitude"
+character (len = *), parameter :: REC_NAME = "time"
+integer :: lon_dimid, lat_dimid, rec_dimid
+
+integer :: start(NDIMS), count(NDIMS)
+real :: lats(NLATS), lons(NLONS)
+integer :: lon_varid, lat_varid
+
+character (len = *), parameter :: WIND_NAME = "Wind_f_inst" 
+character (len = *), parameter :: U_NAME    = "U_f_inst"
+character (len = *), parameter :: V_NAME    = "V_f_inst"
+integer :: wind_varid,u_varid,v_varid
+integer :: dimids(NDIMS)
+
+! We recommend that each variable carry a "units" attribute.
+character (len = *), parameter :: UNITS = "units"
+character (len = *), parameter :: WIND_UNITS = "m s-1"
+character (len = *), parameter :: LAT_UNITS = "degrees_north"
+character (len = *), parameter :: LON_UNITS = "degrees_east"
+real    :: wind_in(NLONS, NLATS), u_out(NLONS, NLATS), v_out(NLONS, NLATS)
+
+! Loop indices
+integer :: lat, lon, rec, i
+character*100   :: infile, outfile1, outfile2  
+
+call getarg(1,infile)
+call getarg(2,outfile1)
+call getarg(3,outfile2)
+
+! READ IN DUMMY VARIABLE
+print *,"Working on infile: ",infile
+ierr = nf90_open(infile,nf90_nowrite, ncid)
+ierr = nf90_inq_varid(ncid,WIND_NAME,wind_varid)
+ierr = nf90_get_var(ncid,wind_varid,wind_in)
+ierr = nf90_inq_varid(ncid,"lat",lat_varid)
+ierr = nf90_inq_varid(ncid,"lon",lon_varid)
+ierr = nf90_get_var(ncid,lat_varid,lats)
+ierr = nf90_get_var(ncid,lon_varid,lons)
+ierr = nf90_close(ncid)
+! WRITE WIND SPEED to U component
+u_out = wind_in
+v_out = 0.0
+! Create the file. 
+ierr = nf90_create(outfile1, nf90_clobber, ncid) 
+! Define the dimensions. The record dimension is defined to have
+! unlimited length - it can grow as needed. In this example it is
+! the time dimension.
+ierr = nf90_def_dim(ncid, LAT_NAME, NLATS, lat_dimid) 
+ierr = nf90_def_dim(ncid, LON_NAME, NLONS, lon_dimid)
+ierr = nf90_def_dim(ncid, REC_NAME, NF90_UNLIMITED, rec_dimid)
+! Define the coordinate variables. We will only define coordinate
+! variables for lat and lon.  Ordinarily we would need to provide
+! an array of dimension IDs for each variable's dimensions, but
+! since coordinate variables only have one dimension, we can
+! simply provide the address of that dimension ID (lat_dimid) and
+! similarly for (lon_dimid).
+ierr = nf90_def_var(ncid, LAT_NAME, NF90_REAL, lat_dimid, lat_varid) 
+ierr = nf90_def_var(ncid, LON_NAME, NF90_REAL, lon_dimid, lon_varid)
+! Assign units attributes to coordinate variables.
+ierr = nf90_put_att(ncid,lat_varid,UNITS,LAT_UNITS)
+ierr = nf90_put_att(ncid,lon_varid,UNITS,LON_UNITS)
+! The dimids array is used to pass the dimids of the dimensions of
+! the netCDF variables. Both of the netCDF variables we are creating
+! share the same four dimensions. In Fortran, the unlimited
+! dimension must come last on the list of dimids.
+dimids = (/ lon_dimid, lat_dimid, rec_dimid /)
+! Define the netCDF variables for the u variable
+ierr = nf90_def_var(ncid, U_NAME, NF90_REAL, dimids, u_varid)
+! Assign units attributes to the netCDF variables.
+ierr = nf90_put_att(ncid, u_varid, UNITS, WIND_UNITS) 
+ierr = nf90_put_att(ncid, u_varid, "missing_value",-9999)
+! End define mode
+ierr = nf90_enddef(ncid)
+! Write the coordinate variable data. This will put the latitudes
+! and longitudes of our data grid into the netCDF file.
+ierr = nf90_put_var(ncid, lat_varid, lats) 
+ierr = nf90_put_var(ncid, lon_varid, lons) 
+! These settings tell netcdf to write one timestep of data. (The
+! setting of start(4) inside the loop below tells netCDF which
+! timestep to write.)
+count = (/ NLONS, NLATS, 1 /)
+start = (/ 1, 1, 1 /)
+do rec = 1, NRECS
+   start(3) = rec
+   ierr =  nf90_put_var(ncid, u_varid, u_out, start = start, &
+                            count = count) 
+end do
+ierr = nf90_close(ncid)
+
+print *, "*** SUCCESS writing example file ", outfile1
+
+! Do the same for V component
+! Create the file. 
+ierr = nf90_create(outfile2, nf90_clobber, ncid)
+! the time dimension.
+ierr = nf90_def_dim(ncid, LAT_NAME, NLATS, lat_dimid)
+ierr = nf90_def_dim(ncid, LON_NAME, NLONS, lon_dimid)
+ierr = nf90_def_dim(ncid, REC_NAME, NF90_UNLIMITED, rec_dimid)
+ierr = nf90_def_var(ncid, LAT_NAME, NF90_REAL, lat_dimid, lat_varid)
+ierr = nf90_def_var(ncid, LON_NAME, NF90_REAL, lon_dimid, lon_varid)
+ierr = nf90_put_att(ncid,lat_varid,UNITS,LAT_UNITS)
+ierr = nf90_put_att(ncid,lon_varid,UNITS,LON_UNITS)
+dimids = (/ lon_dimid, lat_dimid, rec_dimid /)
+ierr = nf90_def_var(ncid, V_NAME, NF90_REAL, dimids, v_varid)
+ierr = nf90_put_att(ncid, v_varid, UNITS, WIND_UNITS)
+ierr = nf90_put_att(ncid, v_varid, "missing_value",-9999)
+ierr = nf90_enddef(ncid)
+ierr = nf90_put_var(ncid, lat_varid, lats)
+ierr = nf90_put_var(ncid, lon_varid, lons)
+count = (/ NLONS, NLATS, 1 /)
+start = (/ 1, 1, 1 /)
+do rec = 1, NRECS
+   start(3) = rec
+   ierr =  nf90_put_var(ncid, v_varid, v_out, start = start, &
+                            count = count)
+end do
+ierr = nf90_close(ncid)
+
+print *, "*** SUCCESS writing example file ", outfile2
+end

HRLDAS/HRLDAS_forcing/run/examples/GLDAS/create_UV_netcdf.perl

@@ -0,0 +1,62 @@
+#!/usr/bin/perl
+
+$filename = "create_UV_netcdf.exe";
+system ("make clean");
+system ("make");
+@nums = ("00","01","02","03","04","05","06","07","08","09", 
+         "10","11","12","13","14","15","16","17","18","19", 
+	 "20","21","22","23","24","25","26","27","28","29", 
+	 "30","31");
+
+@yrs = ("00");
+
+$day_start = 275;
+$day_end = 275;
+
+@hrs = ("00","03","06","09","12","15","18","21");
+
+@noleap_days = (0,31,59,90,120,151,181,212,243,273,304,334,365);
+@leap_days   = (0,31,60,91,121,152,182,213,244,274,305,335,366);
+
+$data_dir = "/glade/work/zhezhang/GLDAS/extracted";
+$results_dir = "/glade/work/zhezhang/GLDAS/extracted";
+
+for $yy (@yrs)
+ {
+
+# This will be the jday time loop
+
+for($julday=$day_start;$julday<=$day_end;$julday++)
+
+ { 
+
+ # This little section finds the text month and day
+
+ @modays = @noleap_days;
+ if($yy == "00" || $yy == "04" || $yy == "08" || $yy == "12") {@modays = @leap_days}
+
+ for($mo=1;$mo<=12;$mo++)
+  {
+    if($julday>$modays[$mo-1] && $julday<=$modays[$mo]) 
+     {
+       $mon = $mo;
+       $day = $julday - $modays[$mo-1];
+     }
+  }
+
+for $hr (@hrs)
+ {
+
+  $file_in  =    "$data_dir/Wind/GLDAS_Wind_f_inst_20$yy$nums[$mon]$nums[$day]$hr";
+  $file1_out =  " $results_dir/U/GLDAS_U_f_inst_20$yy$nums[$mon]$nums[$day]$hr";
+  $file2_out =  " $results_dir/V/GLDAS_V_f_inst_20$yy$nums[$mon]$nums[$day]$hr";
+
+  print ("$file_in \n");
+  print ("$file1_out \n");
+  print ("$file2_out \n");
+  system ("./$filename $file_in $file1_out $file2_out");
+ }
+ }
+ } # End of outer time loop
+
+system("rm $filename");