1

I want to read the variable “tref” for desired latitude and longitude (which I index), 24th M, 4th L, and all S from a NetCDF file:

filename = 'http://iridl.ldeo.columbia.edu/SOURCES/.Models/.NMME/.NCEP-CFSv2/.HINDCAST/.MONTHLY/tref/1+index/.tref/tref/tref/dods';
ncdisp(filename) % so the structure of tref is (X*Y*M*L*S which means lon*lat*ensemble*leadtime*date) 
lon = ncread(filename, 'X');
lat = ncread(filename, 'Y');
date = ncread(filename, 'S');
leadtimes = ncread(filename, 'L');
ensembels = ncread(filename, 'M');
ind_lon = find(lon>=44 & lon<=64);
ind_lat = find(lat>=24 & lat<=40);
tref = ncread(filename, 'tref', [ind_lon ind_lat 24 4 inf]);

I used the following command, but it was not successful. The error is:

Error using horzcat Dimensions of arrays being concatenated are not consistent.

Here are the results of the ncdisp()

Source:
           http://iridl.ldeo.columbia.edu/SOURCES/.Models/.NMME/.NCEP-CFSv2/.HINDCAST/.MONTHLY/tref/1+index/.tref/tref/tref/dods
Format:
           classic
Global Attributes:
           Conventions = 'IRIDL'
Dimensions:
           L = 10
           M = 24
           S = 348
           X = 360
           Y = 181
Variables:
    M   
           Size:       24x1
           Dimensions: M
           Datatype:   single
           Attributes:
                       standard_name = 'realization'
                       long_name     = 'Ensemble Member'
                       pointwidth    = 1
                       gridtype      = 0
                       units         = 'unitless'
    X   
           Size:       360x1
           Dimensions: X
           Datatype:   single
           Attributes:
                       standard_name = 'longitude'
                       pointwidth    = 1
                       gridtype      = 1
                       units         = 'degree_east'
    L   
           Size:       10x1
           Dimensions: L
           Datatype:   single
           Attributes:
                       standard_name = 'forecast_period'
                       long_name     = 'Lead'
                       pointwidth    = 1
                       gridtype      = 0
                       units         = 'months'
    S   
           Size:       348x1
           Dimensions: S
           Datatype:   single
           Attributes:
                       standard_name = 'forecast_reference_time'
                       long_name     = 'Forecast Start Time'
                       pointwidth    = 0
                       calendar      = '360'
                       gridtype      = 0
                       units         = 'months since 1960-01-01'
    Y   
           Size:       181x1
           Dimensions: Y
           Datatype:   single
           Attributes:
                       standard_name = 'latitude'
                       pointwidth    = 1
                       gridtype      = 0
                       units         = 'degree_north'
    tref
           Size:       360x181x24x10x348
           Dimensions: X,Y,M,L,S
           Datatype:   single
           Attributes:
                       pointwidth    = 0
                       standard_name = 'air_temperature'
                       PTVersion     = 2
                       long_name     = 'Reference Temperature'
                       process       = 'Spectral Statistical Interpolation (SSI) analysis from "Final" run.'
                       PDS_TimeRange = 3
                       gribNumBits   = 21
                       center        = 'US Weather Service - National Met. Center'
                       grib_name     = 'TMP'
                       gribcenter    = 7
                       scale_min     = 206.4713
                       gribparam     = 11
                       scale_max     = 315.4577
                       GRIBgridcode  = 3
                       gribleveltype = 105
                       gribfield     = 1
                       units         = 'Kelvin_scale'

Can you help me with this issue, please?

Behzad
  • 15
  • 5
  • 1
    In which line do you get the error? Is it this one: "tref = ncread(filename, 'tref', [ind_lon ind_lat 24 4 inf]);" ? – GrapefruitIsAwesome Nov 08 '22 at 02:17
  • @GrapefruitIsAwesome Yes i got error on this line. – Behzad Nov 08 '22 at 06:18
  • 1
    What are the dimensions of ind_lon & ind_lat? I suspect they are column vectors and you are trying to concatenate them into a row vector – GrapefruitIsAwesome Nov 08 '22 at 12:17
  • 1
    I think the command is not properly written: MatLab expects Fortran type of reading where you first give the start and secondly give the count. Basically, you should have [ind_lon(1) ind_lat(1) 1 1 1] and for the count you give the number of elements how many you want to read: [length(ind_lon) length(ind_lat) 24 4 inf]. Unfortunately, reading this dataset seems to be complicated as I was not able to do even the simplest start = [1 1 1 1 1] and count = [1 1 1 1 1]... – msi_gerva Nov 08 '22 at 17:13
  • Thank you@msi_gerva; Yes you're right it does not even show the results for start = [1 1 1 1 1] and count = [1 1 1 1 1 ]. Do you know what I can do? any suggestion? Thank you. @msi_gerva – Behzad Nov 08 '22 at 23:31
  • 1
    Unfortunately, do not know what could be the issue. I tried to use `ncks` to cut a subset of data and this also fails: `ncks -v tref -d X,44,46 -d Y,29,31 -d S,1,2 -d L,0,1 -d M,1,1 http://iridl.ldeo.columbia.edu/SOURCES/.Models/.NMME/.NCEP-CFSv2/.HINDCAST/.MONTHLY/tref/1+index/.tref/tref/tref/dods test.nc`. The error is: `nco_err_exit(): ERROR Error code is -40. Translation into English with nc_strerror(-40) is "NetCDF: Index exceeds dimension bound"`. Perhaps some netcdf expert can give his/her opinion in this regard... – msi_gerva Nov 09 '22 at 11:06

0 Answers0