/[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.3 - (show annotations) (download)
Thu Sep 22 04:43:34 2005 UTC (18 years, 10 months ago) by enderton
Branch: MAIN
CVS Tags: checkpoint57t_post, checkpoint57v_post, checkpint57u_post
Changes since 1.2: +169 -180 lines
Modification enabling:
(1)  Reading files from directories other than where rdmnc is present.
(2)  Reading field where one iteration spans multiple files.

1 function [S] = rdmnc_mod2(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 % Original Author: Alistair Adcroft
31 % Modifications: Daniel Enderton
32
33 % Initializations
34 file={};
35 filepaths={};
36 files={};
37 varlist={};
38 iters=[];
39
40 % Process function arguments
41 for iarg=1:nargin;
42 arg=varargin{iarg};
43 if ischar(arg)
44 if isempty(dir(char(arg)))
45 varlist{end+1}=arg;
46 else
47 file{end+1}=arg;
48 end
49 else
50 if isempty(iters)
51 iters=arg;
52 else
53 error(['The only allowed numeric argument is iterations',...
54 ' to read in as a vector for the last arguement']);
55 end
56 end
57 end
58
59 % Create list of filenames
60 for eachfile=file
61 filepathtemp=eachfile{:};
62 indecies = find(filepathtemp=='/');
63 if ~isempty(indecies)
64 filepathtemp = filepathtemp(1:indecies(end));
65 else
66 filepathtemp = '';
67 end
68 expandedEachFile=dir(char(eachfile{1}));
69 for i=1:size(expandedEachFile,1);
70 if expandedEachFile(i).isdir==0
71 files{end+1}=expandedEachFile(i).name;
72 filepaths{end+1}=filepathtemp;
73 end
74 end
75 end
76
77
78 % If iterations unspecified, open all the files and make list of all the
79 % iterations that appear, use this iterations list for data extraction.
80 if isempty(iters)
81 iters = [];
82 for ieachfile=1:length(files)
83 eachfile = [filepaths{ieachfile},files{ieachfile}];
84 nc=netcdf(char(eachfile),'read');
85 iters = [iters,nc{'iter'}(:)];
86 close(nc);
87 end
88 iters = unique(iters);
89 end
90
91 % Cycle through files for data extraction.
92 S.attributes=[];
93 for ieachfile=1:length(files)
94 eachfile = [filepaths{ieachfile},files{ieachfile}];
95 nc=netcdf(char(eachfile),'read');
96 S=rdmnc_local(nc,varlist,iters,S);
97 close(nc);
98 end
99
100 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
101 % Local functions %
102 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
103
104 function [A] = read_att(nc);
105 allatt=ncnames(att(nc));
106 if ~isempty(allatt)
107 for attr=allatt;
108 A.(char(attr))=nc.(char(attr))(:);
109 end
110 else
111 A = 'none';
112 end
113
114 function [i0,j0,fn] = findTileOffset(S);
115 fn=0;
116 if isfield(S,'attributes') & isfield(S.attributes,'global')
117 G=S.attributes.global;
118 tn=G.tile_number;
119 snx=G.sNx; sny=G.sNy; nsx=G.nSx; nsy=G.nSy; npx=G.nPx; npy=G.nPy;
120 ntx=nsx*npx;nty=nsy*npy;
121 gbi=mod(tn-1,ntx); gbj=(tn-gbi-1)/ntx;
122 i0=snx*gbi; j0=sny*gbj;
123 if isfield(S.attributes.global,'exch2_myFace')
124 fn=G.exch2_myFace;
125 end
126 else
127 i0=0;j0=0;
128 end
129 %[snx,sny,nsx,nsy,npx,npy,ntx,nty,i0,j0,fn];
130
131 function [S] = rdmnc_local(nc,varlist,iters,S)
132
133 fiter = nc{'iter'}(:); % File iterations present
134 [fii,dii] = ismember(fiter,iters); fii = find(fii); % File iteration index
135 dii = dii(find(dii ~= 0)); % Data interation index
136
137 % No variables specified? Default to all
138 if isempty(varlist), varlist=ncnames(var(nc)); end
139
140 % Attributes for structure
141 if iters>0; S.iters_read_from_file=iters; end
142 S.attributes.global=read_att(nc);
143
144 % Read variable data
145 for ivar=1:size(varlist,2);
146 cvar=char(varlist{ivar});
147 if isempty(nc{cvar})
148 disp(['No such variable ''',cvar,''' in MNC file ',name(nc)]);
149 continue
150 end
151
152 dims=ncnames(dim(nc{cvar})); % Dimensions
153
154 if dims{1}=='T'
155 if isempty(find(fii)), return, end
156 tmpdata=nc{cvar}(find(fii),:);
157 it = length(dims);
158 else
159 tmpdata=nc{cvar}(:);
160 it = 0;
161 end
162
163 tmpdata=squeeze(permute(tmpdata,[9:-1:1]));
164 [ni nj nk nm nn no np]=size(tmpdata);
165 [ni nj nk nm nn no np];
166
167 [i0,j0,fn]=findTileOffset(S);
168 cdim=dims{end}; if cdim(1)~='X'; i0=0; end
169 cdim=dims{end}; if cdim(1)=='Y'; i0=j0; j0=0; end
170 if length(dims)>1;
171 cdim=dims{end-1}; if cdim(1)~='Y'; j0=0; end
172 else
173 j0=0;
174 end
175
176 Sstr = '';
177 for istr = 1:7
178 if istr == it, Sstr = [Sstr,'dii,'];
179 elseif istr == 1, Sstr = [Sstr,'i0+(1:ni),'];
180 elseif istr == 2, Sstr = [Sstr,'j0+(1:nj),'];
181 elseif istr == 3, Sstr = [Sstr,'(1:nk),'];
182 elseif istr == 4, Sstr = [Sstr,'(1:nm),'];
183 elseif istr == 5, Sstr = [Sstr,'(1:nn),'];
184 elseif istr == 6, Sstr = [Sstr,'(1:no),'];
185 elseif istr == 7, Sstr = [Sstr,'(1:np),'];
186 else, error('Can''t handle this many dimensions!');
187 end
188 end
189 eval(['S.(cvar)(',Sstr(1:end-1),')=tmpdata;'])
190 %S.(cvar)(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata;
191 S.attributes.(cvar)=read_att(nc{cvar});
192 end
193
194 if isempty(S)
195 error('Something didn''t work!!!');
196 end

  ViewVC Help
Powered by ViewVC 1.1.22