da_roughem_adj.inc
References to this file elsewhere.
1 subroutine da_roughem_adj(ifreq,gx2,tk,theta,remv,remh, &
2 ADJ_gx2,ADJ_tk,ADJ_remv,ADJ_remh )
3
4 !-----------------------------------------------------------------------
5 ! Purpose: TBD
6 !-----------------------------------------------------------------------
7
8 implicit none
9
10 !-----------------------------------------------------------------
11 ! Input :: ADJ_tk, ADJ_gx2
12 ! Output :: ADJ_remv,ADJ_remh, remv,remh
13 !
14 ! Calculates rough-surface emissivity of ocean surface at SSM/I
15 ! frequencies.
16
17 integer, intent(in) :: ifreq
18 real , intent(in) :: tk, theta, gx2
19 real , intent(in) :: ADJ_remv,ADJ_remh
20 real , intent(inout) :: ADJ_tk,ADJ_gx2
21 real, intent(out) :: remv,remh
22
23 real :: ssw
24
25 real a19v(4),a22v(4),a37v(4),a85v(4)
26 real a19h(4),a22h(4),a37h(4),a85h(4)
27 real f(4)
28 real semv,semh,ADJ_semv,ADJ_semh,remv_save,remh_save
29 real :: tp,g,x1,x2,x3,x4,dtheta
30 real :: ADJ_tp,ADJ_g,ADJ_x1,ADJ_x2,ADJ_x3,ADJ_x4
31
32 data a19v/ -0.111E+01, 0.713E+00, -0.624E-01, 0.212E-01 /
33 data a19h/ 0.812E+00, -0.215E+00, 0.255E-01, 0.305E-02 /
34 data a22v/ -0.134E+01, 0.911E+00, -0.893E-01, 0.463E-01 /
35 data a22h/ 0.958E+00, -0.350E+00, 0.566E-01, -0.262E-01 /
36 data a37v/ -0.162E+01, 0.110E+01, -0.730E-01, 0.298E-01 /
37 data a37h/ 0.947E+00, -0.320E+00, 0.624E-01, -0.300E-01 /
38 data a85v/ -0.145E+01, 0.808E+00, -0.147E-01, -0.252E-01 /
39 data a85h/ 0.717E+00, -0.702E-01, 0.617E-01, -0.243E-01 /
40
41 data f/ 19.35, 22.235, 37.0, 85.5 /
42
43 semv=0.0
44 semh=0.0
45 ADJ_semv=0.0
46 ADJ_semh=0.0
47 remv_save=0.0
48 remh_save=0.0
49 tp=0.0
50 g=0.0
51 x1=0.0
52 x2=0.0
53 x3=0.0
54 x4=0.0
55 dtheta=0.0
56 ADJ_tp=0.0
57 ADJ_g=0.0
58 ADJ_x1=0.0
59 ADJ_x2=0.0
60 ADJ_x3=0.0
61 ADJ_x4=0.0
62
63 tp = tk/t_roughem
64 dtheta = theta-53.0
65 g = 0.5* gx2
66 x1 = g
67 x2 = tp*g
68 x3 = dtheta* g
69 x4 = tp*x3
70
71 if (ifreq .eq. 1) then
72 remv = x1*a19v(1) + x2*a19v(2) + x3*a19v(3) + x4*a19v(4)
73 remh = x1*a19h(1) + x2*a19h(2) + x3*a19h(3) + x4*a19h(4)
74 else if (ifreq .eq. 2) then
75 remv = x1*a22v(1) + x2*a22v(2) + x3*a22v(3) + x4*a22v(4)
76 remh = x1*a22h(1) + x2*a22h(2) + x3*a22h(3) + x4*a22h(4)
77 else if (ifreq .eq. 3) then
78 remv = x1*a37v(1) + x2*a37v(2) + x3*a37v(3) + x4*a37v(4)
79 remh = x1*a37h(1) + x2*a37h(2) + x3*a37h(3) + x4*a37h(4)
80 else if (ifreq .eq. 4) then
81 remv = x1*a85v(1) + x2*a85v(2) + x3*a85v(3) + x4*a85v(4)
82 remh = x1*a85h(1) + x2*a85h(2) + x3*a85h(3) + x4*a85h(4)
83 end if
84
85 ssw=36.5
86 call spemiss(f(ifreq),tk,theta,ssw,semv,semh)
87
88 remv_save = remv
89 remh_save = remh
90
91 ! start
92
93 ADJ_semh = ADJ_remh
94
95 ADJ_semv = ADJ_remv
96
97 remv = remv_save
98 remh = remh_save
99
100 call da_spemiss_adj(f(ifreq),tk,theta,ssw,semv,semh, &
101 ADJ_tk,ADJ_semv,ADJ_semh )
102
103
104 if (ifreq .eq. 1) then
105 ADJ_x1 = ADJ_remh*a19h(1)
106 ADJ_x2 = ADJ_remh*a19h(2)
107 ADJ_x3 = ADJ_remh*a19h(3)
108 ADJ_x4 = ADJ_remh*a19h(4)
109
110 ADJ_x1 = ADJ_remv*a19v(1) + ADJ_x1
111 ADJ_x2 = ADJ_remv*a19v(2) + ADJ_x2
112 ADJ_x3 = ADJ_remv*a19v(3) + ADJ_x3
113 ADJ_x4 = ADJ_remv*a19v(4) + ADJ_x4
114 else if (ifreq .eq. 2) then
115 ADJ_x1 = ADJ_remh*a22h(1)
116 ADJ_x2 = ADJ_remh*a22h(2)
117 ADJ_x3 = ADJ_remh*a22h(3)
118 ADJ_x4 = ADJ_remh*a22h(4)
119
120 ADJ_x1 = ADJ_remv*a22v(1) + ADJ_x1
121 ADJ_x2 = ADJ_remv*a22v(2) + ADJ_x2
122 ADJ_x3 = ADJ_remv*a22v(3) + ADJ_x3
123 ADJ_x4 = ADJ_remv*a22v(4) + ADJ_x4
124 else if (ifreq .eq. 3) then
125 ADJ_x1 = ADJ_remh*a37h(1)
126 ADJ_x2 = ADJ_remh*a37h(2)
127 ADJ_x3 = ADJ_remh*a37h(3)
128 ADJ_x4 = ADJ_remh*a37h(4)
129
130 ADJ_x1 = ADJ_remv*a37v(1) + ADJ_x1
131 ADJ_x2 = ADJ_remv*a37v(2) + ADJ_x2
132 ADJ_x3 = ADJ_remv*a37v(3) + ADJ_x3
133 ADJ_x4 = ADJ_remv*a37v(4) + ADJ_x4
134 else if (ifreq .eq. 4) then
135 ADJ_x1 = ADJ_remh*a85h(1)
136 ADJ_x2 = ADJ_remh*a85h(2)
137 ADJ_x3 = ADJ_remh*a85h(3)
138 ADJ_x4 = ADJ_remh*a85h(4)
139
140 ADJ_x1 = ADJ_remv*a85v(1) + ADJ_x1
141 ADJ_x2 = ADJ_remv*a85v(2) + ADJ_x2
142 ADJ_x3 = ADJ_remv*a85v(3) + ADJ_x3
143 ADJ_x4 = ADJ_remv*a85v(4) + ADJ_x4
144
145 end if
146
147 ADJ_tp = ADJ_x4*x3
148 ADJ_x3 = tp*ADJ_x4 + ADJ_x3
149 ADJ_g = dtheta*ADJ_x3
150 ADJ_tp = ADJ_x2*g + ADJ_tp
151 ADJ_g = tp*ADJ_x2 + ADJ_g
152 ADJ_g = ADJ_x1 + ADJ_g
153 ADJ_gx2= 0.5*ADJ_g + ADJ_gx2
154 ADJ_tk = ADJ_tp/t_roughem + ADJ_tk
155
156 end subroutine da_roughem_adj
157
158