c2--567--1---------2---------3---------4---------5---------6---------712
      Subroutine altitude(p_lvl, t_lvl, z_phys, n)
c
c     -- Altitude/Physical height at specified pressure levels/N. Nath
c        /April '01
c
c     -- Function:  Given atmospheric temperature T(z), at equidistant 
c        log(p) or z levels, where p is the pressure and z is the 
c        pressure-scaled height (the two variables being related by 
c        d(logp) = -dz/H, where H is the constant scale height referred
c        to z = 0 and standard temperature Ts), this routine computes
c        the physical height z_phys at those levels.  The underlying 
c        equation, which is a consequence of the hydrostatic equation 
c        and the perfect gas law, is: (1 + z_phys/R)**-1 = 1 - 
c        Intgrl((T(z')/Ts)*dz'/R; z'=0:z), where R is the radius of the
c        Earth.
c
c     -- Input & Output Variables:
c
            implicit none
c
            real        p_lvl(n), t_lvl(n)
                              !atmospheric pressure (atmos) and 
                              !  temperature (K), at levels 1 (bottom)
                              !  through n (top) (input)

            real        z_phys(n)
                              !physical heights (cm) (output)

            integer     n     !number of levels (input)
c
c     -- Local Variables:
c
            real        r, h, ts
            parameter   (r = 6.378e8, h = 7.996e5, ts = 273.16)
                              !radius of Earth (cm), scale height (cm),
                              !  standard temperature (K)

            real        delz  !layer thickness (cm)

            real        delzn, sum

            integer     i
c
      delz = h*log(p_lvl(1)/p_lvl(n))/(n - 1)
      delzn = delz/r

      z_phys(1) = 0.0
      sum = 0.0
      do i = 1,n-1
         sum = sum + (t_lvl(i) + t_lvl(i+1))/(2.0*ts)
         z_phys(i+1) = r*(1.0/(1.0 - sum*delzn) - 1.0)
      end do
c
      return 
      end