[ncl-talk] Extracting A Time Series of Forcing Data Information from Specific Locations

Dennis Shea shea at ucar.edu
Fri Mar 9 14:40:22 MST 2018 

Apologies: I was in too much of a rush.

The following uses the files you sent.
======================

; vector of names of all of the files and add full file path

diri = "./"
fili = systemfunc("ls "+diri+"NLDAS_FORA0125_H.A2009*.002.grb")
print(fili)   ; make sure the desired files are listed in the correct order

; read in all files

;ListSetType (grb_file, "join")
setfileoption("grb","SingleElementDimensions","Initial_time") ;
initial_time0_hours

grb_file = addfiles(fili,"r")
data = grb_file[:]->DSWRF_110_SFC
data!0 = "time"
printVarSummary(data)

;====
   LAT  = (/ 48.5625/)   ; read from text file
   LON  = (/-98.9375/)
   LAT at units = "degrees_north"
   LON at units = "degrees_east"

; Simple approx:

   print(data(:,{LAT},{LON}) )
;=====
; or better use:

;
https://www.ncl.ucar.edu/Document/Functions/Contributed/linint2_points_Wrap.shtml

points = linint2_points(data&lon_110,data&lat_110,data, False, LON,LAT, 0)
printVarSummary(points)
print("--------------")

+++++++++++++++++++++++++++++
Variable: fili
Type: string
Total Size: 24 bytes
            3 values
Number of Dimensions: 1
Dimensions and sizes:    [3]
Coordinates:
(0)    ./NLDAS_FORA0125_H.A20090316.0100.002.grb
(1)    ./NLDAS_FORA0125_H.A20090316.0200.002.grb
(2)    ./NLDAS_FORA0125_H.A20090316.0300.002.grb

Variable: data
Type: float
Total Size: 1247232 bytes
            311808 values
Number of Dimensions: 3
Dimensions and sizes:    [time | 3] x [lat_110 | 224] x [lon_110 | 464]
Coordinates:
            time: [1833841..1833843]
            lat_110: [25.063..52.938]
            lon_110: [-124.938..-67.063]
Number Of Attributes: 12
  sub_center :    NESDIS Office of Research and Applications
  center :    US National Weather Service - NCEP (WMC)
  long_name :    Downward shortwave radiation flux
  units :    W/m^2
  _FillValue :    1e+20
  level_indicator :    1
  gds_grid_type :    0
  parameter_table_version :    130
  parameter_number :    204
  model :    MESO NAM Model (currently 12 km)
  forecast_time :    0
  forecast_time_units :    hours
------
Variable: data (subsection)       <====  print(data(:,{LAT},{LON}) )
Type: float
Total Size: 12 bytes
            3 values
Number of Dimensions: 1
Dimensions and sizes:    *[time | 3]*
Coordinates:
            time: [1833841..1833843]
Number Of Attributes: 14
  lon_110 :    -98.938
  lat_110 :    48.563
  forecast_time_units :    hours
  forecast_time :    0
  model :    MESO NAM Model (currently 12 km)
  parameter_number :    204
  parameter_table_version :    130
  gds_grid_type :    0
  level_indicator :    1
  _FillValue :    1e+20
  units :    W/m^2
  long_name :    Downward shortwave radiation flux
  center :    US National Weather Service - NCEP (WMC)
  sub_center :    NESDIS Office of Research and Applications
==========

Variable: points
Type: float
Total Size: 12 bytes
            3 values
Number of Dimensions: 2
Dimensions and sizes:    *[time | 3] x [pts | 1]*
Coordinates:
            time: [1833841..1833843]
            pts: [0..0]
Number Of Attributes: 14
  _FillValue :    1e+20
  sub_center :    NESDIS Office of Research and Applications
  center :    US National Weather Service - NCEP (WMC)
  long_name :    Downward shortwave radiation flux
  units :    W/m^2
  level_indicator :    1
  gds_grid_type :    0
  parameter_table_version :    130
  parameter_number :    204
  model :    MESO NAM Model (currently 12 km)
  forecast_time :    0
  forecast_time_units :    hours
  xcoord :    -98.9375
  ycoord :    48.5625










On Fri, Mar 9, 2018 at 1:59 PM, Dennis Shea <shea at ucar.edu> wrote:

> [1]
>
> The following is incorrect. grb_list is a one-dimensional *array* of file
> names containing strings (text).  You have it as a single string.
>
>
>     grb_file = addfiles("grb_list","r")
>
>
> Use
>
>
>     grb_file = addfiles( grb_list ,"r")  ; pass the 1-d array of file names
> =============
> [2]
> You want time series, N'est ce pas? Remove the following:
>
>     ListSetType (grb_file, "join")            ; NO: NOT if you want a time
> serie
> =============
> [3]
>
>
> diri = ncargpath("/depot/phig/apps/LIS_VIC/NLDAS_Forcing/discover/
> nobackup/projects/lis/MET_FORCING/NLDAS2.FORCING/2009/075/")
>
> grb_file = addfiles("grb_list","r")
>
> data      = grb_file[:]->DSWRF_110_SFC
>
> printVarSummary(data)
>
>
> ariable: data
> Type: float
> Total Size: 2494464 bytes
>             623616 values
> Number of Dimensions: 2
> Dimensions and sizes:    [lat_110 | 1344] x [lon_110 | 464]
> Coordinates:
>             lat_110: [25.063..52.938]
>             lon_110: [-124.938..-67.063]
>
>
> =====
>
>    LAT  = (/ 48.5625 /)   ; read from text file
>
>    LAT at units = "degrees_north"
>
>    LON = (/-98.9375/)
>
>    LON&units = "degrees_east"
>
>    npts = dimsizes(LAT)
>
>
> ; Simple approx:
>
>
>    print( data(:,{LAT),{LON}) )
>
> =====
>
>
> or better use:
>
> https://www.ncl.ucar.edu/Document/Functions/Contributed/linint2_points_
> Wrap.shtml
>
>
> points = *linint2_points*(data&lon,data&lat,data, False, LON,LAT, 0)
>
> printVarSummary(points)
>
> print("--------------")
>
>
> Good Luck
>
> On Fri, Mar 9, 2018 at 9:37 AM, Smith, Stuart <smit1770 at purdue.edu> wrote:
>
>> Hello,
>>
>>
>>
>> I would like to compare two sets of forcing data. One is in a .grb file
>> while the other is a .txt. When working the GRIB files (1 hour timesteps),
>> I would like my .ncl script to reads in a time series of data, with the
>> goal to plot specific weather data variables from specific locations to
>> compare the .txt files.
>>
>>
>>
>> However, when running my .ncl script listed below I receive a
>> Segmentation fault error. I’m not sure if this is because I am not defining
>> the latitude and longitude coordinates soon enough, causing an error by
>> reading in too much data. I do not know how to print or plot specific
>> latitude (48.5625) and longitude (-98.9375) points and print their forcing
>> variables. Could you please help improve my .ncl script? I have also
>> uploaded an example of the .grb files I am working with. Thank you for your
>> time.
>>
>>
>>
>>
>>
>> begin
>>
>>
>>
>> ;Go to directory of interest
>>
>>
>>
>> diri = ncargpath("/depot/phig/apps/LIS_VIC/NLDAS_Forcing/discover/
>> nobackup/projects/lis/MET_FORCING/NLDAS2.FORCING/2009/075/")
>>
>>
>>
>>
>>
>> ;list names of all of the files and add ful file path
>>
>> grb_list = systemfunc("ls "+diri+"NLDAS_FORA0125_H.*.002.grb")
>>
>>
>>
>> ;print(grb_list)   ; make sure the desired files are listed in the
>> correct order
>>
>>
>>
>> ; read in all files
>>
>>
>>
>> grb_file = addfiles("grb_list","r")
>>
>> ListSetType (grb_file, "join")
>>
>> data = grb_file[:]->DSWRF_110_SFC
>>
>> data at lat=grb_file[:]->lat_110  ;Lat. of interest 48.5625
>>
>> data at lon=grb_file[:]->lon_110  ;Lon. of interest -98.9375
>>
>>
>>
>> print(data)
>>
>>
>>
>>
>>
>> end
>>
>>
>>
>>
>>
>> Regards,
>>
>>
>>
>> -Stuart
>>
>>
>>
>> _______________________________________________
>> ncl-talk mailing list
>> ncl-talk at ucar.edu
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20180309/7f52fe84/attachment.html>

More information about the ncl-talk mailing list