getnc_s.m

matlacb程序包

  % First find the indices of the data points that are outside the valid
  % range.
    
  if ~isempty(valid_range)
    
    % Check that valid_range is a 2 element vector.
    
    if length(valid_range) == 2
      miss = valid_range;
    else
      error(['The valid_range attribute must be a 2 element vector'])
    end
    
    miss_low = miss(1);
    miss_up = miss(2);
    miss_low_orig = miss_low;
    miss_up_orig = miss_up;
    
    index_miss_low = find (values < miss_low);
    index_miss_up = find (values > miss_up);
  else
    if ~isempty(valid_min)
      miss_low = valid_min;

      % Check that valid_min is a scalar
	
      if length(miss_low) ~= 1
	error(['The valid_min attribute must be a scalar'])
      end
      index_miss_low = find ( values < miss_low );
    end
    
    if ~isempty(valid_max)
      miss_up = valid_max;
  
      % Check that valid_max is a scalar
	
      if length(miss_up) ~= 1
	error(['The valid_max attribute must be a scalar'])
      end
      index_miss_up = find ( values > miss_up );
    end
  end
  
  % Now find the indices of the data points that are the same as the
  % _FillValue.
  
  if ~isempty(fillvalue)
    miss_val = fillvalue;
    
    % Check that _FillValue is a scalar
    
    if length(miss_val) ~= 1
      error(['The _FillValue attribute must be a scalar'])
    end
 
    % Check whether _FillValue is outside the valid range to decide
    % whether to keep going.
    
    keep_going = 1;
    if ~isempty(miss_low_orig)
      if (miss_val < miss_low_orig)
	keep_going = 0;
      end
    end
    if ~isempty(miss_up_orig)
      if (miss_val > miss_up_orig)
	keep_going = 0;
      end
    end
    
    if keep_going
      if var_is_number
	if miss_val < 0
	  miss_val_low = 1.00001*miss_val;
	  miss_val_up = 0.99999*miss_val;
	else
	  miss_val_low = 0.99999*miss_val;
	  miss_val_up = 1.00001*miss_val;
	end
	index_fillvalue = find((miss_val_low <= values) & ...
				(values <= miss_val_up));
      else
	index_fillvalue = find(values == miss_val);
      end
    end
  end

  % Now find the indices of the data points that are the same as the
  % missing_value.
  
  if ~isempty(missing_value)
    miss_val = missing_value;
    
    % Check that missing_value is a scalar
    
    if length(miss_val) ~= 1
      error(['The missing_value attribute must be a scalar'])
    end
 
    % Check whether missing_value is outside the valid range to decide
    % whether to keep going. Also, don't bother doing the search if
    % _FillValue is identical to missing_value.
    
      keep_going = 1;
      if ~isempty(miss_low_orig)
	if (miss_val < miss_low_orig )
	  keep_going = 0;
	end
      end
      if ~isempty(miss_up_orig)
	if (miss_val > miss_up_orig )
	  keep_going = 0;
	end
      end
    if ~isempty(fillvalue)
      if (missing_value == fillvalue)
	keep_going = 0;
      end
    end
    
    if keep_going
      if var_is_number
	if miss_val < 0
	  miss_val_low = 1.00001*miss_val;
	  miss_val_up = 0.99999*miss_val;
	else
	  miss_val_low = 0.99999*miss_val;
	  miss_val_up = 1.00001*miss_val;
	end
	index_missing_value = find((miss_val_low <= values) & ...
				(values <= miss_val_up));
      else
	index_missing_value = find(values == miss_val);
      end
    end
  end
  
  % Combine the arrays of missing value indices into one unordered array.
  % Note that for real numbers the range of the _FillValue and
  % missing_value may intersect both the valid and invalid range and so
  % some indices may appear twice; this does not cause any inaccuracy,
  % although it will result in some inefficiency.
  
  index_miss = [ index_miss_low(:); index_fillvalue(:); ...
		 index_missing_value(:); index_miss_up(:) ];
  len_index_miss = length(index_miss);
  
  if len_index_miss > 0
    if change_miss == 2
      if isnumeric(values)
	values(index_miss) = NaN;
      else
	values(index_miss) = char(0);
      end
    elseif change_miss == 3
      if isnumeric(values)
	values(index_miss) = new_miss;
      else
	if isnumeric(new_miss)
	  values(index_miss) = char(new_miss);
	else
	  values(index_miss) = new_miss;
	end
      end
    elseif change_miss ~= 1
      s = [ 'getnc_s was passed change_miss = ' int2str(change_miss) ];
      error(s)
    end
  end
 case 'java'
  
  % Translate rescale_opts to a 2 element vector.
  
  if length(rescale_opts) == 1
    if rescale_opts == -1
      rescale_opts = [1 1];
    else
      error('rescale_opts must be a 2 element vector or -1')
    end
  end
      
  % Open dataset so that this script will have to do the
  % scale/offset/missing value stuff. This is because the java
  % interface does not deal with missing values properly for character
  % arrays and rescale_opts may force some things to be rescaled and others
  % not. Note that the 2nd argument relates to enhancing
  % the data, i.e., automatically doing scale/offset/missing stuff. The 3rd
  % argument relates to error handling and hopefully is the way to pass a
  % null to the constructor.

  try
    enhance = false;
    ncdJ = ucar.nc2.dataset.NetcdfDataset.openDataset(file, enhance, []);
  catch
    ss = lasterror;
    mess_str = ss.message;
    rcode = -1000000;
    values = error_handle('java', mess_str, rcode, err_opt);
    return
  end
  
  % Get the variable object. Note that this is different according to whether
  % whether varid is numeric or not.
  
  if isnumeric(varid)
    
    % Check that varid is a scalar.
    
    if ndims(varid) > 2
      wrong_size = 1;
    else     
      if all(size(varid) == [1 1])
	wrong_size = 0;
      else
	wrong_size = 1;
      end
    end
    if wrong_size
      mess_str = 'varid must be either a string or a scalar';
      rcode = -1000000;
      values = error_handle(ncdJ, mess_str, rcode, err_opt);
    end
    
    % Check that varid is in the correct range.
    
    variables_list = ncdJ.getVariables();
    nvars = variables_list.size();
    if (varid < 1) || (varid > nvars)
      mess_str = ['Failed trying to access variable number: ' varid];
      rcode = -1000000;
      values = error_handle(ncdJ, mess_str, rcode, err_opt);
      return
    end
    
    % Get the variable object and reset varid to be the name of the variable.
    
    varid_old = varid;
    varJ = variables_list.get(varid_old - 1);
    varid = varJ.getName();
  else
    
    % Try to get the variable object and then check that this succeeded.
    
    varJ = ncdJ.findVariable(varid);
    if isempty(varJ)
      mess_str = ['Failed trying to access variable: ' varid];
      rcode = -1000000;
      values = error_handle(ncdJ, mess_str, rcode, err_opt);
      return
    end
  end
  
  % There are problems translating from netcdf to opendap formats when
  % dealing with character arrays. If varJ.getDataType returns 'char' then an
  % array of char is returned and it can be handled pretty much the same as a
  % numeric array although there are issues with missing values and switching
  % between matlabs numeric and char types. However the opendap data may be
  % seen as an array of strings. In this case varJ.getDataType will be
  % 'String' and things must be handled differently. This is complicated even
  % more by the issue of what is the fastest changing index - hence calls to
  % fliplr and permute.
  
  var_type_str = varJ.getDataType;
  
  % Get the data. We may either get all of the data or only some subset.
  
  if all(bl_corner == -1) && all(tr_corner == -1) && all(stride == -1)
    
    % Get all of variable.
    
    if strcmp(var_type_str, 'String')
      
      % Get the character data as a 1D array. Note all of the fiddling with
      % val_shape_row_orig and val_shape_row_reversed to deal with the
      % difference of fastest changing index between matlab and C.
      
      val_size = varJ.getSize();
      val_rank = varJ.getRank();
      val_shape = varJ.getShape();
      val_shape_row_orig = val_shape(:)';
      val_shape_row_reversed = fliplr(val_shape_row_orig);
      val_int = copyTo1DJavaArray(varJ.read());
      switch class(val_int(1))
       case 'opendap.dap.DString'
	
	% Read the string data allowing for weirdness of the returned
        % value. In particular it may call a null character an empty
        % string. Thus we must read the data one character at a time.
	
	len_val_int = length(val_int);
	values = repmat(' ', len_val_int, 1);
	for ii = 1:len_val_int
	  xx = val_int(ii).getValue();
	  if length(xx) == 0
	    values(ii) = char(0);
	  else
	    values(ii) = char(xx);
	  end
	end
	
	% Turn values into a properly dimensioned array (fastest changing
        % index problem).
	
	if val_rank > 1
	  values = reshape(values, val_shape_row_reversed);
	  values = permute(values, val_rank:-1:1);
	end
	
       otherwise
	error('Cannot handle object')
      end
    else
      values = copyToNDJavaArray(varJ.read());
    end
  else
    
    % Change bl_corner, tr_corner and stride so that they contain the
    % hyperslab description, i.e., replace the -1 values with proper indices
    % and then generate the string readSpec. This contains the fortran-like 
    % information readSpec = '0:7:2,3:8:3';
    
    varRank = varJ.getRank();
    varShape(1:varRank) = varJ.getShape();
    if varRank > 1
      if length(bl_corner) == 1
	bl_corner = repmat(bl_corner(1), varRank, 1);
      end
      if length(tr_corner) == 1
	tr_corner = repmat(tr_corner(1), varRank, 1);
      end
      if length(stride) == 1
	stride = repmat(stride(1), varRank, 1);
      end
    end
    for ii = 1:varRank
      if bl_corner(ii) < 0
	bl_corner(ii) = 1;
	stride(ii) = 1;
	tr_corner(ii) = varShape(ii);
      elseif (bl_corner(ii) == 0) || (bl_corner(ii) > varShape(ii))
	mess_str = '1: hyperslab is badly specified';
	rcode = -1000000;
	values = error_handle(ncdJ, mess_str, rcode, err_opt);
	return
      end
      if tr_corner(ii) == -1
	tr_corner(ii) = varShape(ii);
      elseif (tr_corner(ii) < bl_corner(ii)) || (tr_corner(ii) > varShape(ii))
	mess_str = '2: hyperslab is badly specified';
	rcode = -1000000;
	values = error_handle(ncdJ, mess_str, rcode, err_opt);
	return
      end
    end
    if any(stride < -1) || any(stride == 0)
      mess_str = 'bad stride';
      rcode = -1000000;
      values = error_handle(ncdJ, mess_str, rcode, err_opt);
      return
    end
    stride(find(stride == -1)) = 1;
      
    % Construct the string that specifies the hyperslab and then get the data.
    readSpec = '';
    for ii = 1:varRank
      readSpec = [readSpec num2str(bl_corner(ii) - 1) ':' ...
		  num2str(tr_corner(ii) - 1) ':' ...
		  num2str(stride(ii)) ','];
    end
    readSpec = readSpec(1:(end - 1));
    if strcmp(var_type_str, 'String')
      val_int = copyTo1DJavaArray(varJ.read(readSpec));
      switch class(val_int(1))
       case 'opendap.dap.DString'
	
	% Read the string data allowing for weirdness of the returned
        % value. In particular it may call a null character an empty
        % string. Thus we must read the data one character at a time.
	
	len_val_int = length(val_int);
	values = repmat(' ', len_val_int, 1);
	for ii = 1:len_val_int
	  xx = val_int(ii).getValue();
	  if length(xx) == 0
	    values(ii) = char(0);
	  else
	    values(ii) = char(xx);
	  end
	end
	
	% Turn into properly dimensioned array (fastest changing index
        % problem).
	
	val_rank = varJ.getRank();
	lims_dim = floor((tr_corner(:) - bl_corner(:))./stride(:)) + 1;
	lims_dim = fliplr(lims_dim');
	if val_rank > 1
	  values = reshape(values, lims_dim);
	  values = permute(values, val_rank:-1:1);
	end
	
       otherwise
	error('Cannot handle object')
      end
    else
      values = copyToNDJavaArray(varJ.read(readSpec));
    end
  end
  
  % Make sure that all numeric arrays are returned as double.
  
  if isnumeric(values)
    values = double(values);
  end

  % Possibly squeeze the array.
  
  if squeeze_it
    values = squeeze(values);
  end
  
  % Possibly change the order of the elements according to the value of the
  % variable "order".

  size_val = size(values);
  nd = length(size_val);
  if all(size(order) == 1)
    % order is a scalar
    if order == -2
      values = permute(values, (nd:-1:1));
    end