uWxUtils
uWxUtils.zip contains Delphi Object Pascal source code files containing
all the weather algorithms used in VpLive and VpPressureCalc. The most
important algorithms convert pressure between sensor pressure, station
pressure, altimeter, and sea level reduced pressure. There are three
files: uWxUtils.pas is the main file containing the various functions,
all of which take and return values in metric units, uWxUtilsUS.pas
wraps the functions in uWxUtils.pas with functions that take and return
values in US units, and uWxUtilsVP.pas contains pressure conversion
functions specific to Vantage Pro weather stations. The
uWxUtils.pas file contents are included at the bottom
of this page, with hyperlinked references. This code may be used
by anyone finding it useful, including commercial use.
uWxUtils is freeware.
Visit the Downloads page to get the
latest version.
uWxUtils.pas (with hyperlinked references)
(pas to html formatting via Peter Johnson's
tool
at DelphiDabbler
unit uWxUtils;
interface
type
TWxReal = double;
TWxHumidity = byte;
TSLPAlgorithm = (
paDavisVP,
paUnivie,
paManBar
);
TAltimeterAlgorithm = (
aaASOS,
aaASOS2,
aaMADIS,
aaNOAA,
aaWOB,
aaSMT
);
TVapAlgorithm = (
vaDavisVp,
vaBuck,
vaBuck81,
vaBolton,
vaTetenNWS,
vaTetenMurray,
vaTeten
);
const
DefaultSLPAlgorithm = paManBar;
DefaultAltimeterAlgorithm = aaMADIS;
DefaultVapAlgorithm = vaBolton;
type
TWxUtils = class(TObject)
private
public
class function StationToSensorPressure(pressureHPa: TWxReal;
sensorElevationM: TWxReal; stationElevationM: TWxReal;
currentTempC: TWxReal): TWxReal;
class function StationToAltimeter(pressureHPa: TWxReal; elevationM: TWxReal;
algorithm: TAltimeterAlgorithm = DefaultAltimeterAlgorithm): TWxReal;
class function StationToSeaLevelPressure(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
class function SensorToStationPressure(pressureHPa: TWxReal;
sensorElevationM: TWxReal; stationElevationM: TWxReal;
currentTempC: TWxReal): TWxReal;
class function SeaLevelToStationPressure(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
class function PressureReductionRatio(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
class function ActualVaporPressure(tempC: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
class function SaturationVaporPressure(tempC: TWxReal;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
class function MixingRatio(pressureHPa: TWxReal; tempC: TWxReal; humidity: TWxHumidity): TWxReal;
class function VirtualTempK(pressureHPa: TWxReal; tempC: TWxReal; humidity: TWxHumidity): TWxReal;
class function HumidityCorrection(tempC: TWxReal; elevationM: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
class function DewPoint(tempC: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
class function WindChill(tempC: TWxReal; windSpeedKmph: TWxReal): TWxReal;
class function HeatIndex(tempC: TWxReal; humidity: TWxHumidity): TWxReal;
class function Humidex(tempC: TWxReal; humidity: TWxHumidity): TWxReal;
class function GeopotentialAltitude(geometricAltitudeM: TWxReal): TWxReal;
class function FToC(value: TWxReal): TWxReal;
class function CToF(value: TWxReal): TWxReal;
class function CToK(value: TWxReal): TWxReal;
class function KToC(value: TWxReal): TWxReal;
class function FToR(value: TWxReal): TWxReal;
class function RToF(value: TWxReal): TWxReal;
class function InToHPa(value: TWxReal): TWxReal;
class function HPaToIn(value: TWxReal): TWxReal;
class function FtToM(value: TWxReal): TWxReal;
class function MToFt(value: TWxReal): TWxReal;
class function InToMm(value: TWxReal): TWxReal;
class function MmToIn(value: TWxReal): TWxReal;
class function MToKm(value: TWxReal): TWxReal;
class function KmToM(value: TWxReal): TWxReal;
class function Power(const Base, Exponent: TWxReal): TWxReal;
class function Power10(const exponent: TWxReal): TWxReal;
end;
implementation
const
gravity = 9.80665;
uGC = 8.31432;
moleAir = 0.0289644;
moleWater = 0.01801528;
gasConstantAir = uGC/moleAir;
standardSLP = 1013.25;
standardSlpInHg = 29.921;
standardTempK = 288.15;
earthRadius45 = 6356.766;
standardLapseRate = 0.0065;
standardLapseRateFt = standardLapseRate * 0.3048;
vpLapseRateUS = 0.00275;
manBarLapseRate = 0.0117;
class function TWxUtils.StationToSensorPressure(pressureHPa: TWxReal;
sensorElevationM: TWxReal; stationElevationM: TWxReal;
currentTempC: TWxReal): TWxReal;
begin
Result := InToHPa(HPaToIn(pressureHPa) / Power10(0.00813 * MToFt(sensorElevationM - stationElevationM)
/ FToR(CToF(currentTempC))));
end;
class function TWxUtils.StationToAltimeter(PressureHPa: TWxReal; elevationM: TWxReal;
algorithm: TAltimeterAlgorithm = DefaultAltimeterAlgorithm): TWxReal;
var
geopEl: TWxReal;
k1, k2: TWxReal;
begin
case algorithm of
aaASOS:
begin
Result := InToHPa(Power(Power(HPaToIn(pressureHPa), 0.1903) + (1.313E-5 * MToFt(elevationM)), 5.255));
end;
aaASOS2:
begin
geopEl := GeopotentialAltitude(elevationM);
k1 := standardLapseRate * gasConstantAir / gravity;
k2 := 8.41728638E-5;
Result := Power(Power(pressureHPa, k1) + (k2 * geopEl), 1/k1);
end;
aaMADIS:
begin
k1 := 0.190284;
k2 := 8.4184960528E-5;
Result := Power(Power(pressureHPa - 0.3, k1) + (k2 * elevationM), 1/k1);
end;
aaNOAA:
begin
k1 := 0.190284;
k2 := 8.42288069E-5;
Result := (pressureHPa - 0.3) * Power(1 + (k2 * (elevationM / Power(pressureHPa - 0.3, k1))), 1/k1);
end;
aaWOB:
begin
k1 := standardLapseRate * gasConstantAir / gravity;
k2 := 1.312603E-5;
Result := InToHPa(Power(Power(HPaToIn(pressureHPa), k1) + (k2 * MToFt(elevationM)), 1/k1));
end;
aaSMT:
begin
k1 := 0.190284;
k2 := 4.30899E-5;
geopEl := GeopotentialAltitude(elevationM);
Result := InToHPa((HPaToIn(pressureHPa) - 0.01)
* Power(1 + (k2 * (geopEl / Power(HPaToIn(pressureHPa) - 0.01, k1))), 1/k1));
end;
else Result := pressureHPa;
end;
end;
class function TWxUtils.StationToSeaLevelPressure(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
begin
Result := pressureHPa * PressureReductionRatio(pressureHPa, elevationM, currentTempC,
meanTempC, humidity, algorithm);
end;
class function TWxUtils.SensorToStationPressure(pressureHPa: TWxReal;
sensorElevationM: TWxReal; stationElevationM: TWxReal;
currentTempC: TWxReal): TWxReal;
begin
Result := InToHPa(HPaToIn(pressureHPa) * Power10(0.00813 * MToFt(sensorElevationM - stationElevationM)
/ FToR(CToF(currentTempC))));
end;
class function TWxUtils.SeaLevelToStationPressure(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
begin
Result := pressureHPa / PressureReductionRatio(pressureHPa, elevationM, currentTempC,
meanTempC, humidity, algorithm);
end;
class function TWxUtils.PressureReductionRatio(pressureHPa: TWxReal; elevationM: TWxReal;
currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
var
geopElevationM: TWxReal;
hCorr: TWxReal;
begin
case algorithm of
paUnivie:
begin
geopElevationM := GeopotentialAltitude(elevationM);
Result := Exp(((gravity/gasConstantAir) * geopElevationM)
/ (VirtualTempK(pressureHPa, meanTempC, humidity) + (geopElevationM * standardLapseRate/2)));
end;
paDavisVp:
begin
if (humidity > 0) then begin
hCorr := (9/5) * HumidityCorrection(currentTempC, elevationM, humidity, vaDavisVP);
end else begin
hCorr := 0;
end;
Result := Power(10, (MToFt(elevationM) / (122.8943111 * (CToF(meanTempC) + 460 + (MToFt(elevationM) * vpLapseRateUS/2) + hCorr))));
end;
paManBar:
begin
if (humidity > 0) then begin
hCorr := (9/5) * HumidityCorrection(currentTempC, elevationM, humidity, vaBuck);
end else begin
hCorr := 0;
end;
geopElevationM := GeopotentialAltitude(elevationM);
Result := Exp(geopElevationM * 6.1454E-2 / (CToF(meanTempC) + 459.7 + (geopElevationM * manBarLapseRate / 2) + hCorr));
end;
else Result := 1;
end;
end;
class function TWxUtils.ActualVaporPressure(tempC: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
begin
result := (humidity * SaturationVaporPressure(tempC, algorithm)) / 100;
end;
class function TWxUtils.SaturationVaporPressure(tempC: TWxReal;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
begin
case algorithm of
vaDavisVp: Result := 6.112 * Exp((17.62 * tempC)/(243.12 + tempC));
vaBuck: Result := 6.1121 * Exp((18.678 - (tempC/234.5)) * tempC / (257.14 + tempC));
vaBuck81: Result := 6.1121 * Exp((17.502 * tempC)/(240.97 + tempC));
vaBolton: Result := 6.112 * Exp(17.67 * tempC / (tempC + 243.5));
vaTetenNWS: Result := 6.112 * Power(10,(7.5 * tempC / (tempC + 237.7)));
vaTetenMurray: Result := Power(10, (7.5 * tempC / (237.5 + tempC)) + 0.7858);
vaTeten: Result := 6.1078 * Power(10, (7.5 * tempC / (tempC + 237.3)));
else Result := 0;
end;
end;
class function TWxUtils.MixingRatio(pressureHPa: TWxReal; tempC: TWxReal;
humidity: TWxHumidity): TWxReal;
var
vapPres: TWxReal;
const
k1 = moleWater / moleAir;
begin
vapPres := ActualVaporPressure(tempC, humidity, vaBuck);
Result := 1000 * ((k1 * vapPres) / (pressureHPa - vapPres));
end;
class function TWxUtils.VirtualTempK(pressureHPa: TWxReal; tempC: TWxReal;
humidity: TWxHumidity): TWxReal;
var
vapPres: TWxReal;
const
epsilon = 1 - (moleWater / moleAir);
begin
vapPres := ActualVaporPressure(tempC, humidity, vaBuck);
Result := (CtoK(tempC)) / (1-(epsilon * (vapPres/pressureHPa)));
end;
class function TWxUtils.HumidityCorrection(tempC: TWxReal; elevationM: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
var
vapPress: TWxReal;
begin
vapPress := ActualVaporPressure(tempC, humidity, algorithm);
Result := (vappress * ((2.8322E-9 * Sqr(elevationM)) + (2.225E-5 * elevationM) + 0.10743));
end;
class function TWxUtils.DewPoint(tempC: TWxReal; humidity: TWxHumidity;
algorithm: TVapAlgorithm = DefaultVapAlgorithm): TWxReal;
var
lnVapor: TWxReal;
begin
lnVapor := Ln(ActualVaporPressure(tempc, humidity, algorithm));
case algorithm of
vaDavisVp: Result := ((243.12 * LnVapor) - 440.1) / (19.43 - LnVapor);
else
Result := ((237.7 * LnVapor) - 430.22) / (19.08 - LnVapor);
end;
end;
class function TWxUtils.WindChill(tempC: TWxReal; windSpeedKmph: TWxReal): TWxReal;
var
windPow: TWxReal;
begin
if ((tempC >= 10.0) or (windSpeedKmph <= 4.8)) then begin
Result := tempC;
end else begin
windPow := Power(windSpeedKmph, 0.16);
Result := 13.12 + (0.6215 * tempC) - (11.37 * windPow) + (0.3965 * tempC * windPow);
end;
if (Result > tempC) then Result := tempC;
end;
class function TWxUtils.HeatIndex(tempC: TWxReal; humidity: TWxHumidity): TWxReal;
var
tempF: TWxReal;
tSqrd: TWxReal;
hum: TWxReal;
hSqrd: TWxReal;
begin
tempF := CToF(tempC);
if (tempF < 80) then begin
Result := tempF;
end else begin
tSqrd := tempF * tempF;
hum := humidity;
hSqrd := hum * hum;
Result := 0 - 42.379 + (2.04901523 * tempF) + (10.14333127 * humidity)
- (0.22475541 * tempF * humidity) - (0.00683783 * tSqrd)
- (0.05481717 * hSqrd) + (0.00122874 * tSqrd * humidity)
+ (0.00085282 * tempF * hSqrd) - (0.00000199 * tSqrd * hSqrd);
if ((humidity < 13) and (tempF >= 80) and (tempF <= 112)) then begin
Result := Result - ((13 - humidity)/4) * Sqrt((17 - Abs(tempf - 95))/17);
end else if ((humidity > 85) and (tempF >= 80) and (tempF <= 87)) then begin
Result := Result + ((humidity - 85)/10) * ((87 - tempF)/5);
end;
end;
Result := FToC(Result);
end;
class function TWxUtils.Humidex(tempC: TWxReal; humidity: TWxHumidity): TWxReal;
begin
Result := tempC + ((5/9) * (ActualVaporPressure(tempC, humidity, vaTetenNWS) - 10.0));
end;
class function TWxUtils.GeopotentialAltitude(geometricAltitudeM: TWxReal): TWxReal;
begin
Result := (earthRadius45 * 1000 * geometricAltitudeM) / ((earthRadius45 * 1000) + geometricAltitudeM);
end;
class function TWxUtils.FtoC(value: TWxReal): TWxReal;
begin
Result := ((value - 32) * 5) / 9;
end;
class function TWxUtils.CtoF(value: TWxReal): TWxReal;
begin
Result := ((value * 9) / 5) + 32;
end;
class function TWxUtils.CtoK(value: TWxReal): TWxReal;
begin
Result := 273.15 + value;
end;
class function TWxUtils.KtoC(value: TWxReal): TWxReal;
begin
Result := value - 273.15;
end;
class function TWxUtils.FToR(value: TWxReal): TWxReal;
begin
Result := value + 459.67;
end;
class function TWxUtils.RToF(value: TWxReal): TWxReal;
begin
Result := value - 459.67;
end;
class function TWxUtils.InToHPa(value: TWxReal): TWxReal;
begin
Result := value / 0.02953;
end;
class function TWxUtils.HPaToIn(value: TWxReal): TWxReal;
begin
Result := value * 0.02953;
end;
class function TWxUtils.FtToM(value: TWxReal): TWxReal;
begin
Result := value * 0.3048;
end;
class function TWxUtils.MToFt(value: TWxReal): TWxReal;
begin
Result := value / 0.3048;
end;
class function TWxUtils.InToMm(value: TWxReal): TWxReal;
begin
Result := value * 25.4;
end;
class function TWxUtils.MmToIn(value: TWxReal): TWxReal;
begin
Result := value / 25.4;
end;
class function TWxUtils.MToKm(value: TWxReal): TWxReal;
begin
Result := value * 1.609344;
end;
class function TWxUtils.KmToM(value: TWxReal): TWxReal;
begin
Result := value / 1.609344;
end;
class function TWxUtils.Power(const base, exponent: TWxReal): TWxReal;
begin
if exponent = 0.0 then begin
Result := 1.0;
end else if (base = 0.0) and (exponent > 0.0) then begin
Result := 0.0;
end else begin
Result := Exp(exponent * Ln(base));
end;
end;
class function TWxUtils.Power10(const exponent: TWxReal): TWxReal;
const
ln10 = 2.302585093;
begin
if (exponent = 0.0) then begin
Result := 1.0;
end else begin
Result := Exp(exponent * ln10);
end;
end;
end.