/[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.2 - (show annotations) (download)
Wed Feb 2 10:31:22 2005 UTC (19 years, 4 months ago) by mlosch
Branch: MAIN
CVS Tags: checkpoint57o_post, checkpoint57m_post, checkpoint57s_post, checkpoint57k_post, checkpoint57d_post, checkpoint57g_post, checkpoint57i_post, checkpoint57e_post, checkpoint57g_pre, checkpoint57f_pre, checkpoint57r_post, eckpoint57e_pre, checkpoint57h_done, checkpoint57n_post, checkpoint57p_post, checkpoint57f_post, checkpoint57q_post, checkpoint57j_post, checkpoint57h_pre, checkpoint57l_post, checkpoint57h_post
Changes since 1.1: +15 -7 lines
o make rdmnc also work with Matlab 6.1
  -not pretty, but now it works for both Matlab 6 and 7

1 function [S] = rdmnc(varargin)
2 % S=RDMNC(FILE1,FILE2,...)
3 % S=RDMNC(FILE1,...,ITER)
4 % S=RDMNC(FILE1,...,'VAR1','VAR2',...)
5 % S=RDMNC(FILE1,...,'VAR1','VAR2',...,ITER)
6 %
7 % FILE1 ... is either a single file name (e.g. 'state.nc') or a wild-card
8 % strings expanding to a group of file names (e.g. 'state.*.nc').
9 % There are no assumptions about suffices (e.g. 'state.*' works)
10 % VAR1 ... are model variable names as written in the MNC/netcdf file
11 % ITER is a vector of time levels (consecutive indices referring to snap-shot
12 % in the MNC/netcdf file i.e. 1,2,3,... and not model time)
13 % S is a structure with mutliple fields
14 %
15 % rdmnc() is a rudimentary wrapper for joining and reading netcdf files
16 % produced by MITgcm. It does not give the same flexibility as opening the
17 % netcdf files directly using netcdf(). It is useful for quick loading of
18 % entire model fields which are distributed in multiple netcdf files.
19 % rdmnc() also reads non-MITgcm generated netcdf files.
20 %
21 % e.g. To plot the surface temperature in last time-level written to file
22 % >> S=rdmnd('state.*','XC','YC','RC','T',Inf);
23 % >> imagesc( S.XC, S.YC, S.T(:,:,1)' );
24 %
25 % $Header: /u/gcmpack/MITgcm/utils/matlab/rdmnc.m,v 1.1 2004/03/08 16:45:56 adcroft Exp $
26
27 %Defaults
28 global verbose
29 file={};
30 files={};
31 verbose=0;
32 varlist={};
33 timelevels=[];
34
35 % Process function arguments
36 for iarg=1:nargin;
37 arg=varargin{iarg};
38 if ischar(arg)
39 switch char(arg)
40 case 'verbose',
41 verbose=1;
42 otherwise,
43 if isempty( dir(char(arg)) )
44 varlist{end+1}={arg};
45 else
46 file{end+1}={arg};
47 end
48 end
49 else
50 if isempty(timelevels)
51 timelevels=arg;
52 else
53 error('Specify the time levels in a vector and not as multiple numeric arguments')
54 end
55 end
56 end
57
58 if verbose; disp('Verbose mode enabled'); end
59
60 % Create list of filenames
61 for eachfile=file
62 expandedEachFile=dir(char(eachfile{1}));
63 for i=1:size(expandedEachFile,1);
64 if expandedEachFile(i).isdir==0; files{end+1}=expandedEachFile(i).name; end
65 end
66 end
67
68 % Open each file
69 S.attributes=[];
70 for eachfile={files{1:end}}
71 nc=netcdf(char(eachfile),'read');
72 if ismncfile(nc)
73 S=rdmnc_local(nc,varlist,timelevels,S);
74 else
75 S=rdnetcdf_local(nc,varlist,S);
76 end
77 end
78
79 % -----------------------------------------------------------------------------
80 function [result] = ismncfile(nc);
81 result=~isempty(nc.MITgcm_mnc_ver);
82 %MLresult=~isempty(nc.('MITgcm_mnc_ver'));
83 % -----------------------------------------------------------------------------
84 function [A] = read_att(nc);
85 global verbose
86 allatt=ncnames(att(nc)); if verbose; allatt, end
87 A='none';
88 for attr=allatt;
89 tmp = char(attr);
90 eval(['A.' tmp '= nc.' tmp '(:);'])
91 %ML A.(char(attr))=nc.(char(attr))(:);
92 end
93 % -----------------------------------------------------------------------------
94 function [i0,j0,fn] = findTileOffset(S);
95 global verbose
96 fn=0;
97 if isfield(S,'attributes') & isfield(S.attributes,'global')
98 G=S.attributes.global;
99 tn=G.tile_number;
100 snx=G.sNx; sny=G.sNy; nsx=G.nSx; nsy=G.nSy; npx=G.nPx; npy=G.nPy;
101 ntx=nsx*npx;nty=nsy*npy;
102 gbi=mod(tn-1,ntx); gbj=(tn-gbi-1)/ntx;
103 i0=snx*gbi; j0=sny*gbj;
104 if isfield(S.attributes.global,'exch2_myFace')
105 fn=G.exch2_myFace;
106 end
107 else
108 i0=0;j0=0;
109 end
110 % -----------------------------------------------------------------------------
111 function [S] = rdnetcdf_local(nc,varlist,S)
112 % Read all attributes and variable data from a netcdf file
113 % with special operations for MNC style data
114 global verbose
115
116 % No variables specified? Default to all
117 if isempty(varlist)
118 varlist=ncnames(var(nc)); if verbose; varlist, end
119 end
120
121 % Attributes for structure
122 S.attributes=read_att(nc);
123
124 % Read variable data
125 for ivar=1:size(varlist,2);
126 cvar=char(varlist{ivar});
127 tmpdata=nc{cvar}(:);
128 if isempty(tmpdata)
129 disp(['No such variable ''' cvar ''' in netcdf file' name(nc)])
130 else
131 tmpdata=squeeze(permute(tmpdata,[9:-1:1]));
132 eval(['S.' cvar '=tmpdata;'])
133 eval(['S.attributes.' cvar '=read_att(nc{' cvar '});'])
134 %ML S.(cvar)=tmpdata;
135 %ML S.attributes.(cvar)=read_att(nc{cvar});
136 end
137 end
138 % -----------------------------------------------------------------------------
139 function [S] = rdmnc_local(nc,varlist,timelevels,S)
140 % Read all attributes and variable data from a netcdf file
141 % with special operations for MNC style data
142 global verbose
143
144 nt=size( nc('T'), 1); if verbose; nt, end
145
146 % No variables specified? Default to all
147 if isempty(varlist)
148 varlist=ncnames(var(nc)); if verbose; varlist, end
149 end
150
151 % No iterations specified? Default to all
152 if isempty(timelevels) | isnan(timelevels)
153 timelevels=1:nt;
154 elseif timelevels == Inf
155 timelevels=nt;
156 end
157
158 % Sanity checks
159 if max( find(timelevels > nt) )
160 error('Requested time level is beyond time dimension in netcdf file')
161 end
162
163 % Attributes for structure
164 if timelevels>0; S.timelevels_read_from_file=timelevels; end
165 S.attributes.global=read_att(nc);
166
167 % Read variable data
168 for ivar=1:size(varlist,2);
169 cvar=char(varlist{ivar}); if verbose; cvar, end
170 if isempty(nc{cvar})
171 disp(['No such variable ''' cvar ''' in MNC file ' name(nc)])
172 continue
173 end
174 dims=ncnames(dim(nc{cvar}));
175 if dims{1}=='T'
176 if verbose; disp(['Reading variable ''' cvar ''' with record dimension']); end
177 tmpdata=nc{cvar}(timelevels,:);
178 else
179 if verbose; disp(['Reading variable ''' cvar '''']); end
180 tmpdata=nc{cvar}(:);
181 end
182 if isempty(tmpdata)
183 error(['No such variable ''' cvar ''' in MNC file ' name(nc)])
184 else
185 tmpdata=squeeze(permute(tmpdata,[9:-1:1]));
186 [ni nj nk nm nn no np]=size(tmpdata);
187 if np~=1; error('Wow! This is a very high dimension variable...'); end
188 [i0,j0,fn]=findTileOffset(S);
189 cdim=dims{end}; if cdim(1)~='X'; i0=0; end
190 cdim=dims{end}; if cdim(1)=='Y'; i0=j0; j0=0; end
191 if length(dims)>1;
192 cdim=dims{end-1}; if cdim(1)~='Y'; j0=0; end
193 else
194 j0=0;
195 end
196 eval(['S.' cvar ...
197 '(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata;'])
198 eval(['S.attributes.' cvar ' =read_att(nc{''' cvar '''});'])
199 %ML S.(cvar)(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata;
200 %ML S.attributes.(cvar)=read_att(nc{cvar});
201 end
202 end
203
204 if isempty(S)
205 error('Something didn''t work!!!')
206 end
207 % -----------------------------------------------------------------------------

  ViewVC Help
Powered by ViewVC 1.1.22