da_qc_rad.inc
References to this file elsewhere.
1 subroutine da_qc_rad (ob, iv)
2
3 !---------------------------------------------------------------------------
4 ! Purpose: perform quality control for radiance data.
5 !
6 ! METHOD: seperated QC for each sensor
7 !---------------------------------------------------------------------------
8
9 implicit none
10
11 type (y_type), intent(in) :: ob ! Observation structure.
12 type (ob_type), intent(inout) :: iv ! O-B structure.
13
14 integer :: i, nchan,p,j
15 logical :: amsua, amsub, hirs, msu,airs,hsb
16
17 integer, allocatable :: index(:)
18 integer :: num_tovs_avg
19 integer, allocatable :: excess_count(:)
20 integer, allocatable :: excess_start(:)
21 integer, allocatable :: spare_count(:)
22 integer, allocatable :: spare_start(:)
23 integer :: transfer
24 logical :: copy_found
25
26 #ifdef RTTOV
27 if (trace_use) call da_trace_entry("da_qc_rad")
28
29 allocate (num_tovs_before(iv%num_inst,num_procs))
30 allocate (num_tovs_after(iv%num_inst,num_procs))
31
32 ! Cannot be more total send,receives than combination of processors
33 allocate (tovs_copy_count(iv%num_inst))
34 allocate (tovs_send_pe(iv%num_inst,num_procs*num_procs))
35 allocate (tovs_recv_pe(iv%num_inst,num_procs*num_procs))
36 allocate (tovs_send_start(iv%num_inst,num_procs*num_procs))
37 allocate (tovs_send_count(iv%num_inst,num_procs*num_procs))
38 allocate (tovs_recv_start(iv%num_inst,num_procs*num_procs))
39
40 call da_trace("da_qc_rad", message="allocated tovs redistibution arrays")
41
42 allocate (index(num_procs))
43 allocate (excess_count(num_procs))
44 allocate (excess_start(num_procs))
45 allocate (spare_count(num_procs))
46 allocate (spare_start(num_procs))
47
48 do i = 1, iv%num_inst
49 nchan = iv%instid(i)%nchan
50
51 amsua = trim(inst_name(rtminit_sensor(i))) == 'amsua'
52 amsub = trim(inst_name(rtminit_sensor(i))) == 'amsub'
53 hirs = trim(inst_name(rtminit_sensor(i))) == 'hirs'
54 msu = trim(inst_name(rtminit_sensor(i))) == 'msu'
55 airs = trim(inst_name(rtminit_sensor(i))) == 'airs'
56 hsb = trim(inst_name(rtminit_sensor(i))) == 'hsb'
57
58 if (hirs) then
59 ! 1.0 QC for HIRS
60 call da_qc_hirs(i,nchan,ob,iv)
61 !call da_warning(__FILE__,__LINE__,(/'QC Not implemented for HIRS'/))
62 else if (airs) then
63 call da_qc_airs(i,nchan,ob,iv)
64 !call da_warning(__FILE__,__LINE__,(/'QC Not implemented for AIRS'/))
65 else if ( hsb ) then
66 ! call da_qc_hsb(i,nchan,ob,iv)
67 call da_warning(__FILE__,__LINE__,(/'QC Not implemented for HSB'/))
68 else if (amsua) then
69 call da_qc_amsua(i,nchan,ob,iv)
70 else if ( amsub ) then
71 call da_qc_amsub(i,nchan,ob,iv)
72 else if (msu) then
73 ! call da_qc_msu(i,nchan, ob,iv)
74 call da_warning(__FILE__,__LINE__,(/'QC Not implemented for MSU'/))
75 else
76 write(unit=message(1),fmt='(A,A)') &
77 "Unrecognized instrument",trim(inst_name(rtminit_sensor(i)))
78 call da_error(__FILE__,__LINE__,message(1:1))
79 end if
80
81 ! Report number of observations to other processors via rootproc
82
83 num_tovs_before(i,:) = 0
84 num_tovs_before(i,myproc+1)=iv%instid(i)%num_rad
85 call da_proc_sum_ints(num_tovs_before(i,:))
86
87 #ifdef DM_PARALLEL
88 call wrf_dm_bcast_integer(num_tovs_before(i,:),num_procs)
89 #endif
90
91 num_tovs_after(i,:) = num_tovs_before(i,:)
92
93 if (rootproc .and. print_detail_radiance) then
94 write(unit=message(1),fmt='(A,I1,A)') "Instrument ",i," initial tovs distribution"
95 write(unit=message(2),fmt=*) num_tovs_before(i,:)
96 call da_message(message(1:2))
97 end if
98
99 ! Decide how to reallocate observations
100
101 num_tovs_avg=sum(num_tovs_before(i,:))/num_procs
102
103 call da_trace_int_sort(num_tovs_before(i,:),num_procs,index)
104
105 do p=1,num_procs
106 excess_count(p)=num_tovs_before(i,index(p))-num_tovs_avg
107 spare_count(p)=num_tovs_avg-num_tovs_before(i,index(p))
108 excess_start(p)=num_tovs_avg+1
109 spare_start(p)=num_tovs_before(i,index(p))+1
110 end do
111
112 ! if (rootproc .and. print_detail_radiance) then
113 ! write(unit=stdout,fmt='(A)') "After sort, pe, num,excess, spare, excess_start, spare_start"
114 ! do p=1,num_procs
115 ! write(unit=stdout,fmt='(2I4,5I10)') p,index(p)-1, &
116 ! num_tovs_before(i,index(p)), excess_count(p),spare_count(p), &
117 ! excess_start(p),spare_start(p)
118 ! end do
119 ! end if
120
121 tovs_copy_count(i) = 0
122 tovs_send_start(i,:) = 0
123 tovs_send_count(i,:) = 0
124
125 do
126 copy_found = .false.
127 do p=1,num_procs
128 if (spare_count(p) > tovs_min_transfer) then
129 do j=num_procs,1,-1
130 if (excess_count(j) > tovs_min_transfer ) then
131 copy_found = .true.
132 tovs_copy_count(i)=tovs_copy_count(i)+1
133 tovs_send_pe(i,tovs_copy_count(i)) = index(j)-1
134 tovs_recv_pe(i,tovs_copy_count(i)) = index(p)-1
135 tovs_send_start(i,tovs_copy_count(i)) = excess_start(j)
136 tovs_recv_start(i,tovs_copy_count(i)) = spare_start(p)
137 transfer=min(spare_count(p),excess_count(j))
138 tovs_send_count(i,tovs_copy_count(i)) = transfer
139 num_tovs_after(i,index(p))=num_tovs_after(i,index(p))+transfer
140 num_tovs_after(i,index(j))=num_tovs_after(i,index(j))-transfer
141 spare_count(p)=spare_count(p)-transfer
142 spare_start(p)=spare_start(p)+transfer
143 excess_count(j)=excess_count(j)-transfer
144 excess_start(j)=excess_start(j)+transfer
145 exit
146 end if
147 end do
148 end if
149 end do
150 if (.not. copy_found) exit
151 end do
152
153 if (print_detail_radiance) then
154 ! do p=1,tovs_copy_count(i)
155 ! write (unit=stdout,fmt='(3X,5(A,I5))') "Sending ",tovs_send_count(i,p), &
156 ! " from ",tovs_send_pe(i,p),":",tovs_send_start(i,p)," to ", &
157 ! tovs_recv_pe(i,p),":",tovs_recv_start(i,p)
158 ! end do
159
160 write(unit=message(1),fmt='(A,I1,A)') "Instrument ",i," final tovs distribution"
161 write(unit=message(2),fmt=*) num_tovs_after(i,:)
162 call da_message(message(1:2))
163 end if
164
165 iv % instid(i) % num_rad_glo = sum(num_tovs_after(i,:))
166 end do
167
168 deallocate (index)
169 deallocate (excess_start)
170 deallocate (excess_count)
171 deallocate (spare_start)
172 deallocate (spare_count)
173
174 if (print_detail_radiance) call da_status_rad(iv,__FILE__,__LINE__)
175
176 if (trace_use) call da_trace_exit("da_qc_rad")
177 #endif
178
179 end subroutine da_qc_rad
180
181