da_bias_verif.f90

References to this file elsewhere.
1   PROGRAM da_bias_verif
2 
3   USE RAD_BIAS
4 
5 ! PURPOSE.
6 ! --------
7 ! Apply bias correction to radiance data
8 ! to statistically verify algorithm
9 
10   IMPLICIT NONE
11 
12   TYPE(BIAS) :: tovs
13 
14   REAL(KIND=LONG) :: vmean_dep(JPCHAN), vmean_abs(JPCHAN)
15   REAL(KIND=LONG) :: vec_dep(JPCHAN), vec_abs(JPCHAN)
16   REAL(KIND=LONG) :: vstd_dep(JPCHAN), vstd_abs(JPCHAN)
17   REAL(KIND=LONG) :: vmn(JPCHAN), vstd(JPCHAN), airbias(JPCHAN)
18 
19   REAL :: pred(JPNX)
20 
21   REAL(KIND=LONG) :: xbar(JPNY,JPNX), ybar(JPNY)
22   REAL(KIND=LONG) :: xcov(JPNY,JPNX,JPNX), ycov(JPNY), xycov(JPNY,JPNX)
23 
24   REAL(KIND=LONG) :: coef(JPNY,JPNX), coef0(JPNY)
25   REAL(KIND=LONG) :: reserr(JPNY), rescov(JPNY,JPNY)
26   REAL(KIND=LONG) :: xvec(JPNY,JPNX,JPNX)
27   REAL(KIND=LONG) :: xreser(JPCHAN), xcoef0(JPCHAN), xcoef(JPCHAN,JPNX)
28 
29   REAL(KIND=LONG) :: vmnbd(JPNY,6), vstdbd(JPNY,6)
30   REAL(KIND=LONG) :: vmnbdl(JPNY,6,0:1), vstdbdl(JPNY,6,0:1)
31   REAL(KIND=LONG) :: vmnb0(JPNY,6), vstdb0(JPNY,6)
32   REAL(KIND=LONG) :: vmnbl0(JPNY,6,0:1), vstdbl0(JPNY,6,0:1)
33   REAL(KIND=LONG) :: dvmnbd(JPNY,6), dvstdb(JPNY,6)
34 
35   REAL(KIND=LONG) :: vmnrl(JPCHAN,JPSCAN) = 0.0
36   REAL(KIND=LONG) :: vmnrlb(JPCHAN,JPSCAN,JBAND) = 0.0
37 
38   REAL(KIND=LONG) :: SCORR(JPCHAN)
39 
40   INTEGER :: nsel(10)
41   INTEGER :: nobs(JPCHAN)
42   INTEGER :: nobsy(JPNY),lchan
43 
44   LOGICAL :: LMASK
45 
46   LOGICAL :: lscan = .FALSE.
47 
48   INTEGER :: IR, ibin, I, iband, II, IV, ib, ierr
49   INTEGER :: JS, J, JJ, JCOUNT, JBOX, JMINI, JSCAN, jv, IIV, JJV, sband
50 
51   REAL(KIND=LONG) :: xcorr(JPCHAN)
52   REAL(KIND=LONG) :: coef_year, coef_month, coef_day, coef_time
53 
54   INTEGER :: icount = 0, kscanx, jbandx
55   LOGICAL :: check_limb=.false., check_mask=.false., global
56   REAL    :: FAC = 3.0      ! Number of SD' for QC
57 
58   INTEGER :: nchan,nscan,nband,npred
59 
60   NAMELIST /INPUTS/ global,lscan, check_limb, check_mask, FAC
61 
62 !------------------------------------------------------------------------------
63 !        1.   SETUP.
64 !             -----
65   READ(5,INPUTS,END=100)
66   100 CONTINUE
67   WRITE(6,INPUTS)
68 
69 ! 1.0 read bias correction coefficients
70 !   -------------------------------------------------------------
71      open(UNIT=12,file='bcor.asc',form='formatted')
72      read (12,'(4i6)') nchan,nscan,nband,npred
73      close(12)
74 
75      call read_biascoef(nchan,nscan,nband,npred,global, &
76                       VMNRL(1:nchan,1:nscan),VMNRLB(1:nchan,1:nscan,1:nband), &
77                       xcoef(1:nchan,1:npred),xcoef0(1:nchan), &
78                       nobs(1:nchan),vmean_abs(1:nchan), &
79                       vstd_abs(1:nchan),vmean_dep(1:nchan), &
80                       vstd_dep(1:nchan) )
81 
82 !----------------------------------------------------------------------
83 !        2.   READ IN DATA, Q.C., CALC MEANS AND VARIANCES.
84 !             ---- -- ----- ----- ---- ----- --- ---------
85 
86   READ(UNIT=10,END=300)  tovs%nchan, tovs%npred    ! Read in data
87   REWIND(UNIT=10)
88 
89   allocate(tovs%tb(tovs%nchan))
90   allocate(tovs%omb(tovs%nchan))
91   allocate(tovs%bias(tovs%nchan))
92   allocate(tovs%qc_flag(tovs%nchan))
93   allocate(tovs%cloud_flag(tovs%nchan))
94   allocate(tovs%pred(tovs%npred))
95 
96   300 CONTINUE
97 
98   icount = 0
99 loop2:&
100   DO
101     call da_read_biasprep(tovs,10,ierr)
102     if (ierr == 0) then      ! not end
103          continue
104     elseif (ierr == 1) then  ! end
105          exit
106     else                     ! error
107          stop 'read error in da_veri_bias'
108     endif
109 
110 ! latitude band
111     iband = INT(tovs%lat/30.000001) + 3
112     IF (lscan) THEN
113         JSCAN = tovs%scanpos
114       if (global) then
115         CALL GET_SCORR(JPCHAN,SCORR(1:JPCHAN),tovs%lat,vmnrlb,JSCAN)
116       else
117         SCORR(1:nchan) = vmnrl(1:nchan,JSCAN)
118       end if
119     ELSE
120         JSCAN = KSCANX/2
121         SCORR(1:JPCHAN) = 0.0
122     ENDIF
123 
124     vec_dep(1:nchan) = tovs%omb(1:nchan) - SCORR(1:nchan)
125     vec_abs(1:nchan) =  tovs%tb(1:nchan) - SCORR(1:nchan)
126 
127 ! 3.2 QC:
128       if (tovs%sensor_id == 3) then
129            call qc_amsua(tovs)
130       elseif(tovs%sensor_id == 4) then
131            call qc_amsub(tovs)
132       endif
133 
134 ! 3.3 limb scan check
135 !---------------------
136    if (check_limb) then
137     IF ((tovs%scanpos <= 2) .OR. (tovs%scanpos >= (KSCANX-1))) THEN
138       nsel(2) = nsel(2) + 1
139       CYCLE loop2
140     ENDIF
141    end if
142 
143 ! 3.4 Reject data outside radiosonde mask
144 !------------------------------------------
145     IF (check_mask) THEN
146       CALL MASK(tovs%lat,tovs%lon,LMASK)
147       IF (.NOT. LMASK) THEN
148         nsel(8) = nsel(8) + 1
149         CYCLE loop2
150       ENDIF
151     ENDIF
152 
153 ! 3.5 Reject outliers : facx*sigma, sigma calculated in first pass : loop1
154 !--------------------------------------------------------------------------
155     DO j=1, tovs%nchan
156       IF ( (ABS(vec_dep(j)-vmn(j)) > (vstd(j)*FAC)) ) THEN
157         tovs%qc_flag(j) = -1
158       ENDIF
159     ENDDO
160 
161     do i=1,nchan 
162       airbias(i) = xcoef0(i)
163       do j=1,npred
164         airbias(i) = airbias(i) + tovs%pred(j)*xcoef(i,j)
165       end do
166       vec_dep(i) = vec_dep(i) - airbias(i)
167       vec_abs(i) = vec_abs(i) - airbias(i)
168       if ( tovs%omb(i) == -888888.00 ) vec_dep(i)=tovs%omb(i)
169     end do
170 
171     icount = icount + 1
172     if (nchan == 15) then
173       write(11,'(15i15)') tovs%qc_flag(1:nchan)
174       write(11,'(15f15.3)') tovs%omb(1:nchan)  ! omb no-biascorrection
175       write(11,'(15f15.3)') vec_dep(1:nchan)   ! omb with bias correction
176     else if(nchan==5) then
177       write(11,'(5i15)') tovs%qc_flag(1:nchan)
178       write(11,'(5f15.3)') tovs%omb(1:nchan)  ! omb no-biascorrection
179       write(11,'(5f15.3)') vec_dep(1:nchan)   ! omb with bias correction
180     else if(nchan==19) then
181       write(11,'(19i15)') tovs%qc_flag(1:nchan)
182       write(11,'(19f15.3)') tovs%omb(1:nchan)  ! omb no-biascorrection
183       write(11,'(19f15.3)') vec_dep(1:nchan)   ! omb with bias correction
184     else
185       write(*,*) 'Unknow sensor'
186     end if
187 
188   ENDDO loop2
189 
190   365 CONTINUE
191 
192 !----------------------------------------------------------------------------
193 
194    deallocate(tovs%tb)
195    deallocate(tovs%omb)
196    deallocate(tovs%qc_flag)
197    deallocate(tovs%cloud_flag)
198    deallocate(tovs%pred)
199 
200   END PROGRAM da_bias_verif