Re: problem with NETCDF files produced by NCL

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Mon Mar 03 2014 - 08:37:22 MST

The TRMM V7 uses a new time unit that is not recognized by
by NCL and ncrcat. The error is actually with the UDUNITS package.
UDUNITS support says it will be fixed.

I have a new script that addresses the issue. It will be sent offline.
==============================
Details follow

Dennis,

It does appear that the UDUNITS package doesn't correctly handle a
"naked hour" in a timestamp.

I'll see about fixing this ASAP.

Thanks for sending this in.

> Hello,
> I downloaded 3-hrly TRMM netCDF files from NASA's mirador server.
> A sample file dump for 3B42.20090731.12.7A.nc looks like:
>
> %> ncdump -v time 3B42.20090731.12.7A.nc
>
> netcdf \3B42.20090731.12.7A {
> dimensions:
> time = UNLIMITED ; // (1 currently)
> longitude = 1440 ;
> latitude = 400 ;
> variables:
> double time(time) ;
> time:units = "hours since 2009-7-31 12" ; <========
>
> [SNIP a bunch of variables]
>
> // global attributes:
> :Conventions = "COARDS" ; <=========
> :calendar = "standard" ;
> :comments = "file created by grads using lats4d
> available from http://dao.gsfc.nasa.gov/software/grads/lats4d/" ;
> :model = "geos/das" ;
> :center = "gsfc" ;
> data:
>
> time = 0 ;
> }
> ============
> Each file has its own unique time base units attribute which matches
> the date in the file name.
> ============
> I then used the NCO operator 'ncrcat' to concatenate multiple files.
> The result was a file in which the time variable was not rebased
> correctly. Also, I tried using NCL and it had a problem also.
>
> Rather than go into detail here, please take a look at the
> full description of the issue at the NCO 'help' page. Thread:
>
> ncrcat: time rebase [version "4.3.9"]
>
> http://sourceforge.net/p/nco/discussion/9830/thread/56389da4/
>
> This includes Charlie Zender's response in which he states:
>
> "...A naked integer for the HH hour is not correctly recognized."
>
> "I am uncertain whether this is a problem with UDUnits
> not supporting a format it says it does, or whether TRMM
> data uses a units format that UDUnits does not explicitly support.
> I suspect the latter. If not, please forward this problem to
> Unidata for fixing."
>
> ==============
> I could make a few sample files (~14MB each) available if this
> would be of use.
>
> Regards
> Dennis Shea, NCAR/CGD

On 3/3/14, 1:56 AM, louis Vonder wrote:
> Dear NCL users,
>
> Below appears NCL script used to convert TRMM 3B42 data from HDF to NETCDF.
> This script was proposed by Denis for TRMM 3B42 V6 and was ajusted to work with TRMM 3B42 V7.
> Data are created without problem.
>
> All netcdf files created was merged to a single file by NCO using NCRCAT or with CDO using MERGETIME.
>
> After this operation I tried to compute the annual mean.
> But the file I obtained contain only missing value.
> Even using NCO (ncra) or CDO (timmean), the problem remains the same !!!
>
> I don't know the origin of the problem.
> There is anyway to solved it?
>
>
> Here the header of one file created
>
>
> netcdf \3B42_V7.1998010100 {
> dimensions:
> time = UNLIMITED ; // (1 currently)
> lat = 400 ;
> lon = 1440 ;
> variables:
> double time(time) ;
> time:calendar = "gregorian" ;
> time:units = "hours since 1998-01-01 00:00:0.0" ;
> float lat(lat) ;
> lat:units = "degrees_north" ;
> lat:long_name = "latitude" ;
> float lon(lon) ;
> lon:units = "degrees_east" ;
> lon:long_name = "longitude" ;
> float precip(time, lat, lon) ;
> precip:missing_value = -9999.f ;
> precip:avg_period = "0000-00-00 01:00:00" ;
> precip:delta_t = "0000-00-00 03:00:00" ;
> precip:_FillValue = -9999.f ;
> precip:units = "mm/h" ;
> precip:hdf_name = "precipitation" ;
>
> // global attributes:
> :creation_date = "Sat Sep 21 08:19:32 WAT 2013" ;
> :Conventions = "None" ;
> :source_html = "\n",
> "http://disc.sci.gsfc.nasa.gov/data/datapool/TRMM/01_Data_Products/02_Gridded/06_3-hour_Gpi_Cal_3B_42/index.html\n",
> "" ;
> :title = "TRMM: 3B42 V7" ;
> }
>
>
>
>
> Thank for you for your kind help.
>
>
>
>
>
>
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
>
>
> function parse_3B42_date (fHDF[1]:string)
> ; parse a string with the following structure: 3B42.060101.0.6.HDF
> ; return with a better file name==> 3B42.yyyymmddhh.6.nc
> local perc, per, fhdfc
> begin
> perc = stringtochar(".") ; dimsizes(chr)=2
> per = perc(0) ; period
>
> fhdfc = stringtochar( fHDF )
> iper = ind(fhdfc.eq.per)
>
>
>
> yymmdd = fhdfc(iper(0)+1 : iper(1)-1)
> hh = stringtointeger( chartostring(fhdfc(iper(1)+1 : iper(2)-1)))
>
>
> yyyy = stringtointeger( chartostring(yymmdd(0:3)))
> mm = stringtointeger( chartostring(yymmdd(4:5)))
> dd = stringtointeger( chartostring(yymmdd(6:7)))
>
>
>
> date = new ( 4, "integer", "No_FillValue")
>
>
>
> date(0) = yyyy
> date(1) = mm
> date(2) = dd
> date(3) = hh
>
>
>
>
> return(date)
> end
>
> begin ; begin not required
>
> ;*****************************************************
> ; Create TRMM coordinate variables
> ;*****************************************************
> nlat = 400
> mlon = 1440
> lat = ispan(0,nlat-1,1)*0.25 - 49.875
> lon = ispan(0,mlon-1,1)*0.25 - 179.875
>
> lat!0 = "lat"
> lat&lat = lat
> lat@long_name = "latitude"
> lat@units = "degrees_north"
>
> lon!0 = "lon"
> lon&lon = lon
> lon@long_name = "longitude"
> lon@units = "degrees_east"
> ; miscellaneous
> tunits = "hours since 1998-01-01 00:00:0.0"
> minute = 0
> second = 0.d0
> nline = inttochar(10)
>
> ;*****************************************************
> ; get all 3B42 HDF file names
> ; As returned by "ls", these will not be in chronological order
> ; files were poorly named
> ;*****************************************************
> diri = "./"
> fili = systemfunc("cd "+diri+" ; ls 3B42.*HDF")
> print(fili) ; these may not be in chronological order
> nfil = dimsizes( fili )
> ;*****************************************************
> ; loop over each file and create netCDF
> ;*****************************************************
> do n=0,nfil-1
>
>
> ymdh = parse_3B42_date ( fili(n) )
>
> print(ymdh)
>
>
> f = addfile ( fili(n), "r")
>
>
> PP = f->precipitation
>
>
> P = new((/1, mlon, nlat/), "float")
>
> P(0, :, : ) = PP
>
> P!0 = "time" ; rename to
> P!1 = "lon" ; match model names
> P!2 = "lat" ; convenience only
>
> p = P(time|:, lat|:, lon|:) ; reorder to match model
>
>
>
> ; ADD META DATA
> p@units = "mm/h" ; ?
> p@delta_t = "0000-00-00 03:00:00"
> p@avg_period = "0000-00-00 01:00:00"
> p@_FillValue = -9999.
> p@missing_value = p@_FillValue
>
> ; ADD COORDINATE VARIABLES
> p&lat = lat
> p&lon = lon
>
>
> time = ut_inv_calendar(ymdh(0),ymdh(1),ymdh(2),ymdh(3) \
> ,minute,second,tunits, 0)
> time@calendar = "gregorian"
> time!0 = "time"
> p&time = time
>
>
> ;printVarSummary(p)
>
> ; ADD AUXILARY TIME VARIABLES
> date = ymdh(0)*10000+ ymdh(1)*100 + ymdh(2)
> date!0 = "time"
> date@long_name = "current date"
> date@units = "yyyymmdd"
>
> datesec = ymdh(3)*3600
> datesec!0 = "time"
> datesec@long_name = "current seconds of current date"
> datesec@units = "seconds"
>
> yyyymmddhh = date*100 + ymdh(3)
>
> ;===================================================================
> ; WRITE OUTPUT netCDF: Use efficient method [ie: classic method]
> ;===================================================================
>
>
> ncDir = "./"
> ncName = "3B42_V7."+yyyymmddhh+".nc"
> system ("'rm' -f "+ncDir+ncName)
> fnc = addfile (ncDir+ncName, "c")
>
> ;===================================================================
> ; explicitly declare file definition mode. Improve efficiency.
> ;===================================================================
> setfileoption(fnc, "DefineMode", True)
>
> ;===================================================================
> ; create global attributes of the file
> ;===================================================================
> filAtt = True
>
> filAtt@title = "TRMM: 3B42 V7"
> filAtt@source_html = nline + \
> "http://disc.sci.gsfc.nasa.gov/data/datapool/TRMM/01_Data_Products/02_Gridded/06_3-hour_Gpi_Cal_3B_42/index.html" + nline
> filAtt@Conventions = "None"
> filAtt@creation_date = systemfunc("date")
>
> fileattdef( fnc, filAtt ) ; copy file attributes
> ;===================================================================
> ; predefine the coordinate variables and their dimensionality
> ; note: to get an UNLIMITED record dimension, we set the dimensionality
> ; to -1 and set the unlimited array to True.
> ;===================================================================
> dimNames = (/"time", "lat", "lon"/)
> dimSizes = (/-1 , nlat, mlon/)
> dimUnlim = (/ True, False, False/)
> filedimdef(fnc, dimNames, dimSizes, dimUnlim)
>
>
> ;===================================================================
> ; predefine the the dimensionality of the variables to be written out
> ; filevardef(fnc, "time" ,typeof(time),getvardims(time))
> ;===================================================================
> filevardef(fnc, "time", typeof(time), getvardims(time))
> filevardef(fnc, "lat", typeof(lat), getvardims(lat ))
> filevardef(fnc, "lon", typeof(lon), getvardims(lon ))
> ;filevardef(fnc, "date", typeof(date), getvardims(date))
> ;filevardef(fnc, "datesec", typeof(datesec), getvardims(datesec))
> filevardef(fnc, "precip", typeof(p) , getvardims(p ))
>
> ;===================================================================
> ; Copy attributes associated with each variable to the file
> ; All attributes associated with each variable will be copied.
> ;====================================================================
> filevarattdef(fnc, "time", time) ; copy time attributes
> filevarattdef(fnc, "lat", lat) ; copy lat attributes
> filevarattdef(fnc, "lon", lon) ; copy lon attributes
> ;filevarattdef(fnc, "date", date)
> ;filevarattdef(fnc, "datesec", datesec)
> filevarattdef(fnc, "precip", p) ; copy p attributes
>
> ;===================================================================
> ; explicitly exit file definition mode. **NOT REQUIRED**
> ;===================================================================
> setfileoption(fnc, "DefineMode", False)
>
> ;===================================================================
> ; output only the data values since the dimensionality and such have
> ; been predefined. The "(/", "/)" syntax tells NCL to only output the
> ; data values to the predefined locations on the file.
> ;====================================================================
> fnc->time = (/time/)
> fnc->lat = (/lat /)
> fnc->lon = (/lon /)
> ;fnc->date = (/date/)
> ;fnc->datesec= (/datesec/)
> fnc->precip = (/p/)
>
> delete(P)
> delete(p)
> delete(time)
> delete(PP)
> end do
>
> end ;
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon Mar 3 08:37:28 2014

This archive was generated by hypermail 2.1.8 : Mon Mar 03 2014 - 14:26:18 MST