NCEPLIBS-w3emc 2.12.0
Loading...
Searching...
No Matches
w3ft205.f
Go to the documentation of this file.
1C> @file
2C> @brief Convert (361,91) grid to (45,39) 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 45 by 39 grid. The polar
7C> stereographic map projection is true at 60 deg. n. , The mesh
8C> length is 190.5 km. and the oriention is 60 deg. w. pole
9C> point is at (i,j) = (27,57). new map is awips map 205.
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 lat,lon grid n. hemisphere
17C> 32851 point grid. 360 * 181 one degree grib grid 3 was flipped, greenwish added
18C> to righ side and cut to 361 * 91.
19C> @param[in] INTERP 1 linear interpolation , ne.1 biquadratic
20C> @param[out] APOLA 45*39 grid of northern hemisphere. 1755 point grid is
21C> awips grid type 205
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 45*39 grid orientation
27C> after interpolation. You may use w3fc08() to do this.
28C>
29C> @author Ralph Jones @date 1993-10-19
30 SUBROUTINE w3ft205(ALOLA,APOLA,INTERP)
31C
32 parameter(npts=1755,ii=45,jj=39)
33 parameter(orient=60.0,ipole=27,jpole=57)
34 parameter(xmesh=190.5)
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 xii(1647) = 1.0
83 DO 300 kk = 1,npts
84 angle(kk) = degprd * atan2(xjj(kk),xii(kk))
85 300 CONTINUE
86C
87 DO 400 kk = 1,npts
88 IF (angle(kk).LT.0.0) angle(kk) = angle(kk) + 360.0
89 400 CONTINUE
90C
91 DO 500 kk = 1,npts
92 wlon(kk) = 270.0 + orient - angle(kk)
93 500 CONTINUE
94C
95 DO 600 kk = 1,npts
96 IF (wlon(kk).LT.0.0) wlon(kk) = wlon(kk) + 360.0
97 600 CONTINUE
98C
99 DO 700 kk = 1,npts
100 IF (wlon(kk).GE.360.0) wlon(kk) = wlon(kk) - 360.0
101 700 CONTINUE
102C
103 DO 800 kk = 1,npts
104 w1(kk) = (360.0 - wlon(kk)) / deg + 1.0
105 w2(kk) = xlat(kk) / deg + 1.0
106 800 CONTINUE
107C
108 iswt = 1
109 intrpo = interp
110 GO TO 1000
111C
112C AFTER THE 1ST CALL TO W3FT203 TEST INTERP, IF IT HAS
113C CHANGED RECOMPUTE SOME CONSTANTS
114C
115 900 CONTINUE
116 IF (interp.EQ.intrpo) GO TO 2100
117 intrpo = interp
118C
119 1000 CONTINUE
120 DO 1100 k = 1,npts
121 iv(k) = w1(k)
122 jv(k) = w2(k)
123 xdeli(k) = w1(k) - iv(k)
124 xdelj(k) = w2(k) - jv(k)
125 ip1(k) = iv(k) + 1
126 jy(k,3) = jv(k) + 1
127 jy(k,2) = jv(k)
128 1100 CONTINUE
129C
130 IF (lin) GO TO 1400
131C
132 DO 1200 k = 1,npts
133 ip2(k) = iv(k) + 2
134 im1(k) = iv(k) - 1
135 jy(k,1) = jv(k) - 1
136 jy(k,4) = jv(k) + 2
137 xi2tm(k) = xdeli(k) * (xdeli(k) - 1.0) * .25
138 xj2tm(k) = xdelj(k) * (xdelj(k) - 1.0) * .25
139 1200 CONTINUE
140C
141 1400 CONTINUE
142C
143 IF (lin) GO TO 1700
144C
145 DO 1500 kk = 1,npts
146 IF (jv(kk).LT.2.OR.jv(kk).GT.89) xj2tm(kk) = 0.0
147 1500 CONTINUE
148C
149 1700 CONTINUE
150C
151 IF (.NOT.lin) THEN
152 DO 2000 kk = 1,npts
153 IF (jy(kk,1).LT.1) jy(kk,1) = 1
154 2000 CONTINUE
155 ENDIF
156C
157 2100 CONTINUE
158 IF (lin) THEN
159C
160C LINEAR INTERPOLATION
161C
162 DO 2200 kk = 1,npts
163 eras(kk,2) = (alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
164 & * xdeli(kk) + alola(iv(kk),jy(kk,2))
165 eras(kk,3) = (alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
166 & * xdeli(kk) + alola(iv(kk),jy(kk,3))
167 2200 CONTINUE
168C
169 DO 2300 kk = 1,npts
170 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
171 & * xdelj(kk)
172 2300 CONTINUE
173C
174 ELSE
175C
176C QUADRATIC INTERPOLATION
177C
178 DO 2400 kk = 1,npts
179 eras(kk,1)=(alola(ip1(kk),jy(kk,1))-alola(iv(kk),jy(kk,1)))
180 & * xdeli(kk) + alola(iv(kk),jy(kk,1)) +
181 & ( alola(im1(kk),jy(kk,1)) - alola(iv(kk),jy(kk,1))
182 & - alola(ip1(kk),jy(kk,1))+alola(ip2(kk),jy(kk,1)))
183 & * xi2tm(kk)
184 eras(kk,2)=(alola(ip1(kk),jy(kk,2))-alola(iv(kk),jy(kk,2)))
185 & * xdeli(kk) + alola(iv(kk),jy(kk,2)) +
186 & ( alola(im1(kk),jy(kk,2)) - alola(iv(kk),jy(kk,2))
187 & - alola(ip1(kk),jy(kk,2))+alola(ip2(kk),jy(kk,2)))
188 & * xi2tm(kk)
189 eras(kk,3)=(alola(ip1(kk),jy(kk,3))-alola(iv(kk),jy(kk,3)))
190 & * xdeli(kk) + alola(iv(kk),jy(kk,3)) +
191 & ( alola(im1(kk),jy(kk,3)) - alola(iv(kk),jy(kk,3))
192 & - alola(ip1(kk),jy(kk,3))+alola(ip2(kk),jy(kk,3)))
193 & * xi2tm(kk)
194 eras(kk,4)=(alola(ip1(kk),jy(kk,4))-alola(iv(kk),jy(kk,4)))
195 & * xdeli(kk) + alola(iv(kk),jy(kk,4)) +
196 & ( alola(im1(kk),jy(kk,4)) - alola(iv(kk),jy(kk,4))
197 & - alola(ip1(kk),jy(kk,4))+alola(ip2(kk),jy(kk,4)))
198 & * xi2tm(kk)
199 2400 CONTINUE
200C
201 DO 2500 kk = 1,npts
202 apola(kk) = eras(kk,2) + (eras(kk,3) - eras(kk,2))
203 & * xdelj(kk) + (eras(kk,1) - eras(kk,2)
204 & - eras(kk,3) + eras(kk,4)) * xj2tm(kk)
205 2500 CONTINUE
206C
207C NO POLE POINT
208C
209 ENDIF
210C
211 RETURN
212 END
w3ft205
subroutine w3ft205(alola, apola, interp)
Convert a northern hemisphere 1.0 degree lat.,lon.
Definition w3ft205.f:31