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