Re: Writing netcdf files

From: Almami Johnson <almamij_at_nyahnyahspammersnyahnyah>
Date: Mon Jan 27 2014 - 08:06:03 MST

Dear Dennis and Rick

Thanks for your suggestions. I will try it in few minutes. Sorry that
my question was not clear. What I actually meant is I want to save
each time step (each month of a given year) in a new file. This is the
reason why I thought I needed to loop. Should'nt I create these
arrays?

  paren = new((/nyear,nmon,nlat,nlon/), "float")
  sunshine = new((/nyear,nmon,nlat,nlon/), "float")
  alpha = new((/nyear,nlat,nlon/), "float")
  pot = new((/nyear,nmon,nlat,nlon/), "float")

Thanks again
regards

On 1/27/14, Dennis Shea <shea@ucar.edu> wrote:
> A few comments on your script:
>
> The following
>
> year=(/"1985","1986","1987","1988","1989","1990","1991","1992","1993","1994","1995","1996","1997","1998","1999","2000","2001","2002","2003","2004"/)
> month=(/"01","02","03","04","05","06","07","08","09","10","11","12"/)
>
> could be replaced with
>
> years = tostring(ispan(1985,2004,1)) ; string
> months = sprinti("%0.2i", ispan(1,12,1))
> print(months)
>
> But, given the context of your use, I would suggest
> year = ispan(1985,2004,1) ; integer
> month = ispan(1,12,1) ; integer
> months = sprinti("%0.2i", ispan(1,12,1)) ; string
> or
>
>
> =====
> Remove this ... you do not have a dimension named 'time'
>
> > ; make time and UNLIMITED dimension ; recommended for most
> applications
> > filedimdef(pot_evap,"time",-1,True)
>
> =====
> As mentioned by RickB, your question is not very clear.
>
> You create arrays
>
> > paren = new((/nyear,nmon,nlat,nlon/), "float")
> > sunshine = new((/nyear,nmon,nlat,nlon/), "float")
> > alpha = new((/nyear,nlat,nlon/), "float")
> > pot = new((/nyear,nmon,nlat,nlon/), "float")
>
> Then you do a double 'do loop' and write each time step. Why
>
> I'd suggest you write one file at the end.
>
> =====
>
>
> On 1/27/14, 6:44 AM, brownrig@ucar.edu wrote:
>> Hi,
>>
>> I suspect you may need to explicitly close each NetCDF file after
>> writing to it with a "delete(pot_evap)" at the end of your
>> doubly-nested loop.
>>
>> As for the dimensions, I'm not clear on your intent, but from the
>> loop, were you wanting something more like:
>>
>> pot_evap->pet = pot(i,j,:,:) ;; rather than "pot_eval->pet =
>> pot"
>>
>> Hope that helps...
>> Rick
>>
>>
>>
>> On Mon, 27 Jan 2014 10:58:40 +0000
>> Almami Johnson <almamij@gmail.com> wrote:
>>> Dear NCL users
>>>
>>> I'm trying to calculate monthly potential evapotranspiration based
>>> on
>>> Thornthwaite (1948). I have each monthly term in a file named
>>> something.year.month.nc. I have created this script to output each
>>> month for a given year in a file named output.year.month.nc.
>>> The problem is all the output.year.month.nc are similar which makes
>>> me
>>> think NCL is writing the same thing in all the output files. Could
>>> you
>>> please take a look at this script and tell me what is wrong with it?
>>> Another problem I'm facing is that the ncdump -h of the outputs
>>> looks like this:
>>>
>>> dimensions:
>>> time = UNLIMITED ; // (0 currently)
>>> ncl1 = 20 ;
>>> ncl2 = 12 ;
>>> latitude = 180 ;
>>> longitude = 280 ;
>>> variables:
>>> float pet(ncl1, ncl2, latitude, longitude) ;
>>> pet:units = "mm/month" ;
>>> pet:long_name = "Potential Evap" ;
>>> pet:_FillValue = 9.96921e+36f ;
>>>
>>> My question is how to output without ncl1 and ncl2 dimensions and
>>> have
>>> the time like time = UNLIMITED ; // (1 currently)
>>>
>>> Thanks for your help
>>> Almami
>>>
>>> -------------------------------------------------------------------------------------------------------------------
>>> 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"
>>>
>>> begin
>>>
>>> ;*****************************************************************************************
>>> ; Creating Function to Read at once Precip and Convert Calendar to
>>> Human Readable YYYYMM
>>> ;*****************************************************************************************
>>> undef("getPrecip")
>>> function getPrecip(fname[1]:string, vname[1]:string)
>>> local x
>>>
>>> begin
>>> f = addfile(fname, "r")
>>> x = f->$vname$(0,0,{-10:35},{-30:40}) ; read specified
>>> time period
>>> return(x)
>>> end
>>>
>>> undef("getPrecip1")
>>> function getPrecip1(fname[1]:string, vname[1]:string)
>>> local y
>>>
>>> begin
>>> f = addfile(fname, "r")
>>> y = f->$vname$(0,{-10:35},{-30:40}) ; read specified time
>>> period
>>> return(y)
>>> end
>>>
>>> ;*****************************************************************************************
>>> ; Creating a new variable corresponding to all the 12 months and
>>> getting each
>>> ;*****************************************************************************************
>>>
>>>
>>> reffile = addfile("ENSMEAN_Pre_paren_2004_05.nc", "r")
>>> refparam = reffile->t2m(:,0,{-10:35},{-30:40})
>>> printVarSummary(refparam)
>>>
>>> dim = dimsizes(refparam)
>>> ntime= dim(0)
>>> nlat = dim(1)
>>> nlon = dim(2)
>>>
>>>
>>> year=(/"1985","1986","1987","1988","1989","1990","1991","1992","1993","1994","1995","1996","1997","1998","1999","2000","2001","2002","2003","2004"/)
>>> month=(/"01","02","03","04","05","06","07","08","09","10","11","12"/)
>>> ndays=(/31,28,31,30,31,30,31,31,30,31,30,31/)
>>>
>>> nyear= dimsizes(year)
>>> nmon = dimsizes(month)
>>>
>>> paren = new((/nyear,nmon,nlat,nlon/), "float")
>>> sunshine = new((/nyear,nmon,nlat,nlon/), "float")
>>> alpha = new((/nyear,nlat,nlon/), "float")
>>> pot = new((/nyear,nmon,nlat,nlon/), "float")
>>>
>>>
>>> do i = 0, nyear - 1
>>> do j = 0, nmon - 1
>>>
>>> print("Getting data and preparing for " +year(i)+ " " +month(j)+ "
>>> and for paren alpha and sunshine duration")
>>>
>>> paren(i,j,:,:) =
>>> getPrecip("ENSMEAN_Pre_paren_"+year(i)+"_"+month(j)+".nc","t2m")
>>> printVarSummary(paren(i,j,:,:))
>>>
>>> sunshine(i,j,:,:) =
>>> getPrecip1("ENSMEAN_Pre_sund_"+year(i)+"_"+month(j)+".nc","sund")
>>> printVarSummary(sunshine(i,j,:,:))
>>>
>>> alpha(i,:,:) = getPrecip("ENSMEAN_Pre_alpha_"+year(i)+".nc","t2m")
>>> printVarSummary(alpha(i,:,:))
>>>
>>> print("Now taking potential evapotranspiration for " +year(i)+"
>>> "+month(j)+ " and for paren alpha and sunshine duration")
>>>
>>> pot(i,j,:,:) =
>>> 16*((paren(i,j,:,:))^(alpha(i,:,:)))*((sunshine(i,j,:,:))/12)*((ndays(j))/30)
>>> copy_VarCoords(paren(i,j,:,:),pot(i,j,:,:))
>>> printVarSummary(pot(i,j,:,:))
>>>
>>> ;*******************************************************************************
>>> ; Creating NetCDF files for writing and put time UNLIMITED
>>> ;*******************************************************************************
>>> print("Creating file
>>> ENSMEAN_Pre_pot_evap_"+year(i)+"_"+month(j)+".nc")
>>>
>>> system("/bin/rm -f
>>> ENSMEAN_Pre_pot_evap_"+year(i)+"_"+month(j)+".nc") ; remove any
>>> pre-existing file
>>> pot_evap =
>>> addfile("ENSMEAN_Pre_pot_evap_"+year(i)+"_"+month(j)+".nc","c") ;
>>> open output netCDF file
>>>
>>>
>>> ; make time and UNLIMITED dimension ; recommended for most
>>> applications
>>> filedimdef(pot_evap,"time",-1,True)
>>>
>>> ;************************************************************************************************************
>>> ; Writing down in the NetCDF files
>>> ;************************************************************************************************************
>>>
>>> ; output variables directly
>>> print("Writing down pet in
>>> ENSMEAN_Pre_pot_evap_"+year(i)+"_"+month(j)+".nc")
>>>
>>> pot!2 = "latitude" ; assign named
>>> dimensions
>>> pot!3 = "longitude"
>>>
>>> pot&latitude = refparam&latitude
>>> pot&longitude = refparam&longitude
>>>
>>> pot@long_name = "Potential Evap" ; assign attributes
>>> pot@units = "mm/month"
>>> pot_evap->pet = pot
>>>
>>> end do
>>> 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
Received on Mon Jan 27 08:06:19 2014

This archive was generated by hypermail 2.1.8 : Fri Feb 07 2014 - 16:39:11 MST