da_biasprep.inc
References to this file elsewhere.
1 subroutine da_biasprep(inst,ob,iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: Output information files for bias correction progs
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer , intent(in) :: inst
10 type (y_type) , intent(in) :: ob ! O structure.
11 type (iv_type), intent(in) :: iv ! O-B structure.
12
13 integer :: n,jx,npred,nchan,num_rad,nlevels
14 character(len=80) :: filename
15 character(len=1) :: s1
16 real :: pred(6), q(43)
17 type (bias_type) :: radbias
18
19 num_rad = iv%instid(inst)%info%n2-iv%instid(inst)%info%n1+1
20
21 if (num_rad < 1) return
22
23 if (trace_use) call da_trace_entry("da_biasprep")
24
25 #ifdef DM_PARALLEL
26 write(filename, '(a,i2.2)') 'biasprep_'//trim(iv%instid(inst)%rttovid_string)//'.', myproc
27 #else
28 write(filename, '(a)') 'biasprep_'//trim(iv%instid(inst)%rttovid_string)
29 #endif
30
31 call da_get_unit(biasprep_unit)
32 open(unit=biasprep_unit,FILE=filename,FORM='unformatted')
33
34 !---------------------------------------------------------------------------
35 npred = 4
36 nchan = iv%instid(inst)%nchan
37 nlevels = iv%instid(inst)%nlevels-1
38
39 allocate(radbias%tb(nchan))
40 allocate(radbias%omb(nchan))
41 allocate(radbias%bias(nchan))
42 allocate(radbias%qc_flag(nchan))
43 allocate(radbias%cloud_flag(nchan))
44 allocate(radbias%pred(npred))
45
46 do n=iv%instid(inst)%info%n1,iv%instid(inst)%info%n2
47 if (iv%instid(inst)%info%proc_domain(1,n)) then
48
49 if (rtm_option==rtm_option_rttov) then
50 q(1:43) = iv%instid(inst)%mr(1:43,n)/q2ppmv
51 call da_predictor_rttov(pred(1:npred), npred, iv%instid(inst)%t(1:43,n), &
52 q(1:43), iv%instid(inst)%ts(n))
53 else if (rtm_option==rtm_option_crtm) then
54 #ifdef CRTM
55 call da_predictor_crtm(pred(1:npred), npred, nlevels,iv%instid(inst)%tm(1:nlevels,n), &
56 iv%instid(inst)%qm(1:nlevels,n), iv%instid(inst)%ts(n), &
57 iv%instid(inst)%pm(1:nlevels,n),iv%instid(inst)%pf(0:nlevels,n))
58 #endif
59 end if
60
61 ! transfer information to bias structure
62 radbias%platform_id = iv%instid(inst)%platform_id
63 radbias%satellite_id = iv%instid(inst)%satellite_id
64 radbias%sensor_id = iv%instid(inst)%sensor_id
65
66 read(iv%instid(inst)%info%date_char(n),'(i4,5(a1,i2))') &
67 radbias%year,s1, radbias%month,s1, radbias%day, &
68 s1,radbias%hour, s1,radbias%min, s1,radbias%sec
69
70 radbias%scanline = iv%instid(inst)%scanline(n) ! not available
71 radbias%scanpos = iv%instid(inst)%scanpos(n)
72 radbias%landmask = iv%instid(inst)%landsea_mask(n)
73 radbias%elevation = iv%instid(inst)%info%elv(n)
74 radbias%lat = iv%instid(inst)%info%lat(1,n)
75 radbias%lon = iv%instid(inst)%info%lon(1,n)
76 radbias%surf_flag = iv%instid(inst)%isflg(n)
77 radbias%ps = iv%instid(inst)%ps(n)
78 radbias%t2m = iv%instid(inst)%t2m(n)
79 radbias%q2m = iv%instid(inst)%mr2m(n)/q2ppmv
80 radbias%tsk = iv%instid(inst)%ts(n)
81 radbias%clwp = iv%instid(inst)%clwp(n) ! in mm
82
83 radbias%nchan = nchan
84 radbias%tb(1:nchan) = ob%instid(inst)%tb(1:nchan,n)
85 radbias%omb(1:nchan) = iv%instid(inst)%tb_inv(1:nchan,n)
86 radbias%bias(1:nchan) = 0.0
87
88 radbias%npred = npred
89 radbias%pred(1:npred) = pred(1:npred)
90
91 radbias%qc_flag(1:nchan)= iv%instid(inst)%tb_qc(1:nchan,n)
92 radbias%cloud_flag(1:nchan)= iv%instid(inst)%cloud_flag(1:nchan,n)
93
94 ! set missing data and bad data to missing
95 do jx=1,nchan
96 if (radbias%tb(jx) < 150.0 .or. radbias%tb(jx) > 400.0 ) then
97 radbias%tb(jx) = missing_r
98 radbias%omb(jx) = missing_r
99 end if
100 end do
101
102 !write(unit=biasprep_unit) radbias ! can not compiled with pointer
103
104 call da_write_biasprep(radbias)
105
106 end if
107 end do
108
109 close(unit=biasprep_unit)
110 call da_free_unit(biasprep_unit)
111
112 deallocate(radbias%tb)
113 deallocate(radbias%omb)
114 deallocate(radbias%bias)
115 deallocate(radbias%qc_flag)
116 deallocate(radbias%cloud_flag)
117 deallocate(radbias%pred)
118
119 if (trace_use) call da_trace_exit("da_biasprep")
120
121 end subroutine da_biasprep
122
123