Module atm
!Contains Nick Alley's standard atomosphere routines
contains
!***********************************************************************!
SUBROUTINE stdatmos(h,d)
!  h = geometric altitude, specified by user (m,ft)
!  z = geopotential altitude, returned by subroutine (m,ft)
!  t = temperature, returned by subroutine (K,R)
!  p = pressure, returned by subroutine (N/m^2,lbf/ft^2)
!  d = density, returned by subroutine (kg/m^3,slugs/ft^3)
!  a = speed of sound (ft/s)
!  stdatmos(units,h,z,t,p,d)
USE Precision_Module
IMPLICIT NONE
REAL(prec):: h,z,t,p,d,a
REAL(prec):: hsi,zsi,tsi,psi,dsi,asi
hsi=h*0.3048               ; call statsi(hsi,zsi,tsi,psi,dsi)
z=zsi/0.3048               ; t=tsi*1.8
asi=sqrt(1.4*287.0528*tsi) ; a=asi/0.3048
p=psi*0.02088543           ; d=dsi*0.001940320
RETURN
END SUBROUTINE stdatmos
!***********************************************************************!
SUBROUTINE statsi(h,z,t,p,d)
!     h = geometric altitude, specified by user (m)
!     z = geopotential altitude, returned by subroutine (m)
!     t = temperature, returned by subroutine (K)
!     p = pressure, returned by subroutine (N/m**2)
!     d = density, returned by subroutine (kg/m**3)
USE Precision_Module
IMPLICIT NONE
REAL(prec):: Lt,zsa(0:8),Tsa(0:8),Psa(0:8)
REAL(prec):: h,z,t,p,d,g0,R,Re,ex
INTEGER::i
data zsa/    0.,11000.,20000.,32000.,47000.,52000.,61000.,79000.,9.9e20/
data Tsa/288.15,216.65,216.65,228.65,270.65,270.65,252.65,180.65,180.65/
data g0,R,Re,Psa(0)/9.80665,287.0528,6356766.,101325./
z=Re*h/(Re+h)
do i=1,8
	Lt=-(Tsa(i)-Tsa(i-1))/(zsa(i)-zsa(i-1))
	if(Lt.eq.0.)then
		if(z.le.zsa(i))then
		   t=Tsa(i-1)
		   p=Psa(i-1)*exp(-g0*(z-zsa(i-1))/R/Tsa(i-1))
		   d=p/R/t
		   return
		else
		   Psa(i)=Psa(i-1)*exp(-g0*(zsa(i)-zsa(i-1))/R/Tsa(i-1))
		endif
	else
		ex=g0/R/Lt
		if(z.lt.zsa(i))then
		   t=Tsa(i-1)-Lt*(z-zsa(i-1))
		   p=Psa(i-1)*(t/Tsa(i-1))**ex
		   d=p/R/t
		   return
		else
		   Psa(i)=Psa(i-1)*(Tsa(i)/Tsa(i-1))**ex
		endif
	endif
end do
t=Tsa(8) ; p=0.0 ; d=0.0
RETURN
END SUBROUTINE statsi
!***********************************************************************!
end module atm