1 |
cnh |
1.1 |
code obtained from faulks.lcs.mit.edu:/scratch/adcroft/cubed_sphere |
2 |
|
|
|
3 |
|
|
path('/home/dimitri/matlab/adcroft/bin',path); |
4 |
|
|
tic |
5 |
|
|
[dxg,dyg,dxf,dyf,dxc,dyc,dxv,dyu,Ec,Eu,Ev,Ez,latC,lonC,latG,lonG,... |
6 |
|
|
Q11,Q22,Q12, TUu,TUv,TVu,TVv ]=gengrid_fn(576,4,'conf','c',0,0); |
7 |
|
|
toc |
8 |
|
|
|
9 |
|
|
timings nireas faulks |
10 |
|
|
gengrid_fn(576,2,'conf','c',0,0) 47.9 60.0 |
11 |
|
|
gengrid_fn(576,4,'conf','c',0,0) 163.9 |
12 |
|
|
|
13 |
|
|
current 80S-80N, 1/4-deg isotropic grid: |
14 |
|
|
1440 * 1088 = 1566720 grid points |
15 |
|
|
27.7987 km at Equator, 4.8272 km at 80N |
16 |
|
|
|
17 |
|
|
tripolar grid: |
18 |
|
|
2304 grid boxes around Equator |
19 |
|
|
2*pi*6371/2304 = 17.3742 km grid at Equator |
20 |
|
|
switch to bipolar is at 60N where grid is 17.3742*cos(60) = 8.6871 km |
21 |
|
|
at North Pole grid size is 17.3742/pi = 5.5304 km |
22 |
|
|
smallest distance in Canadian Archipelago is approximately 17.3742/6 = 2.8957 km |
23 |
|
|
worst aspect ratio is approximately 3 |
24 |
|
|
number of grid points in sperical polar region is 2304*960 = 2211840 |
25 |
|
|
number of grid points in bipolar grid is 1152*1152 = 1327104 |
26 |
|
|
total number of grid points is 3538944 |
27 |
|
|
advantages, grid is perfectly orthogonal |
28 |
|
|
there are no singularities |
29 |
|
|
no modifications to sea-ice thermodynamics |
30 |
|
|
|
31 |
|
|
by comparison |
32 |
|
|
the cubed sphere with 2304 grid boxes around Equator |
33 |
|
|
using Alistair's eq 16 has: |
34 |
|
|
576*576*6 = 1990656 grid points (56% of tripolar grid, 127% of 1/4-deg grid) |
35 |
|
|
largest distance is 22.9683 (132% of tripolar grid) |
36 |
|
|
smallest distance is 3.2633 (112% of tripolar grid, 68% of 1/4-deg grid) |
37 |
|
|
worst aspect ratio is 1.6162 (53% of tripolar grid) |
38 |
|
|
but orthogonality is not as good as that of tripolar grid |
39 |
|
|
and sea-ice dynamics need some work |
40 |
|
|
also for the tripolar grid, removing land point will |
41 |
|
|
result in much more dramatic reduction in number of grid |
42 |
|
|
points because the three singularities are over land. |
43 |
|
|
|
44 |
|
|
|
45 |
|
|
|
46 |
|
|
table |
47 |
|
|
first row is min, max, and mean of dxg and dyg |
48 |
|
|
second row is min, max, and mean of dxg./dyg |
49 |
|
|
third row is min, max, and mean of angles, excluding cube corners |
50 |
|
|
|
51 |
|
|
conf, nratio=4 |
52 |
|
|
2.1460 18.0293 15.9182 |
53 |
|
|
0.9172 1.0903 1.0000 |
54 |
|
|
89.9978 98.5344 90.0562 |
55 |
|
|
|
56 |
|
|
conf, nratio=2 |
57 |
|
|
2.1460 18.0293 15.9182 |
58 |
|
|
0.9167 1.0909 1.0000 |
59 |
|
|
89.9978 98.5344 90.0562 |
60 |
|
|
|
61 |
|
|
q=0, nratio=4 |
62 |
|
|
2.2544 17.5336 15.9487 |
63 |
|
|
0.9172 1.0903 1.0006 |
64 |
|
|
89.9977 98.5344 90.0554 |
65 |
|
|
|
66 |
|
|
q=0, nratio=2 |
67 |
|
|
2.2544 17.4069 15.9468 |
68 |
|
|
0.9167 1.0909 1.0006 |
69 |
|
|
89.9977 98.5344 90.0554 |
70 |
|
|
|
71 |
|
|
q=1, nratio=4 (has NaNs) |
72 |
|
|
13.6379 67.1398 16.1959 |
73 |
|
|
0.2031 4.9226 1.0552 |
74 |
|
|
89.9972 98.5343 90.0497 |
75 |
|
|
|
76 |
|
|
q=1, nratio=2 |
77 |
|
|
13.6310 67.1393 16.1925 |
78 |
|
|
0.2031 4.9238 1.0551 |
79 |
|
|
89.9972 98.5343 90.0497 |
80 |
|
|
|
81 |
|
|
q=1/2, nratio=4 |
82 |
|
|
2.8990 21.0171 16.0604 |
83 |
|
|
0.6992 1.4302 1.0137 |
84 |
|
|
89.9974 98.5344 90.0527 |
85 |
|
|
|
86 |
|
|
q=1/2, nratio=2 |
87 |
|
|
2.8958 21.0001 16.0613 |
88 |
|
|
0.6999 1.4287 1.0136 |
89 |
|
|
89.9974 98.5344 90.0527 |
90 |
|
|
|
91 |
|
|
q=7/8, nratio=4 (has NaNs) |
92 |
|
|
13.6379 67.1398 16.1959 |
93 |
|
|
0.2031 4.9226 1.0552 |
94 |
|
|
89.9972 98.5343 90.0497 |
95 |
|
|
|
96 |
|
|
q=7/8, nratio=2 |
97 |
|
|
13.6310 67.1393 16.1925 |
98 |
|
|
0.2031 4.9238 1.0551 |
99 |
|
|
89.9972 98.5343 90.0497 |
100 |
|
|
|
101 |
|
|
q=i3, nratio=4 |
102 |
|
|
13.1144 65.1937 16.1921 |
103 |
|
|
0.2092 4.7807 1.0549 |
104 |
|
|
89.9972 98.5343 90.0497 |
105 |
|
|
|
106 |
|
|
q=i3, nratio=2 |
107 |
|
|
12.4259 62.6094 16.1915 |
108 |
|
|
0.2179 4.5899 1.0546 |
109 |
|
|
89.9973 98.5343 90.0498 |
110 |
|
|
|
111 |
|
|
tan, nratio=4 |
112 |
|
|
3.2633 22.9683 16.0912 |
113 |
|
|
0.6188 1.6162 1.0203 |
114 |
|
|
89.9974 98.5344 90.0519 |
115 |
|
|
|
116 |
|
|
tan, nratio=2 |
117 |
|
|
3.2633 22.9683 16.0912 |
118 |
|
|
0.6188 1.6162 1.0203 |
119 |
|
|
89.9974 98.5344 90.0519 |
120 |
|
|
|
121 |
|
|
tan2, nratio=4 |
122 |
|
|
2.2237 17.7883 15.9327 |
123 |
|
|
0.9174 1.0900 1.0001 |
124 |
|
|
89.9978 98.5344 90.0558 |
125 |
|
|
|
126 |
|
|
tan2, nratio=2 |
127 |
|
|
2.2237 17.7883 15.9327 |
128 |
|
|
0.9169 1.0906 1.0001 |
129 |
|
|
89.9978 98.5344 90.0558 |
130 |
|
|
|
131 |
|
|
new, nratio=4 |
132 |
|
|
2.5095 18.8618 15.9191 |
133 |
|
|
0.8889 1.1250 1.0000 |
134 |
|
|
89.9978 98.5344 90.0561 |
135 |
|
|
|
136 |
|
|
new, nratio=2 |
137 |
|
|
2.8863 20.9484 15.9201 |
138 |
|
|
0.8000 1.2500 1.0002 |
139 |
|
|
89.9978 98.5344 90.0561 |
140 |
|
|
|
141 |
|
|
|
142 |
|
|
|
143 |
|
|
|
144 |
|
|
|
145 |
|
|
|
146 |
|
|
gengrids % Does it all |
147 |
|
|
gengrids_fn(nx,nratio,projection,orientation,plotfigs,writedata) |
148 |
|
|
[dxg,dyg,dxf,dyf,dxc,dyc,dxv,dyu,Ec,Eu,Ev,Ez,latC,lonC,latG,lonG,... |
149 |
|
|
Q11,Q22,Q12, TUu,TUv,TVu,TVv ]=gengrid_fn(24,2,'conf','c',1,0); |
150 |
|
|
|
151 |
|
|
[X,Y,Z,x,y]=read_scubdat(fn,n); % Read SCUB0002.DAT file |
152 |
|
|
[X,Y,Z]=map_xy2xyz(x,y); % Map (x,y) -> (X,Y,Z) |
153 |
|
|
[X,Y,Z]=gmap_xy2xyz(x,y); % Map (x,y) -> (X,Y,Z) |
154 |
|
|
|
155 |
|
|
[lon,lat]=map_xyz2lonlat(X,Y,Z); % Map (X,Y,Z) -> (lon,lat) |
156 |
|
|
[X,Y,Z]=map_lonlat2xyz(lon,lat); % Map (lon,lat) -> (X,Y,Z) |
157 |
|
|
[lat,lon]=calc_geocoords(lx,ly,tile) % Map (x,y) -> (lat,lon) |
158 |
|
|
|
159 |
|
|
[X,Y,Z]=rotate_about_xaxis(lX,lY,lZ,angle); % Rotate about X axis |
160 |
|
|
[X,Y,Z]=rotate_about_yaxis(lX,lY,lZ,angle); % Rotate about Y axis |
161 |
|
|
[X,Y,Z]=rotate_about_zaxis(lX,lY,lZ,angle); % Rotate about Z axis |
162 |
|
|
|
163 |
|
|
V=coord2vector(XG,YG,ZG); % Convert [X,Y,Z] to vector V |
164 |
|
|
|
165 |
|
|
angle=angle_between_coords(X1,Y1,Z1,X2,Y2,Z2); % Angles between coords |
166 |
|
|
angle=angle_between_vectors(V1,V2); % Angles between vectors V1,V2 |
167 |
|
|
excess=excess_of_quad(V1,V2,V3,V4); % Excess of quadrilateral |
168 |
|
|
|
169 |
|
|
q=rescale_coordinate(q,'q=i3'); % Re-scale single coordinate |
170 |
|
|
[dx,dy,E]=calc_fvgrid(lx,ly); % Calculate F.V. grid info |
171 |
|
|
dxg=conf_dx(q); % Calculate conformal dxg |
172 |
|
|
[dxg,dxc,dxf,dxv]=reduce_dx(dx,nratio); % Sum fine-grid dx -> dxg,... |
173 |
|
|
[Ec,Ez,Ev]=reduce_E(dx,nratio); % Sum fine-grid E -> Ec,... |
174 |
|
|
|
175 |
|
|
a=bari(b); % Two point average in I dir. |
176 |
|
|
a=barj(b); % Two point average in J dir. |
177 |
|
|
|
178 |
|
|
write_tiles('LATC',lat) % Write tiled files |
179 |
|
|
displaytiles(lat) % Display tiles in comp. layout |
180 |
|
|
plotcube(X,Y,Z,C) % Display tiles as 3-D sphere |
181 |
|
|
merccube(lonG,latG,C) % Display tiles Mercator proj. |