Re: processing grib files

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Fri Feb 21 2014 - 12:52:06 MST

[1] The following is wrong:
         (/"lat_110","lon_110","TEMP"/)
     As noted by Rick, this should be dimension names not variables
          (/"lat_110","lon_110"/)

[2] More fundamental:

     x = f->FOO ; x has all meta data
     printVarSummary(x)

     X = f->FOO*0.001 ; X will have *no* meta data
     printVarSummary(X)

[3] What you should do to retain meta data:

     x = f->FOO
     x = x*0.001
     x@units = "..."

     Specifically:

     prc = f->A_PCP_110_SFC_acc1h({latS:latN},{lonL:lonR})
     prc = prc*0.001
     prc@units = "..."

     etc

[4] You should use the 'simple' method for writing the netCDF.
     This is far simpler.

      http://www.ncl.ucar.edu/Applications/method_1.shtml

>>> ;**************************************************
>>> ; loops through the directory and converts each
>>> ; file to .nc. variables are converted to correct
>>> ; UEB units.
>>> ;**************************************************
>>> 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 ;latN=42.851414
>>> latS=40.500 ;latS=40.581516
>>> lonL=-112.620 ;lonL=-112.561063
>>> lonR=-110.480 ;lonR=-110.598790
>>>
>>> cmd="ls *.grb"
>>> grbFiles=systemfunc(cmd)
>>> nfiles=dimsizes(grbFiles)
>>>
>>> ; loop through the files and process
>>> do n = 0, nfiles-1
>>> g=addfile(grbFiles(n),"r")
>>> ;print(g)
>>> fn=str_split(grbFiles(n),".")
>>> outfn=fn(0)+"."+fn(1)+"."+fn(2)+"."+fn(3)+".nc"
>>> ;print(outfn)
>>>
>>>
>>> ;define the output file
>>> ncdf_out = addfile("./"+outfn ,"c") ; create output netCDF file
>>> setfileoption(ncdf_out,"DefineMode",True)
>>>
>>> ;setup variable atts
>>> varsize=new((/19,17/),float,"No_FillValue")
>>> latsize=new((/19/),float,"No_FillValue")
>>> lonsize=new((/17/),float,"No_FillValue")
>>>
>>> ;assign file attributes
>>> fAtt = True
>>> fAtt@title = "NLDAS Data - converting .grb to .nc"
>>> fAtt@source_file = grbFiles(n)
>>> fAtt@Conventions = "None"
>>> fAtt@creation_date = systemfunc ("date")
>>> fileattdef(ncdf_out, fAtt )
>>>
>>> ;setup variable dimensions
>>> dimNames = (/"lat_110","lon_110","TEMP"/)
>>> dimSizes =
(/dimsizes(latsize),dimsizes(lonsize),dimsizes(varsize) /)
>>> dimUnlim = (/False, False, False/)
>>> filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
>>>
>>> ; predefine names, type, dimensions
>>>
filevardef(ncdf_out,"lat_110",typeof(latsize),getvardims(latsize))
>>>
filevardef(ncdf_out,"lon_110",typeof(lonsize),getvardims(lonsize))
>>> filevardef(ncdf_out,"TEMP",typeof(varsize),getvardims(varsize))
>>>
>>> ; predefine each variable’s attributes
>>> filevarattdef(ncdf_out,"lat_110" ,g->$"lat_110"$({latS:latN}))
>>> filevarattdef(ncdf_out,"lon_110",g->$"lon_110"$({lonL:llonR}))
>>>
filevarattdef(ncdf_out,"TEMP",g->$"TMP_110_HTGL"$({latS:latN},{lonL:lonR}))
>>> ;setfileoption(ncdf_out,”DefineMode”,False) ; optional
>>>
>>> ;do the conversions and add to ncdf_out
>>> tempc=g->$"TMP_110_HTGL"$({latS:latN},{lonL:lonR})-273.15
>>> ncdf_out->tempc = (/TEMP/)
>>>
>>>
;ncdf_out->$"PRECIP"$=g->A_PCP_110_SFC_acc1h({latS:latN},{lonL:lonR})*0.001
>>>
;ncdf_out->$"SW"$=g->DSWRF_110_SFC({latS:latN},{lonL:lonR})*(1.0/3.6)
>>>
;ncdf_out->$"LW"$=g->DLWRF_110_SFC({latS:latN},{lonL:lonR})*(1.0/3.6)
>>>
>>> ;u=g->U_GRD_110_HTGL({latS:latN},{lonL:lonR})
>>> ;v=g->V_GRD_110_HTGL({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"WIND"$=sqrt((u*u)+(v*v))
>>>
>>> ;press=g->PRES_110_SFC({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"PRESS"$=press
>>>
>>> ;spfh=g->SPF_H_110_HTGL({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"SPFH"$=spfh
>>>
>>> ;esat=611.0*exp((17.3*tempc)/(237.3+tempc))
>>> ;ncdf_out->$"RH"$=(press*spfh)*(esat*(spfh+0.622))
>>>
>>> ; print(ncdf_out)
>>>
>>> end do
>>> end
u1G

On 2/21/14, 11:30 AM, Rick Brownrigg wrote:
> I'm sorry, but I don't follow what you're saying "if I use varsize=new((/19*17/),float,"No_FillValue")…."
>
> I see that TEMP is a 2D variable. Going back to the section of code I cited:
>
>> ;setup variable dimensions
>> dimNames = (/"lat_110","lon_110","TEMP"/)
>> dimSizes = (/dimsizes(latsize),dimsizes(lonsize),dimsizes(varsize) /)
>> dimUnlim = (/False, False, False/)
>> filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
>
> this is defining *dimensions*, not variables. If I understand correctly, you want to define 3 variables and 2 dimensions:
>
> dims:
> lat_110 = 19
> lon_110 = 17
>
> vars:
> lat_110(lat_110)
> lon_110(lon_110)
> TEMP(lat_110, lon_110)
>
> Does this make sense, or am I way off here…
>
> Rick
>
> On Feb 21, 2014, at 11:03 AM, john koudelka <john.koudelka@gmail.com> wrote:
>
>> the TEMP variable is a 2D grid of temperature values. the size of the grid is 19x17 after being clipped.
>>
>> if i use
>> varsize=new((/19*17/),float,"No_FillValue")
>>
>>
>> then i get an error with no explanation at the execution of
>> ; predefine names, type, dimensions
>> filevardef(ncdf_out,"lat_110",typeof(latsize),getvardims(latsize))
>>
>>
>> On Fri, Feb 21, 2014 at 10:43 AM, Rick Brownrigg <brownrig@ucar.edu> wrote:
>> Hi John,
>>
>> It looks like perhaps the problem is in this section of code (?):
>>
>> ;setup variable dimensions
>> dimNames = (/"lat_110","lon_110","TEMP"/)
>> dimSizes = (/dimsizes(latsize),dimsizes(lonsize),dimsizes(varsize) /)
>> dimUnlim = (/False, False, False/)
>> filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
>>
>> As written, its an attempt to create the array dimSizes as (/ 19, 17, (/ 19, 17 /) /) which doesn't make sense, and hence the error. What is the intended size of the dimension TEMP that you are trying to create -- it should be a scalar value?
>>
>> Hope that helps…
>> Rick
>>
>> On Feb 21, 2014, at 10:24 AM, john koudelka <john.koudelka@gmail.com> wrote:
>>
>>> hi all,
>>>
>>> i've run into a problem trying to subset and convert some grib files to netcdf and need some help.
>>>
>>> in process i'm doing some unit conversions and the resulting datasets are missing attributes. i've tried to follow examples for creating and defining a netcdf file as closely as possible, but it seems i'm not doing something correctly here.
>>>
>>> the current error is complaining that the array sizes aren't the same, which they aren't:
>>> "fatal:_NclBuildArray: each element of a literal array must have the same dimension sizes, at least one item doesn't"
>>>
>>> i've searched the archives, but haven't been able to get a good grasp on how people solved this problem. what's the proper way to do this?
>>>
>>>
>>> thanks in advance for your help,
>>> john
>>>
>>> ;**************************************************
>>> ; loops through the directory and converts each
>>> ; file to .nc. variables are converted to correct
>>> ; UEB units.
>>> ;**************************************************
>>> 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 ;latN=42.851414
>>> latS=40.500 ;latS=40.581516
>>> lonL=-112.620 ;lonL=-112.561063
>>> lonR=-110.480 ;lonR=-110.598790
>>>
>>> cmd="ls *.grb"
>>> grbFiles=systemfunc(cmd)
>>> nfiles=dimsizes(grbFiles)
>>>
>>> ; loop through the files and process
>>> do n = 0, nfiles-1
>>> g=addfile(grbFiles(n),"r")
>>> ;print(g)
>>> fn=str_split(grbFiles(n),".")
>>> outfn=fn(0)+"."+fn(1)+"."+fn(2)+"."+fn(3)+".nc"
>>> ;print(outfn)
>>>
>>>
>>> ;define the output file
>>> ncdf_out = addfile("./"+outfn ,"c") ; create output netCDF file
>>> setfileoption(ncdf_out,"DefineMode",True)
>>>
>>> ;setup variable atts
>>> varsize=new((/19,17/),float,"No_FillValue")
>>> latsize=new((/19/),float,"No_FillValue")
>>> lonsize=new((/17/),float,"No_FillValue")
>>>
>>> ;assign file attributes
>>> fAtt = True
>>> fAtt@title = "NLDAS Data - converting .grb to .nc"
>>> fAtt@source_file = grbFiles(n)
>>> fAtt@Conventions = "None"
>>> fAtt@creation_date = systemfunc ("date")
>>> fileattdef(ncdf_out, fAtt )
>>>
>>> ;setup variable dimensions
>>> dimNames = (/"lat_110","lon_110","TEMP"/)
>>> dimSizes = (/dimsizes(latsize),dimsizes(lonsize),dimsizes(varsize) /)
>>> dimUnlim = (/False, False, False/)
>>> filedimdef(ncdf_out,dimNames,dimSizes,dimUnlim)
>>>
>>> ; predefine names, type, dimensions
>>> filevardef(ncdf_out,"lat_110",typeof(latsize),getvardims(latsize))
>>> filevardef(ncdf_out,"lon_110",typeof(lonsize),getvardims(lonsize))
>>> filevardef(ncdf_out,"TEMP",typeof(varsize),getvardims(varsize))
>>>
>>> ; predefine each variable’s attributes
>>> filevarattdef(ncdf_out,"lat_110" ,g->$"lat_110"$({latS:latN}))
>>> filevarattdef(ncdf_out,"lon_110",g->$"lon_110"$({lonL:llonR}))
>>> filevarattdef(ncdf_out,"TEMP",g->$"TMP_110_HTGL"$({latS:latN},{lonL:lonR}))
>>> ;setfileoption(ncdf_out,”DefineMode”,False) ; optional
>>>
>>> ;do the conversions and add to ncdf_out
>>> tempc=g->$"TMP_110_HTGL"$({latS:latN},{lonL:lonR})-273.15
>>> ncdf_out->tempc = (/TEMP/)
>>>
>>> ;ncdf_out->$"PRECIP"$=g->A_PCP_110_SFC_acc1h({latS:latN},{lonL:lonR})*0.001
>>> ;ncdf_out->$"SW"$=g->DSWRF_110_SFC({latS:latN},{lonL:lonR})*(1.0/3.6)
>>> ;ncdf_out->$"LW"$=g->DLWRF_110_SFC({latS:latN},{lonL:lonR})*(1.0/3.6)
>>>
>>> ;u=g->U_GRD_110_HTGL({latS:latN},{lonL:lonR})
>>> ;v=g->V_GRD_110_HTGL({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"WIND"$=sqrt((u*u)+(v*v))
>>>
>>> ;press=g->PRES_110_SFC({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"PRESS"$=press
>>>
>>> ;spfh=g->SPF_H_110_HTGL({latS:latN},{lonL:lonR})
>>> ;ncdf_out->$"SPFH"$=spfh
>>>
>>> ;esat=611.0*exp((17.3*tempc)/(237.3+tempc))
>>> ;ncdf_out->$"RH"$=(press*spfh)*(esat*(spfh+0.622))
>>>
>>> ; print(ncdf_out)
>>>
>>> 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
>
>
>
>
> _______________________________________________
> 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 Feb 21 12:52:14 2014

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