/[MITgcm]/MITgcm/utils/matlab/merccube.m
ViewVC logotype

Annotation of /MITgcm/utils/matlab/merccube.m

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


Revision 1.2 - (hide annotations) (download)
Wed Aug 15 15:47:19 2001 UTC (22 years, 9 months ago) by adcroft
Branch: MAIN
CVS Tags: checkpoint40pre8
Changes since 1.1: +14 -2 lines
Updated comments so we know exactly what this does.

1 adcroft 1.1 function [] = merccube(XX,YY,C)
2 adcroft 1.2 % merccube(x,y,c)
3 adcroft 1.1 %
4 adcroft 1.2 % Plots cubed-sphere data in mercator projection. (x,y) are
5 adcroft 1.1 % coordinates, c is cell-centered scalar to be plotted.
6 adcroft 1.2 % All arrays (x,y,c) should have dimensions of NxNx6 or 6NxN.
7     %
8     % The default plotting mode is shading faceted. Using this or
9     % shading flat, (x,y) should be the coordinates of grid-corners
10     % and can legitimately have dimension (N+1)x(N+1)x6.
11     %
12     % If using shading interp, then (x,y) must be the coordinates of
13     % the cell centers.
14 adcroft 1.1 %
15     % e.g.
16     %
17     % xg=rdmds('XG');
18     % yg=rdmds('YG');
19     % ps=rdmds('Eta.0000000000');
20     % mercube(xg,yg,ps);
21     % colorbar;xlabel('Longitude');ylabel('Latitude');
22 adcroft 1.2 %
23     % xc=rdmds('XC');
24     % yc=rdmds('YC');
25     % mercube(xc,yc,ps);shading interp
26 adcroft 1.1
27     if max(max(max(YY)))-min(min(min(YY))) < 3*pi
28     X=XX*180/pi;
29     Y=YY*180/pi;
30     else
31     X=XX;
32     Y=YY;
33     end
34     Q=C;
35    
36     if ndims(X)==2 & size(X,1)==6*size(X,2)
37     disp('1');
38     [nx ny nt]=size(X);
39     X=permute( reshape(X,[nx/6 6 ny]),[1 3 2]);
40     Y=permute( reshape(Y,[nx/6 6 ny]),[1 3 2]);
41     Q=permute( reshape(Q,[nx/6 6 ny]),[1 3 2]);
42     elseif ndims(X)==3 & size(X,2)==6
43     X=permute( X,[1 3 2]);
44     Y=permute( Y,[1 3 2]);
45     Q=permute( Q,[1 3 2]);
46     elseif ndims(X)==3 & size(X,3)==6
47     [nx ny nt]=size(X);
48     else
49     size(XX)
50     size(YY)
51     size(C)
52     error('Dimensions should be 2 or 3 dimensions: NxNx6, 6NxN or Nx6xN');
53     end
54    
55     if size(X,1)==size(Q,1)
56     whos
57     X(end+1,:,:)=NaN;
58     X(:,end+1,:)=NaN;
59     X(end,:,[1 3 5])=X(1,:,[2 4 6]);
60     X(:,end,[2 4 6])=X(:,1,[3 5 1]);
61     X(:,end,[1 3 5])=squeeze(X(1,end:-1:1,[3 5 1]));
62     X(end,:,[2 4 6])=squeeze(X(end:-1:1,1,[4 6 2]));
63     Y(end+1,:,:)=NaN;
64     Y(:,end+1,:)=NaN;
65     Y(end,:,[1 3 5])=Y(1,:,[2 4 6]);
66     Y(:,end,[2 4 6])=Y(:,1,[3 5 1]);
67     Y(:,end,[1 3 5])=squeeze(Y(1,end:-1:1,[3 5 1]));
68     Y(end,:,[2 4 6])=squeeze(Y(end:-1:1,1,[4 6 2]));
69     end
70     [nx ny nt]=size(X);
71    
72     Q(end+1,:,:)=0;
73     Q(:,end+1,:)=0;
74     Q(end,:,[1 3 5])=Q(1,:,[2 4 6]);
75     Q(:,end,[2 4 6])=Q(:,1,[3 5 1]);
76     Q(:,end,[1 3 5])=squeeze(Q(1,end:-1:1,[3 5 1]));
77     Q(end,:,[2 4 6])=squeeze(Q(end:-1:1,1,[4 6 2]));
78    
79     hnx=ceil(nx/2);
80     hny=ceil(ny/2);
81    
82     for k=1:6;
83     i=1:hnx;
84     x=longitude(X(i,:,k));
85     pcolor(x,Y(i,:,k),Q(i,:,k))
86     axis([-180 180 -90 90])
87     hold on
88     if max(max(max(x)))>180
89     pcolor(x-360,Y(i,:,k),Q(i,:,k))
90     end
91     i=hnx:nx;
92     x=longitude(X(i,:,k));
93     pcolor(x,Y(i,:,k),Q(i,:,k))
94     if max(max(max(x)))>180
95     pcolor(x-360,Y(i,:,k),Q(i,:,k))
96     end
97     end
98     hold off

  ViewVC Help
Powered by ViewVC 1.1.22