--- MITgcm_contrib/enderton/Diagnostics/DiagLoad.m 2005/02/02 16:20:12 1.5 +++ MITgcm_contrib/enderton/Diagnostics/DiagLoad.m 2005/07/05 18:57:48 1.9 @@ -1,7 +1,7 @@ function [data,xax,yax,time,pltslc] = ... DiagLoad(fln,trl,dat,dad,grd,itr,tst,flu,ddf,gdf,avg,slc,pst,... - LoadGridData,DiagDebug,GridSuffix,ZcordFile,Index,Dim,... - Grads,Year0Iter,Months); + LoadGridData,DiagDebug,GridSuffix,ZcordFile,Index,Dim,Vector,Mate,... + Grads,Year0Iter,SecPerYear,Months,FieldName); % Load parameters. @@ -25,6 +25,11 @@ else , mnchandle = 'tave.*' ; end elseif isequal(dat,'Ins') filesuffix = ''; mnchandle = 'state.*'; +elseif isequal(dat(1:2),'Mn') + filesuffix = ''; + if ismember(dat,{'MnA','MnO'}), mnchandle = 'monitor.*'; end + if isequal(dat,'MnI'), mnchandle = 'monitor_sice.*'; end + if isequal(dat,'MnL'), mnchandle = 'monitor_land.*'; end elseif ismember(dat,{'','Non','None'}) filesuffix = ''; mnchandle = ''; else @@ -46,7 +51,8 @@ % Open grads files if ~isequal(Grads,0) if DiagDebug, disp([' Debug -- Loading GRADS field.']); end - [data,xax,yax,zax,months,time,dim] = DiagLoadGradsData(Grads,dad,fln,DiagDebug); + [data,xax,yax,zax,months,time,dim] = ... + DiagLoadGradsData(Grads,dad,fln,DiagDebug); if DiagDebug, disp([' Debug -- ''data'' size after load: ',mat2str(size(data))]); end data = DiagAverage(data,fln,avg,months,ddf,dim); if DiagDebug, disp([' Debug -- ''data'' size after averaging: ',mat2str(size(data))]); end @@ -56,7 +62,12 @@ else error('Loading grads data only set to handle surface plots!'); end - + +% Load monitor data. +elseif isequal(dat(1:2),'Mn') + [data,time] = ... + DiagLoadMonitor(fln,mnchandle,dad,itr,tst,SecPerYear,DiagDebug); + xax = time; yax = NaN; pltslc = 'timfld'; % Load AIM physics data. The data is all stuffed into one file called % 'aimPhy' which currently only comes in MDS format. When loading a field, @@ -69,7 +80,7 @@ % stresses (on atmosphere) must converted to a lat-lon grid. elseif ismember(fln,aimparameters) if isequal(ddf,'MNC') - error('Aim physics reader not set to hande ''MNC'' data format.') + error('Aim physics reader not set to handle ''MNC'' data format.') end data = DiagLoad_Local('aimPhy',dat,dad,grd,itr,ddf,filesuffix,mnchandle); if isequal(fln,'TotP') @@ -94,7 +105,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,flu,... - ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Load coupled fields. Like AIM, the data is all stuffed into one file % called either 'cplFld.*' (MDS) or 'cpl_tave.*' (MNC). Variables within @@ -127,7 +138,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,flu,... - ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Load ice parameters. Note that as the ice variable names appear to be in % constant flux, this part quickly becomes out of date. @@ -142,7 +153,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,flu,... - ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Load generic fields where no modifications are needed. elseif ismember(fln,{'S','T','Temp','aim_RH','phiHyd','Conv'}) @@ -150,7 +161,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Read vertical velocities, eta. Convert from [P] to [hPa] in atmosphere. elseif ismember(fln,{'W','wVel','Eta','ETA'}) @@ -159,16 +170,16 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Load horizontal velocities, convert to lat-lon grid. -elseif ismember(fln,{'U','V','uVel','vVel'}) +elseif ismember(fln,{'U','V','uVel','vVel','fizhi_U','fizhi_V'}) if isequal(dat,'Tav') U = DiagLoad_Local('uVel',dat,dad,grd,itr,ddf,filesuffix,mnchandle); V = DiagLoad_Local('vVel',dat,dad,grd,itr,ddf,filesuffix,mnchandle); else isequal(dat,'Ins') - U = DiagLoad_Local('U',dat,dad,grd,itr,ddf,filesuffix,mnchandle); - V = DiagLoad_Local('V',dat,dad,grd,itr,ddf,filesuffix,mnchandle); + U = DiagLoad_Local([fln(1:end-1),'U'],dat,dad,grd,itr,ddf,filesuffix,mnchandle); + V = DiagLoad_Local([fln(1:end-1),'V'],dat,dad,grd,itr,ddf,filesuffix,mnchandle); end if isequal(ddf,'MNC') U = U(1:end-1,:,:,:); @@ -179,14 +190,14 @@ [XC,XG,YC,YG,Ylat,ZC,ZG,RAC,drC,drF,HFacC,HFacW,HFacS,dxG,dyG,dxC,dyC] = ... DiagLoadGridData(LoadGridData,grd,gdf,flu,GridSuffix,ZcordFile); [U,V]=uvcube2latlon(XC,YC,U,V,XL,YL); - if ismember(fln,{'U','uVel'}) + if ismember(fln,{'U','uVel','fizhi_U'}) data = U; - elseif ismember(fln,{'V','vVel'}) + elseif ismember(fln,{'V','vVel','fizhi_V'}) data = V; end [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Meridional overturning. elseif ismember(fln,{'Bol','Psi','Res'}) @@ -269,7 +280,7 @@ end [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,flu,ddf,gdf,... - avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % Read in variabilities. Here we must make some simple calculations. For @@ -285,7 +296,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); elseif isequal(fln,'Tstd') T = DiagLoad_Local('T' ,dat,dad,grd,itr,ddf,filesuffix,mnchandle); TT = DiagLoad_Local('TT',dat,dad,grd,itr,ddf,filesuffix,mnchandle); @@ -293,7 +304,7 @@ data = DiagAverage(data,fln,avg,months,ddf,Dim); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); elseif isequal(fln,'KEpri') U = DiagLoad_Local('uVel',dat,dad,grd,itr,ddf,filesuffix,mnchandle); V = DiagLoad_Local('vVel',dat,dad,grd,itr,ddf,filesuffix,mnchandle); @@ -312,30 +323,85 @@ data = sqrt(abs((U.*U + V.*V) - (UU + VV))); [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); % If the field is not recognized, try best to try and open it. else disp([' Unknown field, attempting to load.']); data = DiagLoad_Local(fln,dat,dad,grd,itr,ddf,filesuffix,mnchandle); if DiagDebug, disp([' Debug -- ''data'' size after load: ',mat2str(size(data))]); end - if ~isequal(Index,0) + if isequal(Vector,1) | isequal(Vector,2) | isequal(Vector,3) + if isequal(Vector,1) | isequal(Vector,3) + if ~isequal(Index,0) + if isequal(Dim,2), U = squeeze(data(:,:,Index,:)); + elseif isequal(Dim,3), U = squeeze(data(:,:,:,Index,:)); end + end + if ~isequal(Mate,0) + if isequal(Dim,2), V = squeeze(data(:,:,Mate,:)); + elseif isequal(Dim,3), V = squeeze(data(:,:,:,Mate,:)); end + end + elseif isequal(Vector,2) + if ~isequal(Index,0) + if isequal(Dim,2), V = squeeze(data(:,:,Index,:)); + elseif isequal(Dim,3), V = squeeze(data(:,:,:,Index,:)); end + end + if ~isequal(Mate,0) + if isequal(Dim,2), U = squeeze(data(:,:,Mate,:)); + elseif isequal(Dim,3), U = squeeze(data(:,:,:,Mate,:)); end + end + end + U = DiagAverage(U,fln,avg,months,ddf,Dim); + V = DiagAverage(V,fln,avg,months,ddf,Dim); + [XC,XG,YC,YG,Ylat,ZC,ZG,RAC,drC,drF,HFacC,HFacW,HFacS,dxG,dyG,dxC,dyC] = ... + DiagLoadGridData(LoadGridData,grd,gdf,flu,GridSuffix,ZcordFile); + if isequal(Vector,1) | isequal(Vector,2) + [U,V]=uvcube2latlon(XC,YC,U,V,XL,YL); + end + if isequal(Vector,1) + data = U; + elseif isequal(Vector,2) + data = V; + elseif isequal(Vector,3) + if isequal(flu,'A'), delM = - 1000*drF ./ g; + elseif isequal(flu,'O'), delM = drF; end + [psiavg,mskG,ylat]=calc_PsiCube(delM,U.*HFacW,V.*HFacS,dxG,dyG); + data = psiavg(2:end-1,:)'; + if isequal(flu,'A'), data = data./1e6; end + xax = ylat; yax = ZG; pltslc = 'lathgt'; + end + if DiagDebug, disp([' Debug -- ''data'' size after vector manipulation: ',mat2str(size(data))]); end + else +% Sequence here for Vector = 0 + if ~isequal(Index,0) if isequal(Dim,2), data = squeeze(data(:,:,Index,:)); elseif isequal(Dim,3), data = squeeze(data(:,:,:,Index,:)); end + end + if DiagDebug, disp([' Debug -- ''data'' size after indexing: ',mat2str(size(data))]); end + data = DiagAverage(data,fln,avg,months,ddf,Dim); +% Actual Temperature -- When the FieldName says ActT we assume that +% the field that is pointed to in Index is potential temperature + if isequal(FieldName,'ActT') + [XC,XG,YC,YG,Ylat,ZC,ZG,RAC,drC,drF,HFacC,HFacW,HFacS,dxG,dyG,dxC,dyC] = ... + DiagLoadGridData(LoadGridData,grd,gdf,flu,GridSuffix,ZcordFile); + pres = NaN.*zeros(size(data)); + for iz=1:length(ZC), pres(:,:,iz)=ZC(iz); end + temp=data.*(pres./presrefA).^(RdA/cpA); + data=temp; + end + + if DiagDebug, disp([' Debug -- ''data'' size after average: ',mat2str(size(data))]); end end - if DiagDebug, disp([' Debug -- ''data'' size after indexing: ',mat2str(size(data))]); end - data = DiagAverage(data,fln,avg,months,ddf,Dim); - if DiagDebug, disp([' Debug -- ''data'' size after average: ',mat2str(size(data))]); end + if ~isequal(Vector,3) [data,xax,yax,pltslc] = ... DiagSlice(data,fln,trl,dat,dad,grd,itr,tst,... - flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile); + flu,ddf,gdf,avg,slc,pst,LoadGridData,GridSuffix,ZcordFile,Vector,FieldName); if DiagDebug, disp([' Debug -- ''data'' size after slice: ',mat2str(size(data))]); end if DiagDebug, disp([' Debug -- ''xax'' size: ',mat2str(size(xax))]); end if DiagDebug, disp([' Debug -- ''yax'' size: ',mat2str(size(yax))]); end if DiagDebug, disp([' Debug -- ''pltslc'': ',pltslc]); end + end end - %-------------------------------------------------------------------------% % Local functions % %-------------------------------------------------------------------------%