|
@@ -0,0 +1,148 @@
|
|
|
+#!/usr/bin/env python
|
|
|
+
|
|
|
+from aman.com.Weather import Weather
|
|
|
+
|
|
|
+import math
|
|
|
+import scipy.interpolate
|
|
|
+
|
|
|
+class WeatherModel:
|
|
|
+ def __init__(self, gaforId, weather : Weather):
|
|
|
+ self.gafor = gaforId
|
|
|
+ self.weather = weather
|
|
|
+ self.windDirectionModel = None
|
|
|
+ self.windSpeedModel = None
|
|
|
+ self.lastWeatherUpdate = None
|
|
|
+ self.minimumAltitude = 1000000
|
|
|
+ self.maximumAltitude = -1
|
|
|
+
|
|
|
+ # create the density interpolation model
|
|
|
+ # the density model is based on https://aerotoolbox.com/atmcalc/
|
|
|
+ altitudes = [
|
|
|
+ 50000,
|
|
|
+ 45000,
|
|
|
+ 40000,
|
|
|
+ 38000,
|
|
|
+ 36000,
|
|
|
+ 34000,
|
|
|
+ 32000,
|
|
|
+ 30000,
|
|
|
+ 28000,
|
|
|
+ 26000,
|
|
|
+ 24000,
|
|
|
+ 22000,
|
|
|
+ 20000,
|
|
|
+ 18000,
|
|
|
+ 16000,
|
|
|
+ 15000,
|
|
|
+ 14000,
|
|
|
+ 13000,
|
|
|
+ 12000,
|
|
|
+ 11000,
|
|
|
+ 10000,
|
|
|
+ 9000,
|
|
|
+ 8000,
|
|
|
+ 7000,
|
|
|
+ 6000,
|
|
|
+ 5000,
|
|
|
+ 4000,
|
|
|
+ 3000,
|
|
|
+ 2000,
|
|
|
+ 1000,
|
|
|
+ 0
|
|
|
+ ]
|
|
|
+ densities = [
|
|
|
+ 0.18648,
|
|
|
+ 0.23714,
|
|
|
+ 0.24617,
|
|
|
+ 0.33199,
|
|
|
+ 0.36518,
|
|
|
+ 0.39444,
|
|
|
+ 0.42546,
|
|
|
+ 0.45831,
|
|
|
+ 0.402506,
|
|
|
+ 0.432497,
|
|
|
+ 0.464169,
|
|
|
+ 0.60954,
|
|
|
+ 0.65269,
|
|
|
+ 0.69815,
|
|
|
+ 0.74598,
|
|
|
+ 0.77082,
|
|
|
+ 0.79628,
|
|
|
+ 0.82238,
|
|
|
+ 0.84914,
|
|
|
+ 0.87655,
|
|
|
+ 0.90464,
|
|
|
+ 0.93341,
|
|
|
+ 0.96287,
|
|
|
+ 0.99304,
|
|
|
+ 1.02393,
|
|
|
+ 1.05555,
|
|
|
+ 1.08791,
|
|
|
+ 1.12102,
|
|
|
+ 1.1549,
|
|
|
+ 1.18955,
|
|
|
+ 1.225
|
|
|
+ ]
|
|
|
+ self.densityModel = scipy.interpolate.interp1d(altitudes, densities)
|
|
|
+
|
|
|
+ def calculateTAS(self, altitude : int, ias : int):
|
|
|
+ if altitude >= 50000:
|
|
|
+ altitude = 49999
|
|
|
+ if altitude <= 0:
|
|
|
+ altitude = 1
|
|
|
+
|
|
|
+ # calculation based on https://aerotoolbox.com/airspeed-conversions/
|
|
|
+ return ias * math.sqrt(1.225 / self.densityModel(altitude).item())
|
|
|
+
|
|
|
+ def updateWindModel(self):
|
|
|
+ if None == self.lastWeatherUpdate or self.lastWeatherUpdate != self.weather.provider.updateTime:
|
|
|
+ self.lastWeatherUpdate = self.weather.provider.updateTime
|
|
|
+
|
|
|
+ self.minimumAltitude = 1000000
|
|
|
+ self.maximumAltitude = -1
|
|
|
+ self.windDirectionModel = None
|
|
|
+ self.windSpeedModel = None
|
|
|
+
|
|
|
+ if None != self.weather.provider.windData and self.gafor in self.weather.provider.windData:
|
|
|
+ altitudes = []
|
|
|
+ directions = []
|
|
|
+ speeds = []
|
|
|
+
|
|
|
+ # collect the data for the wind model
|
|
|
+ for level in self.weather.provider.windData[self.gafor]:
|
|
|
+ altitudes.append(level[0])
|
|
|
+ directions.append(level[1])
|
|
|
+ speeds.append(level[2])
|
|
|
+
|
|
|
+ # define the thresholds for later boundary checks
|
|
|
+ if self.minimumAltitude > level[0]:
|
|
|
+ self.minimumAltitude = level[0]
|
|
|
+ if self.maximumAltitude < level[0]:
|
|
|
+ self.maximumAltitude = level[0]
|
|
|
+
|
|
|
+ # calculate the models
|
|
|
+ if 1 < len(altitudes):
|
|
|
+ self.windDirectionModel = scipy.interpolate.interp1d(altitudes, directions)
|
|
|
+ self.windSpeedModel = scipy.interpolate.interp1d(altitudes, speeds)
|
|
|
+
|
|
|
+ def calculateGS(self, altitude : int, ias : int, heading : int):
|
|
|
+ self.updateWindModel()
|
|
|
+ tas = self.calculateTAS(altitude, ias)
|
|
|
+
|
|
|
+ # initialize the wind data
|
|
|
+ if None != self.windDirectionModel and None != self.windSpeedModel:
|
|
|
+ direction = 0.0
|
|
|
+ speed = 0.0
|
|
|
+ if None != self.windSpeedModel and None != self.windDirectionModel:
|
|
|
+ if self.maximumAltitude <= altitude:
|
|
|
+ altitude = self.maximumAltitude - 1
|
|
|
+ if self.minimumAltitude >= altitude:
|
|
|
+ altitude = self.minimumAltitude + 1
|
|
|
+ direction = self.windDirectionModel(altitude).item()
|
|
|
+ speed = self.windSpeedModel(altitude).item()
|
|
|
+ else:
|
|
|
+ speed = 0
|
|
|
+ direction = 0
|
|
|
+
|
|
|
+ # calculate the ground speed based on the headwind component
|
|
|
+ return tas + speed * math.cos(math.radians(direction) - math.radians(heading))
|