[ncl-talk] Problem reading data from NCEP DODS server
Mary Haley
haley at ucar.edu
Wed May 18 09:06:54 MDT 2016
Chris,
Yes, the "all(ismissing(sst))" is a way to trap this special case.
I'm not sure what you are asking in your second question. When you subset
"sst" with the given lat/lon range, then the lat/lon coordinate arrays
should also be subsetted automatically.
I just tried to run your script again to verify that you indeed have
coordinate arrays attached to your data, but it appears the file is no
longer there.
You can verify this yourself by doing a "printVarSummary" on "sst" after
you read it off a file.
Here's a little sample script that you should be able to run as-is:
a = addfile("$NCARG_ROOT/lib/ncarg/data/cdf/uv300.nc","r")
u = a->U(1,:,:)
printVarSummary(u)
u_subset = a->U(1,{0:30},{-10:10})
printVarSummary(u_subset)
The output (trimmed) will look like this:
Variable: u
Type: float
. . .
Dimensions and sizes: [lat | 64] x [lon | 128]
Coordinates:
lat: [-87.8638..87.8638]
lon: [-180..177.1875]
Number Of Attributes: 5
time : 7
_FillValue : -999
long_name : Zonal Wind
short_name : U
units : m/s
Variable: u_subset
Type: float
. . .
Dimensions and sizes: [lat | 11] x [lon | 7]
Coordinates:
lat: [1.395307..29.30136]
lon: [-8.4375..8.4375]
Number Of Attributes: 5
time : 7
_FillValue : -999
long_name : Zonal Wind
short_name : U
units : m/s
Note that for "u_subset", the lat/lon coordinates have the new subsetted
range, based on the closest *less than* values to the lat/lon range I gave
it.
--Mary
On Tue, May 17, 2016 at 5:54 PM, Herbster, Christopher G. <herbstec at erau.edu
> wrote:
> Mary,
>
>
>
> Thanks for checking what seems obvious in hindsight. Feeling a bit
> sheepish.
>
>
>
> Can I trap this with an:
>
>
>
> if(all(ismissing(sst))) then
>
> … ; don’t do much
>
> else
>
> … ; make plots
>
> end if
>
>
>
> Boy have I wasted some time thinking this had to do with the longitude
> being 74 – 434 degrees east. I had originally thought that I was somehow
> only getting data over land, where SST would certainly be missing.
>
>
>
> Also, is there a preferred way to have the meta data for the subset data
> have ranges that are indicative of the subset, rather than the whole domain?
>
>
>
> Thanks again for catching my oversight,
>
>
>
> Chris H.
>
> --
>
>
>
> Dr. Christopher G. Herbster
>
> Associate Professor
>
> Director of Science and Technology
>
> for the ERAU Weather Center
>
> Applied Aviation Sciences
>
> Embry-Riddle Aeronautical Univ.
>
> 600 S. Clyde Morris Blvd.
>
> Daytona Beach, FL 32114-3900
>
> 386.226.6444 Office
>
> 386.226.6446 Weather Center
>
> http://wx.erau.edu/
>
>
>
> Schedule at: http://wx.erau.edu/faculty/herbster/Schedules/
>
>
>
> *From:* Mary Haley [mailto:haley at ucar.edu]
> *Sent:* Tuesday, May 17, 2016 5:58 PM
> *To:* Herbster, Christopher G. <herbstec at erau.edu>
> *Cc:* ncl-talk at ucar.edu
> *Subject:* Re: [ncl-talk] Problem reading data from NCEP DODS server
>
>
>
> Chris,
>
>
>
> I think the issue is simply that the first time step of "sst" is all
> missing.
>
>
>
> If you select a different time step, then it looks fine:
>
>
>
> ; sst = f->sst(mytime,:,Slat:Nlat,Wlon:Elon)
>
>
>
> sst = f->sst(:,:,Slat:Nlat,Wlon:Elon) ; NOTE: READING ALL TIME STEPS
>
> print("done reading sst")
>
> print("dimsizes sst = " + str_join(dimsizes(sst),","))
>
>
>
> do nt=0,dimsizes(sst(:,0,0,0))-1
>
> printMinMax(sst(nt,:,:,:),0)
>
> end do
>
>
>
> This is the output:
>
>
>
> (0) dimsizes sst = 65,1,265,481
>
> nt = 0
>
> sea_surface_temperature (c) : min=1.26765e+30 max=1.26765e+30
>
> nt = 1
>
> sea_surface_temperature (c) : min=21.6418 max=32.2542
>
> nt = 2
>
> sea_surface_temperature (c) : min=21.6017 max=32.3401
>
> nt = 3
>
> sea_surface_temperature (c) : min=21.5984 max=32.3429
>
> nt = 4
>
> sea_surface_temperature (c) : min=21.6628 max=32.1836
>
> nt = 5
>
> sea_surface_temperature (c) : min=22.0586 max=32.491
>
> nt = 6
>
> sea_surface_temperature (c) : min=22.2371 max=32.6994
>
>
>
> --Mary
>
>
>
>
>
> On Tue, May 17, 2016 at 10:28 AM, Herbster, Christopher G. <
> herbstec at erau.edu> wrote:
>
> Hi folks,
>
>
>
> I’m trying to retrieve some data for a subset region from the Global Real
> Time Ocean Forecast System (RTOFS), and I can’t figure out why I get all
> missing values when I try to read the data.
>
>
>
> This seems to be the case for me even if I try to read this whole domain.
>
>
>
> Here is some code:
>
> ------
>
> ; read-rtofs-dods.ncl
>
> load "$NCARG_ROOT/lib/ncarg/nclscripts/contrib/cd_string.ncl"
>
> ;
>
>
>
> Dev = "png"
>
>
>
> ;
>
> myDate = 20160516
>
>
>
> ;NOMADS = "http://nomads.ncep.noaa.gov:9090/dods/rtofs/rtofs_global
> "+myDate+"/rtofs_glo_2ds_forecast_daily_prog"
>
>
>
> NOMADS = "http://nomads.ncep.noaa.gov:9090/dods/rtofs/rtofs_global
> "+myDate+"/rtofs_glo_2ds_forecast_3hrly_prog"
>
>
>
> filename = NOMADS ; = url
>
> print("Fetching URL: ")
>
> print(" "+filename)
>
>
>
> exists = isfilepresent(filename)
>
> if(.not.exists) then
>
> print("OPeNDAP isfilepresent test unsuccessful.")
>
> print("Either file doesn't exist, or NCL does not have OPeNDAP
> capabilities on this system")
>
> exit
>
> else
>
> print("OPeNDAP isfilepresent test successful.")
>
> f = addfile(filename,"r")
>
> vars = getfilevarnames(f)
>
> print(vars) ; should be (in any order):
>
>
>
> if(.not.any(ismissing(vars))) then
>
> do i=0,dimsizes(vars)-1
>
> printFileVarSummary (f,vars(i))
>
> end do
>
> end if
>
> end if
>
>
>
> ; Define boundary of data
>
> ; Assign lat/lon boundaries to the data
>
> ; Longitudes between 0 and 74 E need to have 360 added to them
>
> ;
>
> myNlat = 32.0
>
> mySlat = 10.0
>
> myWlon = 260.0
>
> myElon = 300.0
>
>
>
> ;myNlat = 26.0
>
> ;mySlat = 21.0
>
> ;myWlon = 260.0
>
> ;myElon = 300.0
>
>
>
> ; Full domain
>
> ;myNlat = 90.0
>
> ;mySlat = -90.0
>
> ;myWlon = 74.0
>
> ;myElon = 434.0
>
>
>
>
>
>
>
> ; Get lat/lon/time info first
>
> lat = f->lat
>
> print ("Done reading lat")
>
> lon = f->lon
>
> print ("Done reading lon")
>
> time = f->time(0)
>
> print ("Done reading time")
>
> ;
>
> ; Use just one time when "0"
>
> mytime = 0
>
> ;mytime = ":"
>
>
>
> Nlat = closest_val(myNlat, lat)
>
> print ("Nlat = " + Nlat + " Lat = " + lat(Nlat) )
>
>
>
> Slat = closest_val(mySlat, lat)
>
> print ("Slat = " + Slat + " Lat = " + lat(Slat) )
>
>
>
> Wlon = closest_val(myWlon, lon)
>
> print("Wlon = " + Wlon + " Lon = " + lon(Wlon) )
>
>
>
> Elon = closest_val(myElon, lon)
>
> print("Elon = " + Elon + " Lon = " + lon(Elon) )
>
> ;
>
> ;exit
>
> ;
>
> mylat = lat(Slat:Nlat)
>
> mylon = lon(Wlon:Elon)
>
>
>
> printVarSummary (mylat)
>
> printVarSummary (mylon)
>
>
>
> print (mylat)
>
>
>
> ;exit
>
>
>
>
>
> delete(lat)
>
> delete(lon)
>
>
>
> lat = f->lat(Slat:Nlat)
>
> lon = f->lon(Wlon:Elon)
>
> print ("Done reading subset lat/lon")
>
>
>
> printVarSummary(lat)
>
> print (lat)
>
>
>
> printVarSummary(lon)
>
> print (lon)
>
>
>
> ;exit
>
>
>
>
>
> sst = f->sst(mytime,:,Slat:Nlat,Wlon:Elon)
>
> ;sst = f->sst
>
> print ("Done reading sst")
>
> printVarSummary(sst)
>
> print (sst)
>
>
>
> ; exit ; Uncomment to stop after sst data are read
>
>
>
> u_velocity = f->u_velocity(mytime,:,Slat:Nlat,Wlon:Elon)
>
> ;u_velocity = f->u_velocity
>
> print ("Done reading U")
>
> printVarSummary(u_velocity)
>
> print (u_velocity)
>
> v_velocity = f->v_velocity(mytime,:,Slat:Nlat,Wlon:Elon)
>
> ;v_velocity = f->v_velocity
>
> print ("Done reading V")
>
> printVarSummary(v_velocity)
>
> print (v_velocity)
>
> speed = sqrt((u_velocity*u_velocity) + (v_velocity*v_velocity))
>
> print ("Done calculating Speed")
>
> print (speed)
>
>
>
>
>
> exit
>
> ;;;;;;
>
> ;;;;;; END of code snippet
>
>
>
> ------
>
>
>
> If you run this code, I suggest a redirect into a file, as there are a lot
> of lines written to standard out.
>
>
>
> What I don’t understand is that I can read a subset of the lat and lon
> data, but when I try to get any other data they all come up as missing/fill
> values.
>
>
>
> Any help is GREATLY appreciated and I have really chased my tail trying to
> get past this.
>
>
>
> Thanks,
>
> Chris Herbster
>
> --
>
>
>
> Dr. Christopher G. Herbster
>
> Associate Professor
>
> Director of Science and Technology
>
> for the ERAU Weather Center
>
> Applied Aviation Sciences
>
> Embry-Riddle Aeronautical Univ.
>
> 600 S. Clyde Morris Blvd.
>
> Daytona Beach, FL 32114-3900
>
> 386.226.6444 Office
>
> 386.226.6446 Weather Center
>
> http://wx.erau.edu/
>
>
>
> Schedule at: http://wx.erau.edu/faculty/herbster/Schedules/
>
>
>
>
> _______________________________________________
> 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/20160518/186eee57/attachment.html
More information about the ncl-talk
mailing list