[ncl-talk] Problem reading data from NCEP DODS server

[Herbster, Christopher G.]

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/

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20160517/ca0c36fb/attachment.html