Warning, /verification/tutorial_global_oce_latlon/diags_matlab/mit_loadglobal.m is written in an unsupported language. File is not indexed.
view on githubraw file Latest commit add29e06 on 2018-01-31 20:35:05 UTC
051ee7f715 Jean*0001 % m-file: mit_loadglobal.m
0002 % plot various diagnostics of the global 4x4 degree runs
0003
0004 % uses
0005 % mit_loadgrid
0006 % mit_getparm
0007 % mit_oceanmasks
0008 % mit_timesteps
0009 % mit_loadglobaldata
0010 % mit_globalmovie
0011 % mit_plotstreamfunctions
0012 % mit_plotzonalvelocity
0013 % mit_globalmean
0014 % mit_zonalmean
0015 % mit_meridflux
0016 % mit_plotmeandrift
0017
0018 % read in all grid files, etc. This has to be done at the very beginning!
0019 grd = mit_loadgrid('.');
0020 dname = grd.dname;
0021 % add some masks
0022 grd = mit_oceanmasks(grd);
0023 %
0024 if ~exist('meanfields','var')
0025 meanfields = 1; % read in mean fields instead of instantaneous ones,
0026 % include some more diagnostics
0027 % meanfields = 1 is the default
0028 end
0029 %
0030 mit_timesteps
0031
0032 if meanfields
0033 disp('reading averaged fields (*tave.*)')
0034 %
0035 mit_loadglobaldata
0036 %
0037 u=rdmds('uVeltave',timesteps(2:nt)');
0038 v=rdmds('vVeltave',timesteps(2:nt)');
0039 w=rdmds('wVeltave',timesteps(2:nt)');
0040 t=rdmds('Ttave',timesteps(2:nt)');
0041 s=rdmds('Stave',timesteps(2:nt)');
0042 eta=rdmds('ETAtave',timesteps(2:nt)').*repmat(squeeze(grd.cmask(:,:,1)),[1 1 nt-1]);
0043 k=1;
0044 u=cat(4,zeros([nx ny nz 1]),u);
0045 v=cat(4,zeros([nx ny nz 1]),v);
0046 w=cat(4,zeros([nx ny nz 1]),w);
0047 % hack
0048 t=cat(4,rdmds('T',timesteps(1)'),t);
0049 s=cat(4,rdmds('S',timesteps(1)'),s);
0050 eta=cat(3,0.*grd.cmask(:,:,1),eta);
0051 else
0052 % load snap shots
0053 disp('reading snap shot fields')
0054 u=rdmds('U',timesteps');
0055 v=rdmds('V',timesteps');
0056 w=rdmds('W',timesteps');
0057 t=rdmds('T',timesteps');
0058 s=rdmds('S',timesteps');
0059 eta=rdmds('Eta',timesteps').*repmat(squeeze(grd.cmask(:,:,1)),[1 1 nt]);
0060 end
0061 % orientation and masking
0062 if ~strcmp(grd.buoyancy,'OCEANIC');
0063 u = u(:,:,end:-1:1,:);
0064 v = v(:,:,end:-1:1,:);
0065 w = w(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
0066 t = t(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
0067 s = s(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
0068 else
0069 % $$$ u = u;
0070 % $$$ v = v;
0071 w = w.*repmat(grd.cmask,[1 1 1 nt]);
0072 t = t.*repmat(grd.cmask,[1 1 1 nt]);
0073 s = s.*repmat(grd.cmask,[1 1 1 nt]);
0074 end
0075
0076 % transport through Drake Passage
0077 TDP = NaN*ones(nt,1);
0078 kx = 73;
0079 kyg = [5 6 7];
0080 %kyg = [4 5 6 7];
0081 da = grd.dz*grd.dyg(kx,kyg);
0082 for k=kt;
0083 TDP(k) = sum(nansum(squeeze(u(kx,kyg,:,k))'.*da))*1e-6; % in Sv
0084 TDPz(:,k) = nansum(squeeze(u(kx,kyg,:,k)).*da')'*1e-6; % in Sv
0085 end
0086
0087 % set some depth
0088 if ~meanfields
0089 mit_globalmovie
0090 end
0091
0092 % $$$ [x y z]=meshgrid(grd.lonc,grd.grd.latc,-grd.zc);
0093 % $$$ tk = permute(squeeze(t(:,:,:,k)),[2 1 3]);
0094 % $$$ xplane = [0:30:360];
0095 % $$$ yplane = 30;
0096 % $$$ zplane = -[50 500:500:6000];
0097 % $$$ slice(x,y,z,tk,xplane,yplane,zplane);
0098 % $$$ shading flat, colorbar('h')
0099
0100 mit_plotstreamfunctions
0101
0102 mit_plotzonalvelocity
0103
0104 if meanfields
913aa74540 Patr*0105 % plot some mean fields
0106 figure('PaperPosition',[0.31 0.25 10.5 7.88],'PaperOrientation','landscape')
0107 k=kmax;
0108 mySST = t(:,:,1,k);
0109 mySSS = s(:,:,1,k);
0110 myETA = eta(:,:,k);
0111 myVEL = sqrt(u(:,:,1,k).*u(:,:,1,k)+v(:,:,1,k).*v(:,:,1,k));
0112 %
0113 clear sh;
0114 sh(1) = subplot(2,2,1);imagesc(grd.lonc,grd.latc,mySST'.*grd.hfacc(:,:,1)');
0115 title(['Sea Surface Temperature [degC]'])
0116 sh(2) = subplot(2,2,2);imagesc(grd.lonc,grd.latc,mySSS'.*grd.hfacc(:,:,1)');
0117 title(['Sea Surface Salinity [psu]'])
0118 sh(3) = subplot(2,2,3);imagesc(grd.lonc,grd.latc,myETA');
0119 title(['Sea Surface Height [m]'])
0120 sh(4) = subplot(2,2,4);imagesc(grd.lonc,grd.latc,myVEL'.*grd.hfacc(:,:,1)');
0121 title(['Sea Surface Speed [m/s]'])
0122 %
0123 axis(sh,'image'); axis(sh,'xy')
0124 set(gcf,'currentAxes',sh(1));colorbar('h')
0125 set(gcf,'currentAxes',sh(2));colorbar('h')
0126 set(gcf,'currentAxes',sh(3));colorbar('h')
0127 set(gcf,'currentAxes',sh(4));colorbar('h')
0128 end
0129
0130 if meanfields
051ee7f715 Jean*0131 % compute actual heat and E-P fluxes from time averages
0132 % this only makes sense when you load the time averaged fields
0133 tauTheta = mit_getparm('data','tauThetaClimRelax');
0134 tauSalt = mit_getparm('data','tauSaltClimRelax');
0135 rhonil = mit_getparm('data','rhonil');
0136 if isempty(rhonil)
0137 rhonil = 1035;
0138 end
0139 Cp = 3994;
0140 if isempty(tauTheta)
0141 tauTheta = 0;
0142 end
0143 if tauTheta==0
0144 recip_tauTheta = 0.;
0145 else
0146 recip_tauTheta = 1./tauTheta;
0147 end
0148 if isempty(tauSalt)
0149 tauSalt = 0;
0150 end
0151 if tauSalt==0
0152 recip_tauSalt = 0.;
0153 else
0154 recip_tauSalt = 1./tauSalt;
0155 end
0156 figure('PaperPosition',[0.31 0.25 10.5 7.88],'PaperOrientation','landscape')
0157 k=kmax;
0158 qnet_diag = (t(:,:,1,k)-tdmean)*grd.dz(1).*grd.hfacc(:,:,1)*Cp*rhonil*recip_tauTheta;
0159 empr_diag = - (s(:,:,1,k)-sdmean)*grd.dz(1).*grd.hfacc(:,:,1)*recip_tauSalt/35;
0160 clear sh;
0161 sh(1) = subplot(2,2,1);imagesc(grd.lonc,grd.latc,qnet_diag');
0162 title(['Q_{diag} from T-T_{lev} at the surface'])
0163 sh(2) =subplot(2,2,2);imagesc(grd.lonc,grd.latc,qnet'.*grd.cmask(:,:,1)');
0164 title('annual mean Q_{net} (data)')
0165 sh(3) =subplot(2,2,3);imagesc(grd.lonc,grd.latc,empr_diag');
0166 title(['(E-P)_{diag} from S-S_{lev} at the surface'])
0167 sh(4) =subplot(2,2,4);imagesc(grd.lonc,grd.latc,empr'.*grd.cmask(:,:,1)');
0168 title('annual mean E-P (data)')
0169 axis(sh,'image'); axis(sh,'xy')
0170 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
0171 ', ' tuname ' = ' num2str(tim(k))])
0172 set(sh(1:2),'clim',[-1 1]*200);
0173 set(gcf,'currentAxes',sh(1));colorbar('h')
0174 set(gcf,'currentAxes',sh(2));colorbar('h')
0175 set(sh(3:4),'clim',[-1 1]*1.e-7);
0176 set(gcf,'currentAxes',sh(3));colorbar('h')
0177 set(gcf,'currentAxes',sh(4));colorbar('h')
0178
0179 clear sh
0180
0181 if length(tim) > 1
0182 rac3d = repmat(grd.rac,[1 1 nz]).*grd.hfacc;
0183 % horizontal averages of t and s in a hovmoeller-type diagramme
0184 mdt = repmat(NaN,grd.nz,length(kt));
0185 mdtt= mdt;
0186 mdtstd= mdt;
0187 mdt_alt = mdt;
0188 mdtt_alt= mdt;
0189 mdtstd_alt = mdt;
0190 mdt_pac = mdt;
0191 mdtt_pac= mdt;
0192 mdtstd_pac = mdt;
0193 mdt_sou = mdt;
0194 mdtt_sou= mdt;
0195 mdtstd_sou = mdt;
0196 mds = mdt;
0197 mdss= mdt;
0198 mdsstd = mdt;
0199 mds_alt = mdt;
0200 mdss_alt= mdt;
0201 mdsstd_alt = mdt;
0202 mds_pac = mdt;
0203 mdss_pac= mdt;
0204 mdsstd_pac = mdt;
0205 mds_sou = mdt;
0206 mdss_sou= mdt;
0207 mdsstd_sou = mdt;
0208 indopac_hfacc = ...
0209 change(change(grd.pacific_hfacc,'==',NaN,0) ...
0210 +change(grd.indic_hfacc,'==',NaN,0),'==',0,NaN);
0211 indopac_hfacc(:,1:12,:) = NaN;
0212 southern_hfacc = grd.hfacc;
0213 southern_hfacc(:,13:end,:) = NaN;
0214 atl_hfacc = grd.atlantic_hfacc;
0215 atl_hfacc(:,1:12,:) = NaN;
0216 for k=kt;
0217 dt = t(:,:,:,k)-tdatamean; % anomalies
0218 ds = s(:,:,:,k)-sdatamean; % anomalies
0219 [mdt(:,k) mdtt(:,k) mdtstd(:,k)] = mit_globalmean(dt,rac3d);
0220 [mds(:,k) mdss(:,k) mdsstd(:,k)] = mit_globalmean(ds,rac3d);
0221 [mdt_atl(:,k) mdtt_atl(:,k) mdtstd_atl(:,k)] = mit_globalmean(dt,rac3d.*atl_hfacc);
0222 [mds_atl(:,k) mdss_atl(:,k) mdsstd_atl(:,k)] = mit_globalmean(ds,rac3d.*atl_hfacc);
0223 [mdt_pac(:,k) mdtt_pac(:,k) mdtstd_pac(:,k)] = mit_globalmean(dt,rac3d.*indopac_hfacc);
0224 [mds_pac(:,k) mdss_pac(:,k) mdsstd_pac(:,k)] = mit_globalmean(ds,rac3d.*indopac_hfacc);
0225 [mdt_sou(:,k) mdtt_sou(:,k) mdtstd_sou(:,k)] = mit_globalmean(dt,rac3d.*southern_hfacc);
0226 [mds_sou(:,k) mdss_sou(:,k) mdsstd_sou(:,k)] = mit_globalmean(ds,rac3d.*southern_hfacc);
0227 end
0228
0229 figure('PaperPosition',[0.31 0.25 10.5 7.88],'PaperOrientation','landscape')
0230 if length(tim) > 2
0231 clear sh
0232 sh(1) = subplot(2,2,1);
0233 contourf(tim(2:end),-grd.zc,mdt(:,2:end),20);
0234 caxis([-1 1]*max(abs(caxis))); shading flat; colorbar('v')
0235 title('T-T_{lev} horizontally averaged')
0236 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
0237 sh(2) = subplot(2,2,2);
0238 contourf(tim(2:end),-grd.zc,mds(:,2:end),20);
0239 caxis([-1 1]*max(abs(caxis))); shading flat; colorbar('v')
0240 title('S-S_{lev} horizontally averaged')
0241 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
0242 %
0243 sh(3) = subplot(2,2,3);
0244 contourf(tim(2:end),-grd.zc,mdtt(:,2:end),20);
0245 shading flat; colorbar('v')
0246 title('|T-T_{lev}| horizontally averaged');
0247 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
0248 sh(4) = subplot(2,2,4);
0249 contourf(tim(2:end),-grd.zc,mdss(:,2:end),20);
0250 shading flat; colorbar('v')
0251 title('|S-S_{lev}| horizontally averaged')
0252 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
0253 set(sh,'layer','top');
0254 suptitle(['experiment ' dname])
0255 end %if length(tim) > 2
0256
0257 k=kmax;
0258 figure('PaperOrientation','landscape','PaperPosition',[0.25 0.3125 10.5 7.875])
0259 clear sh
0260 sh(1) = subplot(2,4,1);
0261 plot(mdt(:,k),-grd.zc,'b-',mdt(:,k)+mdtstd(:,k),-grd.zc,'r--',mdt(:,k)-mdtstd(:,k),-grd.zc,'r--');
0262 xlabel('T-T_{lev} [degC]'); ylabel('depth [m]'); title('global')
0263 legend('horiz. mean','std',3);
0264 sh(2) = subplot(2,4,2);
0265 plot(mdt_atl(:,k),-grd.zc,'b-',mdt_atl(:,k)+mdtstd_atl(:,k),-grd.zc,'r--',mdt_atl(:,k)-mdtstd_atl(:,k),-grd.zc,'r--');
0266 xlabel('T-T_{lev} [degC]'); %ylabel('depth [m]');
0267 title('atlantic')
0268 sh(3) = subplot(2,4,3);
0269 plot(mdt_pac(:,k),-grd.zc,'b-',mdt_pac(:,k)+mdtstd_pac(:,k),-grd.zc,'r--',mdt_pac(:,k)-mdtstd_pac(:,k),-grd.zc,'r--');
0270 xlabel('T-T_{lev} [degC]'); %ylabel('depth [m]');
0271 title('indo-pacific')
0272 sh(4) = subplot(2,4,4);
0273 plot(mdt_sou(:,k),-grd.zc,'b-',mdt_sou(:,k)+mdtstd_sou(:,k),-grd.zc,'r--',mdt_sou(:,k)-mdtstd_sou(:,k),-grd.zc,'r--');
0274 xlabel('T-T_{lev} [degC]'); ylabel('depth [m]'); title('southern')
0275 %
0276 sh(5) = subplot(2,4,5);
0277 plot(mds(:,k),-grd.zc,'b-',mds(:,k)+mdsstd(:,k),-grd.zc,'r--',mds(:,k)-mdsstd(:,k),-grd.zc,'r--');
0278 xlabel('S-S_{lev} [PSU]'); ylabel('depth [m]');
0279 sh(6) = subplot(2,4,6);
0280 plot(mds_atl(:,k),-grd.zc,'b-',mds_atl(:,k)+mdsstd_atl(:,k),-grd.zc,'r--',mds_atl(:,k)-mdsstd_atl(:,k),-grd.zc,'r--');
0281 xlabel('S-S_{lev} [PSU]'); %ylabel('depth [m]');
0282 sh(7) = subplot(2,4,7);
0283 plot(mds_pac(:,k),-grd.zc,'b-',mds_pac(:,k)+mdsstd_pac(:,k),-grd.zc,'r--',mds_pac(:,k)-mdsstd_pac(:,k),-grd.zc,'r--');
0284 xlabel('S-S_{lev} [PSU]'); %ylabel('depth [m]');
0285 sh(8) = subplot(2,4,8);
0286 plot(mds_sou(:,k),-grd.zc,'b-',mds_sou(:,k)+mdsstd_sou(:,k),-grd.zc,'r--',mds_sou(:,k)-mdsstd_sou(:,k),-grd.zc,'r--');
0287 xlabel('S-S_{lev} [PSU]'); ylabel('depth [m]');
0288 set(sh,'xgrid','on','ygrid','on')
0289 set([sh(4); sh(8)],'YAxisLocation','right')
0290 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
0291 ', ' tuname ' = ' num2str(tim(k))])
0292 end %if length(tim) > 1
0293
0294 % zonal mean of temperature and salinity for different ocean basins
0295 k=kmax;
0296 tdzm_glo = mit_zonalmean(tdatamean,grd.hfacc,grd.dxc);
0297 tdzm_atl = mit_zonalmean(tdatamean,grd.atlantic_hfacc,grd.dxc);
0298 tdzm_pac = mit_zonalmean(tdatamean,grd.pacific_hfacc,grd.dxc);
0299 tdzm_ind = mit_zonalmean(tdatamean,grd.indic_hfacc,grd.dxc);
0300 sdzm_glo = mit_zonalmean(sdatamean,grd.hfacc,grd.dxc);
0301 sdzm_atl = mit_zonalmean(sdatamean,grd.atlantic_hfacc,grd.dxc);
0302 sdzm_pac = mit_zonalmean(sdatamean,grd.pacific_hfacc,grd.dxc);
0303 sdzm_ind = mit_zonalmean(sdatamean,grd.indic_hfacc,grd.dxc);
0304
0305 figure('PaperPosition',[0.25 0.368552 8 10.2629]);
0306 tlev = -5:.5:5;
0307 slev = -1:.1:1;
0308 clear sh clh
0309 clh = cell(8,1);
0310 delay = 1;
0311 for k = kmax
0312 tzm_glo = mit_zonalmean(t(:,:,:,k),grd.hfacc,grd.dxc);
0313 tzm_atl = mit_zonalmean(t(:,:,:,k),grd.atlantic_hfacc,grd.dxc);
0314 tzm_pac = mit_zonalmean(t(:,:,:,k),grd.pacific_hfacc,grd.dxc);
0315 tzm_ind = mit_zonalmean(t(:,:,:,k),grd.indic_hfacc,grd.dxc);
0316 szm_glo = mit_zonalmean(s(:,:,:,k),grd.hfacc,grd.dxc);
0317 szm_atl = mit_zonalmean(s(:,:,:,k),grd.atlantic_hfacc,grd.dxc);
0318 szm_pac = mit_zonalmean(s(:,:,:,k),grd.pacific_hfacc,grd.dxc);
0319 szm_ind = mit_zonalmean(s(:,:,:,k),grd.indic_hfacc,grd.dxc);
0320 caxt = [min(tzm_glo(:)-tdzm_glo(:)) max(tzm_glo(:)-tdzm_glo(:))];
0321 tlev = -5:.5:5;
0322 if max(abs(caxt)) < 1;
0323 tlev = -1:.1:1;
0324 end
0325 if max(abs(caxt)) < .1;
0326 tlev = .1*tlev;
0327 end
0328 caxs = [min(szm_glo(:)-sdzm_glo(:)) max(szm_glo(:)-sdzm_glo(:))];
0329 slev = -1.:.05:1;
0330 if max(abs(caxs)) < 0.2;
0331 slev = -.20:.01:.20;
0332 end
0333 if max(abs(caxs)) < .02;
0334 slev = .1*slev;
0335 end
0336 sh(1) = subplot(4,2,1);
0337 [cs h] = contourf(grd.latc,-grd.zc,(tzm_glo-tdzm_glo)',tlev);
0338 if ~isempty(h);
0339 clh{1} = clabel(cs,h,[0 0]);
0340 end
0341 title('\theta-\theta_{lev} [degC]: global ocean')
0342 sh(3) = subplot(4,2,3);
0343 [cs h] = contourf(grd.latc,-grd.zc,(tzm_atl-tdzm_atl)',tlev);
0344 if ~isempty(h);
0345 clh{3} = clabel(cs,h,[0 0]);
0346 end
0347 title('atlantic ocean')
0348 sh(5) = subplot(4,2,5);
0349 [cs h] = contourf(grd.latc,-grd.zc,(tzm_pac-tdzm_pac)',tlev);
0350 if ~isempty(h);
0351 clh{5} = clabel(cs,h,[0 0]);
0352 end
0353 title('pacific ocean')
0354 sh(7) = subplot(4,2,7);
0355 [cs h] = contourf(grd.latc,-grd.zc,(tzm_ind-tdzm_ind)',tlev);
0356 if ~isempty(h);
0357 clh{7} = clabel(cs,h,[0 0]);
0358 end
0359 title('indian ocean')
0360 set(sh(1:2:end),'clim',[tlev(1) tlev(end)])
0361 colorbar('h')
0362 sh(2) = subplot(4,2,2);
0363 [cs h] = contourf(grd.latc,-grd.zc,(szm_glo-sdzm_glo)',slev);
0364 if ~isempty(h);
0365 clh{2} = clabel(cs,h,[0 0]);
0366 end
0367 title('S-S_{lev} [PSU]: global ocean')
0368 sh(4) = subplot(4,2,4);
0369 [cs h] = contourf(grd.latc,-grd.zc,(szm_atl-sdzm_atl)',slev);
0370 if ~isempty(h);
0371 clh{4} = clabel(cs,h,[0 0]);
0372 end
0373 title('atlantic ocean')
0374 sh(6) = subplot(4,2,6);
0375 [cs h] = contourf(grd.latc,-grd.zc,(szm_pac-sdzm_pac)',slev);
0376 if ~isempty(h);
0377 clh{6} = clabel(cs,h,[0 0]);
0378 end
0379 title('pacific ocean')
0380 sh(8) = subplot(4,2,8);
0381 [cs h] = contourf(grd.latc,-grd.zc,(szm_ind-sdzm_ind)',slev);
0382 if ~isempty(h);
0383 clh{8} = clabel(cs,h,[0 0]);
0384 end
0385 title('indian ocean')
0386 set(sh(2:2:end),'clim',[slev(1) slev(end)])
0387 colorbar('h')
0388 set(sh,'layer','top')
0389 if ~isempty(cat(1,clh{:}))
0390 set(cat(1,clh{:}),'fontsize',8);
0391 end
0392
0393 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
0394 ', ' tuname ' = ' num2str(tim(k)) ', zonal averages'])
0395 drawnow;
0396 pause(delay)
0397 end
0398 clear clh sh
0399
0400 % meridional transport of heat (internal energy) and fresh water
0401 % diagnosed from surface fluxes
0402 petawatts = 1e-15;
0403 sverdrups = 1e-6;
0404 heat_data = -mit_meridflux(qnet,grd.dxc,grd.dyc)*petawatts;
0405 heat_diag = -mit_meridflux(qnet_diag,grd.dxc,grd.dyc)*petawatts;
0406 freshwater_data = -mit_meridflux(empr,grd.dxc,grd.dyc)*sverdrups;
0407 freshwater_diag = -mit_meridflux(empr_diag,grd.dxc,grd.dyc)*sverdrups;
0408 figure
0409 latgf = [2*grd.latc-grd.latg];
0410 clear sh
0411 sh(1) = subplot(2,1,1);
0412 hh = plot(latgf,heat_data,'-',latgf,heat_diag,'--', ...
0413 latgf,heat_data+heat_diag,'-.');
0414 title('heat (internal energy) flux [PW]')
0415 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
0416 ', ' tuname ' = ' num2str(tim(k))])
0417 legend('surface flux','restoring','net',2)
0418 sh(2) = subplot(2,1,2);
0419 fwh = plot(latgf,freshwater_data,'-',latgf,freshwater_diag,'--', ...
0420 latgf,freshwater_data+freshwater_diag,'-.');
0421 title('freshwater flux [Sv]')
0422 legend('surface flux','restoring','net',3)
0423 set(sh,'xlim',[latgf(1) latgf(end)],'xgrid','on','ygrid','on');
0424
0425 end % meanfields
0426
0427 mit_plotmeandrift
0428
0429 fname = [grd.dname '_' sprintf('%u',tim(kmax)) timeunit];
0430 in = findstr(fname,'\');
0431 fname(in) = [];
0432 fn = gcf;
0433 if meanfields
0434 fname = [fname '.ps'];
2b5e163001 Patr*0435 % change required by newer Matlab versions
0436 % f0 = fn-7;
0437 f0=fn.Number-7;
051ee7f715 Jean*0438 else
0439 fname = [fname '_snap.ps'];
2b5e163001 Patr*0440 % change required by newer Matlab versions
0441 % f0 = fn-4;
0442 f0=fn.Number-4;
051ee7f715 Jean*0443 end
0444 print(f0,'-dpsc',fname);
0445 for k = f0+1:fn
0446 print(k,'-dpsc',fname,'-append')
0447 end
