da_write_obs_etkf.inc

References to this file elsewhere.
1 subroutine da_write_obs_etkf(ob, iv, re)
2 
3    !-------------------------------------------------------------------------
4    ! Purpose: Writes out components of iv=O-B structure.
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(inout)  :: re      ! residual vector.
12       
13    integer                       :: n, k, num_obs, ios
14    integer                       :: ounit     ! Output unit           
15    character(len=20)             :: filename
16 
17    if (trace_use) call da_trace_entry("da_write_obs_etkf")
18 
19    !-------------------------------------------------------------------------
20    ! Fix output unit
21    !-------------------------------------------------------------------------
22 
23    call da_get_unit(ounit)
24 
25 #ifdef DM_PARALLEL
26     write(unit=filename, fmt='(a,i3.3)') 'ob.etkf.', myproc
27 #else
28     write(unit=filename, fmt='(a)') 'ob.etkf.000'
29 #endif
30 
31    open (unit=ounit,file=trim(filename),form='formatted',status='replace', &
32       iostat=ios)
33    if (ios /= 0) then
34       call da_error(__FILE__,__LINE__, &
35          (/"Cannot open ETKF observation file"//filename/))
36    end if
37 
38    ! [1] Transfer surface obs:
39 
40    if (iv%info(synop)%nlocal > 0) then
41       num_obs = 0
42       do n = 1, iv%info(synop)%nlocal
43          if (iv%info(synop)%proc_domain(1,n)) num_obs = num_obs + 1
44       end do
45       if (num_obs > 0) then
46          num_obs = 0
47          do n = 1, iv%info(synop)%nlocal
48             if (iv%info(synop)%proc_domain(1,n)) then
49                num_obs = num_obs + 1
50                if ( iv%synop(n)%u%qc >= 0 ) then
51                   write(ounit,'(3f17.7)') ob%synop(n)%u, iv%synop(n)%u%inv, iv%synop(n)%u%error
52                end if
53                if ( iv%synop(n)%v%qc >= 0 ) then
54                   write(ounit,'(3f17.7)') ob%synop(n)%v, iv%synop(n)%v%inv, iv%synop(n)%v%error
55                end if
56                if ( iv%synop(n)%t%qc >= 0 ) then
57                   write(ounit,'(3f17.7)') ob%synop(n)%t, iv%synop(n)%t%inv, iv%synop(n)%t%error
58                end if
59                if ( iv%synop(n)%p%qc >= 0 ) then
60                   write(ounit,'(3f17.7)') ob%synop(n)%p, iv%synop(n)%p%inv, iv%synop(n)%p%error
61                end if
62                if ( iv%synop(n)%q%qc >= 0 ) then
63                   write(ounit,'(3f17.7)') ob%synop(n)%q, iv%synop(n)%q%inv, iv%synop(n)%q%error
64                end if
65             end if
66          end do
67       end if
68    end if
69 
70    ! [2] Transfer metar obs:
71 
72    if (iv%info(metar)%nlocal > 0) then
73       num_obs = 0
74       do n = 1, iv%info(metar)%nlocal
75          if (iv%info(metar)%proc_domain(1,n)) num_obs = num_obs + 1
76       end do
77       if (num_obs > 0) then
78          num_obs = 0
79          do n = 1, iv%info(metar)%nlocal
80             if (iv%info(metar)%proc_domain(1,n)) then
81                num_obs = num_obs + 1
82                if ( iv%metar(n)%u%qc >= 0 ) then
83                   write(ounit,'(3f17.7)') ob%metar(n)%u, iv%metar(n)%u%inv, iv%metar(n)%u%error
84                end if
85                if ( iv%metar(n)%v%qc >= 0 ) then
86                   write(ounit,'(3f17.7)') ob%metar(n)%v, iv%metar(n)%v%inv, iv%metar(n)%v%error
87                end if
88                if ( iv%metar(n)%t%qc >= 0 ) then
89                   write(ounit,'(3f17.7)') ob%metar(n)%t, iv%metar(n)%t%inv, iv%metar(n)%t%error
90                end if
91                if ( iv%metar(n)%p%qc >= 0 ) then
92                   write(ounit,'(3f17.7)') ob%metar(n)%p, iv%metar(n)%p%inv, iv%metar(n)%p%error
93                end if
94                if ( iv%metar(n)%q%qc >= 0 ) then
95                   write(ounit,'(3f17.7)') ob%metar(n)%q, iv%metar(n)%q%inv, iv%metar(n)%q%error
96                end if
97             end if
98          end do
99       end if
100    end if
101 
102    ! [3] Transfer ships obs:
103 
104    if (iv%info(ships)%nlocal > 0) then
105       num_obs = 0
106       do n = 1, iv%info(ships)%nlocal
107          if (iv%info(ships)%proc_domain(1,n)) num_obs = num_obs + 1
108       end do
109       if (num_obs > 0) then
110          num_obs = 0
111          do n = 1, iv%info(ships)%nlocal
112             if (iv%info(ships)%proc_domain(1,n)) then
113                num_obs = num_obs + 1
114                if ( iv%ships(n)%u%qc >= 0 ) then
115                   write(ounit,'(3f17.7)') ob%ships(n)%u, iv%ships(n)%u%inv, iv%ships(n)%u%error
116                end if
117                if ( iv%ships(n)%v%qc >= 0 ) then
118                   write(ounit,'(3f17.7)') ob%ships(n)%v, iv%ships(n)%v%inv, iv%ships(n)%v%error
119                end if
120                if ( iv%ships(n)%t%qc >= 0 ) then
121                   write(ounit,'(3f17.7)') ob%ships(n)%t, iv%ships(n)%t%inv, iv%ships(n)%t%error
122                end if
123                if ( iv%ships(n)%p%qc >= 0 ) then
124                   write(ounit,'(3f17.7)') ob%ships(n)%p, iv%ships(n)%p%inv, iv%ships(n)%p%error
125                end if
126                if ( iv%ships(n)%q%qc >= 0 ) then
127                   write(ounit,'(3f17.7)') ob%ships(n)%q, iv%ships(n)%q%inv, iv%ships(n)%q%error
128                end if
129             end if
130          end do
131       end if
132    end if
133 
134   ! [4.1] Transfer Geo AMVs Obs:
135 
136    if (iv%info(geoamv)%nlocal > 0) then
137       num_obs = 0
138       do n = 1, iv%info(geoamv)%nlocal
139         if (iv%info(geoamv)%proc_domain(1,n)) num_obs = num_obs + 1
140       end do
141       if (num_obs > 0) then
142          num_obs = 0
143          do n = 1, iv%info(geoamv)%nlocal
144             if (iv%info(geoamv)%proc_domain(1,n)) then
145                num_obs = num_obs + 1
146                do k = 1, iv%info(geoamv)%levels(n)
147                   if ( iv%geoamv(n)%u(k)%qc >= 0 ) then
148                      write(ounit,'(3f17.7)') ob%geoamv(n)%u(k), iv%geoamv(n)%u(k)%inv, iv%geoamv(n)%u(k)%error
149                   end if
150                   if ( iv%geoamv(n)%v(k)%qc >= 0 ) then
151                      write(ounit,'(3f17.7)') ob%geoamv(n)%v(k), iv%geoamv(n)%v(k)%inv, iv%geoamv(n)%v(k)%error
152                   end if
153                end do
154             end if
155          end do
156       end if
157    end if
158 
159   ! [4.2] Transfer Polar AMVs Obs:
160 
161    if (iv%info(polaramv)%nlocal > 0) then
162       num_obs = 0
163       do n = 1, iv%info(polaramv)%nlocal
164         if (iv%info(polaramv)%proc_domain(1,n)) num_obs = num_obs + 1
165       end do
166       if (num_obs > 0) then
167          num_obs = 0
168          do n = 1, iv%info(polaramv)%nlocal
169             if (iv%info(polaramv)%proc_domain(1,n)) then
170                num_obs = num_obs + 1
171                do k = 1, iv%info(polaramv)%levels(n)
172                   if ( iv%polaramv(n)%u(k)%qc >= 0 ) then
173                      write(ounit,'(3f17.7)') ob%polaramv(n)%u(k), iv%polaramv(n)%u(k)%inv, iv%polaramv(n)%u(k)%error
174                   end if
175                   if ( iv%polaramv(n)%v(k)%qc >= 0 ) then
176                      write(ounit,'(3f17.7)') ob%polaramv(n)%v(k), iv%polaramv(n)%v(k)%inv, iv%polaramv(n)%v(k)%error
177                   end if
178                end do
179             end if
180          end do
181       end if
182    end if
183 
184    ! [5] Transfer gpspw obs:
185 
186    if (iv%info(gpspw)%nlocal > 0) then
187       num_obs = 0
188       do n = 1, iv%info(gpspw)%nlocal
189          if (iv%info(gpspw)%proc_domain(1,n)) num_obs = num_obs + 1
190       end do
191       if (num_obs > 0) then
192          num_obs = 0
193          do n = 1, iv%info(gpspw)%nlocal
194             if (iv%info(gpspw)%proc_domain(1,n)) then
195                num_obs = num_obs + 1
196                if ( iv%gpspw(n)%tpw%qc >= 0 ) then
197                   write(ounit,'(3f17.7)') ob%gpspw(n)%tpw, iv%gpspw(n)%tpw%inv, iv%gpspw(n)%tpw%error
198                end if
199             end if
200          end do
201       end if
202    end if
203 
204    ! [6] Transfer sonde obs:
205 
206    if (iv%info(sound)%nlocal > 0) then
207       num_obs = 0
208       do n = 1, iv%info(sound)%nlocal
209         if (iv%info(sound)%proc_domain(1,n)) num_obs = num_obs + 1
210       end do
211       if (num_obs > 0) then
212          num_obs = 0
213          do n = 1, iv%info(sound)%nlocal
214             if (iv%info(sound)%proc_domain(1,n)) then
215                num_obs = num_obs + 1
216                do k = 1, iv%info(sound)%levels(n)
217                   if ( iv%sound(n)%u(k)%qc >= 0 ) then
218                      write(ounit,'(3f17.7)') ob%sound(n)%u(k), iv%sound(n)%u(k)%inv, iv%sound(n)%u(k)%error
219                   end if
220                   if ( iv%sound(n)%v(k)%qc >= 0 ) then
221                      write(ounit,'(3f17.7)') ob%sound(n)%v(k), iv%sound(n)%v(k)%inv, iv%sound(n)%v(k)%error
222                   end if
223                   if ( iv%sound(n)%t(k)%qc >= 0 ) then
224                      write(ounit,'(3f17.7)') ob%sound(n)%t(k), iv%sound(n)%t(k)%inv, iv%sound(n)%t(k)%error
225                   end if
226                   if ( iv%sound(n)%q(k)%qc >= 0 ) then
227                      write(ounit,'(3f17.7)') ob%sound(n)%q(k), iv%sound(n)%q(k)%inv, iv%sound(n)%q(k)%error
228                   end if
229                end do
230             end if
231          end do
232       end if
233    end if
234 
235    if (iv%info(sonde_sfc)%nlocal > 0) then
236       num_obs = 0
237       do n = 1, iv%info(sonde_sfc)%nlocal
238         if (iv%info(sonde_sfc)%proc_domain(1,n)) num_obs = num_obs + 1
239       end do
240       if (num_obs > 0) then
241          num_obs = 0
242          do n = 1, iv%info(sonde_sfc)%nlocal
243             if (iv%info(sonde_sfc)%proc_domain(1,n)) then
244                num_obs = num_obs + 1
245                   if ( iv%sonde_sfc(n)%u%qc >= 0 ) then
246                      write(ounit,'(3f17.7)') ob%sonde_sfc(n)%u, iv%sonde_sfc(n)%u%inv, iv%sonde_sfc(n)%u%error
247                   end if
248                   if ( iv%sonde_sfc(n)%v%qc >= 0 ) then
249                      write(ounit,'(3f17.7)') ob%sonde_sfc(n)%v, iv%sonde_sfc(n)%v%inv, iv%sonde_sfc(n)%v%error
250                   end if
251                   if ( iv%sonde_sfc(n)%t%qc >= 0 ) then
252                      write(ounit,'(3f17.7)') ob%sonde_sfc(n)%t, iv%sonde_sfc(n)%t%inv, iv%sonde_sfc(n)%t%error
253                   end if
254                   if ( iv%sonde_sfc(n)%p%qc >= 0 ) then
255                      write(ounit,'(3f17.7)') ob%sonde_sfc(n)%p, iv%sonde_sfc(n)%p%inv, iv%sonde_sfc(n)%p%error
256                   end if
257                   if ( iv%sonde_sfc(n)%q%qc >= 0 ) then
258                      write(ounit,'(3f17.7)') ob%sonde_sfc(n)%q, iv%sonde_sfc(n)%q%inv, iv%sonde_sfc(n)%q%error
259                   end if
260             end if
261          end do
262       end if
263    end if
264 
265   ! [7] Transfer airep obs:
266 
267    if (iv%info(airep)%nlocal > 0) then
268       num_obs = 0
269       do n = 1, iv%info(airep)%nlocal
270         if (iv%info(airep)%proc_domain(1,n)) num_obs = num_obs + 1
271       end do
272       if (num_obs > 0) then
273          num_obs = 0
274          do n = 1, iv%info(airep)%nlocal
275             if (iv%info(airep)%proc_domain(1,n)) then
276                num_obs = num_obs + 1
277                do k = 1, iv%info(airep)%levels(n)
278                   if ( iv%airep(n)%u(k)%qc >= 0 ) then
279                      write(ounit,'(3f17.7)') ob%airep(n)%u(k), iv%airep(n)%u(k)%inv, iv%airep(n)%u(k)%error
280                   end if
281                   if ( iv%airep(n)%v(k)%qc >= 0 ) then
282                      write(ounit,'(3f17.7)') ob%airep(n)%v(k), iv%airep(n)%v(k)%inv, iv%airep(n)%v(k)%error
283                   end if
284                   if ( iv%airep(n)%t(k)%qc >= 0 ) then
285                      write(ounit,'(3f17.7)') ob%airep(n)%t(k), iv%airep(n)%t(k)%inv, iv%airep(n)%t(k)%error
286                   end if
287                end do
288             end if
289          end do
290       end if
291    end if
292 
293    ! [8] Transfer pilot obs:
294 
295    if (iv%info(pilot)%nlocal > 0) then
296       num_obs = 0
297       do n = 1, iv%info(pilot)%nlocal
298         if (iv%info(pilot)%proc_domain(1,n)) num_obs = num_obs + 1
299       end do
300       if (num_obs > 0) then
301          num_obs = 0
302          do n = 1, iv%info(pilot)%nlocal
303             if (iv%info(pilot)%proc_domain(1,n)) then
304                num_obs = num_obs + 1
305                do k = 1, iv%info(pilot)%levels(n)
306                   if ( iv%pilot(n)%u(k)%qc >= 0 ) then
307                      write(ounit,'(3f17.7)') ob%pilot(n)%u(k), iv%pilot(n)%u(k)%inv, iv%pilot(n)%u(k)%error
308                   end if
309                   if ( iv%pilot(n)%v(k)%qc >= 0 ) then
310                      write(ounit,'(3f17.7)') ob%pilot(n)%v(k), iv%pilot(n)%v(k)%inv, iv%pilot(n)%v(k)%error
311                   end if
312                end do
313             end if
314          end do
315       end if
316    end if
317 
318    ! [9] Transfer SSM/I obs:SSMI:
319 
320    if (iv%info(ssmi_rv)%nlocal > 0) then
321       num_obs = 0
322       do n = 1, iv%info(ssmi_rv)%nlocal
323          if (iv%info(ssmi_rv)%proc_domain(1,n)) num_obs = num_obs + 1
324       end do
325       if (num_obs > 0) then
326          num_obs = 0
327          do n = 1, iv%info(ssmi_rv)%nlocal
328             if (iv%info(ssmi_rv)%proc_domain(1,n)) then
329                num_obs = num_obs + 1
330                if ( iv%ssmi_rv(n)%speed%qc >= 0 ) then
331                   write(ounit,'(3f17.7)') ob%ssmi_rv(n)%speed, iv%ssmi_rv(n)%speed%inv, &
332                                           iv%ssmi_rv(n)%speed%error
333                end if
334                if ( iv%ssmi_rv(n)%tpw%qc >= 0 ) then
335                   write(ounit,'(3f17.7)') ob%ssmi_rv(n)%tpw, iv%ssmi_rv(n)%tpw%inv, &
336                                           iv%ssmi_rv(n)%tpw%error
337                end if
338             end if
339          end do
340       end if
341    end if
342 
343 ! SSM/I TB not coded.
344 
345    ! [10] Transfer satem obs:
346 
347    if (iv%info(satem)%nlocal > 0) then
348       num_obs = 0
349       do n = 1, iv%info(satem)%nlocal
350         if (iv%info(satem)%proc_domain(1,n)) num_obs = num_obs + 1
351       end do
352       if (num_obs > 0) then
353          num_obs = 0
354          do n = 1, iv%info(satem)%nlocal
355             if (iv%info(satem)%proc_domain(1,n)) then
356                num_obs = num_obs + 1
357                do k = 1, iv%info(satem)%levels(n)
358                   if ( iv%satem(n)%thickness(k)%qc >= 0 ) then
359                      write(ounit,'(3f17.7)') ob%satem(n)%thickness(k), iv%satem(n)%thickness(k)%inv, &
360                                              iv%satem(n)%thickness(k)%error
361                   end if
362                end do
363             end if
364          end do
365       end if
366    end if
367 
368 !  SSMT1 SSMT2 not coded.
369 
370   ! [11] Transfer scatterometer obs:
371 
372    if (iv%info(qscat)%nlocal > 0) then
373       num_obs = 0
374       do n = 1, iv%info(qscat)%nlocal
375          if (iv%info(qscat)%proc_domain(1,n)) num_obs = num_obs + 1
376       end do
377       if (num_obs > 0) then
378          num_obs = 0
379          do n = 1, iv%info(qscat)%nlocal
380             if (iv%info(qscat)%proc_domain(1,n)) then
381                num_obs = num_obs + 1
382                if ( iv%qscat(n)%u%qc >= 0 ) then
383                   write(ounit,'(3f17.7)') ob%qscat(n)%u, iv%qscat(n)%u%inv, iv%qscat(n)%u%error
384                end if
385                if ( iv%qscat(n)%v%qc >= 0 ) then
386                   write(ounit,'(3f17.7)') ob%qscat(n)%v, iv%qscat(n)%v%inv, iv%qscat(n)%v%error
387                end if
388             end if
389          end do
390       end if
391    end if
392 
393   ! [12] Transfer profiler obs:
394 
395    if (iv%info(profiler)%nlocal > 0) then
396       num_obs = 0
397       do n = 1, iv%info(profiler)%nlocal
398         if (iv%info(profiler)%proc_domain(1,n)) num_obs = num_obs + 1
399       end do
400       if (num_obs > 0) then
401          num_obs = 0
402          do n = 1, iv%info(profiler)%nlocal
403             if (iv%info(profiler)%proc_domain(1,n)) then
404                num_obs = num_obs + 1
405                do k = 1, iv%info(profiler)%levels(n)
406                   if ( iv%profiler(n)%u(k)%qc >= 0 ) then
407                      write(ounit,'(3f17.7)') ob%profiler(n)%u(k), iv%profiler(n)%u(k)%inv, iv%profiler(n)%u(k)%error
408                   end if
409                   if ( iv%profiler(n)%v(k)%qc >= 0 ) then
410                      write(ounit,'(3f17.7)') ob%profiler(n)%v(k), iv%profiler(n)%v(k)%inv, iv%profiler(n)%v(k)%error
411                   end if
412                end do
413             end if
414          end do
415       end if
416    end if
417 
418    ! Transfer Buoy obs:
419 
420    if (iv%info(buoy)%nlocal > 0) then
421       num_obs = 0
422       do n = 1, iv%info(buoy)%nlocal
423          if (iv%info(buoy)%proc_domain(1,n)) num_obs = num_obs + 1
424       end do
425       if (num_obs > 0) then
426          num_obs = 0
427          do n = 1, iv%info(buoy)%nlocal
428             if (iv%info(buoy)%proc_domain(1,n)) then
429                num_obs = num_obs + 1
430                if ( iv%buoy(n)%u%qc >= 0 ) then
431                   write(ounit,'(3f17.7)') ob%buoy(n)%u, iv%buoy(n)%u%inv, iv%buoy(n)%u%error
432                end if
433                if ( iv%buoy(n)%v%qc >= 0 ) then
434                   write(ounit,'(3f17.7)') ob%buoy(n)%v, iv%buoy(n)%v%inv, iv%buoy(n)%v%error
435                end if
436                if ( iv%buoy(n)%t%qc >= 0 ) then
437                   write(ounit,'(3f17.7)') ob%buoy(n)%t, iv%buoy(n)%t%inv, iv%buoy(n)%t%error
438                end if
439                if ( iv%buoy(n)%p%qc >= 0 ) then
440                   write(ounit,'(3f17.7)') ob%buoy(n)%p, iv%buoy(n)%p%inv, iv%buoy(n)%p%error
441                end if
442                if ( iv%buoy(n)%q%qc >= 0 ) then
443                   write(ounit,'(3f17.7)') ob%buoy(n)%q, iv%buoy(n)%q%inv, iv%buoy(n)%q%error
444                end if
445             end if
446          end do
447       end if
448    end if
449 
450    ! Transfer TC bogus obs:
451 
452    if (iv%info(bogus)%nlocal > 0) then
453       num_obs = 0
454       do n = 1, iv%info(bogus)%nlocal
455         if (iv%info(bogus)%proc_domain(1,n)) num_obs = num_obs + 1
456       end do
457       if (num_obs > 0) then
458          num_obs = 0
459          do n = 1, iv%info(bogus)%nlocal
460             if (iv%info(bogus)%proc_domain(1,n)) then
461                num_obs = num_obs + 1
462                do k = 1, iv%info(bogus)%levels(n)
463                   if ( iv%bogus(n)%u(k)%qc >= 0 ) then
464                      write(ounit,'(3f17.7)') ob%bogus(n)%u(k), iv%bogus(n)%u(k)%inv, iv%bogus(n)%u(k)%error
465                   end if
466                   if ( iv%bogus(n)%v(k)%qc >= 0 ) then
467                      write(ounit,'(3f17.7)') ob%bogus(n)%v(k), iv%bogus(n)%v(k)%inv, iv%bogus(n)%v(k)%error
468                   end if
469                   if ( iv%bogus(n)%t(k)%qc >= 0 ) then
470                      write(ounit,'(3f17.7)') ob%bogus(n)%t(k), iv%bogus(n)%t(k)%inv, iv%bogus(n)%t(k)%error
471                   end if
472                   if ( iv%bogus(n)%q(k)%qc >= 0 ) then
473                      write(ounit,'(3f17.7)') ob%bogus(n)%q(k), iv%bogus(n)%q(k)%inv, iv%bogus(n)%q(k)%error
474                   end if
475                end do
476             end if
477          end do
478       end if
479    end if
480 
481    ! Transfer AIRS retrievals:
482 
483    if (iv%info(airsr)%nlocal > 0) then
484       num_obs = 0
485       do n = 1, iv%info(airsr)%nlocal
486         if (iv%info(airsr)%proc_domain(1,n)) num_obs = num_obs + 1
487       end do
488       if (num_obs > 0) then
489          num_obs = 0
490          do n = 1, iv%info(airsr)%nlocal
491             if (iv%info(airsr)%proc_domain(1,n)) then
492                num_obs = num_obs + 1
493                do k = 1, iv%info(airsr)%levels(n)
494                   if ( iv%airsr(n)%t(k)%qc >= 0 ) then
495                      write(ounit,'(3f17.7)') ob%airsr(n)%t(k), iv%airsr(n)%t(k)%inv, iv%airsr(n)%t(k)%error
496                   end if
497                   if ( iv%airsr(n)%q(k)%qc >= 0 ) then
498                      write(ounit,'(3f17.7)') ob%airsr(n)%q(k), iv%airsr(n)%q(k)%inv, iv%airsr(n)%q(k)%error
499                   end if
500                end do
501             end if
502          end do
503       end if
504    end if
505 
506    ! Transfer gpsref obs:
507  
508    if (iv%info(gpsref)%nlocal > 0) then
509       num_obs = 0
510       do n = 1, iv%info(gpsref)%nlocal
511         if (iv%info(gpsref)%proc_domain(1,n)) num_obs = num_obs + 1
512       end do
513       if (num_obs > 0) then
514          num_obs = 0
515          do n = 1, iv%info(gpsref)%nlocal
516             if (iv%info(gpsref)%proc_domain(1,n)) then
517                num_obs = num_obs + 1
518                do k = 1, iv%info(gpsref)%levels(n)
519                   if ( iv%gpsref(n)%ref(k)%qc >= 0 ) then
520                      write(ounit,'(3f17.7)') ob%gpsref(n)%ref(k), iv%gpsref(n)%ref(k)%inv, iv%gpsref(n)%ref(k)%error
521                   end if
522                end do
523             end if
524          end do
525       end if
526    end if
527 
528    close (ounit)
529    call da_free_unit(ounit)
530 
531    if (trace_use) call da_trace_exit("da_write_obs_etkf")
532 
533 end subroutine da_write_obs_etkf
534