Re: Writing netcdf files

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 07:29:33 2014