Re: use cd_calendar when it is necessary to read an variable from many directories

From: Yi-Chih Huang <dscpln_at_nyahnyahspammersnyahnyah>
Date: Fri Nov 08 2013 - 01:06:15 MST

Dave,

    The commands is working. Although I read the NCL documents, yet I
don't fully understand the meaning of the syntax's.

    Anyway, Thank you very much,

           Yi-Chih

On Fri, Nov 8, 2013 at 3:48 AM, Dave Allured - NOAA Affiliate <
dave.allured@noaa.gov> wrote:

> Sorry, I thought of a better solution if you do not need to have the
> numeric time coordinate values [39460..39795] etc. for a later purpose.
> Just compute YYYYMM inside the loop, and the time rebasing problem is
> automatically solved:
>
> 1. After this line:
> xx= new((/nyr,ds(0),ds(1),ds(2)/),float,"No_FillValue")
>
> Add this:
> YYYYMM = new ( (/nyr,ds(0)/), double)
>
> 2. In the read loop, after this line:
> xx(ny,:,:,:)= xs
>
> Add this:
> YYYYMM(ny,:) = cd_calendar(xs&$xs!0$, 1)
>
> 3. After the loop, add this. YYYYMM is a 2-D array, NOT 4-D:
> print(YYYYMM(0:2,0:2))
>
> --Dave
>
> On Thu, Nov 7, 2013 at 11:36 AM, Dave Allured - NOAA Affiliate <
> dave.allured@noaa.gov> wrote:
>
>> Yi-Chih,
>>
>> Yes, you are overwriting the previous time coordinates during array
>> concatenation, xx=xs(ny,:,:,:)= etc. Preserving the right metadata can be
>> tricky when concatenating source arrays, each with their own metadata.
>> Save the time coordinates from each file in a separate array, like this:
>>
>> 1. After this line:
>> xx= new((/nyr,ds(0),ds(1),ds(2)/),float,"No_FillValue")
>>
>> Add this:
>> xtime = new ( (/nyr,ds(0)/), double)
>>
>> 2. In the read loop, after this line:
>> xx(ny,:,:,:)= xs
>>
>> Add this. The critical metadata on the time coordinate variable is also
>> copied:
>> xtime(ny,:) = xs&$xs!0$
>>
>> 3. After the loop, add this. YMDHMS becomes a 2-D array, NOT 4-D:
>> YMDHMS = cd_calendar(xtime, 1)
>> print(YMDHMS(0:2,0:2))
>>
>> Note, this method require that the "units" attribute of the time
>> coordinates be identical between all files. In your case I think this is a
>> good assumption. "hours since 1800-01-01 00:00" is a commonly used fixed
>> value. If this is not the case, then you can add a little more code to
>> perform time rebasing.
>>
>> --Dave
>>
>> On Wed, Nov 6, 2013 at 10:00 PM, Yi-Chih Huang <dscpln@gmail.com> wrote:
>>
>>> Dave,
>>>
>>> Currently my working environment is one variable and one year per
>>> directory. I am new to NCL so that I can't change the data for now
>>> although I have a feeling it would be easier to work on one variable per
>>> file. My input data and the script I am developing are as below.
>>>
>>> Thanks much,
>>>
>>> Yi-Chih
>>>
>>> #####
>>> Copyright (C) 1995-2013 - All Rights Reserved
>>> University Corporation for Atmospheric Research
>>> NCAR Command Language Version 6.1.2
>>> The use of this software is governed by a License Agreement.
>>> See http://www.ncl.ucar.edu/ for more details.
>>>
>>> Variable: f
>>> Type: file
>>> filename: ipvhnl.gdas.2008
>>> path: ipvhnl.gdas.2008.grib2
>>> file global attributes:
>>> dimensions:
>>> initial_time0_hours = 12
>>> lat_0 = 361
>>> lon_0 = 720
>>> variables:
>>> float TMP_P8_L1_GLL0 ( initial_time0_hours, lat_0, lon_0 )
>>> center : US National Weather Service - NCEP (WMC)
>>> production_status : Operational products
>>> long_name : Temperature
>>> units : K
>>> _FillValue : 1e+20
>>> grid_type : Latitude/longitude
>>> parameter_discipline_and_category : Meteorological products,
>>> Temperature
>>> parameter_template_discipline_category_number : ( 8, 0,
>>> 0, 0 )
>>> level_type : Ground or water surface
>>> level : 0
>>> type_of_statistical_processing : Average of N
>>> uninitialized analyses, starting at reference time.
>>> statistical_process_duration : not applicable (intervals of 6
>>> hours)
>>> number_in_average : <ARRAY of 12 elements>
>>> forecast_time : 0
>>> forecast_time_units : hours
>>>
>>> double initial_time0_hours ( initial_time0_hours )
>>> long_name : initial time
>>> units : hours since 1800-01-01 00:00
>>>
>>> double initial_time0_encoded ( initial_time0_hours )
>>> long_name : initial time encoded as double
>>> units : yyyymmddhh.hh_frac
>>>
>>> float lat_0 ( lat_0 )
>>> long_name : latitude
>>> grid_type : Latitude/Longitude
>>> units : degrees_north
>>> Dj : 0.5
>>> Di : 0.5
>>> Lo2 : 359.5
>>> La2 : -90
>>> Lo1 : 0
>>> La1 : 90
>>>
>>> float lon_0 ( lon_0 )
>>> long_name : longitude
>>> grid_type : Latitude/Longitude
>>> units : degrees_east
>>> Dj : 0.5
>>> Di : 0.5
>>> Lo2 : 359.5
>>> La2 : -90
>>> Lo1 : 0
>>> La1 : 90
>>>
>>> string initial_time0 ( initial_time0_hours )
>>> long_name : Initial time of first record
>>> units : mm/dd/yyyy (hh:mm)
>>>
>>>
>>> #####
>>> load "$SysE/lib/ncl/helper_libs.ncl"
>>>
>>> begin
>>> latS = -5.0
>>> latN = 5.0
>>> lonE = 120.0
>>> lonW = 280.0
>>>
>>> yrStrt = 1979
>>> yrLast = 2008
>>>
>>> dir1= "/fs3/SysE_DB/nmm/CFSR/Monthly/Mean/"
>>> dir2=
>>> (/"1979/","1980/","1981/","1982/","1983/","1984/","1985/","1986/","1987/","1988/","1989/",\
>>>
>>> "1990/","1991/","1992/","1993/","1994/","1995/","1996/","1997/","1998/","1999/","2000/",\
>>>
>>> "2001/","2002/","2003/","2004/","2005/","2006/","2007/","2008/"/)
>>> yr= systemfunc("cd "+dir1+" ; ls ") ; input years
>>> print(yr)
>>> varName= (/"sst"/)
>>> nyr= dimsizes(dir2)
>>> print(nyr)
>>>
>>> dir = str_concat(dir1+dir2(0))
>>> fil= systemfunc("cd "+dir+" ; ls sst.nc") ;;;;;;;; read SST
>>> f= addfile (dir+fil, "r")
>>> x= f->$varName(0)$
>>> dims = getfilevardims(f,varName(0))
>>> xs= x(:,{latS:latN},{lonE:lonW}) ; tropical Pacific
>>> ds= dimsizes(xs)
>>> xx= new((/nyr,ds(0),ds(1),ds(2)/),float,"No_FillValue")
>>>
>>> delete (fil)
>>> do ny= 0,nyr-1
>>> ;print(ny)
>>> ;print(dir2(ny))
>>> dir = str_concat(dir1+dir2(ny))
>>> fil= systemfunc("cd "+dir+" ; ls sst.nc")
>>> f= addfile (dir+fil, "r")
>>>
>>> x= f->$varName(0)$
>>> dims = getfilevardims(f,varName(0))
>>> xs= x(:,{latS:latN},{lonE:lonW}) ; tropical Pacific
>>> xx(ny,:,:,:)= xs
>>> end do
>>>
>>> xx!0= "year"
>>> xx!1= "month"
>>> printVarSummary(xx)
>>> ; YMDHMS = cd_calendar(xx, 1)
>>> ;print(YMDHMS)
>>>
>>> ; in = addfile("/fs3/yhuang/precip.mon.mean.nc","r") ;
>>> 197901-201306
>>> ; YMDHMS = cd_calendar(in->time, 0) ; all times
>>> ; index for desired time
>>> ; iYear = ind(YMDHMS(:,0).ge.yrStrt .and. YMDHMS(:,0).le.yrLast)
>>> ; read time and region
>>> precip = in->precip(iYear,{latS:latN},{lonE:lonW})
>>> pcpClm = clmMonTLL( precip )
>>>
>>> xt= precip(0:23,{latS:latN},{lonE:lonW})
>>> ds= dimsizes(xt)
>>> nyr= yrLast -yrStrt +1
>>> pcpm= new((/ds(0),ds(2)/),float)
>>>
>>> yyyymmM = cd_calendar(precip&time,-1) ; verify times are correct
>>>
>>> dimp = dimsizes(precip)
>>> ntim = dimp(0)
>>> nlat = dimp(1)
>>> mlon = dimp(2)
>>>
>>> pcpt = dim_avg_n_Wrap( pcpClm,1 ) ; (time,lat)
>>> do nm= 0,11
>>> pcpm(nm,:) = pcpt(nm,:)
>>> pcpm(nm+12,:)= pcpt(nm,:)
>>> end do
>>>
>>>
>>>
>>> On Thu, Nov 7, 2013 at 1:07 PM, Dave Allured - NOAA Affiliate <
>>> dave.allured@noaa.gov> wrote:
>>>
>>>> Yi-Chih,
>>>>
>>>> It seems that you have concatenated your 30 input directories in an
>>>> unusual way. You created two time dimensions (year and month), where most
>>>> existing NCL examples use a single time dimension ("time") that encodes
>>>> complete dates. The two dimensional case can also be handled, but we need
>>>> to better understand your input data.
>>>>
>>>> If you are still interested in getting cd_calendar to work, please show
>>>> the complete metadata for one of the original input files. Use ncdump -h
>>>> on the OS command line, or else printVarSummary(f) where f is the result of
>>>> addfile in NCL. Showing only xx is not enough, because I think there
>>>> are already mistakes in how this array was created. It looks like you
>>>> may accidentally be overwriting the date/time information in memory, for
>>>> all input files except the last one.
>>>>
>>>> Also please show the complete original NCL script that "reads one
>>>> variable xx from 30 directories for 30 years". I think this is where some
>>>> corrections will be needed.
>>>>
>>>> --Dave
>>>>
>>>>
>>>> On Wed, Nov 6, 2013 at 7:24 PM, Yi-Chih Huang <dscpln@gmail.com> wrote:
>>>>
>>>>> Dennis,
>>>>>
>>>>> Thanks much for the reply. I am trying to make the following
>>>>> commands working when the variable has to be created from 30 directories.
>>>>>
>>>>> YMDHMS = cd_calendar(in->time, 0) ; all times
>>>>> iYear = ind(YMDHMS(:,0).ge.yrStrt .and. YMDHMS(:,0).le.yrLast)
>>>>>
>>>>> The xx is created in the way below. The information month:
>>>>> [39460..39795] is from the months in 2008. If it is not possible to use
>>>>> cd_calendar when the variable has to be created from 30 directories, maybe
>>>>> I should give up cd_calendar and write the script in my own way.
>>>>>
>>>>> xs= x(:,{latS:latN},{lonE:lonW}) ; tropical Pacific
>>>>> xx(ny,:,:,:)= xs
>>>>>
>>>>> Many thanks,
>>>>>
>>>>> Yi-Chih
>>>>>
>>>>>
>>>>> On Thu, Nov 7, 2013 at 10:35 AM, Dennis Shea <shea@ucar.edu> wrote:
>>>>>
>>>>>>
>>>>>>
>>>>>> Hi Yi-Chih
>>>>>>
>>>>>> Please do not send direct emails to me or Mary.
>>>>>> Always send to ncl-talk.
>>>>>>
>>>>>> In particular, I am not paid to answer ncl-talk questions.
>>>>>> I do it voluntarily to help people. However, getting personal
>>>>>> emails makes the communication one-on-one. Now,
>>>>>> your problem is my problem. I just don't have the time,
>>>>>> I have my own job to do. Sending to ncl-talk means
>>>>>> many people see the issue and, hopefully, will respomd.
>>>>>>
>>>>>> The honest truth is that I and others are not sure what
>>>>>> your objective is.
>>>>>>
>>>>>> ---
>>>>>> Naming the dimensions is not important.
>>>>>>
>>>>>> **It is the *units* of the 'time' variable that are important**
>>>>>>
>>>>>>
>>>>>> The documentation states that the 'time' variable
>>>>>>
>>>>>> "Must have a "units" attribute string in the format
>>>>>> "units since basetime", for example, "days since 1971-1-1"
>>>>>>
>>>>>> time is usually a 'one dimensional array', not 2D (year,month)
>>>>>>
>>>>>> Normal netCDF style usage is
>>>>>>
>>>>>> xx(time,lat,lon) ===> xx(360,21,321)
>>>>>>
>>>>>> where the time units would be something like
>>>>>>
>>>>>> "days since ???? ...."
>>>>>>
>>>>>> ===
>>>>>> Perhaps an alternative would be ...
>>>>>>
>>>>>> You created the variable. You must know the years used.
>>>>>> If you know the start and end year
>>>>>>
>>>>>> yrStrt =
>>>>>> yrLast =
>>>>>> year = ispan(yrStrt,yrLast,1)
>>>>>> year!0 = "year"
>>>>>>
>>>>>> month = ispan(1,12,1)
>>>>>> month!0 = "month"
>>>>>>
>>>>>> xx!0 = "year"
>>>>>> xx&year = year
>>>>>>
>>>>>> xx!1 = "month"
>>>>>> xx&month= month
>>>>>>
>>>>>> printVarSummary(xx)
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> ---
>>>>>>
>>>>>> > Variable: xx
>>>>>> > Coordinates:
>>>>>> > month: [39460..39795]
>>>>>> lon: [120..280]
>>>>>>
>>>>>> ==============
>>>>>>
>>>>>> (39795-39460+1) = 336 ..... not 360 so not wure what this represents.
>>>>>>
>>>>>>
>>>>>> On 11/6/13 5:58 PM, Yi-Chih Huang wrote:
>>>>>>
>>>>>>> Dennis,
>>>>>>>
>>>>>>> I gave the names of year and month to the first two dimensions
>>>>>>> of xx.
>>>>>>> But there is no change to YMDHMS. How to give coordinates to xx to
>>>>>>> make cd_calendar
>>>>>>> working?
>>>>>>>
>>>>>>> Thanks much,
>>>>>>>
>>>>>>> Yi-Chih
>>>>>>>
>>>>>>> ######
>>>>>>> xx!0= "year"
>>>>>>> xx!1= "month"
>>>>>>> printVarSummary(xx)
>>>>>>> YMDHMS = cd_calendar(xx, 1)
>>>>>>> print(YMDHMS)
>>>>>>>
>>>>>>> ######
>>>>>>> Variable: xx
>>>>>>> Type: float
>>>>>>> Total Size: 9707040 bytes
>>>>>>> 2426760 values
>>>>>>> Number of Dimensions: 4
>>>>>>> Dimensions and sizes: [year | 30] x [month | 12] x [lat | 21] x
>>>>>>> [lon |
>>>>>>> 321]
>>>>>>> Coordinates:
>>>>>>> month: [39460..39795]
>>>>>>> lat: [-5.. 5]
>>>>>>> lon: [120..280]
>>>>>>> Number Of Attributes: 16
>>>>>>> DataFreq : Monthly
>>>>>>> center : US National Weather Service - NCEP (WMC)
>>>>>>> production_status : Operational products
>>>>>>> long_name : Temperature
>>>>>>> units : K
>>>>>>> grid_type : Latitude/longitude
>>>>>>> parameter_discipline_and_category : Meteorological products,
>>>>>>> Temperature
>>>>>>> parameter_template_discipline_category_number : ( 8, 0, 0,
>>>>>>> 0 )
>>>>>>> level_type : Ground or water surface
>>>>>>> level : 0
>>>>>>> type_of_statistical_processing : Average of N uninitialized
>>>>>>> analyses, starting at reference time.
>>>>>>> statistical_process_duration : not applicable (intervals
>>>>>>> of 6
>>>>>>> hours)
>>>>>>> number_in_average : <ARRAY of 12 elements>
>>>>>>> forecast_time : 0
>>>>>>> forecast_time_units : hours
>>>>>>> _FillValue : 1e+20
>>>>>>>
>>>>>>>
>>>>>>> Variable: YMDHMS
>>>>>>> Type: double
>>>>>>> Total Size: 19414080 bytes
>>>>>>> 2426760 values
>>>>>>> Number of Dimensions: 4
>>>>>>> Dimensions and sizes: [30] x [12] x [21] x [321]
>>>>>>> Coordinates:
>>>>>>> Number Of Attributes: 1
>>>>>>> _FillValue : 9.969209968386869e+36
>>>>>>> (0,0,0,0) 9.969209968386869e+36
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> On Wed, Nov 6, 2013 at 10:40 PM, Dennis Shea <shea@ucar.edu> wrote:
>>>>>>>
>>>>>>> Your 'xx' is a standard variable. It is not a 'time' variable as
>>>>>>>> required by the 'cd_calendar' function:
>>>>>>>> http://www.ncl.ucar.edu/Document/Functions/Built-in/
>>>>>>>> cd_calendar.shtml
>>>>>>>>
>>>>>>>> The documentation states that the 'time' variable
>>>>>>>>
>>>>>>>> "Must have a "units" attribute string in the format
>>>>>>>> "units since basetime", for example, "days since 1971-1-1"
>>>>>>>>
>>>>>>>> Hence, your usage
>>>>>>>> YMDHMS = cd_calendar(xx, 1)
>>>>>>>> is not correct.
>>>>>>>>
>>>>>>>> You should do a printVarSummary(xx)
>>>>>>>>
>>>>>>>> Is there coordinate information for the named dimensions 'year'
>>>>>>>> and 'month'? If there is no information, you must create.
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> On 11/5/13 11:46 PM, Yi-Chih Huang wrote:
>>>>>>>>
>>>>>>>> Hello,
>>>>>>>>>
>>>>>>>>> I read one variable xx from 30 directories for 30 years. The
>>>>>>>>> dimension
>>>>>>>>> of the variable xx is (year,month,lat,lon). I am trying to use
>>>>>>>>> cd_calendar
>>>>>>>>> to get the information of year and month. However, every element
>>>>>>>>> of
>>>>>>>>> YMDHMS
>>>>>>>>> is 9.969209968386869e+36. In addition, the error below occurred
>>>>>>>>> when the
>>>>>>>>> command "iYear = ind(YMDHMS(:,0).ge.yrStrt .and.
>>>>>>>>> YMDHMS(:,0).le.yrLast)" was executed. What is the right way to use
>>>>>>>>> cd_calendar when it is necessary to read an variable from many
>>>>>>>>> directories?
>>>>>>>>>
>>>>>>>>> YMDHMS = cd_calendar(xx, 1)
>>>>>>>>>
>>>>>>>>> (29,11,20,320) 9.969209968386869e+36
>>>>>>>>> fatal:Number of subscripts do not match number of dimensions of
>>>>>>>>> variable,(2) Subscripts used, (4) Subscripts expected
>>>>>>>>> ^Mfatal:["Execute.c":8128]:Execute: Error occurred at or near
>>>>>>>>> line 51 in
>>>>>>>>> file sstTime.ncl
>>>>>>>>>
>>>>>>>>> Thanks much,
>>>>>>>>>
>>>>>>>>> Yi-Chih
>>>>>>>>>
>>>>>>>>
>

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Nov 8 01:06:29 2013

This archive was generated by hypermail 2.1.8 : Mon Nov 11 2013 - 09:45:33 MST