da_read_bufrtovs.inc

References to this file elsewhere.
1 subroutine da_read_bufrtovs (obstype,iv, xp,infile)
2 
3    !---------------------------------------------------------------------------
4    !  Purpose: read in NCEP bufr tovs 1b data to innovation structure
5    !
6    !   METHOD: use F90 sequential data structure to avoid reading file twice  
7    !            so that da_scan_bufrtovs is not necessary any more.
8    !            1. read file radiance data in sequential data structure
9    !            2. do gross QC check
10    !            3. assign sequential data structure to innovation structure
11    !               and deallocate sequential data structure
12    !---------------------------------------------------------------------------
13 
14    use da_control
15 
16    implicit none
17 
18    character(5)      ,  intent (in)  :: obstype
19    character(20)     ,  intent (in)  :: infile
20    type (xpose_type) ,  intent (in)  :: xp
21    type (ob_type)    ,  intent (inout) :: iv
22 
23 #ifdef BUFR
24 
25    integer          :: iost
26    integer(i_kind), allocatable :: nread(:)
27 
28    integer(i_kind),parameter:: n1bhdr=14
29    integer(i_kind),parameter:: maxinfo=12
30    integer(i_kind),parameter:: maxchanl=100
31 
32    logical hirs2,hirs3,hirs4,msu,amsua,amsub,mhs
33    logical outside, outside_all, iuse
34    integer :: inst
35 
36    character(10) date
37    character(8) subset,subfgn
38    character(80) hdr1b
39    integer(i_kind) ihh,i,k,ifov,idd,ireadmg,ireadsb
40    integer(i_kind) iret,idate,im,iy,nchan
41 
42    ! thinning variables
43    integer(i_kind) itt,itx,iobs,iout
44    real(r_kind) terrain,timedif,crit,dist
45    real(r_kind) dlon_earth,dlat_earth
46 
47    real(r_kind) tbmin,tbmax, tbbad
48    real(r_kind) panglr,rato
49    ! real(r_kind) rmask
50    real(r_kind) step,start
51 
52    real(r_double),dimension(maxinfo+maxchanl):: data1b8
53    real(r_double),dimension(n1bhdr):: bfr1bhdr
54 
55    ! Instrument triplet, follow the convension of RTTOV 
56    integer   :: platform_id, satellite_id, sensor_id
57 
58    ! pixel information
59    integer   ::  year,month,day,hour,minute,second  ! observation time
60    real      ::  obs_time
61    ! real      ::  rlat, rlon                         !  lat/lon in degrees   for Anfovs
62    real      ::  satzen, satazi, solzen     !  scan angles for Anfovs
63    integer   ::  landsea_mask
64    real      ::  srf_height
65    ! channels' bright temperature
66    real , allocatable ::   tb_inv(:)                    !  bright temperatures
67    !  end type bright_temperature
68 
69    type (datalink_type), pointer    :: head, p, current
70 
71    integer                        ::  ifgat
72    type(info_type)                ::  info
73    type(model_loc_type)           ::  loc
74 
75    data hdr1b /'SAID FOVN YEAR MNTH DAYS HOUR MINU SECO CLAT CLON SAZA SOZA HOLS LSQL'/
76    !  data hdr1b /'FOVN YEAR MNTH DAYS HOUR MINU SECO CLAT CLON SAZA SOZA HOLS LSQL SLNM BEARAZ'/
77 
78    data tbmin,tbmax,tbbad / 50.0_r_kind, 550.0_r_kind, -9.99e11_r_kind /
79    integer :: num_tovs_local, num_tovs_file, num_tovs_global, num_tovs_selected
80    integer :: num_tovs_thinned, num_tovs_used
81    integer :: lnbufr
82    integer :: n
83 
84    call da_trace_entry("da_read_bufrtovs")
85 
86    ! Initialize variables
87 
88    nchan = 20
89    allocate(nread(1:rtminit_nsensor))
90    nread(1:rtminit_nsensor) = 0
91 
92    ! Set various variables depending on type of data to be read
93 
94    platform_id  = 1                 !! for NOAA series
95 
96    hirs2=     obstype == 'hirs2'
97    hirs3=     obstype == 'hirs3'
98    hirs4=     obstype == 'hirs4'
99    msu=       obstype == 'msu  '
100    amsua=     obstype == 'amsua'
101    amsub=     obstype == 'amsub'
102    mhs=       obstype == 'mhs  '
103 
104    if (hirs2) then
105       sensor_id    =  0
106       step   = 1.80_r_kind
107       start  = -49.5_r_kind
108       nchan=nchan_hirs2
109       subfgn='NC021021'
110       rato=1.1363987_r_kind
111    else if (hirs3) then 
112       sensor_id    =  0
113       step   = 1.80_r_kind
114       start  = -49.5_r_kind
115       nchan=nchan_hirs3
116       subfgn='NC021025'
117    else if (msu) then
118       sensor_id    =  1
119       step   = 9.474_r_kind
120       start  = -47.37_r_kind
121       nchan=nchan_msu
122       subfgn='NC021022'
123       rato=1.1363987_r_kind
124    else if (amsua) then
125       sensor_id    =  3
126       step   = three + one/three
127       start  = -48.33_r_kind
128       nchan=nchan_amsua
129       subfgn='NC021023'
130    else if (amsub)  then
131       sensor_id    =  4
132       step   = 1.1_r_kind
133       start  = -48.95_r_kind
134       nchan=nchan_amsub
135       subfgn='NC021024'
136    end if
137 
138    allocate (tb_inv(nchan))
139 
140    ! 0.0  Open unit to satellite bufr file and read file header
141    !--------------------------------------------------------------
142 
143    call da_get_unit(lnbufr)
144    open(unit=lnbufr,file=trim(infile),form='unformatted', &
145       iostat = iost, status = 'old')
146    if (iost /= 0) then
147       call da_warning(__FILE__,__LINE__, &
148          (/"Cannot open file "//trim(infile)/))
149       call da_free_unit(lnbufr)
150       call da_trace_exit("da_read_bufrtovs")
151       return
152    end if
153 
154    call openbf(lnbufr,'IN',lnbufr)
155    call datelen(10)
156    call readmg(lnbufr,subset,idate,iret)
157    if (subset /= subfgn) then
158       message(1)='The file title does not match the data subset'
159       write(unit=message(2),fmt=*) &
160          'infile=', lnbufr, infile,' subset=', subset, ' subfgn=',subfgn
161       call da_error(__FILE__,__LINE__,message(1:2))
162    end if
163 
164    iy=0
165    im=0
166    idd=0
167    ihh=0
168    write(unit=date,fmt='( i10)') idate
169    read(unit=date,fmt='(i4,3i2)') iy,im,idd,ihh
170    write(unit=stdout,fmt=*) &
171       'Bufr file date is ',iy,im,idd,ihh,infile
172 
173    ! Loop to read bufr file and assign information to a sequential structure
174    !-------------------------------------------------------------------------
175 
176    allocate (head)
177    !  allocate ( head % tb_inv (1:nchan) )
178    nullify  ( head % next )
179    p => head
180 
181    num_tovs_file     = 0    ! number of obs in file
182    num_tovs_global   = 0    ! number of obs within whole domain
183    num_tovs_local    = 0    ! number of obs within tile
184    num_tovs_thinned = 0    ! number of obs rejected by thinning
185    num_tovs_used     = 0    ! number of obs entered into innovation computation
186    num_tovs_selected = 0    ! number of obs limited for debuging
187 
188    if (tovs_start > 1) then
189       write (unit=stdout,fmt='(A,I6)') "   Skipping tovs obs before", tovs_start
190    end if
191 
192    obs: do while (ireadmg(lnbufr,subset,idate)==0 .and. subset==subfgn)
193       do while (ireadsb(lnbufr)==0)
194 
195          ! 1.0     Read header record and data record
196 
197          call ufbint(lnbufr,bfr1bhdr,n1bhdr,1,iret,hdr1b)
198          call ufbrep(lnbufr,data1b8,1,nchan,iret,'TMBR')
199          ! call ufbrep(lnbufr,data1b8,1,1,iret,'BEARAZ')
200 
201          ! check if observation outside range
202 
203          num_tovs_file = num_tovs_file + 1
204 
205          ! 2.0     Extract observation location and other required information
206          !     QC1:  judge if data is in the domain, read next record if not
207          !------------------------------------------------------------------------
208          ! rlat = bfr1bhdr(bufr_lat)
209          ! rlon = bfr1bhdr(bufr_lat)
210          ! if (rlon < 0.0) rlon = rlon+360.0
211 
212          info%lat  =  bfr1bhdr(bufr_lat)
213          info%lon  =  bfr1bhdr(bufr_lon)
214          call da_ll_to_xy (info, loc, xp, outside, outside_all)
215 
216          if (outside_all) cycle
217 
218          !  3.0     Extract other information
219          !------------------------------------------------------
220          !  3.1     Extract satellite id and scan position. 
221    
222          satellite_id = nint(bfr1bhdr(bufr_satellite_id))-191
223          ifov = nint(bfr1bhdr(bufr_ifov))    
224 
225          !  QC2:  limb pixel rejected (not implemented)
226 
227          !  3.2     Extract date information.
228     
229          year   = bfr1bhdr(bufr_year)   
230          month  = bfr1bhdr(bufr_month)  
231          day    = bfr1bhdr(bufr_day)    
232          hour   = bfr1bhdr(bufr_hour)   
233          minute = bfr1bhdr(bufr_minute) 
234          second = bfr1bhdr(bufr_second) 
235 
236          write(unit=info%date_char, fmt='(i4,a,i2.2,a,i2.2,a,i2.2,a,i2.2,a,i2.2)')  &
237             year, '-', month, '-', day, '_', hour, ':', minute, ':', second
238 
239          !  QC3: time consistency check with the background date
240 
241          if (year <= 99) then
242             if (year < 78) then
243                year = year + 2000
244             else
245                year = year + 1900
246             end if
247          end if
248 
249          call da_get_julian_time(year,month,day,hour,minute,obs_time)
250 
251          if (obs_time < time_slots(0) .or.  &
252             obs_time >= time_slots(num_fgat_time)) cycle
253 
254          ! 3.2.1 determine FGAT index ifgat
255    
256          do ifgat=1,num_fgat_time
257             if (obs_time >= time_slots(ifgat-1) .and.  &
258                 obs_time  < time_slots(ifgat)) exit
259          end do
260 
261          ! 3.3 Find wrfvar instrument index from RTTOV instrument triplet
262          !     go to next data if id is not in the lists
263 
264          inst = 0
265          do i = 1, rtminit_nsensor
266             if (platform_id  == rtminit_platform(i) &
267                .and. satellite_id == rtminit_satid(i)    &
268                .and. sensor_id    == rtminit_sensor(i)) then
269                inst = i
270                exit
271             end if
272          end do
273          if (inst == 0) cycle
274 
275          num_tovs_global = num_tovs_global + 1
276 
277          if (num_tovs_global < tovs_start) then
278             cycle
279          end if
280 
281          if (num_tovs_global > tovs_end) then
282             write (unit=stdout,fmt='(A,I6)') "   Skipping radiance obs after", tovs_end
283             exit obs
284          end if
285 
286          num_tovs_selected = num_tovs_selected + 1
287  
288          if (num_tovs_selected > max_tovs_input) then
289             write(unit=message(1),fmt='(A,I10,A)') &
290                "Max number of tovs",max_tovs_input," reached"
291             call da_warning(__FILE__,__LINE__,message(1:1))
292             num_tovs_selected = num_tovs_selected-1
293             exit obs
294          end if
295 
296          if (outside) cycle ! No good for this PE
297          num_tovs_local = num_tovs_local + 1
298 
299          !  Make Thinning
300          !  Map obs to thinning grid
301          !-------------------------------------------------------------------
302          if (thinning) then
303             dlat_earth = info%lat
304             dlon_earth = info%lon
305             if (dlon_earth<zero) dlon_earth = dlon_earth+r360
306             if (dlon_earth>=r360) dlon_earth = dlon_earth-r360
307             dlat_earth = dlat_earth*deg2rad
308             dlon_earth = dlon_earth*deg2rad           
309             timedif = 0. !2.0_r_kind*abs(tdiff)        ! range:  0 to 6
310             terrain = 0.01_r_kind*abs(bfr1bhdr(13))
311             crit = 1. !0.01_r_kind+terrain + timedif !+ 10._r_kind*float(iskip)
312             call map2grids(inst,dlat_earth,dlon_earth,crit,iobs,itx,1,itt,iout,iuse)
313             if (.not. iuse) then
314                num_tovs_thinned=num_tovs_thinned+1
315                cycle
316             end if
317          end if
318 
319          num_tovs_used = num_tovs_used + 1
320          nread(inst) = nread(inst) + 1
321 
322          ! 3.4 extract satellite and solar angle
323    
324          panglr=(start+float(ifov-1)*step)*deg2rad
325          if (hirs2 .or. msu) then
326             satzen = asin(rato*sin(panglr))*rad2deg
327             satzen = abs(satzen)
328          else
329             satzen = bfr1bhdr(bufr_satzen) !*deg2rad   ! local zenith angle
330             satzen = abs(satzen)
331             ! if (amsua .and. ifov .le. 15) satzen=-satzen
332             ! if (amsub .and. ifov .le. 45) satzen=-satzen
333             ! if (hirs3 .and. ifov .le. 28) satzen=-satzen
334          end if
335          satazi = panglr*rad2deg            ! look angle
336          ! if (satazi<0.0) satazi = satazi+360.0
337          solzen = bfr1bhdr(bufr_solzen)              ! solar zenith angle
338 
339          ! 3.5 extract surface information
340    
341          srf_height = bfr1bhdr(bufr_station_height)          ! station height
342          landsea_mask = nint(bfr1bhdr(bufr_landsea_mask))  ! 0:land ; 1:sea (same as RTTOV)
343          ! rmask=one                          ! land
344          ! if (nint(bfr1bhdr(bufr_landsea_mask))==1) rmask=0.0_r_kind   ! reverse the land/sea mask in bufr
345          ! landsea_mask = rmask+.001_r_kind             ! land sea mask
346 
347          info%elv = srf_height
348 
349          ! 3.6 extract channel bright temperature
350    
351          tb_inv(1:nchan) = data1b8(1:nchan)
352          do k = 1, nchan
353             if ( tb_inv(k) < tbmin .or. tb_inv(k) > tbmax) &
354                tb_inv(k) = missing_r
355          end do
356 
357          !  4.0   assign information to sequential radiance structure
358          !--------------------------------------------------------------------------
359          allocate (p % tb_inv (1:nchan))
360          p%info             = info
361          p%loc              = loc
362          p%landsea_mask     = landsea_mask
363          p%scanpos          = ifov
364          p%satzen           = satzen
365          p%satazi           = satazi
366          p%solzen           = solzen
367          p%tb_inv(1:nchan)  = tb_inv(1:nchan)
368          p%sensor_index     = inst
369          p%ifgat            = ifgat
370 
371          allocate (p%next)   ! add next data
372 
373          p => p%next
374          nullify (p%next)
375       end do
376    end do obs
377 
378    iv%total_obs         = iv%total_obs + num_tovs_used
379    iv%total_rad_pixel   = iv%total_rad_pixel + num_tovs_used
380    iv%total_rad_channel = iv%total_rad_channel + num_tovs_used*nchan
381 
382    write(unit=stdout,fmt='(a)') 'num_tovs_file num_tovs_global num_tovs_local num_tovs_used num_tovs_thinned'
383    write(unit=stdout,fmt='(5i10)') num_tovs_file,num_tovs_global, num_tovs_local,num_tovs_used,num_tovs_thinned
384 
385    deallocate(tb_inv)  
386    call closbf(lnbufr)
387    call da_free_unit(lnbufr)
388 
389    !  5.0 allocate innovation radiance structure
390    !----------------------------------------------------------------  
391    
392    do i = 1, iv%num_inst
393       if (nread(i) < 1) cycle
394       iv%instid(i)%num_rad = nread(i)
395       write(UNIT=stdout,FMT='(a,i3,2x,a,3x,i10)') &
396         'Allocating space for radiance innov structure', &
397          i, iv%instid(i)%rttovid_string, iv%instid(i)%num_rad
398       if (rtm_option == rtm_option_crtm) then
399          iv%instid(i)%nlevels = xp%kme-xp%kms+1
400          iv%instid(i)%nchannels=nchan
401       end if
402 
403 #ifdef RTTOV
404       if (rtm_option == rtm_option_rttov) then
405          call rttov_setupchan(1, nchan, coefs(i), &   ! in
406             iv%instid(i)%nfrequencies,iv%instid(i)%nchannels, &
407             iv%instid(i)%nbtout)      ! out
408       end if
409 #endif
410 
411       allocate (iv%instid(i)%info(iv%instid(i)%num_rad))
412       allocate (iv%instid(i)%loc(iv%instid(i)%num_rad))
413       allocate (iv%instid(i)%loc_i(iv%instid(i)%num_rad))
414       allocate (iv%instid(i)%loc_j(iv%instid(i)%num_rad))
415       allocate (iv%instid(i)%loc_k(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
416       allocate (iv%instid(i)%loc_dz(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
417       allocate (iv%instid(i)%loc_dzm(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
418       allocate (iv%instid(i)%zk(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
419       allocate (iv%instid(i)%t(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
420       allocate (iv%instid(i)%mr(iv%instid(i)%nlevels,iv%instid(i)%num_rad))
421       allocate (iv%instid(i)%loc_dx(iv%instid(i)%num_rad))
422       allocate (iv%instid(i)%loc_dy(iv%instid(i)%num_rad))
423       allocate (iv%instid(i)%loc_dxm(iv%instid(i)%num_rad))
424       allocate (iv%instid(i)%loc_dym(iv%instid(i)%num_rad))
425       allocate (iv%instid(i)%tm(xp%kms:xp%kme,iv%instid(i)%num_rad))
426       allocate (iv%instid(i)%qm(xp%kms:xp%kme,iv%instid(i)%num_rad))
427       allocate (iv%instid(i)%qrn(xp%kms:xp%kme,iv%instid(i)%num_rad))
428       allocate (iv%instid(i)%qcw(xp%kms:xp%kme,iv%instid(i)%num_rad))
429       allocate (iv%instid(i)%qci(xp%kms:xp%kme,iv%instid(i)%num_rad))
430       allocate (iv%instid(i)%qsn(xp%kms:xp%kme,iv%instid(i)%num_rad))
431       allocate (iv%instid(i)%qgr(xp%kms:xp%kme,iv%instid(i)%num_rad))
432       allocate (iv%instid(i)%pm(xp%kms:xp%kme,iv%instid(i)%num_rad))
433       allocate (iv%instid(i)%pf(0:xp%kme,iv%instid(i)%num_rad))
434       allocate (iv%instid(i)%u10(iv%instid(i)%num_rad))
435       allocate (iv%instid(i)%v10(iv%instid(i)%num_rad))
436       allocate (iv%instid(i)%t2m(iv%instid(i)%num_rad))
437       allocate (iv%instid(i)%q2m(iv%instid(i)%num_rad))
438       allocate (iv%instid(i)%mr2m(iv%instid(i)%num_rad))
439       allocate (iv%instid(i)%psfc(iv%instid(i)%num_rad))
440       allocate (iv%instid(i)%ts(iv%instid(i)%num_rad))
441       allocate (iv%instid(i)%smois(iv%instid(i)%num_rad))
442       allocate (iv%instid(i)%tslb(iv%instid(i)%num_rad))
443       allocate (iv%instid(i)%snowh(iv%instid(i)%num_rad))
444       allocate (iv%instid(i)%isflg(iv%instid(i)%num_rad))
445       allocate (iv%instid(i)%landsea_mask(iv%instid(i)%num_rad))
446       allocate (iv%instid(i)%elevation(iv%instid(i)%num_rad))
447       allocate (iv%instid(i)%vegfra(iv%instid(i)%num_rad))
448       allocate (iv%instid(i)%vegtyp(iv%instid(i)%num_rad))
449       allocate (iv%instid(i)%soiltyp(iv%instid(i)%num_rad))
450       allocate (iv%instid(i)%clwp(iv%instid(i)%num_rad))
451       allocate (iv%instid(i)%ps(iv%instid(i)%num_rad))
452       allocate (iv%instid(i)%tb_xb(nchan,iv%instid(i)%num_rad))
453       allocate (iv%instid(i)%tb_qc(nchan,iv%instid(i)%num_rad))
454       allocate (iv%instid(i)%tb_inv(nchan,iv%instid(i)%num_rad))
455       allocate (iv%instid(i)%tb_error(nchan,iv%instid(i)%num_rad))
456       allocate (iv%instid(i)%emiss(iv%instid(i)%nchannels,iv%instid(i)%num_rad))
457       allocate(iv%instid(i)%scanpos(iv%instid(i)%num_rad))
458       allocate(iv%instid(i)%scanline(iv%instid(i)%num_rad))
459       allocate(iv%instid(i)%ifgat(iv%instid(i)%num_rad))
460       allocate(iv%instid(i)%cloud_flag(nchan,iv%instid(i)%num_rad))
461       allocate(iv%instid(i)%satzen(iv%instid(i)%num_rad))
462       allocate(iv%instid(i)%satazi(iv%instid(i)%num_rad))
463       allocate(iv%instid(i)%solzen(iv%instid(i)%num_rad))
464       allocate(iv%instid(i)%solazi(iv%instid(i)%num_rad))
465       allocate(iv%instid(i)%proc_domain(iv%instid(i)%num_rad))
466       if (rtm_option == rtm_option_crtm) then
467          allocate(iv%instid(i)%water_coverage(iv%instid(i)%num_rad))
468          allocate(iv%instid(i)%land_coverage(iv%instid(i)%num_rad))
469          allocate(iv%instid(i)%ice_coverage(iv%instid(i)%num_rad))
470          allocate(iv%instid(i)%snow_coverage(iv%instid(i)%num_rad))
471          !allocate(iv%instid(i)%ps_jacobian(nchan,iv%instid(i)%num_rad))
472          !allocate(iv%instid(i)%t_jacobian(nchan,xp%kte,iv%instid(i)%num_rad))
473          !allocate(iv%instid(i)%q_jacobian(nchan,xp%kte,iv%instid(i)%num_rad))
474       end if
475    end do
476    
477    !  6.0 assign sequential structure to innovation structure
478    !-------------------------------------------------------------
479    nread(1:rtminit_nsensor) = 0 
480    p => head
481    ! do while ( associated(p) )
482 
483    do n = 1, num_tovs_used
484       i = p%sensor_index
485       nread(i) = nread(i) + 1
486       iv%instid(i)%info(nread(i))         = p%info
487       iv%instid(i)%loc(nread(i))          = p%loc
488       iv%instid(i)%loc_i(nread(i))        = p%loc%i
489       iv%instid(i)%loc_j(nread(i))        = p%loc%j
490       iv%instid(i)%loc_k(:,nread(i))      = 0
491       iv%instid(i)%loc_dx(nread(i))       = p%loc%dx
492       iv%instid(i)%loc_dy(nread(i))       = p%loc%dy
493       iv%instid(i)%loc_dz(:,nread(i))     = 0.0
494       iv%instid(i)%loc_dxm(nread(i))      = p%loc%dxm
495       iv%instid(i)%loc_dym(nread(i))      = p%loc%dym
496       iv%instid(i)%loc_dzm(:,nread(i))    = 0.0
497       ! z done in da_get_innov_vector_rad
498       iv%instid(i)%t(:,nread(i))          = 0.0
499       iv%instid(i)%mr(:,nread(i))         = 0.0
500       iv%instid(i)%tm(:,nread(i))         = 0.0
501       iv%instid(i)%qm(:,nread(i))         = 0.0
502       iv%instid(i)%qrn(:,nread(i))        = 0.0
503       iv%instid(i)%qcw(:,nread(i))        = 0.0
504       iv%instid(i)%qci(:,nread(i))        = 0.0
505       iv%instid(i)%qsn(:,nread(i))        = 0.0
506       iv%instid(i)%qgr(:,nread(i))        = 0.0
507       iv%instid(i)%pm(:,nread(i))         = 0.0
508       iv%instid(i)%u10(nread(i))          = 0.0
509       iv%instid(i)%v10(nread(i))          = 0.0
510       iv%instid(i)%t2m(nread(i))          = 0.0
511       iv%instid(i)%q2m(nread(i))          = 0.0
512       iv%instid(i)%mr2m(nread(i))         = 0.0
513       iv%instid(i)%psfc(nread(i))         = 0.0
514       iv%instid(i)%ts(nread(i))           = 0.0
515       iv%instid(i)%smois(nread(i))        = 0.0
516       iv%instid(i)%tslb(nread(i))         = 0.0
517       iv%instid(i)%snowh(nread(i))        = 0.0
518       iv%instid(i)%isflg(nread(i))        = 0
519       iv%instid(i)%soiltyp(nread(i))      = 0.0
520       iv%instid(i)%landsea_mask(nread(i)) = p%landsea_mask
521       iv%instid(i)%elevation(nread(i))    = 0.0
522       iv%instid(i)%vegfra(nread(i))       = 0.0
523       iv%instid(i)%vegtyp(nread(i))       = 0.0
524       iv%instid(i)%clwp(nread(i))         = 0.0
525       iv%instid(i)%ps(nread(i))           = 0.0
526       iv%instid(i)%tb_xb(:,nread(i))      = 0
527       iv%instid(i)%tb_inv(:,nread(i))     = p%tb_inv(:) 
528       iv%instid(i)%tb_qc(:,nread(i))      = 0
529       iv%instid(i)%tb_error(:,nread(i))   = 500.
530       iv%instid(i)%emiss(:,nread(i))      = 0.0
531       iv%instid(i)%scanpos(nread(i))      = p%scanpos 
532       ! iv%instid(i)%scanline(nread(i)    = p%scanline
533       iv%instid(i)%scanline(nread(i))     = 0
534       iv%instid(i)%ifgat(nread(i))        = p%ifgat
535       iv%instid(i)%cloud_flag(:,nread(i)) = 1        ! no cloud
536       iv%instid(i)%satzen(nread(i))       = p%satzen  
537       iv%instid(i)%satazi(nread(i))       = p%satazi  
538       iv%instid(i)%solzen(nread(i))       = p%solzen  
539       ! iv%instid(i)%solazi(nread(i))     = p%solazi  
540       iv%instid(i)%solazi(nread(i))       = 0.0
541       iv%instid(i)%proc_domain(nread(i))  = .false.
542 
543       current => p
544       p => p%next
545 
546       ! free current data
547       deallocate ( current % tb_inv )
548       deallocate ( current )
549    end do
550 
551    deallocate ( p )
552 
553    deallocate (nread)
554 
555    ! check if sequential structure has been freed
556    !
557    ! p => head
558    ! do i = 1, num_rad_selected
559    !    write (unit=stdout,fmt=*)  i, p%tb_inv(1:nchan)
560    !    p => p%next
561    ! end do
562 
563    call da_trace_exit("da_read_bufrtovs")
564 #endif
565 
566 end subroutine da_read_bufrtovs
567 
568