da_transform_through_rf.inc
References to this file elsewhere.
1 subroutine da_transform_through_rf(grid,mz, rf_alpha, val,field)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: Control variable transform through the recursive filter.
5 !
6 ! Method: 1) Transform to nondimensional v_hat space.
7 ! 2) Perform recursive filter in non-dimensional space (ds = 1).
8 ! 3) Scale recursive filter output.
9 ! 4) Transform filtered field to dimensional space.
10 ! 5) Optionally scale by background error.
11 !---------------------------------------------------------------------------
12
13 implicit none
14
15 type(domain), intent(inout) :: grid
16 integer, intent(in) :: mz ! Vertical truncation.
17 real, intent(in) :: rf_alpha(mz) ! RF scale parameter.
18 real, intent(in) :: val(jds:jde,mz) ! Error standard deviation.
19 real, intent(inout) :: field(ims:ime,jms:jme,kms:kme) ! Field to be transformed.
20
21 integer :: rf_passes_over_two ! rf_passes / 2
22 integer :: i, j, m, n, pass ! Loop counters.
23 real :: p_x(ims:ime,jms:jme)! sqrt(Grid box area).
24 real :: val_j(grid%xp%jtsy:grid%xp%jtey)
25
26 !----------------------------------------------------------------------
27 ! [1.0]: Initialise:
28 !----------------------------------------------------------------------
29
30 if (trace_use_dull) call da_trace_entry("da_transform_through_rf")
31
32 rf_passes_over_two = rf_passes / 2
33
34 ! [1.1] Define inner product (square root of grid box area):
35 p_x(its:ite,jts:jte) = sqrt(grid%xb%grid_box_area(its:ite,jts:jte))
36
37 grid%xp%v1z(:,:,:) = 0.0
38
39 ! [1.2] Transform to nondimensional v_hat space:
40
41 do m = 1, mz
42 grid%xp%v1z(its:ite,jts:jte,m) = field(its:ite,jts:jte,m) / p_x(its:ite,jts:jte)
43 end do
44
45 !-------------------------------------------------------------------------
46 ! [2.0]: Perform 1D recursive filter in x-direction:
47 !-------------------------------------------------------------------------
48
49 ! [2.1] Apply (i',j',k -> i,j',k') (grid%xp%v1z -> grid%xp%v1x)
50 ! convert from vertical column to x-stripe
51
52 call da_transpose_z2x (grid)
53
54 ! [2.2] Apply 1D filter in x direction:
55
56 n = grid%xp%itex - grid%xp%itsx + 1
57 do m = grid%xp%ktsx, min(grid%xp%ktex,mz)
58 do j = grid%xp%jtsx, grid%xp%jtex
59 do pass = 1, rf_passes_over_two
60 call da_recursive_filter_1d(pass, rf_alpha(m), grid%xp%v1x(grid%xp%itsx:grid%xp%itex,j,m), n)
61 end do
62 end do
63 end do
64
65 !-------------------------------------------------------------------------
66 ! [3.0]: Perform 1D recursive filter in y-direction:
67 !-------------------------------------------------------------------------
68
69 ! [3.1] Apply (i, j' ,k' -> i', j ,k') (grid%xp%v1x -> grid%xp%v1y)
70 ! convert from vertical column to y-stripe
71
72 call da_transpose_x2y (grid)
73
74 ! [3.2] Apply 1D filter in y direction:
75
76 n = grid%xp%jtey - grid%xp%jtsy + 1
77 do m = grid%xp%ktsy, min(grid%xp%ktey,mz)
78 do i = grid%xp%itsy, grid%xp%itey
79 val_j(grid%xp%jtsy:grid%xp%jtey) = grid%xp%v1y(i,grid%xp%jtsy:grid%xp%jtey,m)
80 do pass = 1, rf_passes_over_two
81 call da_recursive_filter_1d(pass, rf_alpha(m), val_j, n)
82 end do
83 grid%xp%v1y(i,grid%xp%jtsy:grid%xp%jtey,m) = val_j(grid%xp%jtsy:grid%xp%jtey)
84 end do
85 end do
86
87 !-------------------------------------------------------------------------
88 ! [4.0]: Perform 1D recursive filter in y-direction:
89 !-------------------------------------------------------------------------
90
91 ! [4.1] Apply (i',j,k' -> i',j',k) (grid%xp%v1y -> grid%xp%v1z)
92 ! convert from y-stripe to vertical column.
93
94 call da_transpose_y2z (grid)
95
96 ! [4.2] Transform filtered field to dimensional space:
97
98 do m = 1, mz
99 field(its:ite,jts:jte,m) = grid%xp%v1z(its:ite,jts:jte,m) * p_x(its:ite,jts:jte)
100 end do
101
102 ! [4.3] Optionally scale by background error:
103
104 ! be_s % val = Gridpoint standard deviation - only required for
105 ! vert_corr = vert_corr_1 as scaling is performed in vertical transform
106 ! for vert_corr = vert_corr_2:
107
108 if (vert_corr == vert_corr_1) then
109 do m = 1, mz
110 do i = its, ite
111 field(i,jts:jte,m) = field(i,jts:jte,m) * val(jts:jte,m)
112 end do
113 end do
114 end if
115
116 if (trace_use_dull) call da_trace_exit("da_transform_through_rf")
117
118 end subroutine da_transform_through_rf
119
120