John,
Try the script attached below.
Wei
---------
;**************************************************
begin
ncfiles=systemfunc("ls *.nc")
if(dimsizes(ncfiles) .gt. 0) then
system("rm *.nc") ;remove any pre-existing file
end if
;specify extent
latN=42.852
latS=40.500
lonL=-112.620
lonR=-110.480
cmd="ls *.grb"
grbFiles=systemfunc(cmd)
nfiles=dimsizes(grbFiles)
;create time var
;should be 'hours since 2010-10-01T00.00
;time should be provided as floating numbers
start=str_split(grbFiles(0),".")
year=str_get_cols(start(1), 1, 4)
month=str_get_cols(start(1), 5, 6)
day=str_get_cols(start(1), 7, 8)
hour=str_get_cols(start(2),0,1)
time = fspan(0,nfiles,1)
time!0 = "time"
time@units = "hours since "+year+"-"+month+"-"+day+"T00.00"
time&time = time
timec=tocharacter(time)
;print(time)
;create and define the output file
ncdf_out = addfile("./test.nc" ,"c")
ncdf_out@title="Temperature Data: degC"
ncdf_out@source="NLDAS Forcing2 - hourly"
ncdf_out@creation_date = systemfunc ("date")
setfileoption(ncdf_out,"DefineMode",True)
;get dimensions and setup file
g=addfile(grbFiles(0),"r")
l=g->lat_110
lat=l({latS:latN})
lat!0="lat"
lon=g->lon_110({lonL:lonR})
lon!0="lon"
nlat=dimsizes(lat)
nlon=dimsizes(lon)
dimNames = (/"lat","lon","time"/)
dimSizes = (/nlat,nlon,dimsizes(time)/)
dimUnlim = (/False,False,True/)
filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
;setup the variable
t=g->TMP_110_HTGL
tc = new((/nlat,nlon,dimsizes(time)/), float)
tc!0="lat"
tc!1="lon"
tc!2="time"
tc&time=time
tc&lat=lat
tc&lon=lon
tc@sub_center=t@sub_center
tc@center=t@center
tc@long_name=t@long_name
tc@units="degC"
tc@level_indicator=t@level_indicator
tc@gds_grid_type=t@gds_grid_type
tc@parameter_table_version=t@parameter_table_version
tc@parameter_number=t@parameter_number
tc@model=t@model
tc@level=t@level
;setup file variables
filevardef(ncdf_out, "lat", typeof(lat), getvardims(lat))
filevardef(ncdf_out, "lon", typeof(lon), getvardims(lon))
filevardef(ncdf_out, "time", typeof(time),getvardims(time))
filevardef(ncdf_out, "TEMP", typeof(tc), getvardims(tc))
ncdf_out->lat = (/lat/)
ncdf_out->lon = (/lon/)
ncdf_out->time = (/time/)
; loop through the files and process
do n = 0, nfiles-1
g=addfile(grbFiles(n),"r")
;do the conversions and add to ncdf_out
; TEMPERATURE
t=g->TMP_110_HTGL
tc=t({latS:latN},{lonL:lonR})-273.15 ;tc(,,n)=t({latS:latN},{lonL:lonR})-273.15
;this spot reserved for regridding to higher resolution
ncdf_out->TEMP(:,:,n)=(/tc/)
end do
print(ncdf_out)
end
huangwei@ucar.edu
VETS/CISL
National Center for Atmospheric Research
P.O. Box 3000 (1850 Table Mesa Dr.)
Boulder, CO 80307-3000 USA
(303) 497-8924
On Mar 14, 2014, at 2:15 PM, john koudelka <john.koudelka@gmail.com> wrote:
> I need some guidance as to how to properly consolidate a time series of temperature grids into a single netcdf file.
>
> Not sure what I'm doing wrong here, printing the ncdf_out at the end shows only 1 item in the time dimension. Doesn't matter if I place the following code in the do loop or after it. I don't know how to associate the grid with the time array and I couldn't understand how this should be done. Also, the model I'm working with needs data to be dimensioned as: y,x,time.
>
> Variable: ncdf_out
> Type: file
> filename: test
> path: ./test.nc
> file global attributes:
> title : Temperature Data: degC
> source : NLDAS Forcing2 - hourly
> creation_date : Fri Mar 14 14:01:28 MDT 2014
> dimensions:
> lat = 19
> lon = 17
> time = 1
> variables:
> float lat ( lat )
>
> float lon ( lon )
>
> float time ( time )
>
> float TEMP ( lat, lon, time )
> _FillValue : 9.96921e+36
>
>
>
> Many thanks in advance for your help,
> john
>
>
>
>
> Here's my code:
> ;**************************************************
> begin
> ncfiles=systemfunc("ls *.nc")
> if(dimsizes(ncfiles) .gt. 0) then
> system("rm *.nc") ;remove any pre-existing file
> end if
>
> ;specify extent
> latN=42.852
> latS=40.500
> lonL=-112.620
> lonR=-110.480
>
> cmd="ls *.grb"
> grbFiles=systemfunc(cmd)
> nfiles=dimsizes(grbFiles)
>
> ;create time var
> ;should be 'hours since 2010-10-01T00.00
> ;time should be provided as floating numbers
> start=str_split(grbFiles(0),".")
> year=str_get_cols(start(1), 1, 4)
> month=str_get_cols(start(1), 5, 6)
> day=str_get_cols(start(1), 7, 8)
> hour=str_get_cols(start(2),0,1)
>
> time = fspan(0,nfiles,1)
> time!0 = "time"
> time@units = "hours since "+year+"-"+month+"-"+day+"T00.00"
> time&time = time
> timec=tocharacter(time)
> ;print(time)
>
>
> ;create and define the output file
> ncdf_out = addfile("./test.nc" ,"c")
> ncdf_out@title="Temperature Data: degC"
> ncdf_out@source="NLDAS Forcing2 - hourly"
> ncdf_out@creation_date = systemfunc ("date")
> setfileoption(ncdf_out,"DefineMode",True)
>
>
> ;get dimensions and setup file
> g=addfile(grbFiles(0),"r")
> l=g->lat_110
> lat=l({latS:latN})
> lat!0="lat"
> lon=g->lon_110({lonL:lonR})
> lon!0="lon"
> nlat=dimsizes(lat)
> nlon=dimsizes(lon)
>
> dimNames = (/"lat","lon","time"/)
> dimSizes = (/nlat,nlon,dimsizes(time)/)
> dimUnlim = (/False,False,False/)
> filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
>
>
> ;setup the variable
> t=g->TMP_110_HTGL
> tc = new((/nlat,nlon,dimsizes(time)/), float)
> tc!0="lat"
> tc!1="lon"
> tc!2="time"
> tc&time=time
> tc&lat=lat
> tc&lon=lon
> tc@sub_center=t@sub_center
> tc@center=t@center
> tc@long_name=t@long_name
> tc@units="degC"
> tc@level_indicator=t@level_indicator
> tc@gds_grid_type=t@gds_grid_type
> tc@parameter_table_version=t@parameter_table_version
> tc@parameter_number=t@parameter_number
> tc@model=t@model
> tc@level=t@level
>
> ;setup file variables
> filevardef(ncdf_out, "lat" ,typeof(lat),getvardims(lat))
> filevardef(ncdf_out, "lon" ,typeof(lon),getvardims(lon))
> filevardef(ncdf_out, "time" ,typeof(time),getvardims(time))
> filevardef(ncdf_out, "TEMP" ,typeof(tc) ,getvardims(tc))
>
>
> ; loop through the files and process
> do n = 0, nfiles-1
> g=addfile(grbFiles(n),"r")
>
> ;do the conversions and add to ncdf_out
> ; TEMPERATURE
> t=g->TMP_110_HTGL
> tc=t({latS:latN},{lonL:lonR})-273.15 ;tc(,,n)=t({latS:latN},{lonL:lonR})-273.15
> ;this spot reserved for regridding to higher resolution
> ;ncdf_out->TEMP=(/tc/)
>
>
> tc!0="lat"
> tc!1="lon"
> tc!2="time"
> tc&lat=lat
> tc&lon=lon
> tc&time=time
> end do
>
>
> ncdf_out->TEMP=(/tc/)
> print(ncdf_out)
> 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 Fri Mar 14 15:41:33 2014
This archive was generated by hypermail 2.1.8 : Fri Mar 21 2014 - 15:49:21 MDT