NCEPLIBS-ip 4.0.0
polfix_mod.f90
1module polfix_mod
2 implicit none
3
4 private
5 public :: polfixs, polfixv
6
7contains
8
9 SUBROUTINE polfixs(NM,NX,KM,RLAT,IB,LO,GO)
10 !$$$ SUBPROGRAM DOCUMENTATION BLOCK
11 !
12 ! SUBPROGRAM: POLFIXS MAKE MULTIPLE POLE SCALAR VALUES CONSISTENT
13 ! PRGMMR: IREDELL ORG: W/NMC23 DATE: 96-04-10
14 !
15 ! ABSTRACT: THIS SUBPROGRAM AVERAGES MULTIPLE POLE SCALAR VALUES
16 ! ON A LATITUDE/LONGITUDE GRID. BITMAPS MAY BE AVERAGED TOO.
17 !
18 ! PROGRAM HISTORY LOG:
19 ! 96-04-10 IREDELL
20 !
21 ! USAGE: CALL POLFIXS(NM,NX,KM,RLAT,IB,LO,GO)
22 !
23 ! INPUT ARGUMENT LIST:
24 ! NO - INTEGER NUMBER OF GRID POINTS
25 ! NX - INTEGER LEADING DIMENSION OF FIELDS
26 ! KM - INTEGER NUMBER OF FIELDS
27 ! RLAT - REAL (NO) LATITUDES IN DEGREES
28 ! IB - INTEGER (KM) BITMAP FLAGS
29 ! LO - LOGICAL*1 (NX,KM) BITMAPS (IF SOME IB(K)=1)
30 ! GO - REAL (NX,KM) FIELDS
31 !
32 ! OUTPUT ARGUMENT LIST:
33 ! LO - LOGICAL*1 (NX,KM) BITMAPS (IF SOME IB(K)=1)
34 ! GO - REAL (NX,KM) FIELDS
35 !
36 ! ATTRIBUTES:
37 ! LANGUAGE: FORTRAN 90
38 !
39 !$$$
40 IMPLICIT NONE
41 !
42 INTEGER, INTENT(IN ) :: NM, NX, KM
43 INTEGER, INTENT(IN ) :: IB(KM)
44 !
45 LOGICAL*1, INTENT(INOUT) :: LO(NX,KM)
46 !
47 REAL, INTENT(IN ) :: RLAT(NM)
48 REAL, INTENT(INOUT) :: GO(NX,KM)
49 !
50 REAL, PARAMETER :: RLATNP=89.9995
51 REAL, PARAMETER :: RLATSP=-rlatnp
52 !
53 INTEGER :: K, N
54 !
55 REAL :: WNP, GNP, TNP, WSP, GSP, TSP
56 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
57 DO k=1,km
58 wnp=0.
59 gnp=0.
60 tnp=0.
61 wsp=0.
62 gsp=0.
63 tsp=0.
64 ! AVERAGE MULTIPLE POLE VALUES
65 !$OMP PARALLEL DO PRIVATE(N) REDUCTION(+:WNP,GNP,TNP,WSP,GSP,TSP) SCHEDULE(STATIC)
66 DO n=1,nm
67 IF(rlat(n).GE.rlatnp) THEN
68 wnp=wnp+1
69 IF(ib(k).EQ.0.OR.lo(n,k)) THEN
70 gnp=gnp+go(n,k)
71 tnp=tnp+1
72 ENDIF
73 ELSEIF(rlat(n).LE.rlatsp) THEN
74 wsp=wsp+1
75 IF(ib(k).EQ.0.OR.lo(n,k)) THEN
76 gsp=gsp+go(n,k)
77 tsp=tsp+1
78 ENDIF
79 ENDIF
80 ENDDO
81 !$OMP END PARALLEL DO
82 ! DISTRIBUTE AVERAGE VALUES BACK TO MULTIPLE POLES
83 IF(wnp.GT.1) THEN
84 IF(tnp.GE.wnp/2) THEN
85 gnp=gnp/tnp
86 ELSE
87 gnp=0.
88 ENDIF
89 !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
90 DO n=1,nm
91 IF(rlat(n).GE.rlatnp) THEN
92 IF(ib(k).NE.0) lo(n,k)=tnp.GE.wnp/2
93 go(n,k)=gnp
94 ENDIF
95 ENDDO
96 !$OMP END PARALLEL DO
97 ENDIF
98 IF(wsp.GT.1) THEN
99 IF(tsp.GE.wsp/2) THEN
100 gsp=gsp/tsp
101 ELSE
102 gsp=0.
103 ENDIF
104 !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
105 DO n=1,nm
106 IF(rlat(n).LE.rlatsp) THEN
107 IF(ib(k).NE.0) lo(n,k)=tsp.GE.wsp/2
108 go(n,k)=gsp
109 ENDIF
110 ENDDO
111 !$OMP END PARALLEL DO
112 ENDIF
113 ENDDO
114 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
115 END SUBROUTINE polfixs
116
117
118 SUBROUTINE polfixv(NM,NX,KM,RLAT,RLON,IB,LO,UO,VO)
119 !$$$ SUBPROGRAM DOCUMENTATION BLOCK
120 !
121 ! SUBPROGRAM: POLFIXV MAKE MULTIPLE POLE VECTOR VALUES CONSISTENT
122 ! PRGMMR: IREDELL ORG: W/NMC23 DATE: 96-04-10
123 !
124 ! ABSTRACT: THIS SUBPROGRAM AVERAGES MULTIPLE POLE VECTOR VALUES
125 ! ON A LATITUDE/LONGITUDE GRID. BITMAPS MAY BE AVERAGED TOO.
126 ! VECTORS ARE ROTATED WITH RESPECT TO THEIR LONGITUDE.
127 !
128 ! PROGRAM HISTORY LOG:
129 ! 96-04-10 IREDELL
130 !
131 ! USAGE: CALL POLFIXV(NM,NX,KM,RLAT,RLON,IB,LO,UO,VO)
132 !
133 ! INPUT ARGUMENT LIST:
134 ! NM - INTEGER NUMBER OF GRID POINTS
135 ! NX - INTEGER LEADING DIMENSION OF FIELDS
136 ! KM - INTEGER NUMBER OF FIELDS
137 ! RLAT - REAL (NM) LATITUDES IN DEGREES
138 ! RLON - REAL (NM) LONGITUDES IN DEGREES
139 ! IB - INTEGER (KM) BITMAP FLAGS
140 ! LO - LOGICAL*1 (NX,KM) BITMAPS (IF SOME IB(K)=1)
141 ! UO - REAL (NX,KM) U-WINDS
142 ! VO - REAL (NX,KM) V-WINDS
143 !
144 ! OUTPUT ARGUMENT LIST:
145 ! LO - LOGICAL*1 (NX,KM) BITMAPS (IF SOME IB(K)=1)
146 ! UO - REAL (NX,KM) U-WINDS
147 ! VO - REAL (NX,KM) V-WINDS
148 !
149 ! ATTRIBUTES:
150 ! LANGUAGE: FORTRAN 90
151 !
152 !$$$
153 IMPLICIT NONE
154 !
155 INTEGER, INTENT(IN ) :: IB(KM), NM, NX, KM
156 !
157 LOGICAL*1, INTENT(INOUT) :: LO(NX,KM)
158 !
159 REAL, INTENT(IN ) :: RLAT(NM), RLON(NM)
160 REAL, INTENT(INOUT) :: UO(NX,KM), VO(NX,KM)
161 !
162 REAL, PARAMETER :: RLATNP=89.9995
163 REAL, PARAMETER :: RLATSP=-rlatnp
164 REAL, PARAMETER :: PI=3.14159265358979
165 REAL, PARAMETER :: DPR=180./pi
166 !
167 INTEGER :: K, N
168 !
169 REAL :: CLON(NM),SLON(NM)
170 REAL :: TNP, UNP, VNP, WNP
171 REAL :: TSP, USP, VSP, WSP
172 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
173 !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
174 DO n=1,nm
175 clon(n)=cos(rlon(n)/dpr)
176 slon(n)=sin(rlon(n)/dpr)
177 ENDDO
178 !$OMP END PARALLEL DO
179 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
180 DO k=1,km
181 wnp=0.
182 unp=0.
183 vnp=0.
184 tnp=0.
185 wsp=0.
186 usp=0.
187 vsp=0.
188 tsp=0.
189 ! AVERAGE MULTIPLE POLE VALUES
190 !$OMP PARALLEL DO PRIVATE(N) REDUCTION(+:WNP,UNP,VNP,TNP,WSP,USP,VSP,TSP) SCHEDULE(STATIC)
191 DO n=1,nm
192 IF(rlat(n).GE.rlatnp) THEN
193 wnp=wnp+1
194 IF(ib(k).EQ.0.OR.lo(n,k)) THEN
195 unp=unp+(clon(n)*uo(n,k)-slon(n)*vo(n,k))
196 vnp=vnp+(slon(n)*uo(n,k)+clon(n)*vo(n,k))
197 tnp=tnp+1
198 ENDIF
199 ELSEIF(rlat(n).LE.rlatsp) THEN
200 wsp=wsp+1
201 IF(ib(k).EQ.0.OR.lo(n,k)) THEN
202 usp=usp+(clon(n)*uo(n,k)+slon(n)*vo(n,k))
203 vsp=vsp+(-slon(n)*uo(n,k)+clon(n)*vo(n,k))
204 tsp=tsp+1
205 ENDIF
206 ENDIF
207 ENDDO
208 !$OMP END PARALLEL DO
209 ! DISTRIBUTE AVERAGE VALUES BACK TO MULTIPLE POLES
210 IF(wnp.GT.1) THEN
211 IF(tnp.GE.wnp/2) THEN
212 unp=unp/tnp
213 vnp=vnp/tnp
214 ELSE
215 unp=0.
216 vnp=0.
217 ENDIF
218 !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
219 DO n=1,nm
220 IF(rlat(n).GE.rlatnp) THEN
221 IF(ib(k).NE.0) lo(n,k)=tnp.GE.wnp/2
222 uo(n,k)=clon(n)*unp+slon(n)*vnp
223 vo(n,k)=-slon(n)*unp+clon(n)*vnp
224 ENDIF
225 ENDDO
226 !$OMP END PARALLEL DO
227 ENDIF
228 IF(wsp.GT.1) THEN
229 IF(tsp.GE.wsp/2) THEN
230 usp=usp/wsp
231 vsp=vsp/wsp
232 ELSE
233 usp=0.
234 vsp=0.
235 ENDIF
236 !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
237 DO n=1,nm
238 IF(rlat(n).LE.rlatsp) THEN
239 IF(ib(k).NE.0) lo(n,k)=tsp.GE.wsp/2
240 uo(n,k)=clon(n)*usp-slon(n)*vsp
241 vo(n,k)=slon(n)*usp+clon(n)*vsp
242 ENDIF
243 ENDDO
244 !$OMP END PARALLEL DO
245 ENDIF
246 ENDDO
247 ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
248 END SUBROUTINE polfixv
249end module polfix_mod