da_read_pseudo_rad.inc
References to this file elsewhere.
1 subroutine da_read_pseudo_rad (iv)
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 type (iv_type), intent (inout) :: iv
19
20 ! Instrument triplet, follow the convension of RTTOV
21 integer :: platform_id, satellite_id, sensor_id
22
23 real, allocatable :: tb_inv(:) ! bright temperatures
24 type (datalink_type) :: p
25
26 logical :: outside, outside_all
27 integer :: i,k,n,nchan,inst, alloc_stat
28 type(info_type) :: info
29 type(model_loc_type) :: loc
30
31 call da_trace_entry("da_read_pseudo_rad")
32
33 ! Initialize variables
34
35 platform_id = pseudo_rad_platid
36 satellite_id = pseudo_rad_satid
37 sensor_id = pseudo_rad_senid
38 if (sensor_id == 0) then
39 nchan=19 !nchan_hirs
40 else if (sensor_id == 1) then
41 nchan=nchan_msu
42 else if (sensor_id == 3) then
43 nchan=nchan_amsua
44 else if (sensor_id == 4) then
45 nchan=nchan_amsub
46 else if (sensor_id == 15) then
47 nchan=nchan_mhs
48 else if (sensor_id == 10) then
49 nchan=nchan_ssmis
50 end if
51
52 inst = 1 ! single instrument
53
54 allocate (tb_inv(nchan))
55
56 info%lat = pseudo_rad_lat ! in degree
57 info%lon = pseudo_rad_lon
58 info%date_char = "0000-00-00_00:00:00"
59 call da_llxy (info, loc, outside, outside_all)
60
61 if (outside) then
62 iv%info(radiance)%nlocal = 0
63 iv%info(radiance)%ntotal = 0
64 iv%info(radiance)%ptotal(1) = 0
65 iv%instid(inst)%num_rad = 0
66 iv%instid(inst)%info%nlocal = 0
67 else
68 iv%info(radiance)%nlocal = 1
69 iv%info(radiance)%ntotal = 1
70 iv%info(radiance)%ptotal(1) = 1
71 iv%instid(inst)%num_rad = 1
72 iv%instid(inst)%info%nlocal = 1
73 end if
74
75 do k = 1, nchan
76 tb_inv(k) = missing_r
77 end do
78 tb_inv(pseudo_rad_ichan)=pseudo_rad_tb
79
80 ! 5.0 allocate innovation radiance structure
81 !----------------------------------------------------------------
82
83 if ( iv%instid(inst)%num_rad > 0 ) then
84 if (rtm_option == rtm_option_crtm) then
85 iv%instid(inst)%nlevels = kme-kms+1
86 iv%instid(inst)%nchannels=nchan
87 end if
88 write(UNIT=stdout,FMT='(a,i3,2x,a,3x,i10)') &
89 'Allocating space for radiance innov structure', &
90 i, iv%instid(inst)%rttovid_string, iv%instid(inst)%num_rad
91
92 #ifdef RTTOV
93 if (rtm_option == rtm_option_rttov) then
94 call rttov_setupchan(1, nchan, coefs(i), & ! in
95 iv%instid(inst)%nfrequencies,iv%instid(inst)%nchannels, &
96 iv%instid(inst)%nbtout) ! out
97 end if
98 #endif
99
100 call da_allocate_rad_iv(i,nchan,iv)
101
102 ! 6.0 assign sequential structure to innovation structure
103 !-------------------------------------------------------------
104
105 n = 1 ! single obs
106
107 allocate (p%tb_inv(1:nchan), stat=alloc_stat)
108 if ( alloc_stat /= 0 ) CALL da_error(__FILE__,__LINE__,(/"error allocating"/))
109
110 p%info = info
111 p%loc = loc
112 p%landsea_mask = 1
113 p%scanpos = 1
114 p%satzen = 0.0
115 p%satazi = 0.0
116 p%solzen = 0.0
117 p%tb_inv(1:nchan) = tb_inv(1:nchan)
118 p%sensor_index = inst
119 p%ifgat = 1
120
121 call da_initialize_rad_iv (inst, n, iv, p)
122
123 iv%instid(inst)%tb_qc(:,n) = qc_bad
124 iv%instid(inst)%tb_error(pseudo_rad_ichan,n) = pseudo_rad_err
125 iv%instid(inst)%tb_qc(pseudo_rad_ichan,n) = qc_good
126
127 deallocate(p%tb_inv)
128 end if
129
130 deallocate(tb_inv)
131
132 call da_trace_exit("da_read_pseudo_rad")
133
134 end subroutine da_read_pseudo_rad
135