timour_matlab/mscripts/meanw1.m

%clear 
clear path 

global Nx Ny Nz
global lat long dz dm mdep
global delt_su su_its t_su delt
global descriptor this_path
global f deltaf Q beta r_expt r_heat H
global time rots it
global g Cp rho_bar alpha
global u v t w
global iterations


param_file_name =   ...
   input(' Please enter the name of the m-file with the parameters for this run : ','s') ;
feval(param_file_name) ;

iterations

itstart = input(' Please enter start iteration : ','s') 
itend = input(' Please enter end iteration : ','s') 


sizeit=size(iterations);
for i=1:sizeit(1)
iter(i)=eval(iterations(i,1:10));
end
nitstart=find(iter==eval(itstart))
nitend=find(iter==eval(itend))

path   = this_path 
cmdstr=['cd ' path ];
eval(cmdstr);
path=pwd

sumw=zeros(Nx,Ny,Nz);
counter=0;

for i=nitstart:nitend
ufilename=(['U.' iterations((i),1:10) ]) ;
vfilename=(['V.' iterations((i),1:10) ]) ;
wfilename=(['W.' iterations((i),1:10) ]) ;
u=rdmds(ufilename,'b');
v=rdmds(vfilename,'b');
w=rdmds(wfilename,'b');
%w=zeros(Nx,Ny,Nz+1);
%dx=dm;
%dy=dm;
%
%for k=Nz:-1:1,
% w(1:Nx-1,1:Ny-1,k)= w(1:Nx-1,1:Ny-1,k+1)...
%                    -u(2:Nx,1:Ny-1,k)*dz/dx +u(1:Nx-1,1:Ny-1,k)*dz/dx ...
%                    -v(1:Nx-1,2:Ny,k)*dz/dy +v(1:Nx-1,1:Ny-1,k)*dz/dy ;
%end
sumw=sumw+w;
counter=counter+1;
end
meanw=sumw/counter;

cmin=-1.e-4;
cmax=1.e-4,
V=[cmin cmax];

title='mean vertical velocity at the 4th level from the top';
        imagesc(lat,long,squeeze(meanw(:,:,4))');shading flat;axis image;caxis(V);colorbar('horizontal');
        set(gca,'ydir','norm')
        text(0,-30,descriptor);
        text(0,140,title);
        text(0,-20,'timestep');text(50,-20,num2str(it))

figure

meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
meanw=smooth3(meanw);
cmin=-1.e-4;
cmax=1.e-4,
V=[cmin cmax];

title='mean vertical velocity at the 4th level from the top';
        imagesc(lat,long,squeeze(meanw(:,:,4))');shading flat;axis image;caxis(V);colorbar('horizontal');
        set(gca,'ydir','norm')
        text(0,-30,descriptor);
        text(0,140,title);
        text(0,-20,'timestep');text(50,-20,num2str(it))

1	%clear
2	clear path
3
4	global Nx Ny Nz
5	global lat long dz dm mdep
6	global delt_su su_its t_su delt
7	global descriptor this_path
8	global f deltaf Q beta r_expt r_heat H
9	global time rots it
10	global g Cp rho_bar alpha
11	global u v t w
12	global iterations
13
14
15	param_file_name = ...
16	input(' Please enter the name of the m-file with the parameters for this run : ','s') ;
17	feval(param_file_name) ;
18
19	iterations
20
21	itstart = input(' Please enter start iteration : ','s')
22	itend = input(' Please enter end iteration : ','s')
23
24
25	sizeit=size(iterations);
26	for i=1:sizeit(1)
27	iter(i)=eval(iterations(i,1:10));
28	end
29	nitstart=find(iter==eval(itstart))
30	nitend=find(iter==eval(itend))
31
32	path = this_path
33	cmdstr=['cd ' path ];
34	eval(cmdstr);
35	path=pwd
36
37	sumw=zeros(Nx,Ny,Nz);
38	counter=0;
39
40	for i=nitstart:nitend
41	ufilename=(['U.' iterations((i),1:10) ]) ;
42	vfilename=(['V.' iterations((i),1:10) ]) ;
43	wfilename=(['W.' iterations((i),1:10) ]) ;
44	u=rdmds(ufilename,'b');
45	v=rdmds(vfilename,'b');
46	w=rdmds(wfilename,'b');
47	%w=zeros(Nx,Ny,Nz+1);
48	%dx=dm;
49	%dy=dm;
50	%
51	%for k=Nz:-1:1,
52	% w(1:Nx-1,1:Ny-1,k)= w(1:Nx-1,1:Ny-1,k+1)...
53	% -u(2:Nx,1:Ny-1,k)dz/dx +u(1:Nx-1,1:Ny-1,k)dz/dx ...
54	% -v(1:Nx-1,2:Ny,k)dz/dy +v(1:Nx-1,1:Ny-1,k)dz/dy ;
55	%end
56	sumw=sumw+w;
57	counter=counter+1;
58	end
59	meanw=sumw/counter;
60
61	cmin=-1.e-4;
62	cmax=1.e-4,
63	V=[cmin cmax];
64
65	title='mean vertical velocity at the 4th level from the top';
66	imagesc(lat,long,squeeze(meanw(:,:,4))');shading flat;axis image;caxis(V);colorbar('horizontal');
67	set(gca,'ydir','norm')
68	text(0,-30,descriptor);
69	text(0,140,title);
70	text(0,-20,'timestep');text(50,-20,num2str(it))
71
72	figure
73
74	meanw=smooth3(meanw);
75	meanw=smooth3(meanw);
76	meanw=smooth3(meanw);
77	meanw=smooth3(meanw);
78	meanw=smooth3(meanw);
79	meanw=smooth3(meanw);
80	meanw=smooth3(meanw);
81	meanw=smooth3(meanw);
82	meanw=smooth3(meanw);
83	meanw=smooth3(meanw);
84	cmin=-1.e-4;
85	cmax=1.e-4,
86	V=[cmin cmax];
87
88	title='mean vertical velocity at the 4th level from the top';
89	imagesc(lat,long,squeeze(meanw(:,:,4))');shading flat;axis image;caxis(V);colorbar('horizontal');
90	set(gca,'ydir','norm')
91	text(0,-30,descriptor);
92	text(0,140,title);
93	text(0,-20,'timestep');text(50,-20,num2str(it))
94