From 42b10e004cc5be1a3892e883b4d8cc3d3d3f7180 Mon Sep 17 00:00:00 2001 From: Patrick Grawehr Date: Sun, 2 Oct 2022 13:01:58 +0200 Subject: [PATCH 1/3] Add a few new WeatherHelper functions --- .../Common/Iot/Device/Common/WeatherHelper.cs | 79 +++++++++++++++++++ src/devices/Common/tests/WeatherTests.cs | 48 +++++++++++ 2 files changed, 127 insertions(+) diff --git a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs index af26066c62..57b323ca84 100644 --- a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs +++ b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs @@ -331,5 +331,84 @@ public static Pressure CalculateBarometricPressure(Pressure measuredPressure, Te } #endregion + + /// + /// Simplified air density (not taking humidity into account) + /// + /// Measured air pressure + /// Measured temperature + /// Approximate standard air density + public static Density CalculateAirDensity(Pressure airPressure, Temperature temperature) + { + double gasConstant = 287.058; // J / (kg * K), for dry air + var result = airPressure.Pascals / (gasConstant * temperature.Kelvins); + return Density.FromKilogramsPerCubicMeter(result); + } + + /// + /// Calculates the air density + /// + /// Measured air pressure + /// Measured temperature + /// Measured relative humidity + /// Approximate standard air density at sea level + /// From https://de.wikipedia.org/wiki/Luftdichte + public static Density CalculateAirDensity(Pressure airPressure, Temperature temperature, RelativeHumidity humidity) + { + double rs = 287.058; + double rd = 461.523; + var pd = CalculateSaturatedVaporPressureOverWater(temperature); + // It's still called "constant" even though it's not constant + double gasConstant = rs / (1 - (humidity.Percent / 100) * (pd.Pascals / airPressure.Pascals) * (1 - (rs / rd))); + var result = airPressure.Pascals / (gasConstant * temperature.Kelvins); + return Density.FromKilogramsPerCubicMeter(result); + } + + /// + /// Calculates the wind chill temperature - this is the perceived temperature in (heavy) winds at cold temperatures. + /// This is only useful at temperatures below about 20°C, above use instead. + /// Not suitable for wind speeds < 5 km/h. + /// + /// The measured air temperature + /// The wind speed (measured at 10m above ground) + /// The perceived temperature. Note that this is not a real temperature, and the skin will never really reach + /// this temperature. This is more an indication on how fast the skin will reach the air temperature. If the skin + /// reaches a temperature of about -5°C, frostbite might occur. + /// From https://de.wikipedia.org/wiki/Windchill + /// The wind speed is less than zero + public static Temperature CalculateWindchill(Temperature temperature, Speed windSpeed) + { + if (windSpeed < Speed.Zero) + { + throw new ArgumentOutOfRangeException(nameof(windSpeed), "The wind speed cannot be negative"); + } + + double va = temperature.DegreesCelsius; + if (windSpeed < Speed.FromKilometersPerHour(1)) + { + // otherwise, the result is complete crap, because the second and third terms of the equation are 0 + windSpeed = Speed.FromKilometersPerHour(1); + } + + double wct = 13.12 + 0.6215 * va + (0.3965 * va - 11.37) * Math.Pow(windSpeed.KilometersPerHour, 0.16); + return Temperature.FromDegreesCelsius(wct); + } + + /// + /// Calculates the wind force on an object. + /// + /// The denisty of the air, calculated using one of the overloads of + /// The speed of the wind + /// Pressure coefficient for the shape of the object. Use 1 for a rectangular object directly facing the wind + /// The Pressure the wind applies on the object + /// From https://de.wikipedia.org/wiki/Winddruck + public static Pressure CalculateWindForce(Density densityOfAir, Speed windSpeed, double pressureCoefficient = 1.0) + { + double v = windSpeed.MetersPerSecond; + double rho = densityOfAir.KilogramsPerCubicMeter; + + double wd = pressureCoefficient * rho / 2 * (v * v); + return Pressure.FromNewtonsPerSquareMeter(wd); + } } } diff --git a/src/devices/Common/tests/WeatherTests.cs b/src/devices/Common/tests/WeatherTests.cs index 59cc246697..16881e8291 100644 --- a/src/devices/Common/tests/WeatherTests.cs +++ b/src/devices/Common/tests/WeatherTests.cs @@ -185,5 +185,53 @@ public void GetRelativeHumidityFromActualAirTemperature(double inTemp, double in Assert.Equal(outHumidityExpected, result.Percent, 3); } + + [Theory] + [InlineData(1013.25, 35, 1.1455)] + [InlineData(1013.25, 0, 1.2922)] + [InlineData(1013.25, -25, 1.4224)] + public void CalculateAirDensitySimple(double inPressure, double inTemp, double expected) + { + var result = WeatherHelper.CalculateAirDensity(Pressure.FromMillibars(inPressure), + Temperature.FromDegreesCelsius(inTemp)); + + Assert.Equal(expected, result.KilogramsPerCubicMeter, 4); + } + + [Theory] + [InlineData(1013.25, 35, 0, 1.1455)] + [InlineData(1013.25, 35, 50, 1.1335)] + [InlineData(1013.25, 35, 100, 1.1214)] + public void CalculateAirDensity(double inPressure, double inTemp, double inHumidity, double expected) + { + var result = WeatherHelper.CalculateAirDensity(Pressure.FromMillibars(inPressure), + Temperature.FromDegreesCelsius(inTemp), RelativeHumidity.FromPercent(inHumidity)); + + Assert.Equal(expected, result.KilogramsPerCubicMeter, 4); + } + + [Theory] + [InlineData(0, 0, 1.75)] // This is to be expected, see WP article + [InlineData(10, 5, 9.8)] + [InlineData(-10, 20, -17.9)] + [InlineData(-15, 40, -27.4)] + public void CalculateWindchill(double temperature, double windSpeed, double expected) + { + var result = WeatherHelper.CalculateWindchill(Temperature.FromDegreesCelsius(temperature), + Speed.FromKilometersPerHour(windSpeed)); + + Assert.Equal(expected, result.DegreesCelsius, 1); + } + + [Theory] + [InlineData(20, 38.5, 1013, 68.83)] + [InlineData(20, 88.9, 1013, 367.0)] + public void CalculateWindForce(double temperature, double windSpeed, double pressure, double expected) + { + Density airDensity = WeatherHelper.CalculateAirDensity(Pressure.FromHectopascals(pressure), + Temperature.FromDegreesCelsius(temperature)); + var density = WeatherHelper.CalculateWindForce(airDensity, Speed.FromKilometersPerHour(windSpeed), 1.0); + Assert.Equal(expected, density.NewtonsPerSquareMeter, 1); + } } } From d7c3df3c8bb2861e1346762ddc334ef4211b7483 Mon Sep 17 00:00:00 2001 From: Patrick Grawehr Date: Tue, 4 Oct 2022 21:41:03 +0200 Subject: [PATCH 2/3] Review findings Improve comments, add constants --- .../Common/Iot/Device/Common/WeatherHelper.cs | 35 +++++++++++++------ src/devices/Common/tests/WeatherTests.cs | 4 ++- 2 files changed, 28 insertions(+), 11 deletions(-) diff --git a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs index 57b323ca84..5db6a156de 100644 --- a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs +++ b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs @@ -12,6 +12,21 @@ namespace Iot.Device.Common /// public static class WeatherHelper { + /// + /// Gas constant of dry Air, J / (kg * K) + /// + public const double SpecificGasConstantOfAir = 287.058; + + /// + /// Gas constant of vapor, J / (kg * K) + /// + public const double SpecificGasConstantOfVapor = 461.523; + + /// + /// Default atmospheric temperature gradient = 0.0065K/m (or 0.65K per 100m) + /// + public const double DefaultTemperatureGradient = 0.0065; + /// /// The mean sea-level pressure (MSLP) is the average atmospheric pressure at mean sea level /// @@ -183,7 +198,7 @@ public static RelativeHumidity GetRelativeHumidityFromActualAirTemperature(Tempe /// The altitude public static Length CalculateAltitude(Pressure pressure, Pressure seaLevelPressure, Temperature airTemperature) { - double meters = ((Math.Pow(seaLevelPressure.Pascals / pressure.Pascals, 1 / 5.255) - 1) * airTemperature.Kelvins) / 0.0065; + double meters = ((Math.Pow(seaLevelPressure.Pascals / pressure.Pascals, 1 / 5.255) - 1) * airTemperature.Kelvins) / DefaultTemperatureGradient; return Length.FromMeters(meters); } @@ -222,7 +237,7 @@ public static Length CalculateAltitude(Pressure pressure) => /// The estimated absolute sea-level pressure /// solved for sea level pressure public static Pressure CalculateSeaLevelPressure(Pressure pressure, Length altitude, Temperature airTemperature) => - Pressure.FromPascals(Math.Pow((((0.0065 * altitude.Meters) / airTemperature.Kelvins) + 1), 5.255) * pressure.Pascals); + Pressure.FromPascals(Math.Pow((((DefaultTemperatureGradient * altitude.Meters) / airTemperature.Kelvins) + 1), 5.255) * pressure.Pascals); /// /// Calculates the approximate absolute pressure from given sea-level pressure, altitude and air temperature. @@ -232,7 +247,7 @@ public static Pressure CalculateSeaLevelPressure(Pressure pressure, Length altit /// The air temperature at the point for which pressure is being calculated /// The estimated absolute pressure at the given altitude public static Pressure CalculatePressure(Pressure seaLevelPressure, Length altitude, Temperature airTemperature) => - Pressure.FromPascals(seaLevelPressure.Pascals / Math.Pow((((0.0065 * altitude.Meters) / airTemperature.Kelvins) + 1), 5.255)); + Pressure.FromPascals(seaLevelPressure.Pascals / Math.Pow((((DefaultTemperatureGradient * altitude.Meters) / airTemperature.Kelvins) + 1), 5.255)); /// /// Calculates the temperature gradient for the given pressure difference @@ -243,7 +258,7 @@ public static Pressure CalculatePressure(Pressure seaLevelPressure, Length altit /// The standard temperature at the given altitude, when the given pressure difference is known /// solved for temperature public static Temperature CalculateTemperature(Pressure pressure, Pressure seaLevelPressure, Length altitude) => - Temperature.FromKelvins((0.0065 * altitude.Meters) / (Math.Pow(seaLevelPressure.Pascals / pressure.Pascals, 1 / 5.255) - 1)); + Temperature.FromKelvins((DefaultTemperatureGradient * altitude.Meters) / (Math.Pow(seaLevelPressure.Pascals / pressure.Pascals, 1 / 5.255) - 1)); /// /// Calculates the barometric pressure from a raw reading, using the reduction formula from the german met service. @@ -302,7 +317,7 @@ public static Pressure CalculateBarometricPressure(Pressure measuredPressure, Te Length measurementAltitude) { double x = (9.80665 / (287.05 * ((measuredTemperature.Kelvins) + 0.12 * vaporPressure.Hectopascals + - (0.0065 * measurementAltitude.Meters) / 2))) * measurementAltitude.Meters; + (DefaultTemperatureGradient * measurementAltitude.Meters) / 2))) * measurementAltitude.Meters; double barometricPressure = measuredPressure.Hectopascals * Math.Exp(x); return Pressure.FromHectopascals(barometricPressure); } @@ -338,10 +353,10 @@ public static Pressure CalculateBarometricPressure(Pressure measuredPressure, Te /// Measured air pressure /// Measured temperature /// Approximate standard air density + /// From https://de.wikipedia.org/wiki/Luftdichte public static Density CalculateAirDensity(Pressure airPressure, Temperature temperature) { - double gasConstant = 287.058; // J / (kg * K), for dry air - var result = airPressure.Pascals / (gasConstant * temperature.Kelvins); + var result = airPressure.Pascals / (SpecificGasConstantOfAir * temperature.Kelvins); return Density.FromKilogramsPerCubicMeter(result); } @@ -355,8 +370,8 @@ public static Density CalculateAirDensity(Pressure airPressure, Temperature temp /// From https://de.wikipedia.org/wiki/Luftdichte public static Density CalculateAirDensity(Pressure airPressure, Temperature temperature, RelativeHumidity humidity) { - double rs = 287.058; - double rd = 461.523; + double rs = SpecificGasConstantOfAir; + double rd = SpecificGasConstantOfVapor; var pd = CalculateSaturatedVaporPressureOverWater(temperature); // It's still called "constant" even though it's not constant double gasConstant = rs / (1 - (humidity.Percent / 100) * (pd.Pascals / airPressure.Pascals) * (1 - (rs / rd))); @@ -386,7 +401,7 @@ public static Temperature CalculateWindchill(Temperature temperature, Speed wind double va = temperature.DegreesCelsius; if (windSpeed < Speed.FromKilometersPerHour(1)) { - // otherwise, the result is complete crap, because the second and third terms of the equation are 0 + // otherwise, the result is unusable, because the second and third terms of the equation are 0, resulting in a constant offset from the input temperature windSpeed = Speed.FromKilometersPerHour(1); } diff --git a/src/devices/Common/tests/WeatherTests.cs b/src/devices/Common/tests/WeatherTests.cs index 16881e8291..6c2b06ca14 100644 --- a/src/devices/Common/tests/WeatherTests.cs +++ b/src/devices/Common/tests/WeatherTests.cs @@ -211,7 +211,9 @@ public void CalculateAirDensity(double inPressure, double inTemp, double inHumid } [Theory] - [InlineData(0, 0, 1.75)] // This is to be expected, see WP article + // When the wind is calm, the windchill temperature is higher than the actual temperature, because the air close to the skin heats up. + // Explained in https://de.wikipedia.org/wiki/Windchill + [InlineData(0, 0, 1.75)] [InlineData(10, 5, 9.8)] [InlineData(-10, 20, -17.9)] [InlineData(-15, 40, -27.4)] From f96fcfcbd5e09e1e423e9d698825682d8ca25b19 Mon Sep 17 00:00:00 2001 From: Patrick Grawehr Date: Thu, 6 Oct 2022 12:24:21 +0200 Subject: [PATCH 3/3] Change members to internal --- src/devices/Common/Iot/Device/Common/WeatherHelper.cs | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs index 5db6a156de..571b355360 100644 --- a/src/devices/Common/Iot/Device/Common/WeatherHelper.cs +++ b/src/devices/Common/Iot/Device/Common/WeatherHelper.cs @@ -15,17 +15,17 @@ public static class WeatherHelper /// /// Gas constant of dry Air, J / (kg * K) /// - public const double SpecificGasConstantOfAir = 287.058; + internal const double SpecificGasConstantOfAir = 287.058; /// /// Gas constant of vapor, J / (kg * K) /// - public const double SpecificGasConstantOfVapor = 461.523; + internal const double SpecificGasConstantOfVapor = 461.523; /// /// Default atmospheric temperature gradient = 0.0065K/m (or 0.65K per 100m) /// - public const double DefaultTemperatureGradient = 0.0065; + internal const double DefaultTemperatureGradient = 0.0065; /// /// The mean sea-level pressure (MSLP) is the average atmospheric pressure at mean sea level