NCEPLIBS-ip  5.1.0
budget_interp_mod.F90
Go to the documentation of this file.
1 
5 
10 module budget_interp_mod
11  use gdswzd_mod
12  use polfix_mod
13  use ip_grids_mod
14  use ip_grid_descriptor_mod
15  use ip_grid_factory_mod
16  implicit none
17 
18  private
19  public :: interpolate_budget
20 
21  interface interpolate_budget
22  module procedure interpolate_budget_scalar
23  module procedure interpolate_budget_vector
24  end interface interpolate_budget
25 
26  ! Smallest positive real value (use for equality comparisons)
27  REAL :: TINYREAL=tiny(1.0)
28 
29 contains
30 
94  SUBROUTINE interpolate_budget_scalar(IPOPT,grid_in,grid_out, &
95  MI,MO,KM,IBI,LI,GI, &
96  NO,RLAT,RLON,IBO,LO,GO,IRET)
97  class(ip_grid), intent(in) :: grid_in, grid_out
98  INTEGER, INTENT(IN ) :: IBI(KM), IPOPT(20)
99  INTEGER, INTENT(IN ) :: KM, MI, MO
100  INTEGER, INTENT( OUT) :: IBO(KM), IRET, NO
101  !
102  LOGICAL*1, INTENT(IN ) :: LI(MI,KM)
103  LOGICAL*1, INTENT( OUT) :: LO(MO,KM)
104  !
105  REAL, INTENT(IN ) :: GI(MI,KM)
106  REAL, INTENT(INOUT) :: RLAT(MO), RLON(MO)
107  REAL, INTENT( OUT) :: GO(MO,KM)
108  !
109  REAL, PARAMETER :: FILL=-9999.
110  !
111  INTEGER :: I1, J1, I2, J2, IB, JB
112  INTEGER :: IX, JX, IXS, JXS
113  INTEGER :: K, KXS, KXT
114  INTEGER :: LB, LSW, MP, MSPIRAL, MX
115  INTEGER :: N, NB, NB1, NB2, NB3, NB4, NV, NX
116  INTEGER :: N11(MO),N21(MO),N12(MO),N22(MO)
117  !
118  REAL :: GB, LAT(1), LON(1)
119  REAL :: PMP, RB2, RLOB(MO), RLAB(MO), WB
120  REAL :: W11(MO), W21(MO), W12(MO), W22(MO)
121  REAL :: WO(MO,KM), XF, YF, XI, YI, XX, YY
122  REAL :: XPTS(MO),YPTS(MO),XPTB(MO),YPTB(MO)
123  REAL :: XXX(1), YYY(1)
124 
125  logical :: to_station_points
126 
127  class(ip_grid), allocatable :: grid_out2
128 
129  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
130  ! COMPUTE NUMBER OF OUTPUT POINTS AND THEIR LATITUDES AND LONGITUDES.
131  ! DO SUBSECTION OF GRID IF KGDSO(1) IS SUBTRACTED FROM 255.
132  iret=0
133 
134  select type(grid_out)
135  type is(ip_station_points_grid)
136  to_station_points = .true.
137  allocate(grid_out2, source = grid_out)
138  CALL gdswzd(grid_out2, 0,mo,fill,xpts,ypts,rlon,rlat,no)
139  IF(no.EQ.0) iret=3
140  CALL gdswzd(grid_in,-1,no,fill,xpts,ypts,rlon,rlat,nv)
141  IF(nv.EQ.0) iret=2
142  class default
143  to_station_points = .false.
144  allocate(grid_out2, source = grid_out)
145  CALL gdswzd(grid_out2, 0,mo,fill,xpts,ypts,rlon,rlat,no)
146  IF(no.EQ.0) iret=3
147  end select
148 
149  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
150  ! SET PARAMETERS
151  IF(ipopt(1).GT.16) iret=32
152  mspiral=max(ipopt(20),1)
153  nb1=ipopt(1)
154  IF(nb1.EQ.-1) nb1=2
155  IF(iret.EQ.0.AND.nb1.LT.0) iret=32
156  lsw=1
157  IF(ipopt(2).EQ.-2) lsw=2
158  IF(ipopt(1).EQ.-1.OR.ipopt(2).EQ.-1) lsw=0
159  IF(iret.EQ.0.AND.lsw.EQ.1.AND.nb1.GT.15) iret=32
160  mp=ipopt(3+ipopt(1))
161  IF(mp.EQ.-1.OR.mp.EQ.0) mp=50
162  IF(mp.LT.0.OR.mp.GT.100) iret=32
163  pmp=mp*0.01
164  IF(iret.EQ.0) THEN
165  nb2=2*nb1+1
166  rb2=1./nb2
167  nb3=nb2*nb2
168  nb4=nb3
169  IF(lsw.EQ.2) THEN
170  rb2=1./(nb1+1)
171  nb4=(nb1+1)**4
172  ELSEIF(lsw.EQ.1) THEN
173  nb4=ipopt(2)
174  DO ib=1,nb1
175  nb4=nb4+8*ib*ipopt(2+ib)
176  ENDDO
177  ENDIF
178  ELSE
179  nb3=0
180  nb4=1
181  ENDIF
182  DO k=1,km
183  DO n=1,no
184  go(n,k)=0.
185  wo(n,k)=0.
186  ENDDO
187  ENDDO
188  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
189  ! LOOP OVER SAMPLE POINTS IN OUTPUT GRID BOX
190  DO nb=1,nb3
191  ! LOCATE INPUT POINTS AND COMPUTE THEIR WEIGHTS
192  jb=(nb-1)/nb2-nb1
193  ib=nb-(jb+nb1)*nb2-nb1-1
194  lb=max(abs(ib),abs(jb))
195  wb=1
196  IF(lsw.EQ.2) THEN
197  wb=(nb1+1-abs(ib))*(nb1+1-abs(jb))
198  ELSEIF(lsw.EQ.1) THEN
199  wb=ipopt(2+lb)
200  ENDIF
201  IF(abs(wb).GT.tinyreal) THEN
202  !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
203  DO n=1,no
204  xptb(n)=xpts(n)+ib*rb2
205  yptb(n)=ypts(n)+jb*rb2
206  ENDDO
207  !$OMP END PARALLEL DO
208  if(to_station_points)then
209  CALL gdswzd(grid_in, 1,no,fill,xptb,yptb,rlob,rlab,nv)
210  CALL gdswzd(grid_in,-1,no,fill,xptb,yptb,rlob,rlab,nv)
211  else
212  CALL gdswzd(grid_out2, 1,no,fill,xptb,yptb,rlob,rlab,nv)
213  CALL gdswzd(grid_in,-1,no,fill,xptb,yptb,rlob,rlab,nv)
214  endif
215  IF(iret.EQ.0.AND.nv.EQ.0.AND.lb.EQ.0) iret=2
216  !$OMP PARALLEL DO PRIVATE(N,XI,YI,I1,I2,J1,J2,XF,YF) SCHEDULE(STATIC)
217  DO n=1,no
218  xi=xptb(n)
219  yi=yptb(n)
220  IF(abs(xi-fill).GT.tinyreal.AND.abs(yi-fill).GT.tinyreal) THEN
221  i1=int(xi)
222  i2=i1+1
223  j1=int(yi)
224  j2=j1+1
225  xf=xi-i1
226  yf=yi-j1
227  n11(n)=grid_in%field_pos(i1, j1)
228  n21(n)=grid_in%field_pos(i2, j1)
229  n12(n)=grid_in%field_pos(i1, j2)
230  n22(n)=grid_in%field_pos(i2, j2)
231  IF(min(n11(n),n21(n),n12(n),n22(n)).GT.0) THEN
232  w11(n)=(1-xf)*(1-yf)
233  w21(n)=xf*(1-yf)
234  w12(n)=(1-xf)*yf
235  w22(n)=xf*yf
236  ELSE
237  n11(n)=0
238  n21(n)=0
239  n12(n)=0
240  n22(n)=0
241  ENDIF
242  ELSE
243  n11(n)=0
244  n21(n)=0
245  n12(n)=0
246  n22(n)=0
247  ENDIF
248  ENDDO
249  !$OMP END PARALLEL DO
250  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
251  ! INTERPOLATE WITH OR WITHOUT BITMAPS
252  !$OMP PARALLEL DO PRIVATE(K,N,GB) SCHEDULE(STATIC)
253  DO k=1,km
254  DO n=1,no
255  IF(n11(n).GT.0) THEN
256  IF(ibi(k).EQ.0) THEN
257  gb=w11(n)*gi(n11(n),k)+w21(n)*gi(n21(n),k) &
258  +w12(n)*gi(n12(n),k)+w22(n)*gi(n22(n),k)
259  go(n,k)=go(n,k)+wb*gb
260  wo(n,k)=wo(n,k)+wb
261  ELSE
262  IF(li(n11(n),k)) THEN
263  go(n,k)=go(n,k)+wb*w11(n)*gi(n11(n),k)
264  wo(n,k)=wo(n,k)+wb*w11(n)
265  ENDIF
266  IF(li(n21(n),k)) THEN
267  go(n,k)=go(n,k)+wb*w21(n)*gi(n21(n),k)
268  wo(n,k)=wo(n,k)+wb*w21(n)
269  ENDIF
270  IF(li(n12(n),k)) THEN
271  go(n,k)=go(n,k)+wb*w12(n)*gi(n12(n),k)
272  wo(n,k)=wo(n,k)+wb*w12(n)
273  ENDIF
274  IF(li(n22(n),k)) THEN
275  go(n,k)=go(n,k)+wb*w22(n)*gi(n22(n),k)
276  wo(n,k)=wo(n,k)+wb*w22(n)
277  ENDIF
278  ENDIF
279  ENDIF
280  ENDDO
281  ENDDO
282  !$OMP END PARALLEL DO
283  ENDIF
284  ENDDO ! sub-grid points
285  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
286  ! COMPUTE OUTPUT BITMAPS AND FIELDS
287  ! KM is often 1 .. do not do OMP PARALLEL DO here
288  km_loop : DO k=1,km
289  ibo(k)=ibi(k)
290  !$OMP PARALLEL DO PRIVATE(N,LAT,LON,XXX,YYY,NV,XX,YY,IXS,JXS,MX,KXS,KXT,IX,JX,NX) SCHEDULE(STATIC)
291  n_loop : DO n=1,no
292  lo(n,k)=wo(n,k).GE.pmp*nb4
293  IF(lo(n,k)) THEN
294  go(n,k)=go(n,k)/wo(n,k)
295  ELSEIF (mspiral.GT.1) THEN
296  lat(1)=rlat(n)
297  lon(1)=rlon(n)
298  CALL gdswzd(grid_in,-1,1,fill,xxx,yyy,lon,lat,nv)
299  xx=xxx(1)
300  yy=yyy(1)
301  IF(nv.EQ.1)THEN
302  i1=nint(xx)
303  j1=nint(yy)
304  ixs=int(sign(1.,xx-i1))
305  jxs=int(sign(1.,yy-j1))
306  spiral_loop : DO mx=2,mspiral**2
307  kxs=int(sqrt(4*mx-2.5))
308  kxt=mx-(kxs**2/4+1)
309  SELECT CASE(mod(kxs,4))
310  CASE(1)
311  ix=i1-ixs*(kxs/4-kxt)
312  jx=j1-jxs*kxs/4
313  CASE(2)
314  ix=i1+ixs*(1+kxs/4)
315  jx=j1-jxs*(kxs/4-kxt)
316  CASE(3)
317  ix=i1+ixs*(1+kxs/4-kxt)
318  jx=j1+jxs*(1+kxs/4)
319  CASE DEFAULT
320  ix=i1-ixs*kxs/4
321  jx=j1+jxs*(kxs/4-kxt)
322  END SELECT
323  nx=grid_in%field_pos(ix, jx)
324  IF(nx.GT.0.)THEN
325  IF(li(nx,k).OR.ibi(k).EQ.0) THEN
326  go(n,k)=gi(nx,k)
327  lo(n,k)=.true.
328  cycle n_loop
329  ENDIF
330  ENDIF
331  ENDDO spiral_loop
332  ibo(k)=1
333  go(n,k)=0.
334  ELSE
335  ibo(k)=1
336  go(n,k)=0.
337  ENDIF
338  ELSE ! no spiral search option
339  ibo(k)=1
340  go(n,k)=0.
341  ENDIF
342  ENDDO n_loop
343  !$OMP END PARALLEL DO
344  ENDDO km_loop
345 
346  select type(grid_out2)
347  type is(ip_equid_cylind_grid)
348  CALL polfixs(no,mo,km,rlat,ibo,lo,go)
349  end select
350  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
351  END SUBROUTINE interpolate_budget_scalar
352 
423  SUBROUTINE interpolate_budget_vector(IPOPT,grid_in,grid_out, &
424  MI,MO,KM,IBI,LI,UI,VI, &
425  NO,RLAT,RLON,CROT,SROT,IBO,LO,UO,VO,IRET)
426  class(ip_grid), intent(in) :: grid_in, grid_out
427  INTEGER, INTENT(IN ) :: IPOPT(20), IBI(KM)
428  INTEGER, INTENT(IN ) :: KM, MI, MO
429  INTEGER, INTENT( OUT) :: IRET, NO, IBO(KM)
430  !
431  LOGICAL*1, INTENT(IN ) :: LI(MI,KM)
432  LOGICAL*1, INTENT( OUT) :: LO(MO,KM)
433  !
434  REAL, INTENT(IN ) :: UI(MI,KM),VI(MI,KM)
435  REAL, INTENT(INOUT) :: RLAT(MO),RLON(MO)
436  REAL, INTENT( OUT) :: UO(MO,KM),VO(MO,KM)
437  REAL, INTENT( OUT) :: CROT(MO),SROT(MO)
438  !
439  REAL, PARAMETER :: FILL=-9999.
440  !
441  INTEGER :: I1,I2,J1,J2,IB,JB,LSW,MP
442  INTEGER :: K,LB,N,NB,NB1,NB2,NB3,NB4,NV
443  INTEGER :: N11(MO),N21(MO),N12(MO),N22(MO)
444  !
445  LOGICAL :: SAME_GRID
446  !
447  REAL :: CM11,SM11,CM12,SM12
448  REAL :: CM21,SM21,CM22,SM22
449  REAL :: PMP,RB2
450  REAL :: C11(MO),C21(MO),C12(MO),C22(MO)
451  REAL :: S11(MO),S21(MO),S12(MO),S22(MO)
452  REAL :: W11(MO),W21(MO),W12(MO),W22(MO)
453  REAL :: UB,VB,WB,UROT,VROT
454  REAL :: U11,V11,U21,V21,U12,V12,U22,V22
455  REAL :: WI1,WJ1,WI2,WJ2
456  REAL :: WO(MO,KM),XI,YI
457  REAL :: XPTS(MO),YPTS(MO)
458  REAL :: XPTB(MO),YPTB(MO),RLOB(MO),RLAB(MO)
459 
460  logical :: to_station_points
461 
462  class(ip_grid), allocatable :: grid_out2
463 
464  ! Save coeffecients between calls and only compute if grids have changed
465  INTEGER, SAVE :: MIX=-1
466  REAL, ALLOCATABLE, SAVE :: CROI(:),SROI(:)
467  REAL, ALLOCATABLE, SAVE :: XPTI(:),YPTI(:),RLOI(:),RLAI(:)
468 
469  class(ip_grid), allocatable, save :: prev_grid_in
470 
471  iret=0
472 
473  ! Negative grid number means interpolate to subgrid
474  ! The type of the subgrid is calculated by 255 +
475  select type(grid_out)
476  type is(ip_station_points_grid)
477  to_station_points = .true.
478  allocate(grid_out2, source = grid_out)
479  CALL gdswzd(grid_out2, 0,mo,fill,xpts,ypts,rlon,rlat,no,crot,srot)
480  IF(no.EQ.0) iret=3
481  CALL gdswzd(grid_in,-1,no,fill,xpts,ypts,rlon,rlat,nv,crot,srot)
482  IF(nv.EQ.0) iret=2
483  class default
484  to_station_points = .false.
485  allocate(grid_out2, source = grid_out)
486  CALL gdswzd(grid_out2, 0,mo,fill,xpts,ypts,rlon,rlat,no,crot,srot)
487  end select
488 
489  if (.not. allocated(prev_grid_in)) then
490  allocate(prev_grid_in, source = grid_in)
491 
492  same_grid = .false.
493  else
494  same_grid = grid_in == prev_grid_in
495 
496  if (.not. same_grid) then
497  deallocate(prev_grid_in)
498  allocate(prev_grid_in, source = grid_in)
499  end if
500  end if
501 
502  IF(.NOT.same_grid) THEN
503  IF(mix.NE.mi) THEN
504  IF(mix.GE.0) DEALLOCATE(xpti,ypti,rloi,rlai,croi,sroi)
505  ALLOCATE(xpti(mi),ypti(mi),rloi(mi),rlai(mi),croi(mi),sroi(mi))
506  mix=mi
507  ENDIF
508  CALL gdswzd(grid_in, 0,mi,fill,xpti,ypti,rloi,rlai,nv,croi,sroi)
509  ENDIF
510 
511  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
512  ! SET PARAMETERS
513  nb1=ipopt(1)
514  IF(nb1.EQ.-1) nb1=2
515  IF(iret.EQ.0.AND.nb1.LT.0) iret=32
516  lsw=1
517  IF(ipopt(2).EQ.-2) lsw=2
518  IF(ipopt(1).EQ.-1.OR.ipopt(2).EQ.-1) lsw=0
519  IF(iret.EQ.0.AND.lsw.EQ.1.AND.nb1.GT.15) iret=32
520  mp=ipopt(3+ipopt(1))
521  IF(mp.EQ.-1.OR.mp.EQ.0) mp=50
522  IF(mp.LT.0.OR.mp.GT.100) iret=32
523  pmp=mp*0.01
524  IF(iret.EQ.0) THEN
525  nb2=2*nb1+1
526  rb2=1./nb2
527  nb3=nb2*nb2
528  nb4=nb3
529  IF(lsw.EQ.2) THEN
530  rb2=1./(nb1+1)
531  nb4=(nb1+1)**4
532  ELSEIF(lsw.EQ.1) THEN
533  nb4=ipopt(2)
534  DO ib=1,nb1
535  nb4=nb4+8*ib*ipopt(2+ib)
536  ENDDO
537  ENDIF
538  ELSE
539  nb3=0
540  nb4=1
541  ENDIF
542  DO k=1,km
543  DO n=1,no
544  uo(n,k)=0
545  vo(n,k)=0
546  wo(n,k)=0.
547  ENDDO
548  ENDDO
549  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
550  ! LOOP OVER SAMPLE POINTS IN OUTPUT GRID BOX
551  DO nb=1,nb3
552  ! LOCATE INPUT POINTS AND COMPUTE THEIR WEIGHTS AND ROTATIONS
553  jb=(nb-1)/nb2-nb1
554  ib=nb-(jb+nb1)*nb2-nb1-1
555  lb=max(abs(ib),abs(jb))
556  wb=1
557  IF(ipopt(2).EQ.-2) THEN
558  wb=(nb1+1-abs(ib))*(nb1+1-abs(jb))
559  ELSEIF(ipopt(2).NE.-1) THEN
560  wb=ipopt(2+lb)
561  ENDIF
562  IF(abs(wb).GT.tinyreal) THEN
563  !$OMP PARALLEL DO PRIVATE(N) SCHEDULE(STATIC)
564  DO n=1,no
565  xptb(n)=xpts(n)+ib*rb2
566  yptb(n)=ypts(n)+jb*rb2
567  ENDDO
568  !$OMP END PARALLEL DO
569  if(to_station_points)then
570  CALL gdswzd(grid_in, 1,no,fill,xptb,yptb,rlob,rlab,nv)
571  CALL gdswzd(grid_in,-1,no,fill,xptb,yptb,rlob,rlab,nv)
572  else
573  CALL gdswzd(grid_out2, 1,no,fill,xptb,yptb,rlob,rlab,nv)
574  CALL gdswzd(grid_in,-1,no,fill,xptb,yptb,rlob,rlab,nv)
575  endif
576  IF(iret.EQ.0.AND.nv.EQ.0.AND.lb.EQ.0) iret=2
577  !$OMP PARALLEL DO PRIVATE(N,XI,YI,I1,I2,WI1,WI2,J1,J2,WJ1,WJ2,CM11,CM21,CM12,CM22,SM11,SM21,SM12,SM22) &
578  !$OMP SCHEDULE(STATIC)
579  DO n=1,no
580  xi=xptb(n)
581  yi=yptb(n)
582  IF(abs(xi-fill).GT.tinyreal.AND.abs(yi-fill).GT.tinyreal) THEN
583  i1=int(xi)
584  i2=i1+1
585  wi2=xi-i1
586  wi1=1-wi2
587  j1=int(yi)
588  j2=j1+1
589  wj2=yi-j1
590  wj1=1-wj2
591  n11(n) = grid_in%field_pos(i1,j1)
592  n21(n) = grid_in%field_pos(i2, j1)
593  n12(n) = grid_in%field_pos(i1, j2)
594  n22(n) = grid_in%field_pos(i2, j2)
595  IF(min(n11(n),n21(n),n12(n),n22(n)).GT.0) THEN
596  w11(n)=wi1*wj1
597  w21(n)=wi2*wj1
598  w12(n)=wi1*wj2
599  w22(n)=wi2*wj2
600  CALL movect(rlai(n11(n)),rloi(n11(n)),rlat(n),rlon(n),cm11,sm11)
601  CALL movect(rlai(n21(n)),rloi(n21(n)),rlat(n),rlon(n),cm21,sm21)
602  CALL movect(rlai(n12(n)),rloi(n12(n)),rlat(n),rlon(n),cm12,sm12)
603  CALL movect(rlai(n22(n)),rloi(n22(n)),rlat(n),rlon(n),cm22,sm22)
604  c11(n)=cm11*croi(n11(n))+sm11*sroi(n11(n))
605  s11(n)=sm11*croi(n11(n))-cm11*sroi(n11(n))
606  c21(n)=cm21*croi(n21(n))+sm21*sroi(n21(n))
607  s21(n)=sm21*croi(n21(n))-cm21*sroi(n21(n))
608  c12(n)=cm12*croi(n12(n))+sm12*sroi(n12(n))
609  s12(n)=sm12*croi(n12(n))-cm12*sroi(n12(n))
610  c22(n)=cm22*croi(n22(n))+sm22*sroi(n22(n))
611  s22(n)=sm22*croi(n22(n))-cm22*sroi(n22(n))
612  ELSE
613  n11(n)=0
614  n21(n)=0
615  n12(n)=0
616  n22(n)=0
617  ENDIF
618  ELSE
619  n11(n)=0
620  n21(n)=0
621  n12(n)=0
622  n22(n)=0
623  ENDIF
624  ENDDO
625  !$OMP END PARALLEL DO
626  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
627  ! INTERPOLATE WITH OR WITHOUT BITMAPS
628  ! KM IS OFTEN 1 .. DO NO PUT OMP PARALLEL DO HERE
629  DO k=1,km
630  !$OMP PARALLEL DO PRIVATE(N,U11,U12,U21,U22,UB,V11,V12,V21,V22,VB) SCHEDULE(STATIC)
631  DO n=1,no
632  IF(n11(n).GT.0) THEN
633  IF(ibi(k).EQ.0) THEN
634  u11=c11(n)*ui(n11(n),k)-s11(n)*vi(n11(n),k)
635  v11=s11(n)*ui(n11(n),k)+c11(n)*vi(n11(n),k)
636  u21=c21(n)*ui(n21(n),k)-s21(n)*vi(n21(n),k)
637  v21=s21(n)*ui(n21(n),k)+c21(n)*vi(n21(n),k)
638  u12=c12(n)*ui(n12(n),k)-s12(n)*vi(n12(n),k)
639  v12=s12(n)*ui(n12(n),k)+c12(n)*vi(n12(n),k)
640  u22=c22(n)*ui(n22(n),k)-s22(n)*vi(n22(n),k)
641  v22=s22(n)*ui(n22(n),k)+c22(n)*vi(n22(n),k)
642  ub=w11(n)*u11+w21(n)*u21+w12(n)*u12+w22(n)*u22
643  vb=w11(n)*v11+w21(n)*v21+w12(n)*v12+w22(n)*v22
644  uo(n,k)=uo(n,k)+wb*ub
645  vo(n,k)=vo(n,k)+wb*vb
646  wo(n,k)=wo(n,k)+wb
647  ELSE
648  IF(li(n11(n),k)) THEN
649  u11=c11(n)*ui(n11(n),k)-s11(n)*vi(n11(n),k)
650  v11=s11(n)*ui(n11(n),k)+c11(n)*vi(n11(n),k)
651  uo(n,k)=uo(n,k)+wb*w11(n)*u11
652  vo(n,k)=vo(n,k)+wb*w11(n)*v11
653  wo(n,k)=wo(n,k)+wb*w11(n)
654  ENDIF
655  IF(li(n21(n),k)) THEN
656  u21=c21(n)*ui(n21(n),k)-s21(n)*vi(n21(n),k)
657  v21=s21(n)*ui(n21(n),k)+c21(n)*vi(n21(n),k)
658  uo(n,k)=uo(n,k)+wb*w21(n)*u21
659  vo(n,k)=vo(n,k)+wb*w21(n)*v21
660  wo(n,k)=wo(n,k)+wb*w21(n)
661  ENDIF
662  IF(li(n12(n),k)) THEN
663  u12=c12(n)*ui(n12(n),k)-s12(n)*vi(n12(n),k)
664  v12=s12(n)*ui(n12(n),k)+c12(n)*vi(n12(n),k)
665  uo(n,k)=uo(n,k)+wb*w12(n)*u12
666  vo(n,k)=vo(n,k)+wb*w12(n)*v12
667  wo(n,k)=wo(n,k)+wb*w12(n)
668  ENDIF
669  IF(li(n22(n),k)) THEN
670  u22=c22(n)*ui(n22(n),k)-s22(n)*vi(n22(n),k)
671  v22=s22(n)*ui(n22(n),k)+c22(n)*vi(n22(n),k)
672  uo(n,k)=uo(n,k)+wb*w22(n)*u22
673  vo(n,k)=vo(n,k)+wb*w22(n)*v22
674  wo(n,k)=wo(n,k)+wb*w22(n)
675  ENDIF
676  ENDIF
677  ENDIF
678  ENDDO
679  !$OMP END PARALLEL DO
680  ENDDO
681  ENDIF
682  ENDDO
683  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
684  ! COMPUTE OUTPUT BITMAPS AND FIELDS
685  ! KM is often 1, do not put OMP PARALLEL here
686  DO k=1,km
687  ibo(k)=ibi(k)
688  !$OMP PARALLEL DO PRIVATE(N,UROT,VROT) SCHEDULE(STATIC)
689  DO n=1,no
690  lo(n,k)=wo(n,k).GE.pmp*nb4
691  IF(lo(n,k)) THEN
692  uo(n,k)=uo(n,k)/wo(n,k)
693  vo(n,k)=vo(n,k)/wo(n,k)
694  urot=crot(n)*uo(n,k)-srot(n)*vo(n,k)
695  vrot=srot(n)*uo(n,k)+crot(n)*vo(n,k)
696  uo(n,k)=urot
697  vo(n,k)=vrot
698  ELSE
699  ibo(k)=1
700  uo(n,k)=0.
701  vo(n,k)=0.
702  ENDIF
703  ENDDO
704  !$OMP END PARALLEL DO
705  ENDDO
706 
707  select type(grid_out2)
708  type is(ip_equid_cylind_grid)
709  CALL polfixv(no,mo,km,rlat,rlon,ibo,lo,uo,vo)
710  end select
711  ! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
712  END SUBROUTINE interpolate_budget_vector
713 
714 end module budget_interp_mod
budget_interp_mod::interpolate_budget
Definition: budget_interp_mod.F90:21
gdswzd_mod::gdswzd
Definition: gdswzd_mod.F90:37
movect
subroutine movect(FLAT, FLON, TLAT, TLON, CROT, SROT)
This subprogram provides the rotation parameters to move a vector along a great circle from one posit...
Definition: movect.F90:26
budget_interp_mod
Budget interpolation routines for scalars and vectors.
Definition: budget_interp_mod.F90:10
budget_interp_mod::interpolate_budget_scalar
subroutine interpolate_budget_scalar(IPOPT, grid_in, grid_out, MI, MO, KM, IBI, LI, GI, NO, RLAT, RLON, IBO, LO, GO, IRET)
Performs budget interpolation from any grid to any grid (or to random station points) for scalar fiel...
Definition: budget_interp_mod.F90:97
budget_interp_mod::interpolate_budget_vector
subroutine interpolate_budget_vector(IPOPT, grid_in, grid_out, MI, MO, KM, IBI, LI, UI, VI, NO, RLAT, RLON, CROT, SROT, IBO, LO, UO, VO, IRET)
This subprogram performs budget interpolation from any grid to any grid (or to random station points)...
Definition: budget_interp_mod.F90:426
gdswzd_mod
Driver module for gdswzd routines.
Definition: gdswzd_mod.F90:25
ip_grid_descriptor_mod
Users derived type grid descriptor objects to abstract away the raw GRIB1 and GRIB2 grid definitions.
Definition: ip_grid_descriptor_mod.F90:15
ip_grid_factory_mod
Routines for creating an ip_grid given a Grib descriptor.
Definition: ip_grid_factory_mod.F90:11
ip_grids_mod
Re-export the individual grids.
Definition: ip_grids_mod.F90:7
polfix_mod
Make multiple pole scalar values consistent.
Definition: polfix_mod.F90:7
polfix_mod::polfixs
subroutine, public polfixs(NM, NX, KM, RLAT, IB, LO, GO)
Make multiple pole scalar values consistent.
Definition: polfix_mod.F90:30
polfix_mod::polfixv
subroutine, public polfixv(NM, NX, KM, RLAT, RLON, IB, LO, UO, VO)
Make multiple pole vector values consistent,.
Definition: polfix_mod.F90:125