da_transform_vptovv.inc
References to this file elsewhere.
1 subroutine da_transform_vptovv(evec, eval, vertical_wgt, &
2 vp, vv, mz, &
3 kds,kde, &
4 ims,ime, jms,jme, kms,kme, &
5 its,ite, jts,jte, kts,kte)
6
7 !---------------------------------------------------------------------------
8 ! Purpose: Transform from fields on vertical levels to fields on vertical
9 ! EOFS.
10 !
11 ! Method: Perform vv(i,j,m) = L^{-1/2} E^{T} vp(i,j,k) transform.
12 !---------------------------------------------------------------------------
13
14 implicit none
15
16 integer, intent(in) :: mz ! # vertical modes.
17 integer, intent(in) :: kds,kde ! domain dims.
18 integer, intent(in) :: ims,ime, jms,jme, kms,kme ! memory dims.
19 integer, intent(in) :: its,ite, jts,jte, kts,kte ! tile dims
20 real, intent(in) :: evec(jts:jte,kds:kde,1:mz) ! Eigenvectors.
21 real, intent(in) :: eval(jts:jte,1:mz) ! Eigenvalues.
22 real, intent(in) :: vertical_wgt(ims:ime,jms:jme,kms:kme) ! Weighting.
23 real, intent(inout) :: vp(ims:ime,jms:jme,kms:kme)! CV in level space.
24 real, intent(out) :: vv(ims:ime,jms:jme,1:mz) ! CV in EOF space.
25
26 integer :: i, j, m ! Loop counters.
27 real :: ETVp ! E(k,m)^{T}*vp(i,j,k)
28
29 ! Do not add trace, as routine used by gen_be
30
31 !-------------------------------------------------------------------
32 ! [1.0] Apply inner-product weighting if vertical_ip > 0:
33 !-------------------------------------------------------------------
34
35 if (vertical_ip /= 0) then
36 vp(its:ite,jts:jte,kts:kte) = vp(its:ite,jts:jte,kts:kte) * &
37 vertical_wgt(its:ite,jts:jte,kts:kte)
38 end if
39
40 !-------------------------------------------------------------------
41 ! [2.0] Perform vv(i,j,m) = L^{-1/2} E^T vp(i,j,k) transform:
42 !-------------------------------------------------------------------
43
44 do m = 1, mz
45 do j = jts, jte
46 do i = its, ite
47 ETVp = sum(evec(j,kts:kte,m) * vp(i,j,kts:kte))
48 vv(i,j,m) = ETVp / eval(j,m)
49 end do
50 end do
51 end do
52
53 end subroutine da_transform_vptovv
54
55