da_sort_rad.inc
References to this file elsewhere.
1 subroutine da_sort_rad (iv)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: sorting radiance innovation to FGAT time bins
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (ob_type), intent (inout) :: iv
10
11 integer :: i,j, n,t, error
12 type(info_type), pointer :: info(:)
13 integer, allocatable :: ifgat(:),landsea_mask(:)
14 integer, allocatable :: loc_i(:), loc_j(:),loc_k(:,:)
15 real, allocatable :: loc_dx(:),loc_dy(:),loc_dz(:,:)
16 real, allocatable :: loc_dxm(:),loc_dym(:),loc_dzm(:,:)
17 real, allocatable :: scanpos(:), satzen(:), satazi(:), solzen(:)
18 real, allocatable :: tb_inv(:,:)
19
20 if (trace_use) call da_trace_entry("da_sort_rad")
21
22 if (num_fgat_time == 1) then
23 do i=1,rtminit_nsensor
24 iv%ob_numb(:)%radiance(i) = 0
25 iv%ob_numb(num_fgat_time)%radiance(i) = iv%instid(i)%num_rad
26 end do
27 if (trace_use) call da_trace_exit("da_sort_rad")
28 return
29 end if
30
31 do i=1,rtminit_nsensor
32 iv%ob_numb(:)%radiance(i) = 0
33 if (iv%instid(i)%num_rad < 1) cycle
34
35 allocate (info (iv%instid(i)%num_rad))
36 allocate (ifgat (iv%instid(i)%num_rad))
37 allocate (landsea_mask(iv%instid(i)%num_rad))
38 allocate (loc_i (iv%instid(i)%num_rad))
39 allocate (loc_j (iv%instid(i)%num_rad))
40 allocate (loc_k (iv%instid(i)%nchan,iv%instid(i)%num_rad))
41 allocate (loc_dx (iv%instid(i)%num_rad))
42 allocate (loc_dy (iv%instid(i)%num_rad))
43 allocate (loc_dz (iv%instid(i)%nchan,iv%instid(i)%num_rad))
44 allocate (loc_dxm (iv%instid(i)%num_rad))
45 allocate (loc_dym (iv%instid(i)%num_rad))
46 allocate (loc_dzm (iv%instid(i)%nchan,iv%instid(i)%num_rad))
47 allocate (scanpos (iv%instid(i)%num_rad))
48 allocate (satzen (iv%instid(i)%num_rad))
49 allocate (satazi (iv%instid(i)%num_rad))
50 allocate (solzen (iv%instid(i)%num_rad))
51 allocate (tb_inv (iv%instid(i)%nchan,iv%instid(i)%num_rad))
52
53 j = 0
54 do t = 1,num_fgat_time
55 do n = 1,iv%instid(i)%num_rad
56 if (iv%instid(i)%ifgat(n) /= t) cycle
57 j = j + 1
58 ifgat(j) = iv%instid(i)%ifgat(n)
59 info(j) = iv%instid(i)%info(n)
60 loc_i(j) = iv%instid(i)%loc_i(n)
61 loc_j(j) = iv%instid(i)%loc_j(n)
62 loc_k(:,j) = iv%instid(i)%loc_k(:,n)
63 loc_dx(j) = iv%instid(i)%loc_dx(n)
64 loc_dy(j) = iv%instid(i)%loc_dy(n)
65 loc_dz(:,j) = iv%instid(i)%loc_dz(:,n)
66 loc_dxm(j) = iv%instid(i)%loc_dxm(n)
67 loc_dym(j) = iv%instid(i)%loc_dym(n)
68 loc_dzm(:,j) = iv%instid(i)%loc_dzm(:,n)
69 landsea_mask(j) = iv%instid(i)%landsea_mask(n)
70 scanpos(j) = iv%instid(i)%scanpos(n)
71 satzen(j) = iv%instid(i)%satzen(n)
72 satazi(j) = iv%instid(i)%satazi(n)
73 solzen(j) = iv%instid(i)%solzen(n)
74
75 tb_inv(1:iv%instid(i)%nchan,j) = iv%instid(i)%tb_inv(1:iv%instid(i)%nchan,n)
76 end do
77 iv%ob_numb(t)%radiance(i) = j
78 end do
79
80 write(unit=stdout,fmt='(a,2x,a,2x,10i7)') &
81 ' FGAT: ',iv%instid(i)%rttovid_string, iv%ob_numb(1:num_fgat_time)%radiance(i)
82
83 do n = 1,iv%instid(i)%num_rad
84 iv%instid(i)%ifgat(n) = ifgat(n)
85 iv%instid(i)%info(n) = info(n)
86 iv%instid(i)%loc_i(j) = loc_i(n)
87 iv%instid(i)%loc_j(j) = loc_j(n)
88 iv%instid(i)%loc_k(:,j) = loc_k(:,n)
89 iv%instid(i)%loc_dx(j) = loc_dx(n)
90 iv%instid(i)%loc_dy(j) = loc_dy(n)
91 iv%instid(i)%loc_dz(:,j) = loc_dz(:,n)
92 iv%instid(i)%loc_dxm(j) = loc_dxm(n)
93 iv%instid(i)%loc_dym(j) = loc_dym(n)
94 iv%instid(i)%loc_dzm(:,j) = loc_dzm(:,n)
95 iv%instid(i)%landsea_mask(n) = landsea_mask(n)
96 iv%instid(i)%scanpos(n) = scanpos(n)
97 iv%instid(i)%satzen(n) = satzen(n)
98 iv%instid(i)%satazi(n) = satazi(n)
99 iv%instid(i)%solzen(n) = solzen(n)
100
101 iv%instid(i)%tb_inv(1:iv%instid(i)%nchan,n) = &
102 tb_inv(1:iv%instid(i)%nchan,n)
103 end do
104
105 deallocate (info)
106 deallocate (ifgat)
107 deallocate (landsea_mask)
108 deallocate (loc_i)
109 deallocate (loc_j)
110 deallocate (loc_k)
111 deallocate (loc_dx)
112 deallocate (loc_dy)
113 deallocate (loc_dz)
114 deallocate (loc_dxm)
115 deallocate (loc_dym)
116 deallocate (loc_dzm)
117 deallocate (scanpos)
118 deallocate (satzen)
119 deallocate (satazi)
120 deallocate (solzen)
121 deallocate (tb_inv)
122 end do
123
124 if (trace_use) call da_trace_exit("da_sort_rad")
125
126 end subroutine da_sort_rad
127
128