Re: help creating a time netcdf file

Dennis and Wei,

Thanks for your responses.

I'm working with a set of grib files, where each file is a step in a time.
So, for example, my directory contains a file for each hour:
NLDAS_FORA0125_H.A20101001.0000.002.grb
NLDAS_FORA0125_H.A20101001.0100.002.grb
NLDAS_FORA0125_H.A20101001.0200.002.grb
. . .
etc.

Each NLDAS file contains a 2d temperature dataset (I know it contains other
variables, temp is all I'm trying to get right now).

I'm trying to loop through the directory, extract, clip, and regrid each
temperature dataset, then place them into a single netcdf file. Ie,
temp.ncwould contain each of the temperature grids, in (y,x,t) format.

I don't need an unlimited dimension in time, since I obtain a file count at
the beginning of the script, subsequently allowing me to declare the size
of the time dimension.

The model I'm working with is only compatible with netcdf 3, and prefers
time located in the third dimension, hence the insistence to setup the file
as y,x,t. From what I've read, there seems to be difficulty in setting it
up like this instead of t, y, x. Is there a function that would permit me
to swap dimensions around? Or is this only a limitation when trying to set
the 3rd dimension as unlimited? Since I know the size of the time
dimension, will I be OK?

I know this is a bit involved, and I'm still learning NCL and netCDF
concepts, so thank you for your patience and thoughtful guidance.

John

On Fri, Mar 14, 2014 at 4:25 PM, Dennis Shea <shea@ucar.edu> wrote:

> In netCDF-3, the unlimited dimension must be the record
> dimension. The record dimension is the leftmost dimension.
> You could try commenting or deleting that or as Wei has suggested
>
> =======> Wei Huang
> Just try to change line:
>
>
> dimUnlim = (/False,False,True/)
>
> to
>
> dimUnlim = (/False,False,False/)
>
> Wei
> ========
>
>
> or
>
> It is my speculation that the GRIB file you are reading
> has only one time step. By default, NCL will not read
> degenerate (singleton) dimensions. Hence, GRIB files with
> only one time step will not have an explicit 'time' dimension.
>
> For illustration ... Try
>
> http://www.ncl.ucar.edu/Document/Tools/ncl_filedump.shtml
>
> Specifically ....
>
> %> ncl_filedump sample_file.grb | less
>
> Then try with the -itime option
>
> %> ncl_filedump -itime sample_file.grb | less
>
> The latter forces NCL to create singleton dimensions.
> ---
>
> Programmatically .... use setfileoption
> http://www.ncl.ucar.edu/Document/Functions/Built-in/setfileoption.shtml
>
>
> ; Force a time dimension for single element time
> setfileoption("grb","SingleElementDimensions","Initial_time")
>
> ; open the ERA-I GRIB files from NCAR-CISL RDA
> g = addfile("sample_file.grb" , "r")
> t =g->TMP....
> printVarSummary(t)
>
> or better
>
> latN = 42.852
> latS = 40.500
> lonL = -112.620
> lonR = -110.480
>
> t =g->TMP....(:,{latS:latN},{lonL:lonR})
> printVarSummary(t)
>
> you can rename dimensions as needed
>
> t!0 = "timE
> t!1 = "lat"
> t!2 = "lon"
>
> We would recommend writing the netCDF the 'simple way'
>
> http://www.ncl.ucar.edu/Applications/method_1.shtml
>
> ---
> system("/bin/rm -f simple.nc") ; remove any pre-existing file
> ncdf = addfile("simple.nc" ,"c") ; open output netCDF file
>
> ; make time and UNLIMITED dimension
> ; recommended for most applications
> filedimdef(ncdf,"time",-1,True)
>
> ; output variables directly
> ncdf->T = t
>
>
>
>
> On 3/14/14, 3:59 PM, john koudelka wrote:
>
>> Wei, thanks for the suggestions.
>>
>>
>> The script errors with:
>> fatal:NetCDF: NC_UNLIMITED in the wrong index
>> fatal:["NclFile.c":474]:FileAddVar: an error occurred while adding a
>> variable to a file, check to make sure data type is supported by the
>> output
>> format
>>
>> And refers to this line:
>> filevardef(ncdf_out, "TEMP" ,typeof(tc) ,getvardims(tc))
>>
>>
>>
>>
>>
>> On Fri, Mar 14, 2014 at 3:41 PM, Wei Huang <huangwei@ucar.edu> wrote:
>>
>> 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
>>
>>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Sat Mar 15 10:05:25 2014