module Atmosphere
( Atmos(..)
, siAtmosphere
, usAtmosphere
, atmosphere
, siAltitudeFromPressure
) where
import Atmosphere.Constants
data Atmos a = Atmos { atmosTemperature :: a
, atmosPressure :: a
, atmosDensity :: a
, atmosSpeedOfSound :: a
, atmosViscosity :: a
, atmosKinematicViscosity :: a
}
siAtmosphere :: (Floating a, Ord a) => a -> Atmos a
siAtmosphere alt_m =
Atmos { atmosTemperature = temp
, atmosPressure = pressure
, atmosDensity = density
, atmosSpeedOfSound = asound
, atmosViscosity = viscosity
, atmosKinematicViscosity = kinematicViscosity
}
where
alt_km = 0.001*alt_m
(sigma, delta, theta) = atmosphere alt_km
temp = _TZERO * theta
pressure = _PZERO * delta
density = _RHOZERO * sigma
asound = _AZERO * sqrt theta
viscosity = metricViscosity theta
kinematicViscosity = viscosity/density
usAtmosphere :: (Floating a, Ord a) => a -> Atmos a
usAtmosphere alt_ft =
Atmos { atmosTemperature = temp
, atmosPressure = pressure
, atmosDensity = density
, atmosSpeedOfSound = asound
, atmosViscosity = viscosity
, atmosKinematicViscosity = kinematicViscosity
}
where
alt_km = 0.001*_FT2METERS*alt_ft
(sigma, delta, theta) = atmosphere alt_km
temp = _KELVIN2RANKINE*_TZERO*theta
pressure = _PZERO*delta/47.88
density = _RHOZERO*sigma/515.379
asound = (_AZERO/_FT2METERS)*sqrt theta
viscosity=(1.0/_PSF2NSM)*metricViscosity theta
kinematicViscosity = viscosity/density
siAltitudeFromPressure :: (Floating a, Ord a) => a -> a
siAltitudeFromPressure pressureIn = 1000*alt
where
alt = _REARTH / (_REARTH/h 1)
deltaIn = pressureIn / _PZERO
(htabI, tbase, ptabI, tgradI) = getI htpgTable
where
getI [htab'] = htab'
getI (htab0:htab1@(_,_,delta',_):htabs)
| deltaIn < delta' = getI (htab1:htabs)
| otherwise = htab0
getI [] = error "something went wrong"
h
| 0.0 == tgradI = htabI tbase / _GMR * (log (deltaIn / ptabI))
| otherwise = htabI + tbase/tgradI*((deltaIn/ptabI)**(tgradI/_GMR) 1)
metricViscosity :: (Floating a, Ord a) => a -> a
metricViscosity theta = _BETAVISC*sqrt(t*t*t)/(t+_SUTH)
where
t = theta * _TZERO
atmosphere :: (Floating a, Ord a) => a -> (a,a,a)
atmosphere alt = (sigma, delta, theta)
where
h = alt*_REARTH/(alt+_REARTH)
(htabI, tbase, ptabI, tgradI) = getI htpgTable
where
getI [htab'] = htab'
getI (htab0:htab1@(h',_,_,_):htabs)
| h > h' = getI (htab1:htabs)
| otherwise = htab0
getI [] = error "something went wrong"
deltah = h htabI
tlocal = tbase + tgradI*deltah
theta = tlocal/_TZERO
delta
| 0.0 == tgradI = ptabI*exp(_GMR*deltah/tbase)
| otherwise = ptabI*(tbase/tlocal)**(_GMR/tgradI)
sigma = delta/theta