getnc_s.m

读取Network Common Data Form (netCDF)数据

function values = getnc_s(file, varid, corner, end_point, stride, order, ...
      change_miss, new_miss, squeeze_it, rescale_opts)
% GETNC_S returns a hyperslab of values from a netCDF variable.
%
%  function values = getnc_s(file, varid, corner, end_point, stride, order, ...
%        change_miss, new_miss, squeeze_it, rescale_opts)
%
% DESCRIPTION:
% getnc_s is an non-interactive function that gets a hyperslab
% of values from a netCDF variable; its arguments specify what data
% points are to be retrieved.  Thus, getnc_s is suitable for use
% in a matlab script or function file.  In practice its common use is to be
% called by getnc.
%  
% INPUT:
%  file is the name of a netCDF file with or without the .cdf or .nc extent.
%  varid may be an integer or a string.  If it is an integer then it
%    must be the menu number of the n dimensional variable as used
%    by a call to inqnc or getnc.  If it is a string then it should
%    be the name of the variable.
%  corner is a vector of length n specifying the hyperslab corner
%    with the lowest index values (the bottom left-hand corner in a
%    2-space).  The corners refer to the dimensions in the same
%    order that these dimensions are listed in the relevant questions
%    in getnc.m and in the inqnc.m description of the variable.  A
%    negative element means that all values in that direction will be
%    returned.  If a negative scalar is used this means that all of the
%    elements in the array will be returned.
%  end_point is a vector of length n specifying the hyperslab corner
%    with the highest index values (the top right-hand corner in a
%    2-space).  The corners refer to the dimensions in the same order
%    that these dimensions are listed in the relevant questions in
%    getnc.m and in the inqnc.m description of the variable.  An
%    element in the end_point vector wil be ignored if the corresponding
%    element in the corner vector is negative.
%  stride is a vector of length n specifying the interval between
%    accessed values of the hyperslab (sub-sampling) in each of the n
%    dimensions.  A value of 1 accesses adjacent values in the given
%    dimension; a value of 2 accesses every other value; and so on. If
%    no sub-sampling is required in any direction then it is allowable
%    to just pass the scalar 1 (or -1 to be consistent with the corner
%    and end_point notation).
%  order is a vector of length n specifying the order of the dimensions in
%    the returned array.  order = -1 or [1 2 3 .. n] for an n dimensional
%    netCDF variable will return an array with the dimensions in the same
%    order  as described by a call to inqnc(file) from within matlab or
%    'ncdump  -h' from the command line.  Putting order = -2 will reverse
%    this order.  More general permutations are given re-arranging the
%    numbers 1 to n in the vector.
%  change_miss == 1 causes missing values to be returned unchanged.
%    change_miss == 2 causes missing values to be changed to NaN.
%    change_miss == 3 causes missing values to be changed to new_miss
%    (after rescaling if that is necessary).
%    change_miss < 0 produces the default (missing values to be changed to
%    NaN).
%  new_miss is the value given to missing data if change_miss = 3.
%  squeeze_it specifies whether the returned array should be squeezed.
%    That is, when squeeze_it is non-zero then the squeeze function will
%    be applied to the returned array to eliminate singleton array
%    dimensions.  This is the default.  Note also that a 1-d array is
%    returned as a column vector.
%  rescale_opts is a 2 element vector specifying whether or not rescaling is
%    carried out on retrieved variables or attributes. Only use this option
%    if you are sure that you know what you are doing.
%    If rescale_opts(1) == 1 then a variable read in by getnc.m will be
%    rescaled by 'scale_factor' and  'add_offset' if these are attributes of
%    the variable; this is the default. If rescale_opts(1) == 0 then this
%    rescaling will not be done.
%    If rescale_opts(2) == 1 then the attributes '_FillValue', 'valid_range',
%    'valid_min' and 'valid_max' read in by getnc.m (and used to find the
%    missing values of the relevant variable) will be rescaled by
%    'scale_factor' and 'add_offset'; this is the default. If rescale_opts(2)
%    == 0 then this rescaling will not be done.
%
% OUTPUT:
%  values is a scalar, vector or array of values that is read in
%     from the NetCDF file
%
% NOTES:
% 1) In order for getnc to work non-interactively it is only strictly
% necessary to pass the first 2 input arguments to getnc - sensible
% defaults are available for the rest.
% These are:
% corner, end_point = [-1 ... -1], => all elements retrieved
% stride = 1, => all elements retrieved
% order = [1 2 3 .. n] for an n dimensional netCDF variable
% change_miss = 2, => missing values replaced by NaNs
% new_miss = 0;
% squeeze_it = 1; => singleton dimensions will be removed
% 2) It is not acceptable to pass only 3 input arguments since there is
% no default in the case of the corner points being specified but the
% end points not.
% 3) By default the order of the dimensions of a returned array will be the
% same as they appear in the relevant call to 'inqnc' (from matlab) or
% 'ncdump -h' (from the command line).  (This is the opposite to what
% happened in an earlier version of getnc.)  This actually involves getnc
% re-arranging the returned array because the netCDF utilities follow the C
% convention for data storage and matlab follows the fortran convention.
%    To be more explicit, suppose that we use inqnc to examine a netCDF file.
% Then 2 lines in the output might read:
%
% The 3 dimensions are  1) month = 12  2) lat = 90  3) lon = 180.
%     ---  Information about airtemp(month lat lon )  ---
%
% Likewise using 'ncdump -h' from the command line would have the
% corresponding line:
%
%        short airtemp(month, lat, lon);
%
% The simplest possible call to getnc will give, as you would expect,
%
% >> airtemp = getnc('oberhuber', 'airtemp');
% >> size(airtemp) = 12    90   180
%
% Since the netCDF file follows the C convention this means that in the
% actual storage of the data lon is the fastest changing index, followed by
% lat, and then month.  However matlab (and fortran) use the opposite
% convention and so month is the fastest changing index, followed by lat, and
% then lon.  getnc actually used the permute function to reverse the storage
% order.  If efficiency is a concern (because of using very large arrays or
% a large number of small arrays) then passing order == -2 to getnc will
% produce the fastest response by returning the data in its 'natural' order
% (a 180x90x12 array in our example)
%
% 4) If the values are returned in a one-dimensional array then this will be a
% column vector.  This choice provides consistency if there is an
% unlimited dimension.  That is, if the length of the unlimited
% dimension is n then an m x n array will be returned even for n = 1.)
%
% 5) A strange 'feature' of matlab 5 is that it will not tolerate a singleton
% dimension as the final dimension.  Thus, if you chose to have only one
% element in the final dimension this dimension will be 'squeezed' whether
% you want it to be or not - this seems to be unavoidable.
%
% EXAMPLES:
% 1) Get all the elements of the variable, note the order of the dimensions:
% >> airtemp = getnc('oberhuber', 'airtemp');
% >> size(airtemp)
% ans =
%     12    90   180
%
% 2) Get a subsample of the variable, note the stride:
% >> airtemp = getnc('oberhuber', 'airtemp', [-1 1 3], [-1 46 6], [1 5 1]);
% >> size(airtemp)
% ans =
%     12    10     4
%
% 3) Get all the elements of the variable, but with dimensions permuted:
% >> airtemp = getnc('oberhuber', 'airtemp', -1, -1, -1, [2 3 1]);
% >> size(airtemp)
% ans =
%     90   180    12
% 
% 4) Get all the elements of the variable, but with missing values
%    replaced with 1000.  Note that the corner, end_point, stride and
%    order vectors have been replaced by -1 to indicate that the whole
%    range is required in the default order:
% >> airtemp = getnc('oberhuber', 'airtemp', -1, -1, -1, -1, 3, 1000); 
% >> size(airtemp)
% ans =
%     12    90   180
%
% 5) Get a subsample of the variable, a singleton dimension is squeezed:
% >> airtemp = getnc('oberhuber', 'airtemp', [-1 7 -1], [-1 7 -1]);   
% >> size(airtemp)                                                         
% ans =
%     12   180
% 
% 6) Get a subsample of the variable, a singleton dimension is not squeezed:
% >> airtemp = getnc('oberhuber','airtemp',[-1 7 -1],[-1 7 -1],-1,-1,-1,-1,0);
% >> size(airtemp)                                                            
% ans =
%     12     1   180
%
%
% CALLER:   general purpose
% CALLEE:   fill_att.m, check_st.m, y_rescal.m, ncmex.mex
%
% AUTHOR:   J. V. Mansbridge, CSIRO
%---------------------------------------------------------------------

%     Copyright (C), J.V. Mansbridge, 1992
%     Commonwealth Scientific and Industrial Research Organisation
%     Revision $Revision: 1.12 $
%     Author   $Author: man133 $
%     Date     $Date: 2004/11/16 04:50:30 $
%     RCSfile  $RCSfile: getnc_s.m,v $
% 
%--------------------------------------------------------------------

% Written by Jim Mansbridge december 10 1991

% In November 1998 some code was added to deal better with byte type data. Note
% that any values greater than 127 will have 256 subtracted from them. This is
% because on some machines (an SGI running irix6.2 is an example) values are
% returned in the range 0 to 255. Note that in the fix the values less than 128
% are unaltered and so we do not have to know whether the particular problem has
% occurred or not; for machines where there is no problem no values will be
% altered. This is applied to byte type attributes (like _FillValue) as well as
% the variable values.

% Check that there are the correct number of arguments. If the 6th, 7th
% or 8th arguments are missing then they are set here.  If the 3rd, 4th
% or 5th arguments are missing then their defaults will be set later
% when we find out the dimensions of the variable.  It produces an error
% if the corner is set but the end_point is not.

if (nargin == 2) | (nargin == 4) | (nargin == 5)
  order = -1;
  change_miss = 2;
  new_miss = 0;
  squeeze_it = 1;
elseif nargin == 6
  change_miss = 2;
  new_miss = 0;
  squeeze_it = 1;
elseif nargin == 7
  new_miss = 0;
  squeeze_it = 1;
elseif nargin == 8
  squeeze_it = 1;
elseif (nargin ~= 9) & (nargin ~= 10)
  s = [ ' number of input arguments = ' int2str(nargin) ];
  disp(s)
  help getnc_s
  return
end

% I make ncmex calls to find the integers that specify the attribute
% types

nc_byte = ncmex('parameter', 'nc_byte'); %1
nc_char = ncmex('parameter', 'nc_char'); %2
nc_short = ncmex('parameter', 'nc_short'); %3
nc_long = ncmex('parameter', 'nc_long'); %4
nc_float = ncmex('parameter', 'nc_float'); %5
nc_double = ncmex('parameter', 'nc_double'); %6

% Find out whether values should be automatically rescaled or not.

if nargin == 10
  [rescale_var, rescale_att] = y_rescal(rescale_opts);
else
  [rescale_var, rescale_att] = y_rescal;
end

% Set the value of imap.  Note that this is used simply as a
% placeholder in calls to vargetg - its value is never used.

imap = 0;

% Set some constants.

blank = abs(' ');

% Check that the values of the first 2 input arguments are acceptable.

if ~isstr(file)
  error(' FILE is not a string');
end

if ~isstr(varid)
  size_var = size(varid);
  if size_var(1) ~= 1 | size_var(2) ~= 1
    error('ERROR: varid must be a scalar or a string')
  end
  if varid < 0
    error('ERROR: varid is less than zero');
  end
end

% Check that the file is accessible.  If it is then its full name will
% be stored in the variable cdf.  The file may have the extent .cdf or
% .nc and be in the current directory or the common data set (whose
% path can be found by a call to pos_cds.m.  If a compressed form
% of the file is in the current directory then the procedure gives an
% error message and exits.  This is unlike the interactive version ,
% getnc.m.  If the netcdf file is not accessible then the m file is
% exited with an error message.

cdf_list = { '' '.nc' '.cdf'};
ilim = length(cdf_list);
for i = 1:ilim 
  cdf = [ file cdf_list{i} ];
  err = check_nc(cdf);

  if (err == 0) | (err == 4)
    break;
  elseif err == 1
    if i == ilim
      error([ file ' could not be found' ]);
    end
  elseif err == 2
    path_name = pos_cds;
    cdf = [ path_name cdf ];
    break;
  elseif err == 3
    error([ 'exiting because ' cdf ' is in compressed form' ]);
  end
end

% Open the netcdf file.

[cdfid, rcode] = ncmex('OPEN', cdf, 'NC_NOWRITE');
if rcode == -1
  error(['** ERROR ** ncopen: rcode = ' num2str(rcode)])
end

% Suppress all error messages from netCDF 

[rcode] = ncmex('setopts', 0);

% Collect information about the cdf file.

[ndimens, nvars, ngatts, recdim, rcode] =  ncmex('ncinquire', cdfid);
if rcode == -1
  error(['** ERROR ** ncinquire: rcode = ' num2str(rcode)])
end

% Determine the netcdf id number for the required variable.  If varid
% is a string then an appropriate call to ncmex is used to convert it
% to the relevant integer.  Note the ugly way that the letters varid
% have 3 meanings in the one line of code.  If varid is a number then
% it is decremented.  This is done because inqnc & getnc count the
% variables from 1 to nvars whereas the calls to the ncmex routines
% use c-type counting, from 0 to nvars - 1.

if isstr(varid)
  varid = ncmex('varid', cdfid, varid);
  if rcode == -1
    error(['** ERROR ** ncvarid: rcode = ' num2str(rcode)])
  end
else
  varid = varid - 1;
end

% Check the value of varid

if varid < 0 | varid >= nvars
   error([ 'getnc_s was passed varid = ' int2str(varid) ])
end

% Find out info about the variable, in particular find nvdims, the number of
% dimensions that the variable has.

[varnam, vartypv, nvdims, vdims, nvatts, rcode] = ...
    ncmex('ncvarinq', cdfid, varid);
if rcode == -1