/[MITgcm]/MITgcm_contrib/gael/matlab_class/gcmfaces_IO/ncload.m
ViewVC logotype

Annotation of /MITgcm_contrib/gael/matlab_class/gcmfaces_IO/ncload.m

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.1 - (hide annotations) (download)
Sun Jan 24 17:05:13 2016 UTC (9 years, 5 months ago) by gforget
Branch: MAIN
CVS Tags: checkpoint65x, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint66f, checkpoint66e, checkpoint66d, checkpoint66c, checkpoint66b, checkpoint66a, checkpoint66o, HEAD
- move routines from MITprof/profiles_misc/ to gcmfaces/gcmfaces_IO/
  to allow for nctiles IO without having to rely on MITprof
- these wrapper routines were originally provided by F. Roquet
  that work both with the native implememtation of ncetcf in recent
  Matlab versions and the third party mex codes used earlier

1 gforget 1.1 function [] = ncload(fileIn, varargin);
2    
3     % ncload -- Load NetCDF variables.
4     % ncload('fileIn', 'var1', 'var2', ...) loads the
5     % given variables of 'fileIn' into the Matlab
6     % workspace of the "caller" of this routine. If no names
7     % are given, all variables are loaded.
8    
9     global useNativeMatlabNetcdf;
10     if isempty(useNativeMatlabNetcdf); useNativeMatlabNetcdf = ~isempty(which('netcdf.open')); end;
11    
12     f = ncopen(fileIn, 'nowrite');
13     if isempty(f), return, end
14     vars=ncvars(f);
15     if isempty(varargin); varargin = vars; end;
16    
17     if (useNativeMatlabNetcdf);
18    
19     for i = 1:length(varargin)
20     if sum(strcmp(vars,varargin{i}))>0;
21     %get variable
22     varid = netcdf.inqVarID(f,varargin{i});
23     aa=netcdf.getVar(f,varid);
24     %inverse the order of dimensions
25     bb=length(size(aa)); aa=permute(aa,[bb:-1:1]);
26     %replace missing value with NaN
27     [atts]=ncatts(f,varid);
28     if strcmp(atts,'missing_value')&isreal(aa);
29     spval = double(netcdf.getAtt(f,varid,'missing_value'));
30     elseif strcmp(atts,'_FillValue')&isreal(aa);
31     spval = double(netcdf.getAtt(f,varid,'_FillValue'));
32     else;
33     spval=[];
34     end;
35     if ~isempty(spval); aa(aa==spval)=NaN; end;
36     else;
37     aa=[];
38     end;
39     assignin('caller', varargin{i}, aa)
40     end
41    
42    
43    
44    
45     else;%try to use old mex stuff
46    
47    
48     for i = 1:length(varargin)
49     if ~isstr(varargin{i}), varargin{i} = inputname(i+1); end
50     oldfld = f{varargin{i}}(:);
51     spval = f{varargin{i}}.missing_value(:);
52     if isempty(spval);
53     spval = f{varargin{i}}.FillValue_(:);
54     end
55     fld = oldfld;
56     if ~isempty(spval)&~ischar(fld);
57     replace = find(oldfld == spval);
58     nreplace = length(replace);
59     if nreplace>0
60     fld(replace) = NaN*ones(1,nreplace);
61     end %if
62     end %if
63     %NaN-substitution messes up backward compatibility so I comment it out
64     % assignin('caller', varargin{i}, fld);
65     %and revert to oldfld
66     assignin('caller', varargin{i}, oldfld)
67     end
68    
69     end
70    
71     ncclose(f);
72    

  ViewVC Help
Powered by ViewVC 1.1.22