Skip to content

Commit

Permalink
ENH: Transform lat/lon into properties
Browse files Browse the repository at this point in the history
  • Loading branch information
Gui-FernandesBR committed Oct 9, 2022
1 parent fb1a95f commit ff665d7
Showing 1 changed file with 79 additions and 113 deletions.
192 changes: 79 additions & 113 deletions rocketpy/Flight.py
Original file line number Diff line number Diff line change
Expand Up @@ -640,8 +640,10 @@ def __init__(
self.calculatedRailButtonForces = False
self.calculatedGeodesicCoordinates = False
self.calculatedPressureSignals = False
self.latitude = Function(0)
self.longitude = Function(0)
self._drift = Function(0)
self.bearing = Function(0)
self._latitude = Function(0)
self._longitude = Function(0)
# Initialize solver monitors
self.functionEvaluations = []
self.functionEvaluationsPerTimeStep = []
Expand Down Expand Up @@ -2860,6 +2862,81 @@ def maxRailButton2ShearForce(self):
else:
return np.max(self.railButton2ShearForce[:outOfRailTimeIndex])

@cached_property
def drift(self):
"""Rocket horizontal distance to tha launch point, in meters, at each
time step.
Returns
-------
drift: Function
Rocket horizontal distance to tha launch point, in meters, at each
time step.
"""
return (self.x**2 + self.y**2) ** 0.5

@cached_property
def bearing(self):
"""Rocket bearing compass, in degrees, at each time step.
Returns
-------
bearing: Function
Rocket bearing, in degrees, at each time step.
"""
return np.arctan2(self.y, self.x) * 180 / np.pi

@cached_property
def latitude(self):
"""Rocket latitude coordinate, in degrees, at each time step.
Returns
-------
latitude: Function
Rocket latitude coordinate, in degrees, at each time step.
"""
lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians

# Applies the haversine equation to find final lat/lon coordinates
latitude = np.rad2deg(
math.asin(
math.sin(lat1) * math.cos(self.drift / self.env.earthRadius)
+ math.cos(lat1)
* math.sin(self.drift / self.env.earthRadius)
* math.cos(self.bearing)
)
)
latitude = Function(latitude, "Time (s)", "Latitude (°)", "linear")

return latitude

@cached_property
def longitude(self):
"""Rocket longitude coordinate, in degrees, at each time step.
Returns
-------
longitude: Function
Rocket longitude coordinate, in degrees, at each time step.
"""
lat1 = np.deg2rad(self.env.lat) # Launch lat point converted to radians
lon1 = np.deg2rad(self.env.lon) # Launch lon point converted to radians

# Applies the haversine equation to find final lat/lon coordinates
longitude = np.rad2deg(
lon1
+ math.atan2(
math.sin(self.bearing)
* math.sin(self.drift / self.env.earthRadius)
* math.cos(lat1),
math.cos(self.drift / self.env.earthRadius)
- math.sin(lat1) * math.sin(self.latitude),
)
)
longitude = Function(longitude, "Time (s)", "Longitude (°)", "linear")

return longitude

@cached_property
def __calculate_rail_button_forces(self):
"""Calculate the forces applied to the rail buttons while rocket is still
Expand Down Expand Up @@ -2909,113 +2986,6 @@ def __calculate_rail_button_forces(self):

return F11, F12, F21, F22

def __calculate_geodesic_coordinates(self):
"""Calculate the geodesic coordinates (i.e. latitude and longitude) of
the rocket during the whole flight.
Applies the Haversine equation considering a WGS84 ellipsoid.
Returns
-------
None
"""
# Converts x and y positions to lat and lon
# We are currently considering the earth as a sphere.
bearing = []
distance = [
(self.x[i][1] ** 2 + self.y[i][1] ** 2) ** 0.5 for i in range(len(self.x))
]
for i in range(len(self.x)):
# Check if the point is over the grid (i.e. if x*y == 0)
if self.x[i][1] == 0:
if self.y[i][1] < 0:
bearing.append(3.14159265359)
else:
bearing.append(0)
continue
if self.y[i][1] == 0:
if self.x[i][1] < 0:
bearing.append(3 * 3.14159265359 / 2)
elif self.x[i][1] > 0:
bearing.append(3.14159265359 / 2)
else:
continue
continue

# Calculate bearing as the azimuth considering different quadrants
if self.x[i][1] * self.y[i][1] > 0 and self.x[i][1] > 0:
bearing.append(math.atan(abs(self.x[i][1]) / abs(self.y[i][1])))
elif self.x[i][1] * self.y[i][1] < 0 and self.x[i][1] > 0:
bearing.append(
3.14159265359 / 2 + math.atan(abs(self.y[i][1]) / abs(self.x[i][1]))
)
elif self.x[i][1] * self.y[i][1] > 0 and self.x[i][1] < 0:
bearing.append(
3.14159265359 + math.atan(abs(self.x[i][1]) / abs(self.y[i][1]))
)
elif self.x[i][1] * self.y[i][1] < 0 and self.x[i][1] < 0:
bearing.append(
3 * 3.14159265359 / 2
+ math.atan(abs(self.y[i][1]) / abs(self.x[i][1]))
)

# Store values of distance and bearing using appropriate units
# self.distance = distance # Must be in meters
# self.bearing = bearing # Must be in radians

lat1 = (
3.14159265359 * self.env.lat / 180
) # Launch lat point converted to radians
lon1 = (
3.14159265359 * self.env.lon / 180
) # Launch long point converted to radians

R = self.env.earthRadius
lat2 = []
lon2 = []
# Applies the haversine equation to find final lat/lon coordinates
for i in range(len(self.x)):
# Please notice that for distances lower than 1 centimeter the difference on latitude or longitude too small
if abs(self.y[i][1]) < 1e-2:
lat2.append(self.env.lat)
else:
lat2.append(
(180 / 3.14159265359)
* math.asin(
math.sin(lat1) * math.cos(distance[i] / R)
+ math.cos(lat1)
* math.sin(distance[i] / R)
* math.cos(bearing[i])
)
)
if abs(self.x[i][1]) < 1e-2:
lon2.append(self.env.lon)
else:
lon2.append(
(180 / 3.14159265359)
* (
lon1
+ math.atan2(
math.sin(bearing[i])
* math.sin(distance[i] / R)
* math.cos(lat1),
math.cos(distance[i] / R)
- math.sin(lat1) * math.sin(lat2[i]),
)
)
)

latitude = [[self.solution[i][0], lat2[i]] for i in range(len(self.solution))]
longitude = [[self.solution[i][0], lon2[i]] for i in range(len(self.solution))]

# Store final values of lat/lon as a function of time
# TODO: Must be converted to properties
self.latitude = Function(latitude, "Time (s)", "Latitude (°)", "linear")
self.longitude = Function(longitude, "Time (s)", "Longitude (°)", "linear")

self.calculatedGeodesicCoordinates = True

return None

def __calculate_pressure_signal(self):
# Transform parachute sensor feed into functions
for parachute in self.rocket.parachutes:
Expand Down Expand Up @@ -3054,10 +3024,6 @@ def postProcess(self, interpolation="spline", extrapolation="natural"):
None
"""

# Post process other quantities # TODO: Implement boolean to check if already done
if self.calculatedGeodesicCoordinates is not True:
self.__calculate_geodesic_coordinates()

if self.calculatedPressureSignals is not True:
self.__calculate_pressure_signal()

Expand Down

0 comments on commit ff665d7

Please sign in to comment.