var/da/da_ssmi/da_tb_adj.inc

   1 subroutine da_tb_adj(ifreq,theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld,            &
   2 !                  tbv,tbh,                                                  &
   3                   ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,                   &
   4                   ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh            )
   5
   6    !-----------------------------------------------------------------------
   7    ! Purpose: TBD
   8    ! Output : ADJ_p0,  ADJ_ta,   ADJ_gamma, ADJ_sst, ADJ_wv, ADJ_hwv, ADJ_u
   9    !          ADJ_alw, ADJ_zcld
  10    ! Input  : ADJ_tbv, ADJ_tbh,  tbv,  tbh
  11    !-----------------------------------------------------------------------
  12
  13    implicit none
  14
  15    integer, intent(in   ) :: ifreq
  16    real   , intent(in   ) :: theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld
  17    real   , intent(inout) :: ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,     &
  18                              ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld
  19    real   , intent(in   ) :: ADJ_tbv,ADJ_tbh
  20 !   real   , intent(in   ) :: tbv,tbh
  21
  22    real :: freq(4),ebiasv(4),ebiash(4),cf1(4),cf2(4),cg(4)
  23
  24    real :: f,costhet,gx2,tbup,tbdn,tauatm,sigma,remv,remh
  25    real :: effangv,effangh,tbdnv,foam,foam_save,emissv,emissh
  26    real :: dum
  27    real :: refv,refh,semv,semh,scatv,scath,tbdnh
  28    real :: ADJ_gx2,ADJ_tbup,ADJ_tbdn,ADJ_tauatm,ADJ_sigma
  29    real :: ADJ_tremv,ADJ_remh,ADJ_effangv,ADJ_effangh
  30    real :: ADJ_tbdnh,ADJ_dum,ADJ_foam,ADJ_emissv
  31    real :: ADJ_emissh,ADJ_refv,ADJ_refh,ADJ_semv
  32    real :: ADJ_semh,ADJ_scatv,ADJ_scath,ADJ_remv,ADJ_tbdnv
  33    real :: ADJ_theta
  34
  35    real :: fem
  36
  37    data freq/19.35,22.235,37.0,85.5/
  38
  39    ! empirical bias corrections for surface emissivity
  40
  41    data ebiasv/0.0095, 0.005,-0.014, -0.0010/
  42    data ebiash/0.004,   0.0,-0.023, 0.043/
  43
  44
  45    data cf1/0.0015,0.004,0.0027,0.005/
  46    data cg/4.50e-3, 5.200e-3, 5.5e-3, 7.2e-3 /
  47
  48    data cf2/0.0023,0.000,0.0002,-0.006/
  49
  50    ! 'foam' emissivity
  51    data fem /1.0/
  52
  53    if (trace_use) call da_trace_entry("da_tb_adj")
  54
  55    f=0.0;costhet=0.0;gx2=0.0;tbup=0.0;tbdn=0.0;tauatm=0.0
  56    sigma=0.0;remv=0.0;remh=0.0;effangv=0.0;effangh=0.0
  57    tbdnv=0.0;foam=0.0;foam_save=0.0;emissv=0.0;emissh=0.0
  58    dum=0.0;refv=0.0;refh=0.0;semv=0.0;semh=0.0;scatv=0.0
  59    scath=0.0;tbdnh=0.0;ADJ_gx2=0.0;ADJ_tbup=0.0;ADJ_tbdn=0.0
  60    ADJ_tauatm=0.0;ADJ_sigma=0.0;ADJ_tremv=0.0;ADJ_remh=0.0
  61    ADJ_effangv=0.0;ADJ_effangh=0.0;ADJ_tbdnh=0.0;ADJ_dum=0.0
  62    ADJ_foam=0.0;ADJ_emissv=0.0;ADJ_emissh=0.0;ADJ_refv=0.0
  63    ADJ_refh=0.0;ADJ_semv=0.0;ADJ_semh=0.0;ADJ_scatv=0.0
  64    ADJ_scath=0.0;ADJ_remv=0.0;ADJ_tbdnv=0.0
  65    ADJ_theta =0.0
  66
  67
  68    ! write (unit=stdout,fmt=*) 'ifreq',ifreq,theta,p0,ta,gamma,sst,wv,hwv,u,alw,zcld,          &
  69    !              tbv,tbh,                                                  &
  70    !              ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,                   &
  71    !               ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
  72
  73    f = freq(ifreq)
  74    costhet = cos(theta*0.017453)
  75
  76    ! effective surface slope variance
  77
  78    gx2 = cg(ifreq)*    u
  79
  80    ! get upwelling atmospheric brightness temperature
  81
  82    call tbatmos(ifreq,theta,p0,wv,hwv,ta,gamma,alw,zcld, tbup,tbdn,tauatm)
  83
  84    ! convert transmittance to optical depth
  85
  86    sigma = -alog(tauatm)*costhet
  87
  88    ! get rough surface emissivity
  89
  90    call roughem(ifreq,gx2,sst,theta,remv,remh)
  91
  92    ! get effective zenith angles for scattered radiation at surface
  93
  94    call effang(ifreq,theta,gx2,sigma,effangv,effangh)
  95
  96    ! get effective sky brightness temperatures for scattered radiation
  97
  98    call tbatmos(ifreq,effangv,p0,wv,hwv,ta,gamma,alw,zcld, dum,tbdnv,dum)
  99
 100    call tbatmos(ifreq,effangh,p0,wv,hwv,ta,gamma,alw,zcld, dum,tbdnh,dum)
 101
 102    ! compute 'foam' coverage
 103
 104    foam = cf1(ifreq)*    u
 105
 106    if (u .gt. 5.0) then
 107       foam_save = foam
 108       foam =     foam + cf2(ifreq)*(   u-5.0)
 109    end if
 110
 111    ! compute surface emissivities and reflectivity
 112
 113    emissv =     foam*fem + (1.0 - foam)*(remv + ebiasv(ifreq))
 114    emissh =     foam*fem + (1.0 - foam)*(remh + ebiash(ifreq))
 115    refv =   1.0 - emissv
 116    refh =   1.0 - emissh
 117
 118    ! compute surface emission term
 119
 120    semv = sst*emissv
 121    semh = sst*emissh
 122
 123    ! compute surface scattering term
 124
 125    scatv = refv*tbdnv
 126    scath = refh*tbdnh
 127
 128    ! combine to get space-observed brightness temperature
 129
 130    ! tbv =     tbup + tauatm*(semv + scatv)
 131    ! tbh =     tbup + tauatm*(semh + scath)
 132
 133
 134    ! start
 135    ! write (unit=stdout,fmt=*) 'ifreq 1',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,  &
 136    !                 ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
 137
 138
 139    ADJ_tbup   = ADJ_tbh                    !!! first
 140    ADJ_tauatm = ADJ_tbh*(semh + scath)     !!! first
 141    ADJ_semh   = tauatm*ADJ_tbh             !!! first
 142    ADJ_scath  = tauatm*ADJ_tbh             !!! first
 143
 144    ADJ_tbup   = ADJ_tbv                  + ADJ_tbup
 145    ADJ_tauatm = ADJ_tbv*(semv + scatv)   + ADJ_tauatm
 146    ADJ_semv   = tauatm*ADJ_tbv             !!! first
 147    ADJ_scatv  = tauatm*ADJ_tbv             !!! first
 148
 149    ADJ_refh   = ADJ_scath*tbdnh            !!! first
 150    ADJ_tbdnh  = refh*ADJ_scath             !!! first
 151    ADJ_refv   = ADJ_scatv*tbdnv            !!! first
 152    ADJ_tbdnv  = refv*ADJ_scatv             !!! first
 153    ADJ_sst    = ADJ_semh*emissh          + ADJ_sst
 154    ADJ_emissh = sst*ADJ_semh               !!! first
 155    ADJ_sst    = ADJ_semv*emissv          + ADJ_sst
 156    ADJ_emissv = sst*ADJ_semv               !!! first
 157
 158    ADJ_emissh = - ADJ_refh               + ADJ_emissh
 159    ADJ_emissv = - ADJ_refv               + ADJ_emissv
 160
 161    ADJ_foam   =   ADJ_emissh*fem                      !!! first
 162    ADJ_foam   = - ADJ_emissh*(remh + ebiash(ifreq)) + ADJ_foam
 163    ADJ_remh   =   (1.0 - foam)*ADJ_emissh             !!! first
 164    ADJ_foam   =   ADJ_emissv*fem                    + ADJ_foam
 165    ADJ_foam   = - ADJ_emissv*(remv + ebiasv(ifreq)) + ADJ_foam
 166    ADJ_remv   =   (1.0 - foam)*ADJ_emissv             !!! first
 167
 168    if (u .gt. 5.0) then
 169      ADJ_u = cf2(ifreq)*ADJ_foam  + ADJ_u
 170      foam=foam_save
 171    end if
 172    ADJ_u = cf1(ifreq)*ADJ_foam    + ADJ_u
 173
 174    ADJ_dum = 0.0
 175    dum     = 0.0
 176    ADJ_effangh = 0.0
 177    call da_tbatmos_adj(ifreq,effangh,p0,wv,hwv,ta,gamma,alw,    &
 178                        zcld,dum,tbdnh,dum,                      &
 179                        ADJ_effangh,ADJ_p0,ADJ_wv,ADJ_hwv,       &
 180                        ADJ_ta,ADJ_gamma,ADJ_alw,ADJ_zcld,       &
 181                        ADJ_dum,ADJ_tbdnh,ADJ_dum                )
 182    dum     = 0.0
 183    ADJ_dum = 0.0
 184
 185    ! write (unit=stdout,fmt=*) 'ifreq 2',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,  &
 186    !                 ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
 187
 188    ADJ_effangv = 0.0
 189    call da_tbatmos_adj(ifreq,effangv,p0,wv,hwv,ta,gamma,alw,    &
 190                        zcld,dum,tbdnv,dum,                      &
 191                        ADJ_effangv,ADJ_p0,ADJ_wv,ADJ_hwv,       &
 192                        ADJ_ta,ADJ_gamma,ADJ_alw,ADJ_zcld,       &
 193                        ADJ_dum,ADJ_tbdnv,ADJ_dum                )
 194
 195    ADJ_gx2=0.0
 196    ADJ_sigma=0.0
 197    ! write (unit=stdout,fmt=*) 'ifreq 3',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,  &
 198    !                 ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
 199
 200    call da_effang_adj(ifreq,theta,gx2,sigma,effangv,effangh,    &
 201                       ADJ_gx2,ADJ_sigma,ADJ_effangv,ADJ_effangh )
 202
 203    call da_roughem_adj(ifreq,gx2,sst,theta,remv,remh,         &
 204                        ADJ_gx2,ADJ_sst,ADJ_remv,ADJ_remh      )
 205
 206    ADJ_tauatm = - costhet*ADJ_sigma/tauatm + ADJ_tauatm
 207
 208
 209    ! write (unit=stdout,fmt=*) 'ifreq 4',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,  &
 210    !              ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
 211
 212    call da_tbatmos_adj(ifreq,theta,p0,wv,hwv,ta,gamma,alw,zcld, &
 213                        tbup,tbdn,tauatm,                        &
 214                        ADJ_theta,ADJ_p0,ADJ_wv,ADJ_hwv,ADJ_ta,ADJ_gamma,  &
 215                        ADJ_alw,ADJ_zcld,ADJ_tbup,ADJ_tbdn,      &
 216                        ADJ_tauatm                               )
 217
 218    ADJ_theta=0.0   ! first
 219
 220    ADJ_u = cg(ifreq)*ADJ_gx2 + ADJ_u
 221
 222    ! write (unit=stdout,fmt=*) 'ifreq 5',ADJ_p0,ADJ_ta,ADJ_gamma,ADJ_sst,ADJ_wv,  &
 223    !              ADJ_hwv,ADJ_u,ADJ_alw,ADJ_zcld,ADJ_tbv,ADJ_tbh
 224    ! end
 225
 226    if (trace_use) call da_trace_exit("da_tb_adj")
 227
 228 end subroutine da_tb_adj
 229
 230