/[MITgcm]/MITgcm/verification/tutorial_global_oce_latlon/diags_matlab/mit_loadglobal.m
ViewVC logotype

Contents of /MITgcm/verification/tutorial_global_oce_latlon/diags_matlab/mit_loadglobal.m

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


Revision 1.3 - (show annotations) (download)
Sat Aug 12 20:25:13 2006 UTC (17 years, 8 months ago) by jmc
Branch: MAIN
CVS Tags: checkpoint64y, checkpoint64x, checkpoint64z, checkpoint64q, checkpoint64p, checkpoint64s, checkpoint64r, checkpoint64u, checkpoint64t, checkpoint64w, checkpoint64v, checkpoint64i, checkpoint64h, checkpoint64k, checkpoint64j, checkpoint64m, checkpoint64l, checkpoint64o, checkpoint64n, checkpoint64a, checkpoint64c, checkpoint64b, checkpoint64e, checkpoint64d, checkpoint64g, checkpoint64f, checkpoint58u_post, checkpoint58w_post, checkpoint63p, checkpoint63q, checkpoint63r, checkpoint63s, checkpoint63l, checkpoint63m, checkpoint63n, checkpoint63o, checkpoint63h, checkpoint63i, checkpoint63j, checkpoint63k, checkpoint63d, checkpoint63e, checkpoint63f, checkpoint63g, checkpoint63a, checkpoint63b, checkpoint63c, checkpoint64, checkpoint65, checkpoint60, checkpoint61, checkpoint62, checkpoint63, checkpoint58r_post, checkpoint58n_post, checkpoint58x_post, checkpoint58t_post, checkpoint65x, checkpoint65y, checkpoint65r, checkpoint65s, checkpoint65p, checkpoint65q, checkpoint65v, checkpoint65w, checkpoint65t, checkpoint65u, checkpoint65j, checkpoint65k, checkpoint65h, checkpoint65i, checkpoint65n, checkpoint65o, checkpoint65l, checkpoint65m, checkpoint65b, checkpoint65c, checkpoint65a, checkpoint65f, checkpoint65g, checkpoint65d, checkpoint65e, checkpoint58q_post, checkpoint59q, checkpoint59p, checkpoint59r, checkpoint59e, checkpoint59d, checkpoint59g, checkpoint59f, checkpoint59a, checkpoint59c, checkpoint59b, checkpoint59m, checkpoint59l, checkpoint59o, checkpoint59n, checkpoint59i, checkpoint59h, checkpoint59k, checkpoint59j, checkpoint59, checkpoint58o_post, checkpoint62c, checkpoint62b, checkpoint62a, checkpoint62g, checkpoint62f, checkpoint62e, checkpoint62d, checkpoint62k, checkpoint62j, checkpoint62i, checkpoint62h, checkpoint62o, checkpoint62n, checkpoint62m, checkpoint62l, checkpoint62s, checkpoint62r, checkpoint62q, checkpoint62p, checkpoint62w, checkpoint62v, checkpoint62u, checkpoint62t, checkpoint62z, checkpoint62y, checkpoint62x, checkpoint58y_post, checkpoint58v_post, checkpoint58s_post, checkpoint61f, checkpoint61g, checkpoint61d, checkpoint61e, checkpoint61b, checkpoint61c, checkpoint58p_post, checkpoint61a, checkpoint61n, checkpoint61o, checkpoint61l, checkpoint61m, checkpoint61j, checkpoint61k, checkpoint61h, checkpoint61i, checkpoint61v, checkpoint61w, checkpoint61t, checkpoint61u, checkpoint61r, checkpoint61s, checkpoint61p, checkpoint61q, checkpoint61z, checkpoint61x, checkpoint61y
Changes since 1.2: +1 -1 lines
accidentally removed ; put them back.

1 % m-file: mit_loadglobal.m
2 % plot various diagnostics of the global 4x4 degree runs
3
4 % uses
5 % mit_loadgrid
6 % mit_getparm
7 % mit_oceanmasks
8 % mit_timesteps
9 % mit_loadglobaldata
10 % mit_globalmovie
11 % mit_plotstreamfunctions
12 % mit_plotzonalvelocity
13 % mit_globalmean
14 % mit_zonalmean
15 % mit_meridflux
16 % mit_plotmeandrift
17
18 % $Header: $
19 % $Name: $
20
21 % read in all grid files, etc. This has to be done at the very beginning!
22 grd = mit_loadgrid('.');
23 dname = grd.dname;
24 % add some masks
25 grd = mit_oceanmasks(grd);
26 %
27 if ~exist('meanfields','var')
28 meanfields = 1; % read in mean fields instead of instantaneous ones,
29 % include some more diagnostics
30 % meanfields = 1 is the default
31 end
32 %
33 mit_timesteps
34
35 if meanfields
36 disp('reading averaged fields (*tave.*)')
37 %
38 mit_loadglobaldata
39 %
40 u=rdmds('uVeltave',timesteps(2:nt)');
41 v=rdmds('vVeltave',timesteps(2:nt)');
42 w=rdmds('wVeltave',timesteps(2:nt)');
43 t=rdmds('Ttave',timesteps(2:nt)');
44 s=rdmds('Stave',timesteps(2:nt)');
45 eta=rdmds('ETAtave',timesteps(2:nt)').*repmat(squeeze(grd.cmask(:,:,1)),[1 1 nt-1]);
46 k=1;
47 u=cat(4,zeros([nx ny nz 1]),u);
48 v=cat(4,zeros([nx ny nz 1]),v);
49 w=cat(4,zeros([nx ny nz 1]),w);
50 % hack
51 t=cat(4,rdmds('T',timesteps(1)'),t);
52 s=cat(4,rdmds('S',timesteps(1)'),s);
53 eta=cat(3,0.*grd.cmask(:,:,1),eta);
54 else
55 % load snap shots
56 disp('reading snap shot fields')
57 u=rdmds('U',timesteps');
58 v=rdmds('V',timesteps');
59 w=rdmds('W',timesteps');
60 t=rdmds('T',timesteps');
61 s=rdmds('S',timesteps');
62 eta=rdmds('Eta',timesteps').*repmat(squeeze(grd.cmask(:,:,1)),[1 1 nt]);
63 end
64 % orientation and masking
65 if ~strcmp(grd.buoyancy,'OCEANIC');
66 u = u(:,:,end:-1:1,:);
67 v = v(:,:,end:-1:1,:);
68 w = w(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
69 t = t(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
70 s = s(:,:,end:-1:1,:).*repmat(grd.cmask,[1 1 1 nt]);
71 else
72 % $$$ u = u;
73 % $$$ v = v;
74 w = w.*repmat(grd.cmask,[1 1 1 nt]);
75 t = t.*repmat(grd.cmask,[1 1 1 nt]);
76 s = s.*repmat(grd.cmask,[1 1 1 nt]);
77 end
78
79 % transport through Drake Passage
80 TDP = NaN*ones(nt,1);
81 kx = 73;
82 kyg = [5 6 7];
83 %kyg = [4 5 6 7];
84 da = grd.dz*grd.dyg(kx,kyg);
85 for k=kt;
86 TDP(k) = sum(nansum(squeeze(u(kx,kyg,:,k))'.*da))*1e-6; % in Sv
87 TDPz(:,k) = nansum(squeeze(u(kx,kyg,:,k)).*da')'*1e-6; % in Sv
88 end
89
90 % set some depth
91 if ~meanfields
92 mit_globalmovie
93 end
94
95 % $$$ [x y z]=meshgrid(grd.lonc,grd.grd.latc,-grd.zc);
96 % $$$ tk = permute(squeeze(t(:,:,:,k)),[2 1 3]);
97 % $$$ xplane = [0:30:360];
98 % $$$ yplane = 30;
99 % $$$ zplane = -[50 500:500:6000];
100 % $$$ slice(x,y,z,tk,xplane,yplane,zplane);
101 % $$$ shading flat, colorbar('h')
102
103 mit_plotstreamfunctions
104
105 mit_plotzonalvelocity
106
107 if meanfields
108 % compute actual heat and E-P fluxes from time averages
109 % this only makes sense when you load the time averaged fields
110 tauTheta = mit_getparm('data','tauThetaClimRelax');
111 tauSalt = mit_getparm('data','tauSaltClimRelax');
112 rhonil = mit_getparm('data','rhonil');
113 if isempty(rhonil)
114 rhonil = 1035;
115 end
116 Cp = 3994;
117 if isempty(tauTheta)
118 tauTheta = 0;
119 end
120 if tauTheta==0
121 recip_tauTheta = 0.;
122 else
123 recip_tauTheta = 1./tauTheta;
124 end
125 if isempty(tauSalt)
126 tauSalt = 0;
127 end
128 if tauSalt==0
129 recip_tauSalt = 0.;
130 else
131 recip_tauSalt = 1./tauSalt;
132 end
133 figure('PaperPosition',[0.31 0.25 10.5 7.88],'PaperOrientation','landscape')
134 k=kmax;
135 qnet_diag = (t(:,:,1,k)-tdmean)*grd.dz(1).*grd.hfacc(:,:,1)*Cp*rhonil*recip_tauTheta;
136 empr_diag = - (s(:,:,1,k)-sdmean)*grd.dz(1).*grd.hfacc(:,:,1)*recip_tauSalt/35;
137 clear sh;
138 sh(1) = subplot(2,2,1);imagesc(grd.lonc,grd.latc,qnet_diag');
139 title(['Q_{diag} from T-T_{lev} at the surface'])
140 sh(2) =subplot(2,2,2);imagesc(grd.lonc,grd.latc,qnet'.*grd.cmask(:,:,1)');
141 title('annual mean Q_{net} (data)')
142 sh(3) =subplot(2,2,3);imagesc(grd.lonc,grd.latc,empr_diag');
143 title(['(E-P)_{diag} from S-S_{lev} at the surface'])
144 sh(4) =subplot(2,2,4);imagesc(grd.lonc,grd.latc,empr'.*grd.cmask(:,:,1)');
145 title('annual mean E-P (data)')
146 axis(sh,'image'); axis(sh,'xy')
147 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
148 ', ' tuname ' = ' num2str(tim(k))])
149 set(sh(1:2),'clim',[-1 1]*200);
150 set(gcf,'currentAxes',sh(1));colorbar('h')
151 set(gcf,'currentAxes',sh(2));colorbar('h')
152 set(sh(3:4),'clim',[-1 1]*1.e-7);
153 set(gcf,'currentAxes',sh(3));colorbar('h')
154 set(gcf,'currentAxes',sh(4));colorbar('h')
155
156 clear sh
157
158 if length(tim) > 1
159 rac3d = repmat(grd.rac,[1 1 nz]).*grd.hfacc;
160 % horizontal averages of t and s in a hovmoeller-type diagramme
161 mdt = repmat(NaN,grd.nz,length(kt));
162 mdtt= mdt;
163 mdtstd= mdt;
164 mdt_alt = mdt;
165 mdtt_alt= mdt;
166 mdtstd_alt = mdt;
167 mdt_pac = mdt;
168 mdtt_pac= mdt;
169 mdtstd_pac = mdt;
170 mdt_sou = mdt;
171 mdtt_sou= mdt;
172 mdtstd_sou = mdt;
173 mds = mdt;
174 mdss= mdt;
175 mdsstd = mdt;
176 mds_alt = mdt;
177 mdss_alt= mdt;
178 mdsstd_alt = mdt;
179 mds_pac = mdt;
180 mdss_pac= mdt;
181 mdsstd_pac = mdt;
182 mds_sou = mdt;
183 mdss_sou= mdt;
184 mdsstd_sou = mdt;
185 indopac_hfacc = ...
186 change(change(grd.pacific_hfacc,'==',NaN,0) ...
187 +change(grd.indic_hfacc,'==',NaN,0),'==',0,NaN);
188 indopac_hfacc(:,1:12,:) = NaN;
189 southern_hfacc = grd.hfacc;
190 southern_hfacc(:,13:end,:) = NaN;
191 atl_hfacc = grd.atlantic_hfacc;
192 atl_hfacc(:,1:12,:) = NaN;
193 for k=kt;
194 dt = t(:,:,:,k)-tdatamean; % anomalies
195 ds = s(:,:,:,k)-sdatamean; % anomalies
196 [mdt(:,k) mdtt(:,k) mdtstd(:,k)] = mit_globalmean(dt,rac3d);
197 [mds(:,k) mdss(:,k) mdsstd(:,k)] = mit_globalmean(ds,rac3d);
198 [mdt_atl(:,k) mdtt_atl(:,k) mdtstd_atl(:,k)] = mit_globalmean(dt,rac3d.*atl_hfacc);
199 [mds_atl(:,k) mdss_atl(:,k) mdsstd_atl(:,k)] = mit_globalmean(ds,rac3d.*atl_hfacc);
200 [mdt_pac(:,k) mdtt_pac(:,k) mdtstd_pac(:,k)] = mit_globalmean(dt,rac3d.*indopac_hfacc);
201 [mds_pac(:,k) mdss_pac(:,k) mdsstd_pac(:,k)] = mit_globalmean(ds,rac3d.*indopac_hfacc);
202 [mdt_sou(:,k) mdtt_sou(:,k) mdtstd_sou(:,k)] = mit_globalmean(dt,rac3d.*southern_hfacc);
203 [mds_sou(:,k) mdss_sou(:,k) mdsstd_sou(:,k)] = mit_globalmean(ds,rac3d.*southern_hfacc);
204 end
205
206 figure('PaperPosition',[0.31 0.25 10.5 7.88],'PaperOrientation','landscape')
207 if length(tim) > 2
208 clear sh
209 sh(1) = subplot(2,2,1);
210 contourf(tim(2:end),-grd.zc,mdt(:,2:end),20);
211 caxis([-1 1]*max(abs(caxis))); shading flat; colorbar('v')
212 title('T-T_{lev} horizontally averaged')
213 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
214 sh(2) = subplot(2,2,2);
215 contourf(tim(2:end),-grd.zc,mds(:,2:end),20);
216 caxis([-1 1]*max(abs(caxis))); shading flat; colorbar('v')
217 title('S-S_{lev} horizontally averaged')
218 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
219 %
220 sh(3) = subplot(2,2,3);
221 contourf(tim(2:end),-grd.zc,mdtt(:,2:end),20);
222 shading flat; colorbar('v')
223 title('|T-T_{lev}| horizontally averaged');
224 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
225 sh(4) = subplot(2,2,4);
226 contourf(tim(2:end),-grd.zc,mdss(:,2:end),20);
227 shading flat; colorbar('v')
228 title('|S-S_{lev}| horizontally averaged')
229 xlabel(['Time [' timeunit ']']); ylabel('Depth [m]')
230 set(sh,'layer','top');
231 suptitle(['experiment ' dname])
232 end %if length(tim) > 2
233
234 k=kmax;
235 figure('PaperOrientation','landscape','PaperPosition',[0.25 0.3125 10.5 7.875])
236 clear sh
237 sh(1) = subplot(2,4,1);
238 plot(mdt(:,k),-grd.zc,'b-',mdt(:,k)+mdtstd(:,k),-grd.zc,'r--',mdt(:,k)-mdtstd(:,k),-grd.zc,'r--');
239 xlabel('T-T_{lev} [degC]'); ylabel('depth [m]'); title('global')
240 legend('horiz. mean','std',3);
241 sh(2) = subplot(2,4,2);
242 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--');
243 xlabel('T-T_{lev} [degC]'); %ylabel('depth [m]');
244 title('atlantic')
245 sh(3) = subplot(2,4,3);
246 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--');
247 xlabel('T-T_{lev} [degC]'); %ylabel('depth [m]');
248 title('indo-pacific')
249 sh(4) = subplot(2,4,4);
250 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--');
251 xlabel('T-T_{lev} [degC]'); ylabel('depth [m]'); title('southern')
252 %
253 sh(5) = subplot(2,4,5);
254 plot(mds(:,k),-grd.zc,'b-',mds(:,k)+mdsstd(:,k),-grd.zc,'r--',mds(:,k)-mdsstd(:,k),-grd.zc,'r--');
255 xlabel('S-S_{lev} [PSU]'); ylabel('depth [m]');
256 sh(6) = subplot(2,4,6);
257 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--');
258 xlabel('S-S_{lev} [PSU]'); %ylabel('depth [m]');
259 sh(7) = subplot(2,4,7);
260 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--');
261 xlabel('S-S_{lev} [PSU]'); %ylabel('depth [m]');
262 sh(8) = subplot(2,4,8);
263 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--');
264 xlabel('S-S_{lev} [PSU]'); ylabel('depth [m]');
265 set(sh,'xgrid','on','ygrid','on')
266 set([sh(4); sh(8)],'YAxisLocation','right')
267 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
268 ', ' tuname ' = ' num2str(tim(k))])
269 end %if length(tim) > 1
270
271 % zonal mean of temperature and salinity for different ocean basins
272 k=kmax;
273 tdzm_glo = mit_zonalmean(tdatamean,grd.hfacc,grd.dxc);
274 tdzm_atl = mit_zonalmean(tdatamean,grd.atlantic_hfacc,grd.dxc);
275 tdzm_pac = mit_zonalmean(tdatamean,grd.pacific_hfacc,grd.dxc);
276 tdzm_ind = mit_zonalmean(tdatamean,grd.indic_hfacc,grd.dxc);
277 sdzm_glo = mit_zonalmean(sdatamean,grd.hfacc,grd.dxc);
278 sdzm_atl = mit_zonalmean(sdatamean,grd.atlantic_hfacc,grd.dxc);
279 sdzm_pac = mit_zonalmean(sdatamean,grd.pacific_hfacc,grd.dxc);
280 sdzm_ind = mit_zonalmean(sdatamean,grd.indic_hfacc,grd.dxc);
281
282 figure('PaperPosition',[0.25 0.368552 8 10.2629]);
283 tlev = -5:.5:5;
284 slev = -1:.1:1;
285 clear sh clh
286 clh = cell(8,1);
287 delay = 1;
288 for k = kmax
289 tzm_glo = mit_zonalmean(t(:,:,:,k),grd.hfacc,grd.dxc);
290 tzm_atl = mit_zonalmean(t(:,:,:,k),grd.atlantic_hfacc,grd.dxc);
291 tzm_pac = mit_zonalmean(t(:,:,:,k),grd.pacific_hfacc,grd.dxc);
292 tzm_ind = mit_zonalmean(t(:,:,:,k),grd.indic_hfacc,grd.dxc);
293 szm_glo = mit_zonalmean(s(:,:,:,k),grd.hfacc,grd.dxc);
294 szm_atl = mit_zonalmean(s(:,:,:,k),grd.atlantic_hfacc,grd.dxc);
295 szm_pac = mit_zonalmean(s(:,:,:,k),grd.pacific_hfacc,grd.dxc);
296 szm_ind = mit_zonalmean(s(:,:,:,k),grd.indic_hfacc,grd.dxc);
297 caxt = [min(tzm_glo(:)-tdzm_glo(:)) max(tzm_glo(:)-tdzm_glo(:))];
298 tlev = -5:.5:5;
299 if max(abs(caxt)) < 1;
300 tlev = -1:.1:1;
301 end
302 if max(abs(caxt)) < .1;
303 tlev = .1*tlev;
304 end
305 caxs = [min(szm_glo(:)-sdzm_glo(:)) max(szm_glo(:)-sdzm_glo(:))];
306 slev = -1.:.05:1;
307 if max(abs(caxs)) < 0.2;
308 slev = -.20:.01:.20;
309 end
310 if max(abs(caxs)) < .02;
311 slev = .1*slev;
312 end
313 sh(1) = subplot(4,2,1);
314 [cs h] = contourf(grd.latc,-grd.zc,(tzm_glo-tdzm_glo)',tlev);
315 if ~isempty(h);
316 clh{1} = clabel(cs,h,[0 0]);
317 end
318 title('\theta-\theta_{lev} [degC]: global ocean')
319 sh(3) = subplot(4,2,3);
320 [cs h] = contourf(grd.latc,-grd.zc,(tzm_atl-tdzm_atl)',tlev);
321 if ~isempty(h);
322 clh{3} = clabel(cs,h,[0 0]);
323 end
324 title('atlantic ocean')
325 sh(5) = subplot(4,2,5);
326 [cs h] = contourf(grd.latc,-grd.zc,(tzm_pac-tdzm_pac)',tlev);
327 if ~isempty(h);
328 clh{5} = clabel(cs,h,[0 0]);
329 end
330 title('pacific ocean')
331 sh(7) = subplot(4,2,7);
332 [cs h] = contourf(grd.latc,-grd.zc,(tzm_ind-tdzm_ind)',tlev);
333 if ~isempty(h);
334 clh{7} = clabel(cs,h,[0 0]);
335 end
336 title('indian ocean')
337 set(sh(1:2:end),'clim',[tlev(1) tlev(end)])
338 colorbar('h')
339 sh(2) = subplot(4,2,2);
340 [cs h] = contourf(grd.latc,-grd.zc,(szm_glo-sdzm_glo)',slev);
341 if ~isempty(h);
342 clh{2} = clabel(cs,h,[0 0]);
343 end
344 title('S-S_{lev} [PSU]: global ocean')
345 sh(4) = subplot(4,2,4);
346 [cs h] = contourf(grd.latc,-grd.zc,(szm_atl-sdzm_atl)',slev);
347 if ~isempty(h);
348 clh{4} = clabel(cs,h,[0 0]);
349 end
350 title('atlantic ocean')
351 sh(6) = subplot(4,2,6);
352 [cs h] = contourf(grd.latc,-grd.zc,(szm_pac-sdzm_pac)',slev);
353 if ~isempty(h);
354 clh{6} = clabel(cs,h,[0 0]);
355 end
356 title('pacific ocean')
357 sh(8) = subplot(4,2,8);
358 [cs h] = contourf(grd.latc,-grd.zc,(szm_ind-sdzm_ind)',slev);
359 if ~isempty(h);
360 clh{8} = clabel(cs,h,[0 0]);
361 end
362 title('indian ocean')
363 set(sh(2:2:end),'clim',[slev(1) slev(end)])
364 colorbar('h')
365 set(sh,'layer','top')
366 if ~isempty(cat(1,clh{:}))
367 set(cat(1,clh{:}),'fontsize',8);
368 end
369
370 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
371 ', ' tuname ' = ' num2str(tim(k)) ', zonal averages'])
372 drawnow;
373 pause(delay)
374 end
375 clear clh sh
376
377 % meridional transport of heat (internal energy) and fresh water
378 % diagnosed from surface fluxes
379 petawatts = 1e-15;
380 sverdrups = 1e-6;
381 heat_data = -mit_meridflux(qnet,grd.dxc,grd.dyc)*petawatts;
382 heat_diag = -mit_meridflux(qnet_diag,grd.dxc,grd.dyc)*petawatts;
383 freshwater_data = -mit_meridflux(empr,grd.dxc,grd.dyc)*sverdrups;
384 freshwater_diag = -mit_meridflux(empr_diag,grd.dxc,grd.dyc)*sverdrups;
385 figure
386 latgf = [2*grd.latc-grd.latg];
387 clear sh
388 sh(1) = subplot(2,1,1);
389 hh = plot(latgf,heat_data,'-',latgf,heat_diag,'--', ...
390 latgf,heat_data+heat_diag,'-.');
391 title('heat (internal energy) flux [PW]')
392 suptitle(['experiment ' dname ', timestep = ' num2str(timesteps(k)) ...
393 ', ' tuname ' = ' num2str(tim(k))])
394 legend('surface flux','restoring','net',2)
395 sh(2) = subplot(2,1,2);
396 fwh = plot(latgf,freshwater_data,'-',latgf,freshwater_diag,'--', ...
397 latgf,freshwater_data+freshwater_diag,'-.');
398 title('freshwater flux [Sv]')
399 legend('surface flux','restoring','net',3)
400 set(sh,'xlim',[latgf(1) latgf(end)],'xgrid','on','ygrid','on');
401
402 end % meanfields
403
404 mit_plotmeandrift
405
406 fname = [grd.dname '_' sprintf('%u',tim(kmax)) timeunit];
407 in = findstr(fname,'\');
408 fname(in) = [];
409 fn = gcf;
410 if meanfields
411 fname = [fname '.ps'];
412 f0 = fn-7;
413 else
414 fname = [fname '_snap.ps'];
415 f0 = fn-4;
416 end
417 print(f0,'-dpsc',fname);
418 for k = f0+1:fn
419 print(k,'-dpsc',fname,'-append')
420 end

  ViewVC Help
Powered by ViewVC 1.1.22