/[MITgcm]/MITgcm_contrib/timour_matlab/mscripts/radmean.m
ViewVC logotype

Contents of /MITgcm_contrib/timour_matlab/mscripts/radmean.m

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


Revision 1.1 - (show annotations) (download)
Wed Sep 3 21:22:22 2003 UTC (21 years, 10 months ago) by edhill
Branch: MAIN
CVS Tags: HEAD
initial checkin of Timour's MatLAB scripts

1 % create circular sampling rings
2 %t=ones(Nx,Ny,Nz);
3 Mp1=Mp;
4 Dp1=Dp;
5
6 for kk=1:4
7 pv1=smooth3(pv1);
8 Mp1=smooth3(Mp1);
9 Dp1=smooth3(Dp1);
10 end
11
12
13 Nn=Nx/2;
14
15 rad=[1:1:(Nn-1)];
16 ring=ones(Nx,Ny,Nz);
17 weight=zeros(Nn-1,1);
18 Dring=zeros(Nx,Ny,Nz);
19 Mring=zeros(Nx,Ny,Nz);
20 radialmean=zeros(1,Nn-1);
21 XCenter = (Nx+1)/2 ;
22 YCenter = (Ny+1)/2 ;
23
24
25 for k=1:Nn-1
26 for i=1:Nx
27 for j=1:Ny
28 radius(i,j) = sqrt((i-XCenter)^2+(j-YCenter)^2);
29 end
30 end
31 radius(find( round(radius) > (k+1) ))=radius(find (round(radius) > (k+1) ))*0;
32 radius(find( round(radius) < (k) ))=radius(find (round(radius) < (k) ))*0;
33 radius(find(radius>0))=1.0;
34 ring(:,:,k)=radius;
35 weight(k)=length(find(ring(:,:,k) > 0));
36 end
37
38
39 for k=1:Nz
40 for i=1:Nn-1
41 Tring(:,:,k)=meantheta(:,:,k).*ring(:,:,i);
42 Dring(:,:,k)=Dp1(:,:,k).*ring(:,:,i);
43 Mring(:,:,k)=Mp1(:,:,k).*ring(:,:,i);
44 meanDp(i,k)=Nx*Ny*mean(mean(squeeze(Dring(:,:,k))))/weight(i);
45 meanMp(i,k)=Nx*Ny*mean(mean(squeeze(Mring(:,:,k))))/weight(i);
46 meanT(i,k)=Nx*Ny*mean(mean(squeeze(Tring(:,:,k))))/weight(i);
47 end
48 end
49
50 rad=dx:dx:(Nn-1)*dx;
51
52 for k=1:Nz
53 radd(:,k)=rad(:);
54 end
55
56 c=-0.05:0.005:0.005;
57 figure
58 A=meanDp(:,5:26).*radd(:,5:26);
59
60 i=find(A>0.00);
61 A(i)=A(i)/2;
62 %contourf(flipud(A'),c)
63 imagesc(flipud((A(1:Nn-5,:)')));
64 set(gca,'DataAspectRatio',[2,2,2])
65 set(gca,'ydir','norm')
66 xlabel('r ','FontSize',20)
67 ylabel('z ','Rotation',[0],'FontSize',20)
68 set(gca,'XtickLabel','|||')
69 set(gca,'YtickLabel','|||')
70 % title('D_p, azimuthal average ','Fontsize',20);
71 text(20,Nz-5,'D_p, azimuthal average ','Fontsize',20);
72
73 hold on
74 contour(flipud(meanT(1:Nn-5,5:26)')/1000,10);
75
76 %subplot(2,1,2)
77 %contourf(flipud((meanMp(:,5:26).*radd(:,5:26))'),-c);colorbar;
78
79 %--------TEST---------------------
80
81 av=6.28*mean(mean(meanDp(:,5:26).*radd(:,5:26)))*(Nn-1)*Nz*dx*dz
82 av1=mean(mean(mean(Dp1(:,:,5:26))))*Nx*Ny*(26-5)*dx*dy*dz
83 % vol=Nz*dz*3.14*((Nn-1)*dx)*((Nn-1)*dx)
84
85 avm=6.28*mean(mean(meanMp(:,5:26).*radd(:,5:26)))*(Nn-1)*Nz*dx*dz
86 avm1=mean(mean(mean(Mp1(:,:,5:26))))*Nx*Ny*(26-5)*dx*dy*dz

  ViewVC Help
Powered by ViewVC 1.1.22