/[MITgcm]/MITgcm_contrib/enderton/Diagnostics/DiagUtility/rdmnc_mod.m
ViewVC logotype

Contents of /MITgcm_contrib/enderton/Diagnostics/DiagUtility/rdmnc_mod.m

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


Revision 1.1 - (show annotations) (download)
Mon Jan 31 15:43:29 2005 UTC (20 years, 5 months ago) by enderton
Branch: MAIN
CVS Tags: HEAD
 o Initial check in.

1 function [S] = rdmnc_mod(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 filepaths={};
31 files={};
32 verbose=0;
33 varlist={};
34 timelevels=[];
35
36 % Process function arguments
37 for iarg=1:nargin;
38 arg=varargin{iarg};
39 if ischar(arg)
40 switch char(arg)
41 case 'verbose',
42 verbose=1;
43 otherwise,
44 if isempty( dir(char(arg)) )
45 varlist{end+1}={arg};
46 else
47 file{end+1}={arg};
48 end
49 end
50 else
51 if isempty(timelevels)
52 timelevels=arg;
53 else
54 error('Specify the time levels in a vector and not as multiple numeric arguments')
55 end
56 end
57 end
58
59 if verbose; disp('Verbose mode enabled'); end
60
61 % Create list of filenames
62 for eachfile=file
63 filepathtemp=eachfile{:}{:};
64 indecies = find(filepathtemp=='/');
65 if ~isempty(indecies)
66 filepathtemp = filepathtemp(1:indecies(end));
67 else
68 filepathtemp = '';
69 end
70 expandedEachFile=dir(char(eachfile{1}));
71 for i=1:size(expandedEachFile,1);
72 if expandedEachFile(i).isdir==0
73 files{end+1}=expandedEachFile(i).name;
74 filepaths{end+1}=filepathtemp;
75 end
76 end
77 end
78
79 % Open each file
80 S.attributes=[];
81 for ieachfile=1:length(files)
82 eachfile = [filepaths{ieachfile},files{ieachfile}];
83 nc=netcdf(char(eachfile),'read');
84 if ismncfile(nc)
85 S=rdmnc_local(nc,varlist,timelevels,S);
86 else
87 S=rdnetcdf_local(nc,varlist,S);
88 end
89 close(nc);
90 end
91
92 % -----------------------------------------------------------------------------
93 function [result] = ismncfile(nc);
94 result=~isempty(nc.('MITgcm_mnc_ver'));
95 % -----------------------------------------------------------------------------
96 function [A] = read_att(nc);
97 global verbose
98 allatt=ncnames(att(nc)); if verbose; allatt, end
99 A='none';
100 for attr=allatt;
101 A.(char(attr))=nc.(char(attr))(:);
102 end
103 % -----------------------------------------------------------------------------
104 function [i0,j0,fn] = findTileOffset(S);
105 global verbose
106 fn=0;
107 if isfield(S,'attributes') & isfield(S.attributes,'global')
108 G=S.attributes.global;
109 tn=G.tile_number;
110 snx=G.sNx; sny=G.sNy; nsx=G.nSx; nsy=G.nSy; npx=G.nPx; npy=G.nPy;
111 ntx=nsx*npx;nty=nsy*npy;
112 gbi=mod(tn-1,ntx); gbj=(tn-gbi-1)/ntx;
113 i0=snx*gbi; j0=sny*gbj;
114 if isfield(S.attributes.global,'exch2_myFace')
115 fn=G.exch2_myFace;
116 end
117 else
118 i0=0;j0=0;
119 end
120 [snx,sny,nsx,nsy,npx,npy,ntx,nty,i0,j0,fn];
121 % -----------------------------------------------------------------------------
122 function [S] = rdnetcdf_local(nc,varlist,S)
123 % Read all attributes and variable data from a netcdf file
124 % with special operations for MNC style data
125 global verbose
126
127 % No variables specified? Default to all
128 if isempty(varlist)
129 varlist=ncnames(var(nc)); if verbose; varlist, end
130 end
131
132 % Attributes for structure
133 S.attributes=read_att(nc);
134
135 % Read variable data
136 for ivar=1:size(varlist,2);
137 cvar=char(varlist{ivar});
138 tmpdata=nc{cvar}(:);
139 if isempty(tmpdata)
140 disp(['No such variable ''' cvar ''' in netcdf file' name(nc)])
141 else
142 tmpdata=squeeze(permute(tmpdata,[9:-1:1]));
143 S.(cvar)=tmpdata;
144 S.attributes.(cvar)=read_att(nc{cvar});
145 end
146 end
147 % -----------------------------------------------------------------------------
148 function [S] = rdmnc_local(nc,varlist,timelevels,S)
149 % Read all attributes and variable data from a netcdf file
150 % with special operations for MNC style data
151 global verbose
152
153 nt=size( nc('T'), 1); if verbose; nt, end
154
155 % No variables specified? Default to all
156 if isempty(varlist)
157 varlist=ncnames(var(nc)); if verbose; varlist, end
158 end
159
160 % No iterations specified? Default to all
161 if isempty(timelevels) | isnan(timelevels)
162 timelevels=1:nt;
163 elseif timelevels == Inf
164 timelevels=nt;
165 end
166
167 % Sanity checks
168 if max( find(timelevels > nt) )
169 error('Requested time level is beyond time dimension in netcdf file')
170 end
171
172 % Attributes for structure
173 if timelevels>0; S.timelevels_read_from_file=timelevels; end
174 S.attributes.global=read_att(nc);
175
176 % Read variable data
177 for ivar=1:size(varlist,2);
178 cvar=char(varlist{ivar}); if verbose; cvar, end
179 if isempty(nc{cvar})
180 disp(['No such variable ''' cvar ''' in MNC file ' name(nc)])
181 continue
182 end
183 dims=ncnames(dim(nc{cvar}));
184 if dims{1}=='T'
185 if verbose; disp(['Reading variable ''' cvar ''' with record dimension']); end
186 tmpdata=nc{cvar}(timelevels,:);
187 else
188 if verbose; disp(['Reading variable ''' cvar '''']); end
189 tmpdata=nc{cvar}(:);
190 end
191 if isempty(tmpdata)
192 error(['No such variable ''' cvar ''' in MNC file ' name(nc)])
193 else
194 tmpdata=squeeze(permute(tmpdata,[9:-1:1]));
195 [ni nj nk nm nn no np]=size(tmpdata);
196 [ni nj nk nm nn no np];
197 if np~=1; error('Wow! This is a very high dimension variable...'); end
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 S.(cvar)(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata;
207 S.attributes.(cvar)=read_att(nc{cvar});
208 end
209 end
210
211 if isempty(S)
212 error('Something didn''t work!!!')
213 end
214 % -----------------------------------------------------------------------------

  ViewVC Help
Powered by ViewVC 1.1.22