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;
    // The following source code may be freely used, including for commercial purposes
    // Steve Hatchett, SoftWx, Inc.
    // http://www.softwx.com/

    {**********************************************
    This file contains functions for performing various weather related calculations.

    Notes about pressure
      Sensor Pressure           raw pressure indicated by the barometer instrument
      Station Pressure          Sensor Pressure adjusted for any difference between sensor elevation and official station elevation
      Field Pressure     (QFE)  Usually the same as Station Pressure
      Altimeter Setting  (QNH)  Station Pressure adjusted for elevation (assumes standard atmosphere)
      Sea Level Pressure (QFF)  Station Pressure adjusted for elevation, temperature and humidity

    Notes about input parameters:
      currentTemp -   current instantaneous station temperature
      meanTemp -      average of current temp and the temperature 12 hours in
                      the past. If the 12 hour temp is not known, simply pass
                      the same value as currentTemp for the mean temp.
      humidity -      Value should be 0 to 100. For the pressure conversion
                      functions, pass a value of zero if you do not want to
                      the algorithm to include the humidity correction factor
                      in the calculation. If you provide a humidity value
                      > 0, then humidity effect will be included in the
                      calculation.
     elevation -     This should be the geometric altitude of the station
                      (this is the elevation provided by surveys and normally
                      used by people when they speak of elevation). Some
                      algorithms will convert the elevation internally into
                      a geopotential altitude.
      sensorElevation - This should be the geometric altitude of the actual
                      barometric sensor (which could be different than the
                      official station elevation).

    Notes about Sensor Pressure vs. Station Pressure:
      SensorToStationPressure and StationToSensorPressure functions are based
      on an ASOS algorithm. It corrects for a difference in elevation between
      the official station location and the location of the barometetric sensor.
      It turns out that if the elevation difference is under 30 ft, then the
      algorithm will give the same result (a 0 to .01 inHg adjustment) regardless
      of temperature. In that case, the difference can be covered using a simple
      fixed offset. If the difference is 30 ft or greater, there is some effect
      from temperature, though it is small. For example, at a 100ft difference,
      the adjustment will be .13 inHg at -30F and .10 at 100F. The bottom line
      is that while ASOS stations may do this calculation, it is likely unneeded
      for home weather stations, and the station pressure and the sensor pressure
      can be treated as equivalent.
    **********************************************}
    interface

    type
    // by changing this you can change the precision of the floating point calculations
      TWxReal = double;
      TWxHumidity = byte; // humidities expressed as integer percentage, i.e. 88 = 88%

    // formulas/algorithms
      // Sea Level Pressure reduction algorithms
      TSLPAlgorithm = (
        paDavisVP,  // algorithm closely approximates SLP calculation used inside Davis Vantage Pro weather equipment console (http://www.davisnet.com/weather/)
        paUnivie,   // http://www.univie.ac.at/IMG-Wien/daquamap/Parametergencom.html
        paManBar    // from Manual of Barometry (1963)
        );

      // Altimeter algorithms
      TAltimeterAlgorithm = (
        aaASOS,     // formula described in the ASOS training docs
        aaASOS2,    // metric formula that was likely used to derive the aaASOS formula
        aaMADIS,    // apparently the formula used by the MADIS system
        aaNOAA,     // essentially the same as aaSMT with any result differences caused by unit conversion rounding error and geometric vs. geopotential elevation
        aaWOB,      // Weather Observation Handbook (algorithm similar to aaASOS & aaASOS2 - main differences being precision of constants used)
        aaSMT       // Smithsonian Meteorological Tables (1963)
        );

      TVapAlgorithm = (
        vaDavisVp,  // algorithm closely approximates calculation used by Davis Vantage Pro weather stations and software
        vaBuck,     // this and the remaining algorithms described at http://cires.colorado.edu/~voemel/vp.html
        vaBuck81,
        vaBolton,
        vaTetenNWS,
        vaTetenMurray,
        vaTeten
        );

    const
      DefaultSLPAlgorithm = paManBar;
      DefaultAltimeterAlgorithm = aaMADIS;
      DefaultVapAlgorithm = vaBolton;

    type
      TWxUtils = class(TObject)
      private
      public
        // basic pressure functions - all take metric params (hectoPascals(hPa), celsius, meters)
        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 AltimeterToStationPressure(pressureHPa: TWxReal; elevationM: TWxReal;
     //                       algorithm: TAltimeterAlgorithm = DefaultAltimeterAlgorithm): TWxReal;
        class function SeaLevelToStationPressure(pressureHPa: TWxReal; elevationM: TWxReal;
                            currentTempC: TWxReal; meanTempC: TWxReal; humidity: TWxHumidity;
                            algorithm: TSLPAlgorithm = DefaultSLPAlgorithm): TWxReal;
        // low level pressure related functions
        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;

        // temperature related functions
        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;
        // simplified algorithm for geopotential altitude from US Standard Atmosphere 1976
        // assumes latitude 45.5 degrees
        class function GeopotentialAltitude(geometricAltitudeM: TWxReal): TWxReal;

        // basic unit conversion functions
        class function FToC(value: TWxReal): TWxReal;       // Fahrenheit to Celsius
        class function CToF(value: TWxReal): TWxReal;       // Celsius to Fahrenheit
        class function CToK(value: TWxReal): TWxReal;       // Celsius to Kelvin
        class function KToC(value: TWxReal): TWxReal;       // Kelvin to Celsius
        class function FToR(value: TWxReal): TWxReal;       // Fahrenheit to Rankine
        class function RToF(value: TWxReal): TWxReal;       // Rankine to Fahrenheit
        class function InToHPa(value: TWxReal): TWxReal;    // Inches Mercury to HectoPascals/millibars
        class function HPaToIn(value: TWxReal): TWxReal;    // HectoPascals/millibars to Inches Mercury
        class function FtToM(value: TWxReal): TWxReal;      // Feet to Meters
        class function MToFt(value: TWxReal): TWxReal;      // Meters to Feet
        class function InToMm(value: TWxReal): TWxReal;     // Inches to Millimeters
        class function MmToIn(value: TWxReal): TWxReal;     // Millimeters to Inches
        class function MToKm(value: TWxReal): TWxReal;      // Miles to Kilometers
        class function KmToM(value: TWxReal): TWxReal;      // Kilometers to Miles

        // general math functions
        class function Power(const Base, Exponent: TWxReal): TWxReal;
        class function Power10(const exponent: TWxReal): TWxReal;
      end;

    implementation

    const
    // U.S. Standard Atmosphere (1976) constants
      gravity = 9.80665;          // g at sea level at latitude 45.5 degrees in m/sec^2
      uGC = 8.31432;              // universal gas constant in J/mole-K
      moleAir = 0.0289644;        // mean molecular mass of air in kg/mole
      moleWater = 0.01801528;     // molecular weight of water in kg/mole
      gasConstantAir = uGC/moleAir; // (287.053) gas constant for air in J/kgK
      standardSLP = 1013.25;      // standard sea level pressure in hPa
      standardSlpInHg = 29.921;   // standard sea level pressure in inHg
      standardTempK = 288.15;     // standard sea level temperature in Kelvin
      earthRadius45 = 6356.766;   // radius of the earth at latitude 45.5 degrees in km

      standardLapseRate = 0.0065; // standard lapse rate (6.5C/1000m i.e. 6.5K/1000m)

      standardLapseRateFt = standardLapseRate * 0.3048; // (0.0019812) standard lapse rate per foot (1.98C/1000ft)
      vpLapseRateUS = 0.00275;    // lapse rate used by Davis VantagePro (2.75F/1000ft)
      manBarLapseRate = 0.0117;   // lapse rate from Manual of Barometry (11.7F/1000m, which = 6.5C/1000m)

    class function TWxUtils.StationToSensorPressure(pressureHPa: TWxReal;
                            sensorElevationM: TWxReal; stationElevationM: TWxReal;
                            currentTempC: TWxReal): TWxReal;
    begin
      // from ASOS formula specified in US units
      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:
          // see ASOS training at http://www.nwstc.noaa.gov
          // see also http://wahiduddin.net/calc/density_altitude.htm
          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; // approx. 0.190263
            k2 := 8.41728638E-5; // (standardLapseRate / standardTempK) * (Power(standardSLP,  k1)
            Result := Power(Power(pressureHPa, k1) + (k2 * geopEl), 1/k1);
          end;
        aaMADIS:
          // from MADIS API by NOAA Forecast Systems Lab, see http://madis.noaa.gov/madis_api.html
          begin
            k1 := 0.190284; // discrepency with calculated k1 probably because Smithsonian used less precise gas constant and gravity values
            k2 := 8.4184960528E-5; // (standardLapseRate / standardTempK) * (Power(standardSLP, k1)
            Result := Power(Power(pressureHPa - 0.3, k1) + (k2 * elevationM), 1/k1);
          end;
        aaNOAA:
          // see http://www.srh.noaa.gov/elp/wxcalc/formulas/altimeterSetting.html
          begin
            k1 := 0.190284; // discrepency with k1 probably because Smithsonian used less precise gas constant and gravity values
            k2 := 8.42288069E-5; // (standardLapseRate / 288) * (Power(standardSLP, k1SMT);
            Result := (pressureHPa - 0.3) * Power(1 + (k2 * (elevationM / Power(pressureHPa - 0.3, k1))), 1/k1);
          end;
        aaWOB:
          // see http://www.wxqa.com/archive/obsman.pdf
          begin
            k1 := standardLapseRate * gasConstantAir / gravity; // approx. 0.190263
            k2 := 1.312603E-5; //(standardLapseRateFt / standardTempK) * Power(standardSlpInHg, k1);
            Result := InToHPa(Power(Power(HPaToIn(pressureHPa), k1) + (k2 * MToFt(elevationM)), 1/k1));
          end;
        aaSMT:
          // see WMO Instruments and Observing Methods Report No.19 at http://www.wmo.int/pages/prog/www/IMOP/publications/IOM-19-Synoptic-AWS.pdf
          begin
            k1 := 0.190284; // discrepency with calculated value probably because Smithsonian used less precise gas constant and gravity values
            k2 := 4.30899E-5; // (standardLapseRate / 288) * (Power(standardSlpInHg, k1SMT));
            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; // unknown algorithm
      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;
    // see ASOS training at http://www.nwstc.noaa.gov
    begin
      // from US units ASOS formula
      Result := InToHPa(HPaToIn(pressureHPa) * Power10(0.00813 * MToFt(sensorElevationM - stationElevationM)
                                           / FToR(CToF(currentTempC))));
    end;
    {
    // still to do
    class function TWxUtils.AltimeterToStationPressure(pressureHPa: TWxReal; elevationM: TWxReal;
                            algorithm: TAltimeterAlgorithm = DefaultAltimeterAlgorithm): TWxReal;
    begin
    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:
          //  see http://www.univie.ac.at/IMG-Wien/daquamap/Parametergencom.html
          begin
            geopElevationM := GeopotentialAltitude(elevationM);
            Result := Exp(((gravity/gasConstantAir) * geopElevationM)
               / (VirtualTempK(pressureHPa, meanTempC, humidity) + (geopElevationM * standardLapseRate/2)));
          end;
        paDavisVp:
          // see http://www.exploratorium.edu/weather/barometer.html
          begin
            if (humidity > 0) then begin
              hCorr := (9/5) * HumidityCorrection(currentTempC, elevationM, humidity, vaDavisVP);
            end else begin
              hCorr := 0;
            end;
            // In the case of DavisVp, take the constant values literally.
            Result := Power(10, (MToFt(elevationM) / (122.8943111 * (CToF(meanTempC) + 460 + (MToFt(elevationM) * vpLapseRateUS/2) + hCorr))));
          end;
        paManBar:
          // see WMO Instruments and Observing Methods Report No.19 at http://www.wmo.int/pages/prog/www/IMOP/publications/IOM-19-Synoptic-AWS.pdf
          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; // unknown algorithm
      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
      // see http://cires.colorado.edu/~voemel/vp.html   comparison of vapor pressure algorithms
      // see (for DavisVP) http://www.exploratorium.edu/weather/dewpoint.html
      case algorithm of
        vaDavisVp: Result := 6.112 * Exp((17.62 * tempC)/(243.12 + tempC)); // Davis Calculations Doc
        vaBuck: Result := 6.1121 * Exp((18.678 - (tempC/234.5)) * tempC / (257.14 + tempC)); // Buck(1996)
        vaBuck81: Result := 6.1121 * Exp((17.502 * tempC)/(240.97 + tempC)); // Buck(1981)
        vaBolton: Result := 6.112 * Exp(17.67 * tempC / (tempC + 243.5)); // Bolton(1980)
        vaTetenNWS: Result := 6.112 * Power(10,(7.5 * tempC / (tempC + 237.7))); // Magnus Teten see www.srh.weather.gov/elp/wxcalc/formulas/vaporPressure.html
        vaTetenMurray: Result := Power(10, (7.5 * tempC / (237.5 + tempC)) + 0.7858); // Magnus Teten (Murray 1967)
        vaTeten: Result := 6.1078 * Power(10, (7.5 * tempC / (tempC + 237.3))); // Magnus Teten see www.vivoscuola.it/US/RSIGPP3202/umidita/attivita/relhumONA.htm
        else Result := 0; // unknown algorithm
      end;
    end;

    class function TWxUtils.MixingRatio(pressureHPa: TWxReal; tempC: TWxReal;
                        humidity: TWxHumidity): TWxReal;
    var
      vapPres: TWxReal;
    const
      k1 = moleWater / moleAir; // 0.62198
    begin
      // see http://www.wxqa.com/archive/obsman.pdf
      // see also http://www.vivoscuola.it/US/RSIGPP3202/umidita/attivita/relhumONA.htm
      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); // 0.37802
    begin
      // see http://www.univie.ac.at/IMG-Wien/daquamap/Parametergencom.html
      // see also http://www.vivoscuola.it/US/RSIGPP3202/umidita/attivita/relhumONA.htm
      // see also http://wahiduddin.net/calc/density_altitude.htm
      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;
	// see American Meteorological Society Journal
	// see http://www.msc.ec.gc.ca/education/windchill/science_equations_e.cfm
	// see http://www.weather.gov/os/windchill/index.shtml
    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;
    // see http://www.hpc.ncep.noaa.gov/heat_index/hi_equation.html
    var
      tempF: TWxReal;
      tSqrd: TWxReal;
      hum: TWxReal;
      hSqrd: TWxReal;
    begin
      tempF := CToF(tempC);
      if (tempF < 80) then begin
        // heat index algorithm is only valid for temps above 80F
        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);
        // Rothfusz adjustments
        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;

    //************** Conversion Functions ******************

    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;

    //************* General Math Functions ****************

    // raise base to the given power (i.e. base**exponent
    class function TWxUtils.Power(const base, exponent: TWxReal): TWxReal;
    begin
      if exponent = 0.0 then begin
        Result := 1.0; // n**0 = 1
      end else if (base = 0.0) and (exponent > 0.0) then begin
        Result := 0.0; // 0**n = 0, n > 0
      end else begin
        Result := Exp(exponent * Ln(base));
      end;
    end;

    class function TWxUtils.Power10(const exponent: TWxReal): TWxReal;
    const
      ln10 = 2.302585093; // Ln(10);
    begin
      if (exponent = 0.0) then begin
        Result := 1.0;
      end else begin
        Result := Exp(exponent * ln10);
      end;
    end;
    end.