UFO
ufo_satcolumn_mod.F90
Go to the documentation of this file.
2 contains
3 
4 subroutine simulate_column_ob(nlayers_obs, nlayers_model, avgkernel_obs, &
5  prsl_obs, prsl_model, temp_model, zi_model, &
6  profile_model, hofx, &
7  troplev_obs, airmass_tot, airmass_trop)
8  use kinds
11  implicit none
12  integer, intent(in ) :: nlayers_obs, nlayers_model
13  real(kind_real), intent(in ), dimension(nlayers_obs) :: avgkernel_obs
14  real(kind_real), intent(in ), dimension(nlayers_obs) :: prsl_obs
15  real(kind_real), intent(in ), dimension(nlayers_model) :: prsl_model
16  real(kind_real), intent(in ), dimension(nlayers_model) :: profile_model
17  real(kind_real), intent(in ), dimension(nlayers_model) :: temp_model
18  real(kind_real), intent(in ), dimension(nlayers_model+1) :: zi_model
19  real(kind_real), intent( out) :: hofx
20  real(kind_real), intent(in ), optional :: airmass_tot, airmass_trop
21  integer, intent(in ), optional :: troplev_obs
22  real(kind_real) :: airmass_ratio, lnp_ob, wf, dz
23  real(kind_real), dimension(nlayers_obs) :: avgkernel_use
24  real(kind_real), dimension(nlayers_model) :: lnp_model, profile_model_use
25  real(kind_real), dimension(nlayers_obs) :: profile_obslayers
26  integer :: k, wi
27  logical :: troposphere
28 
29  troposphere = .false.
30  ! determine if we are computing a total column or just tropospheric column
31  if (present(airmass_tot) .and. present(airmass_trop) .and. present(troplev_obs)) then
32  troposphere = .true.
33  end if
34 
35  if (troposphere) then
36  ! compute air mass ratio and apply it to the averaging kernel
37  airmass_ratio = airmass_tot / airmass_trop
38  avgkernel_use = avgkernel_obs * airmass_ratio
39 
40  ! set all layers to zero above tropopause layer
41  avgkernel_use(troplev_obs+1:nlayers_obs) = zero
42  else
43  avgkernel_use = avgkernel_obs
44  end if
45 
46  ! need to convert from kg/kg to molec/m3 for each layer's T and P
47  profile_model_use = profile_model * ((avogadro*prsl_model)/(temp_model*gas_constant))
48 
49  ! need to convert from molec/m3 to molec/cm2 to sum later
50  do k=1,nlayers_model
51  dz = (zi_model(k+1) - zi_model(k))
52  profile_model_use(k) = profile_model_use(k) * dz ! convert to molec/m2
53  profile_model_use(k) = profile_model_use(k) / 1.0e4_kind_real ! to molec/cm2
54  end do
55 
56  ! need to interpolate model profile layers to observation layers
57  lnp_model = log(prsl_model)
58  do k=1,nlayers_obs
59  lnp_ob = log(prsl_obs(k))
60  call vert_interp_weights(nlayers_model, lnp_ob, lnp_model, wi, wf)
61  call vert_interp_apply(nlayers_model, profile_model_use, &
62  profile_obslayers(k), wi, wf)
63  end do
64 
65  ! compute hofx as column integral
66  hofx = zero
67  do k=1,nlayers_obs
68  hofx = hofx + (avgkernel_use(k) * profile_obslayers(k))
69  end do
70 
71 end subroutine simulate_column_ob
72 
73 end module satcolumn_mod
satcolumn_mod
Definition: ufo_satcolumn_mod.F90:1
ufo_constants_mod::gas_constant
real(kind_real), parameter, public gas_constant
Definition: ufo_constants_mod.F90:17
vert_interp_mod
Fortran module to perform linear interpolation.
Definition: vert_interp.F90:8
vert_interp_mod::vert_interp_weights
subroutine vert_interp_weights(nlev, obl, vec, wi, wf)
Definition: vert_interp.F90:22
ufo_constants_mod
Definition: ufo_constants_mod.F90:2
satcolumn_mod::simulate_column_ob
subroutine simulate_column_ob(nlayers_obs, nlayers_model, avgkernel_obs, prsl_obs, prsl_model, temp_model, zi_model, profile_model, hofx, troplev_obs, airmass_tot, airmass_trop)
Definition: ufo_satcolumn_mod.F90:8
ufo_constants_mod::avogadro
real(kind_real), parameter, public avogadro
Definition: ufo_constants_mod.F90:16
vert_interp_mod::vert_interp_apply
subroutine vert_interp_apply(nlev, fvec, f, wi, wf)
Definition: vert_interp.F90:73
ufo_constants_mod::zero
real(kind_real), parameter, public zero
Definition: ufo_constants_mod.F90:24