| 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_contrib/enderton/Diagnostics/DiagUtility/rdmnc.m,v 1.1 2005/01/31 15:43:29 enderton 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 |
% ----------------------------------------------------------------------------- |
| 83 |
function [A] = read_att(nc); |
| 84 |
global verbose |
| 85 |
allatt=ncnames(att(nc)); if verbose; allatt, end |
| 86 |
A='none'; |
| 87 |
for attr=allatt; |
| 88 |
A.(char(attr))=nc.(char(attr))(:); |
| 89 |
end |
| 90 |
% ----------------------------------------------------------------------------- |
| 91 |
function [i0,j0,fn] = findTileOffset(S); |
| 92 |
global verbose |
| 93 |
fn=0; |
| 94 |
if isfield(S,'attributes') & isfield(S.attributes,'global') |
| 95 |
G=S.attributes.global; |
| 96 |
tn=G.tile_number; |
| 97 |
snx=G.sNx; sny=G.sNy; nsx=G.nSx; nsy=G.nSy; npx=G.nPx; npy=G.nPy; |
| 98 |
ntx=nsx*npx;nty=nsy*npy; |
| 99 |
gbi=mod(tn-1,ntx); gbj=(tn-gbi-1)/ntx; |
| 100 |
i0=snx*gbi; j0=sny*gbj; |
| 101 |
if isfield(S.attributes.global,'exch2_myFace') |
| 102 |
fn=G.exch2_myFace; |
| 103 |
end |
| 104 |
else |
| 105 |
i0=0;j0=0; |
| 106 |
end |
| 107 |
[snx,sny,nsx,nsy,npx,npy,ntx,nty,i0,j0,fn]; |
| 108 |
% ----------------------------------------------------------------------------- |
| 109 |
function [S] = rdnetcdf_local(nc,varlist,S) |
| 110 |
% Read all attributes and variable data from a netcdf file |
| 111 |
% with special operations for MNC style data |
| 112 |
global verbose |
| 113 |
|
| 114 |
% No variables specified? Default to all |
| 115 |
if isempty(varlist) |
| 116 |
varlist=ncnames(var(nc)); if verbose; varlist, end |
| 117 |
end |
| 118 |
|
| 119 |
% Attributes for structure |
| 120 |
S.attributes=read_att(nc); |
| 121 |
|
| 122 |
% Read variable data |
| 123 |
for ivar=1:size(varlist,2); |
| 124 |
cvar=char(varlist{ivar}); |
| 125 |
tmpdata=nc{cvar}(:); |
| 126 |
if isempty(tmpdata) |
| 127 |
disp(['No such variable ''' cvar ''' in netcdf file' name(nc)]) |
| 128 |
else |
| 129 |
tmpdata=squeeze(permute(tmpdata,[9:-1:1])); |
| 130 |
S.(cvar)=tmpdata; |
| 131 |
S.attributes.(cvar)=read_att(nc{cvar}); |
| 132 |
end |
| 133 |
end |
| 134 |
% ----------------------------------------------------------------------------- |
| 135 |
function [S] = rdmnc_local(nc,varlist,timelevels,S) |
| 136 |
% Read all attributes and variable data from a netcdf file |
| 137 |
% with special operations for MNC style data |
| 138 |
global verbose |
| 139 |
|
| 140 |
nt=size( nc('T'), 1); if verbose; nt, end |
| 141 |
|
| 142 |
% No variables specified? Default to all |
| 143 |
if isempty(varlist) |
| 144 |
varlist=ncnames(var(nc)); if verbose; varlist, end |
| 145 |
end |
| 146 |
|
| 147 |
% No iterations specified? Default to all |
| 148 |
if isempty(timelevels) | isnan(timelevels) |
| 149 |
timelevels=1:nt; |
| 150 |
elseif timelevels == Inf |
| 151 |
timelevels=nt; |
| 152 |
end |
| 153 |
|
| 154 |
% Sanity checks |
| 155 |
if max( find(timelevels > nt) ) |
| 156 |
error('Requested time level is beyond time dimension in netcdf file') |
| 157 |
end |
| 158 |
|
| 159 |
% Attributes for structure |
| 160 |
if timelevels>0; S.timelevels_read_from_file=timelevels; end |
| 161 |
S.attributes.global=read_att(nc); |
| 162 |
|
| 163 |
% Read variable data |
| 164 |
for ivar=1:size(varlist,2); |
| 165 |
cvar=char(varlist{ivar}); if verbose; cvar, end |
| 166 |
if isempty(nc{cvar}) |
| 167 |
disp(['No such variable ''' cvar ''' in MNC file ' name(nc)]) |
| 168 |
continue |
| 169 |
end |
| 170 |
dims=ncnames(dim(nc{cvar})); |
| 171 |
if dims{1}=='T' |
| 172 |
if verbose; disp(['Reading variable ''' cvar ''' with record dimension']); end |
| 173 |
tmpdata=nc{cvar}(timelevels,:); |
| 174 |
else |
| 175 |
if verbose; disp(['Reading variable ''' cvar '''']); end |
| 176 |
tmpdata=nc{cvar}(:); |
| 177 |
end |
| 178 |
if isempty(tmpdata) |
| 179 |
error(['No such variable ''' cvar ''' in MNC file ' name(nc)]) |
| 180 |
else |
| 181 |
tmpdata=squeeze(permute(tmpdata,[9:-1:1])); |
| 182 |
[ni nj nk nm nn no np]=size(tmpdata); |
| 183 |
[ni nj nk nm nn no np]; |
| 184 |
if np~=1; error('Wow! This is a very high dimension variable...'); end |
| 185 |
[i0,j0,fn]=findTileOffset(S); |
| 186 |
cdim=dims{end}; if cdim(1)~='X'; i0=0; end |
| 187 |
cdim=dims{end}; if cdim(1)=='Y'; i0=j0; j0=0; end |
| 188 |
if length(dims)>1; |
| 189 |
cdim=dims{end-1}; if cdim(1)~='Y'; j0=0; end |
| 190 |
else |
| 191 |
j0=0; |
| 192 |
end |
| 193 |
S.(cvar)(i0+(1:ni),j0+(1:nj),(1:nk),(1:nm),(1:nn),(1:no),(1:np))=tmpdata; |
| 194 |
S.attributes.(cvar)=read_att(nc{cvar}); |
| 195 |
end |
| 196 |
end |
| 197 |
|
| 198 |
if isempty(S) |
| 199 |
error('Something didn''t work!!!') |
| 200 |
end |
| 201 |
% ----------------------------------------------------------------------------- |
Parent Directory
| 
Revision Log
| 
Revision Graph