NCEPLIBS-w3emc 2.12.0
Loading...
Searching...
No Matches
w3ft06v.f
Go to the documentation of this file.
1C> @file
2C> @brief Convert (145,37) grid to (65,65) s. hemi. grid.
3C> @author Ralph Jones @date 1985-04-10
4
5C> Convert a southern hemisphere 2.5 degree lat.,lon. 145 by
6C> 37 grid to a polar stereographic 65 by 65 grid. The polar
7C> stereographic map projection is true at 60 deg. s.; The mesh
8C> length is 381 km. and the oriention is 260 deg. w.
9C>
10C> ### Program History Log:
11C> Date | Programmer | Comment
12C> -----|------------|--------
13C> 1985-04-10 | Ralph Jones | Vectorized version of w3ft05.
14C> 1989-10-21 | Ralph Jones | Changes to increase speed.
15C> 1991-07-24 | Ralph Jones | Change to cray cft77 fortran.
16C> 1993-05-31 | Ralph Jones | Recompile so linear interpolation works.
17C>
18C> @param[in] ALOLA - 145*37 gid 2.5 lat,lon grid s. hemishere. 5365 point
19C> grid is o.n.84 type 30 or 1e hex.
20C> @param[in] INTERP - 1 linear interpolation , ne.1 biquadratic.
21C> @param[out] APOLA - 65*65 grid of northern hemi. 4225 point grid is o.n. 84
22C> type 28 or 1c hex.
23C>
24C> @remark
25C> - 1. W1 and w2 are used to store sets of constants which are
26C> reusable for repeated calls to the subroutine.
27C> - 2. Wind components are not rotated to the 65*65 grid orientation
28C> after interpolation. You may use w3fc10 to do this.
29C> - 3. The grid points values on the equator have been extrapolated
30C> outward to all the grid points outside the equator on the 65*65
31C> grid (about 1100 points).
32C>
33C> @author Ralph Jones @date 1985-04-10
34 SUBROUTINE w3ft06v(ALOLA,APOLA,INTERP)
35C
36 REAL R2(4225), WLON(4225)
37 REAL XLAT(4225), XI(65,65), XJ(65,65)
38 REAL XII(4225), XJJ(4225), ANGLE(4225)
39 REAL ALOLA(145,37), APOLA(4225), ERAS(4225,4)
40 REAL W1(4225), W2(4225)
41 REAL XDELI(4225), XDELJ(4225)
42 REAL XI2TM(4225), XJ2TM(4225)
43C
44 INTEGER IV(4225), JV(4225), JY(4225,4)
45 INTEGER IM1(4225), IP1(4225), IP2(4225)
46C
47 LOGICAL LIN
48C
49 SAVE
50C
51 equivalence(xi(1,1),xii(1)),(xj(1,1),xjj(1))
52C
53 DATA degprd/57.2957795/
54 DATA earthr/6371.2/
55 DATA intrpo/99/
56 DATA iswt /0/
57C
58 lin = .false.
59 IF (interp.EQ.1) lin = .true.
60 IF (iswt.EQ.1) GO TO 900
61C
62 orient = 260.0
63 deg = 2.5
64 xmesh = 381.0
65 gi2 = (1.86603 * earthr) / xmesh
66 gi2 = gi2 * gi2
67C
68C NEXT 32 LINES OF CODE PUTS SUBROUTINE W3FB03 IN LINE
69C
70 DO 100 j = 1,65
71 xj1 = j - 33
72 DO 100 i = 1,65
73 xi(i,j) = i - 33
74 xj(i,j) = xj1
75 100 CONTINUE
76C
77 DO 200 kk = 1,4225
78 r2(kk) = xjj(kk) * xjj(kk) + xii(kk) * xii(kk)
79 xlat(kk) = -degprd *
80 & asin((gi2 - r2(kk)) / (gi2 + r2(kk)))
81 200 CONTINUE
82C
83 xii(2113) = 1.0
84 DO 300 kk = 1,4225
85 angle(kk) = degprd * atan2(xjj(kk),xii(kk))
86 300 CONTINUE
87C
88 DO 400 kk = 1,4225
89 IF (angle(kk).LT.0.0) angle(kk) = angle(kk) + 360.0
90 400 CONTINUE
91C
92 DO 500 kk = 1,4225
93 wlon(kk) = angle(kk) + orient - 270.0
94 500 CONTINUE
95C
96 DO 600 kk = 1,4225
97 IF (wlon(kk).GE.360.0) wlon(kk) = wlon(kk) - 360.0
98 600 CONTINUE
99C
100 DO 700 kk = 1,4225
101 IF (wlon(kk).LT.0.0) wlon(kk) = wlon(kk) + 360.0
102 700 CONTINUE
103C
104 xlat(2113) = -90.0
105 wlon(2113) = 0.0
106C
107 DO 800 kk = 1,4225
108 w1(kk) = (360.0 - wlon(kk)) / deg + 1.0
109 w2(kk) = (xlat(kk) + 90.0) / deg + 1.0
110 800 CONTINUE
111C
112 iswt = 1
113 intrpo = interp
114 GO TO 1000
115C
116C AFTER THE 1ST CALL TO W3FT05 TEST INTERP, IF IT HAS
117C CHANGED RECOMPUTE SOME CONSTANTS
118C
119 900 CONTINUE
120 IF (interp.EQ.intrpo) GO TO 2100
121 intrpo = interp
122C
123 1000 CONTINUE
124 DO 1100 k = 1,4225
125 iv(k) = w1(k)
126 jv(k) = w2(k)
127 xdeli(k) = w1(k) - iv(k)
128 xdelj(k) = w2(k) - jv(k)
129 ip1(k) = iv(k) + 1
130 jy(k,3) = jv(k) + 1
131 jy(k,2) = jv(k)
132 1100 CONTINUE
133C
134 IF (lin) GO TO 1400
135C
136 DO 1200 k = 1,4225
137 ip2(k) = iv(k) + 2
138 im1(k) = iv(k) - 1
139 jy(k,1) = jv(k) - 1
140 jy(k,4) = jv(k) + 2
141 xi2tm(k) = xdeli(k) * (xdeli(k) - 1.0) * .25
142 xj2tm(k) = xdelj(k) * (xdelj(k) - 1.0) * .25
143 1200 CONTINUE
144C
145 DO 1300 kk = 1,4225
146 IF (iv(kk).EQ.1) THEN
147 ip2(kk) = 3
148 im1(kk) = 144
149 ELSE IF (iv(kk).EQ.144) THEN
150 ip2(kk) = 2
151 im1(kk) = 143
152 ENDIF
153 1300 CONTINUE
154C
155 1400 CONTINUE
156C
157 IF (lin) GO TO 1700
158C
159 DO 1500 kk = 1,4225
160 IF (jv(kk).LT.2.OR.jv(kk).GT.35) xj2tm(kk) = 0.0
161 1500 CONTINUE
162C
163 DO 1600 kk = 1,4225
164 IF (ip2(kk).LT.1) ip2(kk) = 1
165 IF (im1(kk).LT.1) im1(kk) = 1
166 IF (ip2(kk).GT.145) ip2(kk) = 145
167 IF (im1(kk).GT.145) im1(kk) = 145
168 1600 CONTINUE
169C
170 1700 CONTINUE
171 DO 1800 kk = 1,4225
172 IF (iv(kk).LT.1) iv(kk) = 1
173 IF (ip1(kk).LT.1) ip1(kk) = 1
174 IF (iv(kk).GT.145) iv(kk) = 145
175 IF (ip1(kk).GT.145) ip1(kk) = 145
176 1800 CONTINUE
177C
178C LINEAR INTERPOLATION
179C
180 DO 1900 kk = 1,4225
181 IF (jy(kk,2).LT.1) jy(kk,2) = 1
182 IF (jy(kk,2).GT.37) jy(kk,2) = 37
183 IF (jy(kk,3).LT.1) jy(kk,3) = 1
184 IF (jy(kk,3).GT.37) jy(kk,3) = 37
185 1900 CONTINUE
186C
187 IF (.NOT.lin) THEN
188 DO 2000 kk = 1,4225
189 IF (jy(kk,1).LT.1) jy(kk,1) = 1
190 IF (jy(kk,1).GT.37) jy(kk,1) = 37
191 IF (jy(kk,4).LT.1) jy(kk,4) = 1
192 IF (jy(kk,4).GT.37) jy(kk,4) = 37
193 2000 CONTINUE
194 ENDIF
195C
196 2100 CONTINUE
197 IF (lin) THEN
198C
199C LINEAR INTERPOLATION
200C
201 DO 2200 kk = 1,4225
202 eras(kk,2) = (alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
203 & * xdeli(kk) + alola(iv(kk),jy(kk,2))
204 eras(kk,3) = (alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
205 & * xdeli(kk) + alola(iv(kk),jy(kk,3))
206 2200 CONTINUE
207C
208 DO 2300 kk = 1,4225
209 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
210 & * xdelj(kk)
211 2300 CONTINUE
212C
213 ELSE
214C
215C QUADRATIC INTERPOLATION
216C
217 DO 2400 kk = 1,4225
218 eras(kk,1)=(alola(ip1(kk),jy(kk,1))-alola(iv(kk),jy(kk,1)))
219 & * xdeli(kk) + alola(iv(kk),jy(kk,1)) +
220 & ( alola(im1(kk),jy(kk,1)) - alola(iv(kk),jy(kk,1))
221 & - alola(ip1(kk),jy(kk,1))+alola(ip2(kk),jy(kk,1)))
222 & * xi2tm(kk)
223 eras(kk,2)=(alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
224 & * xdeli(kk) + alola(iv(kk),jy(kk,2)) +
225 & ( alola(im1(kk),jy(kk,2)) - alola(iv(kk),jy(kk,2))
226 & - alola(ip1(kk),jy(kk,2))+alola(ip2(kk),jy(kk,2)))
227 & * xi2tm(kk)
228 eras(kk,3)=(alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
229 & * xdeli(kk) + alola(iv(kk),jy(kk,3)) +
230 & ( alola(im1(kk),jy(kk,3)) - alola(iv(kk),jy(kk,3))
231 & - alola(ip1(kk),jy(kk,3))+alola(ip2(kk),jy(kk,3)))
232 & * xi2tm(kk)
233 eras(kk,4)=(alola(ip1(kk),jy(kk,4))-alola(iv(kk),jy(kk,4)))
234 & * xdeli(kk) + alola(iv(kk),jy(kk,4)) +
235 & ( alola(im1(kk),jy(kk,4)) - alola(iv(kk),jy(kk,4))
236 & - alola(ip1(kk),jy(kk,4))+alola(ip2(kk),jy(kk,4)))
237 & * xi2tm(kk)
238 2400 CONTINUE
239C
240 DO 2500 kk = 1,4225
241 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
242 & * xdelj(kk) + (eras(kk,1) - eras(kk,2)
243 & - eras(kk,3) + eras(kk,4)) * xj2tm(kk)
244 2500 CONTINUE
245C
246 ENDIF
247C
248C SET POLE POINT , WMO STANDARD FOR U OR V
249C
250 apola(2113) = alola(73,1)
251C
252 RETURN
253 END
w3ft06v
subroutine w3ft06v(alola, apola, interp)
Convert a southern hemisphere 2.5 degree lat.,lon.
Definition w3ft06v.f:35