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
99 LOGICAL*1,
INTENT(IN ) :: LI(MI,KM)
100 LOGICAL*1,
INTENT( OUT) :: LO(MO,KM)
102 REAL,
INTENT(IN ) :: GI(MI,KM)
103 REAL,
INTENT(INOUT) :: RLAT(MO), RLON(MO)
104 REAL,
INTENT( OUT) :: GO(MO,KM)
106 REAL,
PARAMETER :: FILL=-9999.
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)
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)
122 logical :: to_station_points
124 class(ip_grid),
allocatable :: grid_out2
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)
137 CALL gdswzd(grid_in,-1,no,fill,xpts,ypts,rlon,rlat,nv)
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)
148 IF(ipopt(1).GT.16) iret=32
149 mspiral=max(ipopt(20),1)
152 IF(iret.EQ.0.AND.nb1.LT.0) iret=32
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
158 IF(mp.EQ.-1.OR.mp.EQ.0) mp=50
159 IF(mp.LT.0.OR.mp.GT.100) iret=32
169 ELSEIF(lsw.EQ.1)
THEN
172 nb4=nb4+8*ib*ipopt(2+ib)
190 ib=nb-(jb+nb1)*nb2-nb1-1
191 lb=max(abs(ib),abs(jb))
194 wb=(nb1+1-abs(ib))*(nb1+1-abs(jb))
195 ELSEIF(lsw.EQ.1)
THEN
201 xptb(n)=xpts(n)+ib*rb2
202 yptb(n)=ypts(n)+jb*rb2
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)
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)
212 IF(iret.EQ.0.AND.nv.EQ.0.AND.lb.EQ.0) iret=2
217 IF(xi.NE.fill.AND.yi.NE.fill)
THEN
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
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
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)
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)
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)
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)
289 lo(n,k)=wo(n,k).GE.pmp*nb4
291 go(n,k)=go(n,k)/wo(n,k)
292 ELSEIF (mspiral.GT.1)
THEN
295 CALL gdswzd(grid_in,-1,1,fill,xxx,yyy,lon,lat,nv)
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))
306 SELECT CASE(mod(kxs,4))
308 ix=i1-ixs*(kxs/4-kxt)
312 jx=j1-jxs*(kxs/4-kxt)
314 ix=i1+ixs*(1+kxs/4-kxt)
318 jx=j1+jxs*(kxs/4-kxt)
320 nx=grid_in%field_pos(ix, jx)
322 IF(li(nx,k).OR.ibi(k).EQ.0)
THEN
343 select type(grid_out2)
344 type is(ip_equid_cylind_grid)
345 CALL polfixs(no,mo,km,rlat,ibo,lo,go)
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)
428 LOGICAL*1,
INTENT(IN ) :: LI(MI,KM)
429 LOGICAL*1,
INTENT( OUT) :: LO(MO,KM)
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)
436 REAL,
PARAMETER :: FILL=-9999.
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)
444 REAL :: CM11,SM11,CM12,SM12
445 REAL :: CM21,SM21,CM22,SM22
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)
457 logical :: to_station_points
459 class(ip_grid),
allocatable :: grid_out2
462 INTEGER,
SAVE :: MIX=-1
463 REAL,
ALLOCATABLE,
SAVE :: CROI(:),SROI(:)
464 REAL,
ALLOCATABLE,
SAVE :: XPTI(:),YPTI(:),RLOI(:),RLAI(:)
466 class(ip_grid),
allocatable,
save :: prev_grid_in
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)
478 CALL gdswzd(grid_in,-1,no,fill,xpts,ypts,rlon,rlat,nv,crot,srot)
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)
486 if (.not.
allocated(prev_grid_in))
then
487 allocate(prev_grid_in, source = grid_in)
491 same_grid = grid_in == prev_grid_in
493 if (.not. same_grid)
then
494 deallocate(prev_grid_in)
495 allocate(prev_grid_in, source = grid_in)
499 IF(.NOT.same_grid)
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))
505 CALL gdswzd(grid_in, 0,mi,fill,xpti,ypti,rloi,rlai,nv,croi,sroi)
512 IF(iret.EQ.0.AND.nb1.LT.0) iret=32
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
518 IF(mp.EQ.-1.OR.mp.EQ.0) mp=50
519 IF(mp.LT.0.OR.mp.GT.100) iret=32
529 ELSEIF(lsw.EQ.1)
THEN
532 nb4=nb4+8*ib*ipopt(2+ib)
551 ib=nb-(jb+nb1)*nb2-nb1-1
552 lb=max(abs(ib),abs(jb))
554 IF(ipopt(2).EQ.-2)
THEN
555 wb=(nb1+1-abs(ib))*(nb1+1-abs(jb))
556 ELSEIF(ipopt(2).NE.-1)
THEN
562 xptb(n)=xpts(n)+ib*rb2
563 yptb(n)=ypts(n)+jb*rb2
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)
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)
573 IF(iret.EQ.0.AND.nv.EQ.0.AND.lb.EQ.0) iret=2
579 IF(xi.NE.fill.AND.yi.NE.fill)
THEN
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
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))
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
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)
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)
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)
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)
687 lo(n,k)=wo(n,k).GE.pmp*nb4
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)
704 select type(grid_out2)
705 type is(ip_equid_cylind_grid)
706 CALL polfixv(no,mo,km,rlat,rlon,ibo,lo,uo,vo)
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...
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.
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.
Make multiple pole scalar values consistent.
subroutine, public polfixs(NM, NX, KM, RLAT, IB, LO, GO)
Make multiple pole scalar values consistent.
subroutine, public polfixv(NM, NX, KM, RLAT, RLON, IB, LO, UO, VO)
Make multiple pole vector values consistent,.