da_read_radar.inc
References to this file elsewhere.
1 subroutine da_read_radar (ob, xp, filename)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: Read the Radar observation file
5 !----------------------------------------------------------------------------------------!
6
7 implicit none
8
9 type (xpose_type), intent(in) :: xp ! Domain decomposition vars.
10 type (ob_type), intent (inout) :: ob
11
12 character(*), intent(in) :: filename
13
14 character (LEN = 120) :: char_total_Radar
15 character (LEN = 120) :: char_ned
16
17 integer :: i, j, n, iost, nlevels, fm
18
19 type (Radar_multi_level_type) :: platform
20
21 logical :: outside
22
23 integer :: total_obs, num_Radar, total_Radar
24 integer ::n_dup, ndup, iunit
25
26 if (trace_use) call da_trace_entry("da_read_radar")
27
28 total_obs = ob%total_obs
29 num_Radar = ob%num_Radar
30
31 ! 1. open file
32
33 call da_get_unit(iunit)
34 open(unit = iunit, &
35 FILE = trim(filename), &
36 FORM = 'FORMATTED', &
37 ACCESS = 'SEQUENTIAL', &
38 iostat = iost, &
39 STATUS = 'OLD')
40
41 if (iost /= 0) then
42 ! Missing file does not matter
43 call da_warning(__FILE__,__LINE__, &
44 (/"Cannot open radar file "//trim(filename)/))
45 return
46 end if
47
48 ! 2. read total radar
49
50 ! 2.1 read first line
51
52 read (unit=iunit, fmt = '(A)', iostat = iost) char_total_Radar
53
54 ! 2.2 process error
55
56 if (iost /= 0) then
57 call da_error(__FILE__,__LINE__, &
58 (/"Cannot read radar file"/))
59 end if
60
61 ! 2.3 total radar number
62
63 read (unit=char_total_Radar (15:17),fmt='(I3)', iostat = iost) total_Radar
64
65 if (print_detail_radar) write (unit=stdout,fmt='(/,A,I3,/)') &
66 ' TOTAL RADAR: ', total_Radar
67
68 ! 2.4 skip one line
69
70 read (unit=iunit, fmt = '(A)', iostat = iost)
71
72 ! 3. read radar data
73
74 do n = 1, total_Radar
75
76 ! 3.1 skip one blank line
77
78 read (unit=iunit, fmt = '(A)', iostat = iost)
79
80 ! 3.2 read header
81
82 read (unit=iunit, fmt = '(A)', iostat = iost) char_ned
83
84 ! 3.3 read header information
85
86 read (unit=char_ned (1:5),fmt='(A5)', iostat = iost) platform % stn % platform
87
88 if (print_detail_radar) write (unit=stdout,fmt='(A)') 'RADAR Observations information'
89
90 read (unit=char_ned (8:19),fmt='(A12)', iostat = iost) platform % stn % name
91
92 if (print_detail_radar) write (unit=stdout,fmt='(A,A5,A,A12)') &
93 ' Expect: ',platform % stn % platform, &
94 'data at station:', platform % stn % name
95
96 read (unit=char_ned(20:27),fmt='(F8.3)', iostat = iost) platform % stn % lon
97
98 read (unit=char_ned (30:37),fmt='(F8.3)', iostat = iost) platform % stn % lat
99
100 read (unit=char_ned (40:47),fmt='(F8.1)', iostat = iost) platform % stn % elv
101
102 if (print_detail_radar) write (unit=stdout,fmt='(A,2(F8.3,2X),F8.1)') &
103 'The station longitude, latutude, and altitude are: ', &
104 platform % stn % lon, &
105 platform % stn % lat, platform % stn % elv
106
107 read (unit=char_ned (50:68),fmt='(A19)', iostat = iost) platform % stn % date_char
108
109 if (print_detail_radar) write (unit=stdout,fmt='(A,A19)') &
110 'The observation time for this data is ', &
111 platform % stn % date_char
112
113 read (unit=char_ned (69:74),fmt='(I6)', iostat = iost) platform % stn % numObs
114
115 if (print_detail_radar) write (unit=stdout,fmt='(A,I6)') &
116 'Total number of Super-observations is ', &
117 platform % stn % numObs
118
119 read (unit=char_ned (75:80),fmt='(I6)', iostat = iost) platform % stn % levels
120
121 if (print_detail_radar) write (unit=stdout,fmt='(A,I6)') &
122 'Vertical layers for each Super-observation is ', &
123 platform % stn % levels
124
125 ! 3.4 skip two lines
126
127 read (unit=iunit, fmt = '(A)', iostat = iost)
128 read (unit=iunit, fmt = '(A)', iostat = iost)
129
130 ! 3.5 loop over records
131
132 reports: do j = 1, platform % stn % numObs
133
134 ! 3.5.1 read station general info
135
136 read (unit = iunit, iostat = iost, &
137 fmt = '(A12,3X,A19,2X,2(F12.3,2X),F8.1,2X,I6)') &
138 platform % info % platform, &
139 platform % info % date_char, &
140 platform % info % lat, &
141 platform % info % lon, &
142 platform % info % elv, &
143 platform % info % levels
144
145 read(platform % info % platform (4:6), '(I3)') fm
146
147 total_obs = total_obs + 1
148
149
150 ! 3.5.2 read each level
151
152 do i = 1, platform % info % levels
153 ! height
154 platform%each(i) = Radar_each_level_type(missing_r, missing, -1.0, &
155 field_type(missing_r, missing, missing_r), & ! rv
156 field_type(missing_r, missing, missing_r)) ! rf
157
158 read (unit = iunit, fmt = '(3X, F12.1, 2(F12.3,I4,F12.3,2X))') &
159 platform % each (i) % height, &
160 platform % each (i) % rv % inv, &
161 platform % each (i) % rv % qc, &
162 platform % each (i) % rv % error, &
163 platform % each (i) % rf % inv, &
164 platform % each (i) % rf % qc, &
165 platform % each (i) % rf % error
166
167 if (platform % each (i) % rv % error == 0.0) then
168 platform % each (i) % rv % error = 1.0
169 end if
170
171 if (platform % each (i) % rf % error == 0.0) then
172 platform % each (i) % rf % error = 1.0
173 end if
174
175 if (platform % each (i) % rv % inv == missing_r .or. &
176 platform % each (i) % rv % error == missing_r) then
177 platform % each (i) % rv % qc = missing_data
178 end if
179
180 if (platform % each (i) % rf % inv == missing_r .or. &
181 platform % each (i) % rf % error == missing_r) then
182 platform % each (i) % rf % qc = missing_data
183 end if
184 end do
185
186 call da_ll_to_xy (platform%info, platform%loc, xp, outside)
187
188 if (outside) then
189 cycle reports
190 end if
191
192 nlevels = platform%info%levels
193
194 if (nlevels > max_ob_levels) then
195 write(unit=message(1),fmt='(A,2I8)') &
196 ' Radar=> nlevels > max_ob_levels:',nlevels, max_ob_levels
197 call da_warning(__FILE__,__LINE__,message(1:1))
198 nlevels = max_ob_levels
199 platform%info%levels = nlevels
200 else if (nlevels < 1) then
201 cycle reports
202 end if
203
204 ! Loop over duplicating obs for global
205 n_dup = 1
206 if (global .and. &
207 (platform%loc%i == xp%ids .or. platform%loc%i == xp%ide)) n_dup= 2
208 do ndup = 1, n_dup
209 select case (fm)
210
211 case (128)
212
213 if (ndup==1) ob%num_Radar_glo = ob%num_Radar_glo + 1
214
215 num_Radar = num_Radar + 1
216
217 if (num_Radar > max_Radar) then
218 write(unit=message(1),fmt='(A,I6,A,I6)') &
219 "Radar #= ",ob % num_Radar, " > max_Radar = ", max_Radar
220 call da_error(__FILE__,__LINE__,message(1:1))
221 end if
222 platform % loc%obs_global_index = ob%num_Radar_glo
223
224 ob % Radar (num_Radar) % stn_loc % lon = platform % stn % lon
225 ob % Radar (num_Radar) % stn_loc % lat = platform % stn % lat
226 ob % Radar (num_Radar) % stn_loc % elv = platform % stn % elv
227
228 ob % Radar (num_Radar) % info = platform % info
229 ob % Radar (num_Radar) % loc = platform % loc
230
231 allocate (ob % Radar (num_Radar) % height (1:nlevels))
232 allocate (ob % Radar (num_Radar) % height_qc(1:nlevels))
233 allocate (ob % Radar (num_Radar) % zk (1:nlevels))
234 allocate (ob % Radar (num_Radar) % rv (1:nlevels))
235 allocate (ob % Radar (num_Radar) % rf (1:nlevels))
236
237 do i = 1, nlevels
238 ob % Radar (num_Radar) % height(i) = platform % each(i) % height
239 ob % Radar (num_Radar) % height_qc(i) = platform % each(i) % height_qc
240 ob % Radar (num_Radar) % zk(i) = platform % each(i) % zk
241 ob % Radar (num_Radar) % rv(i) = platform % each(i) % rv
242 ob % Radar (num_Radar) % rf(i) = platform % each(i) % rf
243 end do
244
245 case default;
246 write(unit=message(1), fmt='(a)') 'Unsaved obs found:'
247 write(unit=message(2), fmt='(2a)') &
248 'platform % info % platform=', platform % info % platform
249 write(unit=message(3), fmt='(a, i3)') &
250 'platform % info % levels=', platform % info % levels
251 call da_warning(__FILE__,__LINE__,message(1:3))
252 end select
253
254 if (global .and. ndup == 1) then
255 if (platform%loc % i >= xp % ide) then
256 platform%loc%i = xp % ids
257 platform%loc%proc_domain = .false.
258 else if (platform%loc % i <= xp % ids) then
259 platform%loc%i = xp % ide
260 platform%loc%proc_domain = .false.
261 end if
262 end if
263 end do ! loop over duplicate
264 end do reports
265 end do
266
267 ob % num_Radar = num_Radar
268
269 ! 4. print out
270
271 if (print_detail_radar) then
272 write(unit=stdout, fmt='(a)') ' '
273 write(unit=stdout, fmt='(5x,a,i6,a)') &
274 ' Read: ', ob % num_Radar, ' Radar reports '
275 end if
276
277 if (ob%ob_numb(ob%current_ob_time)%Radar /= ob%num_Radar) then
278 write(unit=message(1),fmt='(a, i6, 2x, a, i6)') &
279 'Radar obs mismatch ob%ob_numb(ob%current_ob_time)%Radar=', &
280 ob%ob_numb(ob%current_ob_time)%Radar, &
281 'ob%num_Radar=', ob%num_Radar
282 call da_error(__FILE__,__LINE__,message(1:1))
283 end if
284
285 close(iunit)
286 call da_free_unit(iunit)
287
288 if (trace_use) call da_trace_exit("da_read_radar")
289
290
291 end subroutine da_read_radar
292
293