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 (iv_type), intent (inout) :: iv
10
11 integer :: i,j, n,t, error
12 integer, allocatable :: ifgat(:),landsea_mask(:)
13 integer, allocatable :: loc_i(:,:), loc_j(:,:),loc_k(:,:)
14 real, allocatable :: loc_dx(:,:),loc_dy(:,:),loc_dz(:,:)
15 real, allocatable :: loc_dxm(:,:),loc_dym(:,:),loc_dzm(:,:)
16 character (len = 40), allocatable :: name(:)
17 character (len = 12), allocatable :: platform(:)
18 character (len = 5), allocatable :: id(:)
19 character (len = 19), allocatable :: date_char(:)
20 integer, allocatable :: levels(:)
21 real, allocatable :: lat(:,:)
22 real, allocatable :: lon(:,:)
23 real, allocatable :: elv(:)
24 real, allocatable :: pstar(:)
25 real, allocatable :: scanpos(:), satzen(:), satazi(:), solzen(:)
26 real, allocatable :: tb_inv(:,:)
27
28 if (trace_use) call da_trace_entry("da_sort_rad")
29
30 if (num_fgat_time == 1) then
31 do i=1,rtminit_nsensor
32 iv%instid(i)%info%plocal(:) = 0
33 iv%instid(i)%info%plocal(num_fgat_time) = iv%instid(i)%num_rad
34 end do
35 if (trace_use) call da_trace_exit("da_sort_rad")
36 return
37 end if
38
39 do i=1,rtminit_nsensor
40 iv%instid(i)%info%plocal(:) = 0
41 if (iv%instid(i)%num_rad < 1) cycle
42
43 allocate (ifgat (iv%instid(i)%num_rad))
44 allocate (landsea_mask(iv%instid(i)%num_rad))
45 allocate (loc_i (max_ob_levels, iv%instid(i)%num_rad))
46 allocate (loc_j (max_ob_levels, iv%instid(i)%num_rad))
47 allocate (loc_k (iv%instid(i)%nchan,iv%instid(i)%num_rad))
48 allocate (loc_dx (max_ob_levels, iv%instid(i)%num_rad))
49 allocate (loc_dy (max_ob_levels, iv%instid(i)%num_rad))
50 allocate (loc_dz (iv%instid(i)%nchan,iv%instid(i)%num_rad))
51 allocate (loc_dxm (max_ob_levels, iv%instid(i)%num_rad))
52 allocate (loc_dym (max_ob_levels, iv%instid(i)%num_rad))
53 allocate (loc_dzm (iv%instid(i)%nchan,iv%instid(i)%num_rad))
54 allocate (name (iv%instid(i)%num_rad))
55 allocate (platform (iv%instid(i)%num_rad))
56 allocate (id (iv%instid(i)%num_rad))
57 allocate (date_char (iv%instid(i)%num_rad))
58 allocate (levels (iv%instid(i)%num_rad))
59 allocate (lat (iv%instid(i)%nchan,iv%instid(i)%num_rad))
60 allocate (lon (iv%instid(i)%nchan,iv%instid(i)%num_rad))
61 allocate (elv (iv%instid(i)%num_rad))
62 allocate (pstar (iv%instid(i)%num_rad))
63 allocate (scanpos (iv%instid(i)%num_rad))
64 allocate (satzen (iv%instid(i)%num_rad))
65 allocate (satazi (iv%instid(i)%num_rad))
66 allocate (solzen (iv%instid(i)%num_rad))
67 allocate (tb_inv (iv%instid(i)%nchan,iv%instid(i)%num_rad))
68
69 j = 0
70 do t = 1,num_fgat_time
71 do n = 1,iv%instid(i)%num_rad
72 if (iv%instid(i)%ifgat(n) /= t) cycle
73 j = j + 1
74 ifgat(j) = iv%instid(i)%ifgat(n)
75 name(j) = iv%instid(i)%info%name(n)
76 platform(j) = iv%instid(i)%info%platform(n)
77 id(j) = iv%instid(i)%info%id(n)
78 date_char(j) = iv%instid(i)%info%date_char(n)
79 levels(j) = iv%instid(i)%info%levels(n)
80 elv(j) = iv%instid(i)%info%elv(n)
81 pstar(j) = iv%instid(i)%info%pstar(n)
82
83 lat (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%lat(:,n)
84 lon (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%lon(:,n)
85 loc_i (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%i(:,n)
86 loc_j (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%j(:,n)
87 loc_k (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%k(:,n)
88 loc_dx (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dx(:,n)
89 loc_dy (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dy(:,n)
90 loc_dxm(1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dxm(:,n)
91 loc_dym(1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dym(:,n)
92 loc_dz (1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dz(:,n)
93 loc_dzm(1:iv%instid(i)%info%levels(n),j) = iv%instid(i)%info%dzm(:,n)
94
95 landsea_mask(j) = iv%instid(i)%landsea_mask(n)
96 scanpos(j) = iv%instid(i)%scanpos(n)
97 satzen(j) = iv%instid(i)%satzen(n)
98 satazi(j) = iv%instid(i)%satazi(n)
99 solzen(j) = iv%instid(i)%solzen(n)
100
101 tb_inv(1:iv%instid(i)%nchan,j) = iv%instid(i)%tb_inv(1:iv%instid(i)%nchan,n)
102 end do
103 iv%instid(i)%info%plocal(t) = j
104 write (0,*) __FILE__,__LINE__,"i,t,iv%instid(i)%info%plocal(t)",i,t,iv%instid(i)%info%plocal(t)
105 end do
106
107 write(unit=stdout,fmt='(a,2x,a,2x,10i7)') &
108 ' FGAT: ',iv%instid(i)%rttovid_string, iv%info(i)%plocal(1:num_fgat_time)
109
110 do n = 1,iv%instid(i)%num_rad
111 iv%instid(i)%ifgat(n) = ifgat(n)
112
113 iv%instid(i)%info%name(n) = name(n)
114 iv%instid(i)%info%platform(n) = platform(n)
115 iv%instid(i)%info%id(n) = id(n)
116 iv%instid(i)%info%date_char(n) = date_char(n)
117 iv%instid(i)%info%levels(n) = levels(n)
118 iv%instid(i)%info%lat(:,n) = lat(:,n)
119 iv%instid(i)%info%lon(:,n) = lon(:,n)
120 iv%instid(i)%info%elv(n) = elv(n)
121 iv%instid(i)%info%pstar(n) = pstar(n)
122
123 iv%instid(i)%info%i(:,j) = loc_i(:,n)
124 iv%instid(i)%info%j(:,j) = loc_j(:,n)
125 iv%instid(i)%info%k(:,j) = loc_k(:,n)
126 iv%instid(i)%info%dx(:,j) = loc_dx(:,n)
127 iv%instid(i)%info%dy(:,j) = loc_dy(:,n)
128 iv%instid(i)%info%dz(:,j) = loc_dz(:,n)
129 iv%instid(i)%info%dxm(:,j) = loc_dxm(:,n)
130 iv%instid(i)%info%dym(:,j) = loc_dym(:,n)
131 iv%instid(i)%info%dzm(:,j) = loc_dzm(:,n)
132 iv%instid(i)%landsea_mask(n) = landsea_mask(n)
133 iv%instid(i)%scanpos(n) = scanpos(n)
134 iv%instid(i)%satzen(n) = satzen(n)
135 iv%instid(i)%satazi(n) = satazi(n)
136 iv%instid(i)%solzen(n) = solzen(n)
137
138 iv%instid(i)%tb_inv(1:iv%instid(i)%nchan,n) = tb_inv(1:iv%instid(i)%nchan,n)
139 end do
140
141 deallocate (ifgat)
142 deallocate (landsea_mask)
143 deallocate (name)
144 deallocate (platform)
145 deallocate (id)
146 deallocate (date_char)
147 deallocate (levels)
148 deallocate (lat)
149 deallocate (lon)
150 deallocate (elv)
151 deallocate (pstar)
152 deallocate (loc_i)
153 deallocate (loc_j)
154 deallocate (loc_k)
155 deallocate (loc_dx)
156 deallocate (loc_dy)
157 deallocate (loc_dz)
158 deallocate (loc_dxm)
159 deallocate (loc_dym)
160 deallocate (loc_dzm)
161 deallocate (scanpos)
162 deallocate (satzen)
163 deallocate (satazi)
164 deallocate (solzen)
165 deallocate (tb_inv)
166 end do
167
168 if (trace_use) call da_trace_exit("da_sort_rad")
169
170 end subroutine da_sort_rad
171
172