da_write_oa_rad_ascii.inc
References to this file elsewhere.
1 subroutine da_write_oa_rad_ascii ( ob, iv, re )
2
3 !---------------------------------------------------------------------------
4 ! Purpose: write out OMB and OMA vector structure for radiance data.
5 !---------------------------------------------------------------------------
6
7 implicit none
8
9 type (y_type), intent(in) :: ob ! Observation structure.
10 type (iv_type), intent(in) :: iv ! O-B structure.
11 type (y_type), intent(in) :: re ! O-A structure.
12
13 #ifdef RTTOV
14 integer :: n ! Loop counter.
15 integer :: i, k ! Index dimension.
16 integer :: nlevelss ! Number of obs levels.
17
18 integer :: ios, oma_rad_unit
19 character(len=filename_len) :: filename
20 integer :: ndomain
21
22 if (trace_use) call da_trace_entry("da_write_oa_rad_ascii")
23
24 do i = 1, iv%num_inst
25 if (iv%instid(i)%num_rad < 1) cycle
26
27 ! count number of obs within the proc_domain
28 !---------------------------------------------
29 ndomain = 0
30 do n =1,iv%instid(i)%num_rad
31 if (iv%instid(i)%info%proc_domain(1,n)) then
32 ndomain = ndomain + 1
33 end if
34 end do
35
36 #ifdef DM_PARALLEL
37 write(unit=filename, fmt='(a,i2.2)') 'oma_'//trim(iv%instid(i)%rttovid_string)//'.', myproc
38 #else
39 write(unit=filename, fmt='(a)') 'oma_'//trim(iv%instid(i)%rttovid_string)
40 #endif
41 call da_get_unit(oma_rad_unit)
42 open(unit=oma_rad_unit,file=trim(filename),form='formatted',iostat=ios)
43 if (ios /= 0) then
44 call da_error(__FILE__,__LINE__, &
45 (/"Cannot open oma radiance file"//filename/))
46 end if
47 write(unit=oma_rad_unit,fmt='(a,a,i7,a,i5,a)') trim(iv%instid(i)%rttovid_string), &
48 ' number-of-pixels : ', ndomain, &
49 ' channel-number-of-each-pixel : ', iv%instid(i)%nchan, &
50 ' index-of-channels : '
51 write(unit=oma_rad_unit,fmt='(10i5)') iv%instid(i)%ichan
52 write(unit=oma_rad_unit,fmt= *) ' pixel-info : i date scanpos landsea_mask elv lat lon satzen satazi '
53 write(unit=oma_rad_unit,fmt= *) ' xb-surf-info : i t2m mr2m(ppmv) u10 v10 ps ts smois tslb snowh isflg &
54 & soiltyp vegtyp vegfra elev clwp'
55 ndomain = 0
56 do n=1,iv%instid(i)%num_rad
57 if (iv%instid(i)%info%proc_domain(1,n)) then
58 ndomain=ndomain+1
59 write(unit=oma_rad_unit,fmt='(a,i7,2x,a,2i3,f6.0,4f8.2)') 'INFO : ', ndomain, &
60 iv%instid(i)%info%date_char(n), &
61 iv%instid(i)%scanpos(n), &
62 iv%instid(i)%landsea_mask(n), &
63 iv%instid(i)%info%elv(n), &
64 iv%instid(i)%info%lat(1,n), &
65 iv%instid(i)%info%lon(1,n), &
66 iv%instid(i)%satzen(n), &
67 iv%instid(i)%satazi(n)
68 if (iv%instid(i)%isflg(n)==0) then
69 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') ' SEA : ', n, &
70 iv%instid(i)%t2m(n), &
71 iv%instid(i)%mr2m(n), &
72 iv%instid(i)%u10(n), &
73 iv%instid(i)%v10(n), &
74 iv%instid(i)%ps(n), &
75 iv%instid(i)%ts(n), &
76 iv%instid(i)%smois(n), &
77 iv%instid(i)%tslb(n), &
78 iv%instid(i)%snowh(n), &
79 iv%instid(i)%isflg(n), &
80 nint(iv%instid(i)%soiltyp(n)), &
81 nint(iv%instid(i)%vegtyp(n)), &
82 iv%instid(i)%vegfra(n), &
83 iv%instid(i)%elevation(n), &
84 iv%instid(i)%clwp(n)
85 else if (iv%instid(i)%isflg(n)==1) then
86 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') ' ICE : ', n, &
87 iv%instid(i)%t2m(n), &
88 iv%instid(i)%mr2m(n), &
89 iv%instid(i)%u10(n), &
90 iv%instid(i)%v10(n), &
91 iv%instid(i)%ps(n), &
92 iv%instid(i)%ts(n), &
93 iv%instid(i)%smois(n), &
94 iv%instid(i)%tslb(n), &
95 iv%instid(i)%snowh(n), &
96 iv%instid(i)%isflg(n), &
97 nint(iv%instid(i)%soiltyp(n)), &
98 nint(iv%instid(i)%vegtyp(n)), &
99 iv%instid(i)%vegfra(n), &
100 iv%instid(i)%elevation(n), &
101 iv%instid(i)%clwp(n)
102 else if (iv%instid(i)%isflg(n)==2) then
103 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'LAND : ', n, &
104 iv%instid(i)%t2m(n), &
105 iv%instid(i)%mr2m(n), &
106 iv%instid(i)%u10(n), &
107 iv%instid(i)%v10(n), &
108 iv%instid(i)%ps(n), &
109 iv%instid(i)%ts(n), &
110 iv%instid(i)%smois(n), &
111 iv%instid(i)%tslb(n), &
112 iv%instid(i)%snowh(n), &
113 iv%instid(i)%isflg(n), &
114 nint(iv%instid(i)%soiltyp(n)), &
115 nint(iv%instid(i)%vegtyp(n)), &
116 iv%instid(i)%vegfra(n), &
117 iv%instid(i)%elevation(n), &
118 iv%instid(i)%clwp(n)
119 else if (iv%instid(i)%isflg(n)==3) then
120 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'SNOW : ', n, &
121 iv%instid(i)%t2m(n), &
122 iv%instid(i)%mr2m(n), &
123 iv%instid(i)%u10(n), &
124 iv%instid(i)%v10(n), &
125 iv%instid(i)%ps(n), &
126 iv%instid(i)%ts(n), &
127 iv%instid(i)%smois(n), &
128 iv%instid(i)%tslb(n), &
129 iv%instid(i)%snowh(n), &
130 iv%instid(i)%isflg(n), &
131 nint(iv%instid(i)%soiltyp(n)), &
132 nint(iv%instid(i)%vegtyp(n)), &
133 iv%instid(i)%vegfra(n), &
134 iv%instid(i)%elevation(n), &
135 iv%instid(i)%clwp(n)
136 else if (iv%instid(i)%isflg(n)==4) then
137 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'MSEA : ', n, &
138 iv%instid(i)%t2m(n), &
139 iv%instid(i)%mr2m(n), &
140 iv%instid(i)%u10(n), &
141 iv%instid(i)%v10(n), &
142 iv%instid(i)%ps(n), &
143 iv%instid(i)%ts(n), &
144 iv%instid(i)%smois(n), &
145 iv%instid(i)%tslb(n), &
146 iv%instid(i)%snowh(n), &
147 iv%instid(i)%isflg(n), &
148 nint(iv%instid(i)%soiltyp(n)), &
149 nint(iv%instid(i)%vegtyp(n)), &
150 iv%instid(i)%vegfra(n), &
151 iv%instid(i)%elevation(n), &
152 iv%instid(i)%clwp(n)
153 else if (iv%instid(i)%isflg(n)==5) then
154 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'MICE : ', n, &
155 iv%instid(i)%t2m(n), &
156 iv%instid(i)%mr2m(n), &
157 iv%instid(i)%u10(n), &
158 iv%instid(i)%v10(n), &
159 iv%instid(i)%ps(n), &
160 iv%instid(i)%ts(n), &
161 iv%instid(i)%smois(n), &
162 iv%instid(i)%tslb(n), &
163 iv%instid(i)%snowh(n), &
164 iv%instid(i)%isflg(n), &
165 nint(iv%instid(i)%soiltyp(n)), &
166 nint(iv%instid(i)%vegtyp(n)), &
167 iv%instid(i)%vegfra(n), &
168 iv%instid(i)%elevation(n), &
169 iv%instid(i)%clwp(n)
170 else if (iv%instid(i)%isflg(n)==6) then
171 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'MLND : ', n, &
172 iv%instid(i)%t2m(n), &
173 iv%instid(i)%mr2m(n), &
174 iv%instid(i)%u10(n), &
175 iv%instid(i)%v10(n), &
176 iv%instid(i)%ps(n), &
177 iv%instid(i)%ts(n), &
178 iv%instid(i)%smois(n), &
179 iv%instid(i)%tslb(n), &
180 iv%instid(i)%snowh(n), &
181 iv%instid(i)%isflg(n), &
182 nint(iv%instid(i)%soiltyp(n)), &
183 nint(iv%instid(i)%vegtyp(n)), &
184 iv%instid(i)%vegfra(n), &
185 iv%instid(i)%elevation(n), &
186 iv%instid(i)%clwp(n)
187 else if (iv%instid(i)%isflg(n)==7) then
188 write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3)') 'MSNO : ', n, &
189 iv%instid(i)%t2m(n), &
190 iv%instid(i)%mr2m(n), &
191 iv%instid(i)%u10(n), &
192 iv%instid(i)%v10(n), &
193 iv%instid(i)%ps(n), &
194 iv%instid(i)%ts(n), &
195 iv%instid(i)%smois(n), &
196 iv%instid(i)%tslb(n), &
197 iv%instid(i)%snowh(n), &
198 iv%instid(i)%isflg(n), &
199 nint(iv%instid(i)%soiltyp(n)), &
200 nint(iv%instid(i)%vegtyp(n)), &
201 iv%instid(i)%vegfra(n), &
202 iv%instid(i)%elevation(n), &
203 iv%instid(i)%clwp(n)
204 end if
205
206 write(unit=oma_rad_unit,fmt='(a)') 'OBS : '
207 write(unit=oma_rad_unit,fmt='(10f11.2)') ob%instid(i)%tb(:,n)
208 write(unit=oma_rad_unit,fmt='(a)') 'BAK : '
209 write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%tb_xb(:,n)
210 write(unit=oma_rad_unit,fmt='(a)') 'IVBC : '
211 write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%tb_inv(:,n)
212 write(unit=oma_rad_unit,fmt='(a)') 'OMA : '
213 write(unit=oma_rad_unit,fmt='(10f11.2)') re%instid(i)%tb(:,n)
214 write(unit=oma_rad_unit,fmt='(a)') 'EMS : '
215 write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%emiss(1:iv%instid(i)%nchan,n)
216 write(unit=oma_rad_unit,fmt='(a)') 'ERR : '
217 write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%tb_error(:,n)
218 write(unit=oma_rad_unit,fmt='(a)') 'QC : '
219 write(unit=oma_rad_unit,fmt='(10i11)') iv%instid(i)%tb_qc(:,n)
220
221 if (write_profile) then
222 nlevelss = iv%instid(i)%nlevels
223 write(unit=oma_rad_unit,fmt=*) &
224 'RTM_level pres(mb) T(k) Q(ppmv) WRF_level pres(mb) T(k) q(g/kg) clw(g/kg) rain(g/kg)'
225 do k=kts,kte
226 if (k <= nlevelss) then
227 write(unit=oma_rad_unit,fmt='(i3,f10.2,f8.2,e11.4,i3,f10.2,f8.2,3e11.4)') &
228 k, & ! RTTOV levels
229 coefs(i) % ref_prfl_p(k) , &
230 iv%instid(i)%t(k,n) , &
231 iv%instid(i)%mr(k,n), &
232 k, & ! WRF model levels
233 iv%instid(i)%pm(k,n) , &
234 iv%instid(i)%tm(k,n) , &
235 iv%instid(i)%qm(k,n)*1000 , &
236 iv%instid(i)%qcw(k,n)*1000.0, &
237 iv%instid(i)%qrn(k,n)*1000.0
238
239 else
240 write(unit=oma_rad_unit,fmt='(32x,i3,f10.2,f8.2,3e11.4)') k, &
241 iv%instid(i)%pm(k,n) , &
242 iv%instid(i)%tm(k,n) , &
243 iv%instid(i)%qm(k,n)*1000 , &
244 iv%instid(i)%qcw(k,n)*1000.0, &
245 iv%instid(i)%qrn(k,n)*1000.0
246 ! iv%instid(i)%qci(k,n)*1000.0, &
247 ! iv%instid(i)%qsn(k,n)*1000.0, &
248 ! iv%instid(i)%qgr(k,n)*10000.0
249 end if
250 end do ! end loop profile
251 end if ! end if write_profile
252 end if ! end if proc_domain
253 end do ! end do pixels
254 close(unit=oma_rad_unit)
255 call da_free_unit(oma_rad_unit)
256 end do !! end do instruments
257
258 if (trace_use) call da_trace_exit("da_write_oa_rad_ascii")
259 #endif
260
261 end subroutine da_write_oa_rad_ascii
262
263