NCEPLIBS-w3emc 2.12.0
Loading...
Searching...
No Matches
w3ft214.f
Go to the documentation of this file.
1C> @file
2C> @brief Convert (361,91) grid to (97,69) n. hemi. grid
3C> @author Ralph Jones @date 1993-10-19
4
5C> Convert a northern hemisphere 1.0 degree lat.,lon. 361 by
6C> 91 grid to a polar stereographic 97 by 69 grid. The polar
7C> stereographic map projection is true at 60 deg. n. , The mesh
8C> length is 47.625 km. and the oriention is 150 deg. w.
9C> awips grid 214 regional - alaska - double resolution
10C>
11C> ### Program History Log:
12C> Date | Programmer | Comment
13C> -----|------------|--------
14C> 1993-10-19 | Ralph Jones | Initial.
15C>
16C> @param[in] ALOLA 361*91 grid 1.0 deg. lat,lon grid n. hemi.
17C> 32851 point grid. 360 * 181 one degree grib grid 3 was flipped, greenwish added
18C> to right side and cut to 361 * 91.
19C> @param[in] INTERP 1 linear interpolation , ne.1 biquadratic
20C> @param[out] APOLA 97*69 grid of northern hemisphere. 6693 point grid is
21C> awips grid type 214
22C>
23C> @note
24C> - 1. W1 and w2 are used to store sets of constants which are
25C> reusable for repeated calls to the subroutine.
26C> - 2. Wind components are not rotated to the 97*69 grid orientation
27C> after interpolation. you may use w3fc08 to do this.
28C>
29C> @author Ralph Jones @date 1993-10-19
30 SUBROUTINE w3ft214(ALOLA,APOLA,INTERP)
31C
32 parameter(npts=6693,ii=97,jj=69)
33 parameter(orient=150.0,ipole=49,jpole=101)
34 parameter(xmesh=47.625)
35C
36 REAL R2(NPTS), WLON(NPTS)
37 REAL XLAT(NPTS), XI(II,JJ), XJ(II,JJ)
38 REAL XII(NPTS), XJJ(NPTS), ANGLE(NPTS)
39 REAL ALOLA(361,91), APOLA(NPTS), ERAS(NPTS,4)
40 REAL W1(NPTS), W2(NPTS)
41 REAL XDELI(NPTS), XDELJ(NPTS)
42 REAL XI2TM(NPTS), XJ2TM(NPTS)
43C
44 INTEGER IV(NPTS), JV(NPTS), JY(NPTS,4)
45 INTEGER IM1(NPTS), IP1(NPTS), IP2(NPTS)
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.
60C
61 IF (iswt.EQ.1) GO TO 900
62C
63 deg = 1.0
64 gi2 = (1.86603 * earthr) / xmesh
65 gi2 = gi2 * gi2
66C
67C NEXT 32 LINES OF CODE PUTS SUBROUTINE W3FB05 IN LINE
68C
69 DO 100 j = 1,jj
70 xj1 = j - jpole
71 DO 100 i = 1,ii
72 xi(i,j) = i - ipole
73 xj(i,j) = xj1
74 100 CONTINUE
75C
76 DO 200 kk = 1,npts
77 r2(kk) = xjj(kk) * xjj(kk) + xii(kk) * xii(kk)
78 xlat(kk) = degprd *
79 & asin((gi2 - r2(kk)) / (gi2 + r2(kk)))
80 200 CONTINUE
81C
82 DO 300 kk = 1,npts
83 angle(kk) = degprd * atan2(xjj(kk),xii(kk))
84 300 CONTINUE
85C
86 DO 400 kk = 1,npts
87 IF (angle(kk).LT.0.0) angle(kk) = angle(kk) + 360.0
88 400 CONTINUE
89C
90 DO 500 kk = 1,npts
91 wlon(kk) = 270.0 + orient - angle(kk)
92 500 CONTINUE
93C
94 DO 600 kk = 1,npts
95 IF (wlon(kk).LT.0.0) wlon(kk) = wlon(kk) + 360.0
96 600 CONTINUE
97C
98 DO 700 kk = 1,npts
99 IF (wlon(kk).GE.360.0) wlon(kk) = wlon(kk) - 360.0
100 700 CONTINUE
101C
102 DO 800 kk = 1,npts
103 w1(kk) = (360.0 - wlon(kk)) / deg + 1.0
104 w2(kk) = xlat(kk) / deg + 1.0
105 800 CONTINUE
106C
107 iswt = 1
108 intrpo = interp
109 GO TO 1000
110C
111C AFTER THE 1ST CALL TO W3FT214 TEST INTERP, IF IT HAS
112C CHANGED RECOMPUTE SOME CONSTANTS
113C
114 900 CONTINUE
115 IF (interp.EQ.intrpo) GO TO 2100
116 intrpo = interp
117C
118 1000 CONTINUE
119 DO 1100 k = 1,npts
120 iv(k) = w1(k)
121 jv(k) = w2(k)
122 xdeli(k) = w1(k) - iv(k)
123 xdelj(k) = w2(k) - jv(k)
124 ip1(k) = iv(k) + 1
125 jy(k,3) = jv(k) + 1
126 jy(k,2) = jv(k)
127 1100 CONTINUE
128C
129 IF (lin) GO TO 1400
130C
131 DO 1200 k = 1,npts
132 ip2(k) = iv(k) + 2
133 im1(k) = iv(k) - 1
134 jy(k,1) = jv(k) - 1
135 jy(k,4) = jv(k) + 2
136 xi2tm(k) = xdeli(k) * (xdeli(k) - 1.0) * .25
137 xj2tm(k) = xdelj(k) * (xdelj(k) - 1.0) * .25
138 1200 CONTINUE
139C
140 DO 1300 kk = 1,npts
141 IF (iv(kk).EQ.1) THEN
142 ip2(kk) = 3
143 im1(kk) = 360
144 ELSE IF (iv(kk).EQ.360) THEN
145 ip2(kk) = 2
146 im1(kk) = 359
147 ENDIF
148 1300 CONTINUE
149C
150 1400 CONTINUE
151C
152 IF (lin) GO TO 1700
153C
154 DO 1500 kk = 1,npts
155 IF (jv(kk).LT.2.OR.jv(kk).GT.89) xj2tm(kk) = 0.0
156 1500 CONTINUE
157C
158 DO 1600 kk = 1,npts
159 IF (ip2(kk).LT.1) ip2(kk) = 1
160 IF (im1(kk).LT.1) im1(kk) = 1
161 IF (ip2(kk).GT.361) ip2(kk) = 361
162 IF (im1(kk).GT.361) im1(kk) = 361
163 1600 CONTINUE
164C
165 1700 CONTINUE
166 DO 1800 kk = 1,npts
167 IF (iv(kk).LT.1) iv(kk) = 1
168 IF (ip1(kk).LT.1) ip1(kk) = 1
169 IF (iv(kk).GT.361) iv(kk) = 361
170 IF (ip1(kk).GT.361) ip1(kk) = 361
171 1800 CONTINUE
172C
173C LINEAR INTERPOLATION
174C
175 DO 1900 kk = 1,npts
176 IF (jy(kk,2).LT.1) jy(kk,2) = 1
177 IF (jy(kk,2).GT.91) jy(kk,2) = 91
178 IF (jy(kk,3).LT.1) jy(kk,3) = 1
179 IF (jy(kk,3).GT.91) jy(kk,3) = 91
180 1900 CONTINUE
181C
182 IF (.NOT.lin) THEN
183 DO 2000 kk = 1,npts
184 IF (jy(kk,1).LT.1) jy(kk,1) = 1
185 IF (jy(kk,1).GT.91) jy(kk,1) = 91
186 IF (jy(kk,4).LT.1) jy(kk,4) = 1
187 IF (jy(kk,4).GT.91) jy(kk,4) = 91
188 2000 CONTINUE
189 ENDIF
190C
191 2100 CONTINUE
192 IF (lin) THEN
193C
194C LINEAR INTERPOLATION
195C
196 DO 2200 kk = 1,npts
197 eras(kk,2) = (alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
198 & * xdeli(kk) + alola(iv(kk),jy(kk,2))
199 eras(kk,3) = (alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
200 & * xdeli(kk) + alola(iv(kk),jy(kk,3))
201 2200 CONTINUE
202C
203 DO 2300 kk = 1,npts
204 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
205 & * xdelj(kk)
206 2300 CONTINUE
207C
208 ELSE
209C
210C QUADRATIC INTERPOLATION
211C
212 DO 2400 kk = 1,npts
213 eras(kk,1)=(alola(ip1(kk),jy(kk,1))-alola(iv(kk),jy(kk,1)))
214 & * xdeli(kk) + alola(iv(kk),jy(kk,1)) +
215 & ( alola(im1(kk),jy(kk,1)) - alola(iv(kk),jy(kk,1))
216 & - alola(ip1(kk),jy(kk,1))+alola(ip2(kk),jy(kk,1)))
217 & * xi2tm(kk)
218 eras(kk,2)=(alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
219 & * xdeli(kk) + alola(iv(kk),jy(kk,2)) +
220 & ( alola(im1(kk),jy(kk,2)) - alola(iv(kk),jy(kk,2))
221 & - alola(ip1(kk),jy(kk,2))+alola(ip2(kk),jy(kk,2)))
222 & * xi2tm(kk)
223 eras(kk,3)=(alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
224 & * xdeli(kk) + alola(iv(kk),jy(kk,3)) +
225 & ( alola(im1(kk),jy(kk,3)) - alola(iv(kk),jy(kk,3))
226 & - alola(ip1(kk),jy(kk,3))+alola(ip2(kk),jy(kk,3)))
227 & * xi2tm(kk)
228 eras(kk,4)=(alola(ip1(kk),jy(kk,4))-alola(iv(kk),jy(kk,4)))
229 & * xdeli(kk) + alola(iv(kk),jy(kk,4)) +
230 & ( alola(im1(kk),jy(kk,4)) - alola(iv(kk),jy(kk,4))
231 & - alola(ip1(kk),jy(kk,4))+alola(ip2(kk),jy(kk,4)))
232 & * xi2tm(kk)
233 2400 CONTINUE
234C
235 DO 2500 kk = 1,npts
236 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
237 & * xdelj(kk) + (eras(kk,1) - eras(kk,2)
238 & - eras(kk,3) + eras(kk,4)) * xj2tm(kk)
239 2500 CONTINUE
240C
241 ENDIF
242C
243 RETURN
244 END
w3ft214
subroutine w3ft214(alola, apola, interp)
Convert a northern hemisphere 1.0 degree lat.,lon.
Definition w3ft214.f:31