da_roughem_tl.inc
References to this file elsewhere.
1 subroutine da_roughem_tl(ifreq,gx2,tk,theta,remv,remh, TGL_gx2,TGL_tk,TGL_remv,TGL_remh )
2
3 !----------------------------------------------------------------
4 ! Purpose: Calculates rough-surface emissivity of ocean surface at SSM/I
5 ! frequencies.
6 ! Input : TGL_tk, TGL_gx2
7 ! Output : TGL_remv,TGL_remh, remv,remh
8 !----------------------------------------------------------------
9
10 implicit none
11
12 integer, intent(in) :: ifreq
13 real, intent(in) :: tk, theta, gx2, TGL_tk,TGL_gx2
14 real, intent(out) :: TGL_remv,TGL_remh, remv,remh
15
16 real :: a19v(4),a22v(4),a37v(4),a85v(4)
17 real :: a19h(4),a22h(4),a37h(4),a85h(4)
18 real :: f(4)
19 real :: semv,semh,TGL_semv,TGL_semh,ssw
20 real :: tp,g,x1,x2,x3,x4,dtheta
21 real :: TGL_tp,TGL_g,TGL_x1,TGL_x2,TGL_x3,TGL_x4
22
23 data a19v/ -0.111E+01, 0.713E+00, -0.624E-01, 0.212E-01 /
24 data a19h/ 0.812E+00, -0.215E+00, 0.255E-01, 0.305E-02 /
25 data a22v/ -0.134E+01, 0.911E+00, -0.893E-01, 0.463E-01 /
26 data a22h/ 0.958E+00, -0.350E+00, 0.566E-01, -0.262E-01 /
27 data a37v/ -0.162E+01, 0.110E+01, -0.730E-01, 0.298E-01 /
28 data a37h/ 0.947E+00, -0.320E+00, 0.624E-01, -0.300E-01 /
29 data a85v/ -0.145E+01, 0.808E+00, -0.147E-01, -0.252E-01 /
30 data a85h/ 0.717E+00, -0.702E-01, 0.617E-01, -0.243E-01 /
31
32 data f/ 19.35, 22.235, 37.0, 85.5 /
33
34 if (trace_use) call da_trace_entry("da_roughem_tl")
35
36 tp = tk/t_roughem
37 TGL_tp = TGL_tk/t_roughem
38 dtheta = theta-53.0
39 g = 0.5* gx2
40 TGL_g = 0.5*TGL_gx2
41 x1 = g
42 TGL_x1 = TGL_g
43 x2 = tp*g
44 TGL_x2 = TGL_tp*g + tp*TGL_g
45 x3 = dtheta* g
46 TGL_x3 = dtheta*TGL_g
47 x4 = tp*x3
48 TGL_x4 = TGL_tp*x3 + tp*TGL_x3
49
50 if (ifreq == 1) then
51 remv = x1*a19v(1) + x2*a19v(2) + x3*a19v(3) + x4*a19v(4)
52 TGL_remv = TGL_x1*a19v(1) + TGL_x2*a19v(2) + TGL_x3*a19v(3) + TGL_x4*a19v(4)
53 remh = x1*a19h(1) + x2*a19h(2) + x3*a19h(3) + x4*a19h(4)
54 TGL_remh = TGL_x1*a19h(1) + TGL_x2*a19h(2) + TGL_x3*a19h(3) + TGL_x4*a19h(4)
55 else if (ifreq == 2) then
56 remv = x1*a22v(1) + x2*a22v(2) + x3*a22v(3) + x4*a22v(4)
57 TGL_remv = TGL_x1*a22v(1) + TGL_x2*a22v(2) + TGL_x3*a22v(3) + TGL_x4*a22v(4)
58 remh = x1*a22h(1) + x2*a22h(2) + x3*a22h(3) + x4*a22h(4)
59 TGL_remh = TGL_x1*a22h(1) + TGL_x2*a22h(2) + TGL_x3*a22h(3) + TGL_x4*a22h(4)
60 else if (ifreq == 3) then
61 remv = x1*a37v(1) + x2*a37v(2) + x3*a37v(3) + x4*a37v(4)
62 TGL_remv = TGL_x1*a37v(1) + TGL_x2*a37v(2) + TGL_x3*a37v(3) + TGL_x4*a37v(4)
63 remh = x1*a37h(1) + x2*a37h(2) + x3*a37h(3) + x4*a37h(4)
64 TGL_remh = TGL_x1*a37h(1) + TGL_x2*a37h(2) + TGL_x3*a37h(3) + TGL_x4*a37h(4)
65 else if (ifreq == 4) then
66 remv = x1*a85v(1) + x2*a85v(2) + x3*a85v(3) + x4*a85v(4)
67 TGL_remv = TGL_x1*a85v(1) + TGL_x2*a85v(2) + TGL_x3*a85v(3) + TGL_x4*a85v(4)
68 remh = x1*a85h(1) + x2*a85h(2) + x3*a85h(3) + x4*a85h(4)
69 TGL_remh = TGL_x1*a85h(1) + TGL_x2*a85h(2) + TGL_x3*a85h(3) + TGL_x4*a85h(4)
70 end if
71
72 ssw=36.5
73 call da_spemiss_tl(f(ifreq),tk,theta,ssw,semv,semh, TGL_tk,TGL_semv,TGL_semh)
74
75 TGL_remv = TGL_remv + TGL_semv
76 remv = remv + semv
77 TGL_remh = TGL_remh + TGL_semh
78 remh = remh + semh
79
80 if (trace_use) call da_trace_exit("da_roughem_tl")
81
82 end subroutine da_roughem_tl
83
84