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

From: Yi-Chih Huang <dscpln_at_nyahnyahspammersnyahnyah>
Date: Wed Nov 06 2013 - 22:00:02 MST

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
>>>>>>
>>>>>>
>>>>>>
>>>>
>>
>> _______________________________________________
>> 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 Wed Nov 6 22:00:17 2013

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