var/da/da_radiance/da_biasprep.inc

   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)
  17    type (bias_type)   :: radbias
  18    real,allocatable      :: q(:), temp(:), hum(:), pf(:)
  19
  20    num_rad = iv%instid(inst)%info%n2-iv%instid(inst)%info%n1+1
  21
  22    if (num_rad < 1) return
  23
  24    if (trace_use) call da_trace_entry("da_biasprep")
  25
  26 #ifdef DM_PARALLEL
  27    write(filename, '(a,i4.4)') 'biasprep_'//trim(iv%instid(inst)%rttovid_string)//'.', myproc
  28 #else
  29    write(filename, '(a)') 'biasprep_'//trim(iv%instid(inst)%rttovid_string)
  30 #endif
  31
  32    call da_get_unit(biasprep_unit)
  33    open(unit=biasprep_unit,FILE=filename,FORM='unformatted')
  34
  35    !---------------------------------------------------------------------------
  36    npred = 4
  37    nchan = iv%instid(inst)%nchan
  38    nlevels = iv%instid(inst)%nlevels-1
  39    allocate(q(nlevels))
  40    allocate(temp(nlevels))
  41    allocate(hum(nlevels))
  42    allocate(pf(0:nlevels))
  43
  44    allocate(radbias%tb(nchan))
  45    allocate(radbias%omb(nchan))
  46    allocate(radbias%bias(nchan))
  47    allocate(radbias%qc_flag(nchan))
  48    allocate(radbias%cloud_flag(nchan))
  49    allocate(radbias%pred(npred))
  50
  51    do n=iv%instid(inst)%info%n1,iv%instid(inst)%info%n2
  52       if (iv%instid(inst)%info%proc_domain(1,n)) then
  53
  54         if (rtm_option==rtm_option_rttov) then
  55 #ifdef RTTOV
  56          q(1:nlevels) = iv%instid(inst)%mr(1:nlevels,n)/q2ppmv
  57          call da_predictor_rttov(inst, pred(1:npred), npred, nlevels, &
  58             iv%instid(inst)%t(1:nlevels,n), q(1:nlevels), iv%instid(inst)%ts(n))
  59 #endif
  60         else if (rtm_option==rtm_option_crtm) then
  61 #ifdef CRTM
  62 ! FIX? problems with IBM AIX COMPILER
  63          temp(1:nlevels) = iv%instid(inst)%tm(1:nlevels,n)
  64          hum(1:nlevels) = iv%instid(inst)%qm(1:nlevels,n)
  65          pf(0:nlevels) = iv%instid(inst)%pf(0:nlevels,n)
  66          call da_predictor_crtm(pred(1:npred), npred, nlevels,temp, &
  67             hum, iv%instid(inst)%ts(n), pf)
  68 #endif
  69         end if
  70
  71          ! transfer information to bias structure
  72          radbias%platform_id  = iv%instid(inst)%platform_id
  73          radbias%satellite_id = iv%instid(inst)%satellite_id
  74          radbias%sensor_id    = iv%instid(inst)%sensor_id
  75
  76          read(iv%instid(inst)%info%date_char(n),'(i4,5(a1,i2))') &
  77                                    radbias%year,s1, radbias%month,s1, radbias%day, &
  78                                    s1,radbias%hour, s1,radbias%min, s1,radbias%sec
  79
  80          radbias%scanline     = iv%instid(inst)%scanline(n)    ! not available
  81          radbias%scanpos      = iv%instid(inst)%scanpos(n)
  82          radbias%landmask     = iv%instid(inst)%landsea_mask(n)
  83          radbias%elevation    = iv%instid(inst)%info%elv(n)
  84          radbias%lat          = iv%instid(inst)%info%lat(1,n)
  85          radbias%lon          = iv%instid(inst)%info%lon(1,n)
  86          radbias%surf_flag    = iv%instid(inst)%isflg(n)
  87          radbias%ps           = iv%instid(inst)%ps(n)
  88          radbias%t2m          = iv%instid(inst)%t2m(n)
  89          radbias%q2m          = iv%instid(inst)%mr2m(n)/q2ppmv
  90          radbias%tsk          = iv%instid(inst)%ts(n)
  91          radbias%clwp         = iv%instid(inst)%clwp(n)  ! in mm
  92
  93          radbias%nchan        = nchan
  94          radbias%tb(1:nchan)  = ob%instid(inst)%tb(1:nchan,n)
  95          radbias%omb(1:nchan) = ob%instid(inst)%tb(1:nchan,n)-iv%instid(inst)%tb_xb(1:nchan,n)
  96          radbias%bias(1:nchan) = 0.0
  97
  98          radbias%npred         = npred
  99          radbias%pred(1:npred) = pred(1:npred)
 100
 101          radbias%qc_flag(1:nchan)= iv%instid(inst)%tb_qc(1:nchan,n)
 102          radbias%cloud_flag(1:nchan)= iv%instid(inst)%cloud_flag(1:nchan,n)
 103
 104          ! set missing data and bad data to missing
 105          do jx=1,nchan
 106             if (radbias%tb(jx) < 150.0 .or. radbias%tb(jx) > 400.0 ) then
 107                radbias%tb(jx)   = missing_r
 108                radbias%omb(jx)  = missing_r
 109             end if
 110          end do
 111
 112          !write(unit=biasprep_unit) radbias ! can not compiled with pointer
 113
 114          call da_write_biasprep(radbias)
 115
 116       end if
 117    end do
 118
 119    close(unit=biasprep_unit)
 120    call da_free_unit(biasprep_unit)
 121
 122    deallocate(q)
 123    deallocate(temp)
 124    deallocate(hum)
 125    deallocate(pf)
 126
 127    deallocate(radbias%tb)
 128    deallocate(radbias%omb)
 129    deallocate(radbias%bias)
 130    deallocate(radbias%qc_flag)
 131    deallocate(radbias%cloud_flag)
 132    deallocate(radbias%pred)
 133
 134    if (trace_use) call da_trace_exit("da_biasprep")
 135
 136 end subroutine da_biasprep
 137
 138