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
Received on Mon Mar 3 01:56:38 2014
This archive was generated by hypermail 2.1.8 : Mon Mar 03 2014 - 14:26:18 MST