Hi All,
I'm in uncharted waters here, playing with ncdf routines for the first time,
trying to read some Numerical Weather mode data from an OPeNDAP system.
My call to the NCDF_VARGET routine is tripping up on variables defined as :-
FLOAT = Array[1088, 746, 70]
Of that 3D array, I want [*,*,0], the bottom layer of the atmosphere.
My fledging code is given in full below, but here's the problem call :-
zw_id = ncdf_varid(file_id, 'zonal_wnd')
ncdf_varget, file_id, zw_id, zonal_wind, COUNT = [391,391,1], OFFSET = [428,191,0],STRIDE=[1,1,1]
According to the IDL Help, COUNT is a 1-based vector with an element for each dimension of the data to be written. OFFSET seems to be a normal 0 based vector,
and STRIDE is entirely optional.
If I don't provide any of the COUNT, OFFSET or STRIDE keywords, then the entire file is read in properly - no probs there.
If I leave off just the STRIDE keyword, there's no change. Still fails.
In short, I haven't been able to find any combination of the COUNT, OFFSET and STRIDE keywords that works. Remove them all, and I get the entire file, but that's very slow, and simply more data then I need by a factor of ~70!
My wife would say it's obviously because I'm doing something stupid. You may agree. I'm hoping someone out there can point out just what stupid thing it is that I'm doing...
Andrew
NB : This codes relies on wget in IDL v8.5.1 to find the latest dataset online.
If you don't have v8.5.1 you can either spawn a call to wget on your system, or
just just drop the base url into a browser to see what the latest year, month and day are.
PRO ADC_OPENDAP
CLOSE,/ALL
DEVICE,DECOMP=0
base_url = ' http://opendap.bom.gov.au:8080/thredds/dodsC/nmoc_catalogs/A CCESS-R/operational/model_level/latest/'
text = WGET(base_url,/STRING_ARRAY)
first_file_line = text[20]
Access_pos = strpos(first_file_line,'latest/ACCESS-R') + 16
latest_model_run = STRMID(first_file_line,access_pos,10)
print,'latest_model_run = ',latest_model_run
; year =
; month = 7
; day = 3
; hour = 0
var_names = ['lat','lon','lvl','zonal_wnd','merid_wnd','air_temp','verti cal_wnd','pressure','dewpt','relhum']
for forecast_hour = 9,9 Do Begin
hour_file = 'ACCESS-R_' + latest_model_run + '_' + $ String(forecast_hour,form='(i3.3)') + '_model.nc'
this_url = base_url + hour_file
file_id = ncdf_open(this_url)
info = ncdf_inquire(file_id)
Nvars = info.nvars
var_ids = ncdf_varidsinq(file_id)
t = systime(1)
; lat,lon and lvl are single dimension vectors , and are rad in AOK
lat_id = ncdf_varid(file_id, 'lat')
ncdf_varget, file_id, lat_id, lat, COUNT = [391], STRIDE=[1], OFFSET = [191]
help,lat
; lon_id = ncdf_varid(file_id, 'lon')
; ncdf_varget, file_id, lon_id, lon, COUNT = [391], STRIDE=[1], OFFSET = [428]
; lvl_id = ncdf_varid(file_id, 'lvl')
; ncdf_varget, file_id, lvl_id, lvl;
;
; Here's the (first) problem child. The dimension of 70 is or various ;atmospheric heights.
;The error message:-
;% Insufficient number of indices in OFFSET array (<INT Array[3]>).
;ZONAL_WIND FLOAT = Array[1088, 746, 70]
zw_id = ncdf_varid(file_id, 'zonal_wnd')
ncdf_varget, file_id, zw_id, zonal_wind, COUNT = [391,391,1], OFFSET = [428,191,0],STRIDE=[1,1,1]
print,'zw ok'
help,zonal_wind
stop
airtemp_id = ncdf_varid(file_id, 'air_temp')
ncdf_varget, file_id, airtemp_id, air_temp
vert_id = ncdf_varid(file_id, 'vertical_wnd')
ncdf_varget, file_id, vert_id, vertical_wind
merid_id = ncdf_varid(file_id, 'merid_wnd')
ncdf_varget, file_id, merid_id, merid_wind
pressure_id = ncdf_varid(file_id, 'pressure')
ncdf_varget, file_id, pressure_id, pressure
dewpt_id = ncdf_varid(file_id, 'dewpt')
ncdf_varget, file_id, dewpt_id, dewpt
relhum_id = ncdf_varid(file_id, 'relhum')
ncdf_varget, file_id, relhum_id, relhum
print,'time taken = ',systime(1) - t,' seconds'
help,lat,lon,relhum
ncdf_close, file_id
END
Window,xs=1190,ys=820
LOADCT,39
!P.Position=[0.1,0.1,0.9,0.9]
map_set,-24,125,0,limit=[-65,65,16.95,185]
!ORDER = 1
TV,CONGRID(BYTSCL(relhum),1190*.8,820*.8),0.1,0.1,/NORM
END
|