timour_matlab/mscripts/asym2.m

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) ;

% compute the E-W asymmetry of the thermocline and its depth

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

sumps=zeros(Nx,Ny);
counter=0;

for i=nitstart:nitend
psfilename=(['PS.' iterations((i),1:10) ]) ;
ps=rdmds(psfilename,'b');
sumps=sumps+ps;
counter=counter+1;
end

meanps=sumps/counter;


% h=zeros(Nx,Ny);

hh=zeros(Nx,1);

%t0=input('Enter temperature : ')
%meantheta(:,:,Nz)=0.;
%for i=1:Nx,
%for j=1:Ny,
%kk=find(meantheta(i,j,:)<t0);
%if kk(1)>1
%h(i,j)=(kk(1)-1)*dz+dz*(meantheta(i,j,kk(1)-1)-t0)/(meantheta(i,j,kk(1)-1)-meantheta(i,j,kk(1)));
%else
%h(i,j)=0;
%end
%end
%end
% hmax=max(max(h))
% [I,J]=find(h==hmax)
% NNN=J


%
NNN=Ny/2+1
%
hh=meanps(:,NNN);
hmax=max(hh);
I=find(hh==hmax)

ii=find(hh>0);
ii1=size(ii);
istart=ii(1)
iend=ii(ii1(1))
x0=(I-istart)/(iend-istart)


plot(hh);
1	clear path
2
3	global Nx Ny Nz
4	global lat long dz dm mdep
5	global delt_su su_its t_su delt
6	global descriptor this_path
7	global f deltaf Q beta r_expt r_heat H
8	global time rots it
9	global g Cp rho_bar alpha
10	global u v t w
11	global iterations
12
13
14	param_file_name = ...
15	input(' Please enter the name of the m-file with the parameters for this run : ','s') ;
16	feval(param_file_name) ;
17
18	% compute the E-W asymmetry of the thermocline and its depth
19
20	iterations
21
22	itstart = input(' Please enter start iteration : ','s')
23	itend = input(' Please enter end iteration : ','s')
24
25
26	sizeit=size(iterations);
27	for i=1:sizeit(1)
28	iter(i)=eval(iterations(i,1:10));
29	end
30	nitstart=find(iter==eval(itstart))
31	nitend=find(iter==eval(itend))
32
33	path = this_path
34	cmdstr=['cd ' path ];
35	eval(cmdstr);
36	path=pwd
37
38	sumps=zeros(Nx,Ny);
39	counter=0;
40
41	for i=nitstart:nitend
42	psfilename=(['PS.' iterations((i),1:10) ]) ;
43	ps=rdmds(psfilename,'b');
44	sumps=sumps+ps;
45	counter=counter+1;
46	end
47
48	meanps=sumps/counter;
49
50
51
52
53	% h=zeros(Nx,Ny);
54
55	hh=zeros(Nx,1);
56
57	%t0=input('Enter temperature : ')
58	%meantheta(:,:,Nz)=0.;
59	%for i=1:Nx,
60	%for j=1:Ny,
61	%kk=find(meantheta(i,j,:)<t0);
62	%if kk(1)>1
63	%h(i,j)=(kk(1)-1)dz+dz(meantheta(i,j,kk(1)-1)-t0)/(meantheta(i,j,kk(1)-1)-meantheta(i,j,kk(1)));
64	%else
65	%h(i,j)=0;
66	%end
67	%end
68	%end
69	% hmax=max(max(h))
70	% [I,J]=find(h==hmax)
71	% NNN=J
72
73
74	%
75	NNN=Ny/2+1
76	%
77	hh=meanps(:,NNN);
78	hmax=max(hh);
79	I=find(hh==hmax)
80
81	ii=find(hh>0);
82	ii1=size(ii);
83	istart=ii(1)
84	iend=ii(ii1(1))
85	x0=(I-istart)/(iend-istart)
86
87
88	plot(hh);