anchor.diff:               lhgt      = 0.0_r8
                
                ! Convert input water vapor mixing ratio from PPB to unitless
  -             thish2omxr = H2O_profiles(idx_lon,idx_lat,n) / 1.0E9
  +             thish2omxr = lh2o(n) / 1.0E9
   
                ! Calculate mean molecular weight of wet air 
                Mwet = (1.0_r8 - thish2omxr)*Mdry + thish2omxr*Mh2o
  
positive.code:              ! Calculate layer thickness using the following
               ! dz = -(R*T/g) * dlnP
               n1                   = n - 1
               detlnp = log(lpre(n1)) - log(lpre(n))
               lhgt = Rwet * ltmp(n) * detlnp / gplanet ! meter
               
               rho                = layer_press(n) / ltmp(n) / Rstar 
               aircolumn          = rho*lhgt*Navogadro*1.0E-4 ! # air/cm^2
               
               ! -------------------------------------------------------------
               climatology(ixtrack,itimes,n) = aircolumn * Gas_profiles(idx_lon,idx_lat,n) /1.0E9 ! [GAS]/cm^2
               
            END DO
            
            !  Set non-physical entries to zero.
  
negative.code:     ! ---------------
      ! Local variables
      ! ---------------
      INTEGER (KIND=i4) :: itimes, ixtrack, spix, epix, idx_lat, idx_lon, ilevel, n, n1, one, status
      REAL    (KIND=r8)                      :: rho, lhgt, aircolumn, clima_psurf
      REAL    (KIND=r8), DIMENSION (0:CmETA) :: lpre, level_press, clima_lpre
      REAL    (KIND=r8), DIMENSION (1:CmETA) :: ltmp, lh2o, lgas, re_tmp, re_h2o, re_gas, &
           layer_press, clima_layer_press, clima_local_heights
      REAL    (KIND=r8) :: thish2omxr, Mwet, Rwet, detlnp
  
      ! -----------------------
  