/[MITgcm]/MITgcm/utils/matlab/rdmnc.m
ViewVC logotype

Contents of /MITgcm/utils/matlab/rdmnc.m

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


Revision 1.13 - (show annotations) (download)
Mon Feb 19 03:43:36 2007 UTC (17 years, 3 months ago) by jmc
Branch: MAIN
Changes since 1.12: +32 -10 lines
clean-up after replacing the squeeze(); put back some debuging print ;
modified to work with singleton dimension in X ;

1 function [S] = rdmnc(varargin)
2
3 % Usage:
4 % S=RDMNC(FILE1,FILE2,...)
5 % S=RDMNC(FILE1,...,ITER)
6 % S=RDMNC(FILE1,...,'VAR1','VAR2',...)
7 % S=RDMNC(FILE1,...,'VAR1','VAR2',...,ITER)
8 %
9 % Input:
10 % FILE1 Either a single file name (e.g. 'state.nc') or a wild-card
11 % strings expanding to a group of file names (e.g. 'state.*.nc').
12 % There are no assumptions about suffices (e.g. 'state.*' works).
13 % VAR1 Model variable names as written in the MNC/netcdf file.
14 % ITER Vector of iterations in the MNC/netcdf files, not model time.
15 %
16 % Output:
17 % S Structure with fields corresponding to 'VAR1', 'VAR2', ...
18 %
19 % Description:
20 % This function is a rudimentary wrapper for joining and reading netcdf
21 % files produced by MITgcm. It does not give the same flexibility as
22 % opening the netcdf files directly using netcdf(), but is useful for
23 % quick loading of entire model fields which are distributed in multiple
24 % netcdf files.
25 %
26 % Example:
27 % >> S=rdmnd('state.*','XC','YC','T');
28 % >> imagesc( S.XC, S.YC, S.T(:,:,1)' );
29 %
30 % Author: Alistair Adcroft
31 % Modifications: Daniel Enderton
32
33 % $Header: /u/gcmpack/MITgcm/utils/matlab/rdmnc.m,v 1.12 2007/02/17 23:49:43 jmc Exp $
34 % $Name: $
35
36 % Initializations
37 dBug=0;
38 file={};
39 filepaths={};
40 files={};
41 varlist={};
42 iters=[];
43
44 % Process function arguments
45 for iarg=1:nargin;
46 arg=varargin{iarg};
47 if ischar(arg)
48 if isempty(dir(char(arg)))
49 varlist{end+1}=arg;
50 else
51 file{end+1}=arg;
52 end
53 else
54 if isempty(iters)
55 iters=arg;
56 else
57 error(['The only allowed numeric argument is iterations',...
58 ' to read in as a vector for the last arguement']);
59 end
60 end
61 end
62
63 % Create list of filenames
64 for eachfile=file
65 filepathtemp=eachfile{:};
66 indecies = find(filepathtemp=='/');
67 if ~isempty(indecies)
68 filepathtemp = filepathtemp(1:indecies(end));
69 else
70 filepathtemp = '';
71 end
72 expandedEachFile=dir(char(eachfile{1}));
73 for i=1:size(expandedEachFile,1);
74 if expandedEachFile(i).isdir==0
75 files{end+1}=expandedEachFile(i).name;
76 filepaths{end+1}=filepathtemp;
77 end
78 end
79 end
80
81
82 % If iterations unspecified, open all the files and make list of all the
83 % iterations that appear, use this iterations list for data extraction.
84 if isempty(iters)
85 iters = [];
86 for ieachfile=1:length(files)
87 eachfile = [filepaths{ieachfile},files{ieachfile}];
88 nc=netcdf(char(eachfile),'read');
89 nciters = nc{'iter'}(:);
90 if isempty(nciters), nciters = nc{'T'}(:); end
91 iters = [iters,nciters];
92 close(nc);
93 end
94 iters = unique(iters);
95 end
96
97 % Cycle through files for data extraction.
98 S.attributes=[];
99 for ieachfile=1:length(files)
100 eachfile = [filepaths{ieachfile},files{ieachfile}];
101 if dBug > 0, fprintf([' open: ',eachfile]); end
102 nc=netcdf(char(eachfile),'read');
103 S=rdmnc_local(nc,varlist,iters,S,dBug);
104 close(nc);
105 end
106
107
108 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
109 % Local functions %
110 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
111
112 function [A] = read_att(nc);
113 allatt=ncnames(att(nc));
114 if ~isempty(allatt)
115 for attr=allatt;
116 A.(char(attr))=nc.(char(attr))(:);
117 end
118 else
119 A = 'none';
120 end
121
122 function [i0,j0,fn] = findTileOffset(S);
123 fn=0;
124 if isfield(S,'attributes') & isfield(S.attributes,'global')
125 G=S.attributes.global;
126 tn=G.tile_number;
127 snx=G.sNx; sny=G.sNy; nsx=G.nSx; nsy=G.nSy; npx=G.nPx; npy=G.nPy;
128 ntx=nsx*npx;nty=nsy*npy;
129 gbi=mod(tn-1,ntx); gbj=(tn-gbi-1)/ntx;
130 i0=snx*gbi; j0=sny*gbj;
131 if isfield(S.attributes.global,'exch2_myFace')
132 fn=G.exch2_myFace;
133 end
134 else
135 i0=0;j0=0;
136 end
137 %[snx,sny,nsx,nsy,npx,npy,ntx,nty,i0,j0,fn];
138
139 function [S] = rdmnc_local(nc,varlist,iters,S,dBug)
140
141 fiter = nc{'iter'}(:); % File iterations present
142 if isempty(fiter), fiter = nc{'T'}(:); end
143 if isinf(iters); iters = fiter(end); end
144 if isnan(iters); iters = fiter; end
145 [fii,dii] = ismember(fiter,iters); fii = find(fii); % File iteration index
146 dii = dii(find(dii ~= 0)); % Data interation index
147 if dBug > 0,
148 fprintf(' ; fii='); fprintf(' %i',fii);
149 fprintf(' ; dii='); fprintf(' %i',dii); fprintf(' \n');
150 end
151
152 % No variables specified? Default to all
153 if isempty(varlist), varlist=ncnames(var(nc)); end
154
155 % Attributes for structure
156 if iters>0; S.iters_read_from_file=iters; end
157 S.attributes.global=read_att(nc);
158
159 % Read variable data
160 for ivar=1:size(varlist,2)
161
162 cvar=char(varlist{ivar});
163 if isempty(nc{cvar})
164 disp(['No such variable ''',cvar,''' in MNC file ',name(nc)]);
165 continue
166 end
167
168 dims = ncnames(dim(nc{cvar})); % Dimensions
169 if dims{1} == 'T'
170
171 if isempty(find(fii)), disp('Iters not found'); return, end
172
173 tmpdata = nc{cvar}(fii,:);
174 it = length(dims);
175 %- if only 1 time record, 1rst dim get lost; add it back:
176 % if size(nc{cvar},1) == 1,
177 if length(fii) == 1,
178 tmpdata=reshape(tmpdata,[1 size(tmpdata)]);
179 end
180 else
181 tmpdata = nc{cvar}(:);
182 it = 0;
183 end
184
185 nDims=length(size(nc{cvar}));
186 if dBug > 1,
187 fprintf([' var:',cvar,': nDims=%i ('],nDims);fprintf(' %i',size(nc{cvar}));
188 fprintf('):%i,nD=%i,it=%i ;',length(size(tmpdata)),length(dims),it);
189 end
190 if length(dims) > 1,
191 tmpdata=permute(tmpdata,[nDims:-1:1]);
192 end
193 if dBug > 1,
194 % fprintf('(tmpdata:');fprintf(' %i',size(tmpdata)); fprintf(')');
195 end
196 [ni nj nk nm nn no np]=size(tmpdata);
197
198 [i0,j0,fn]=findTileOffset(S);
199 cdim=dims{end}; if cdim(1)~='X'; i0=0; end
200 cdim=dims{end}; if cdim(1)=='Y'; i0=j0; j0=0; end
201 if length(dims)>1;
202 cdim=dims{end-1}; if cdim(1)~='Y'; j0=0; end
203 else
204 j0=0;
205 end
206 if dBug > 1,
207 fprintf(' i0,ni= %i %i; j,nj= %i %i; nk=%i :',i0,ni,j0,nj,nk);
208 end
209
210 Sstr = '';
211 for istr = 1:max(nDims,length(dims)),
212 if istr == it, Sstr = [Sstr,'dii,'];
213 elseif istr == 1, Sstr = [Sstr,'i0+(1:ni),'];
214 elseif istr == 2, Sstr = [Sstr,'j0+(1:nj),'];
215 elseif istr == 3, Sstr = [Sstr,'(1:nk),'];
216 elseif istr == 4, Sstr = [Sstr,'(1:nm),'];
217 elseif istr == 5, Sstr = [Sstr,'(1:nn),'];
218 elseif istr == 6, Sstr = [Sstr,'(1:no),'];
219 elseif istr == 7, Sstr = [Sstr,'(1:np),'];
220 else, error('Can''t handle this many dimensions!');
221 end
222 end
223 eval(['S.(cvar)(',Sstr(1:end-1),')=tmpdata;'])
224 %S.(cvar)(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata;
225 if dBug > 1, fprintf(' %i',size(S.(cvar))); fprintf('\n'); end
226
227 S.attributes.(cvar)=read_att(nc{cvar});
228 end
229
230 if isempty(S)
231 error('Something didn''t work!!!');
232 end

  ViewVC Help
Powered by ViewVC 1.1.22