I'm looking for some help checking that my math/reasoning is correct, with regards to calculating the y value for a corresponding and known x value. I am not confident with maths and my formulas have come from researching on the internet and as a result I am aware I may have gone awry without realising.
I am writing an application which allows the user to analyse line of sight between two locations. In order for this to be effective over longer ranges, I wish to simulate the curvature of the Earth (Assuming that the Earth is a perfect circle). This curve will then have Google Elevation API data placed on top of it. Here is an example of the output
The green line represents the line of sight between the start and end location, relative to their heights above the ground. The orange line is my current simulation of the Earth's curvature and the blue line is the imposed elevation data which has been manipulated to rise and fall with appropriate values depending on the y value of the Earth's curvature.
My code uses a class (called Circle) for handling calculations. The code for the circle class, is as follows -
# The circle class accepts two arguments, the radius of the circle and the length of arc between the locations
# which are being assessed. This generates a representation of the Earth's curvature which will be represented on the
# final graph. This is important to represent as the curvature of the Earth makes a big difference to whether
# or not, line of sight exists between locations.
import numpy as np
class Circle:
def __init__(self, radius_of_circle, length_of_arc):
self.radius = radius_of_circle
self.arc_length = length_of_arc
# Define Circle class methods
# Calculate the central angle, in degrees, by using the arc_length
# Gives angle in degrees at centre of the circle between the two points (beginning and end points of arc_length)
# Returns floating point
def calc_degrees(self):
return (self.arc_length / (np.pi * self.calc_diameter())) * 360
# Calculate the central angle in radians, between two points on the circle
# Returns floating point
def calc_radians(self): # Where theta is the angle between both points at the centre of the circle
return np.radians(self.calc_degrees()) # Convert degrees to radians to work with chord_length formula
# Returns the chord lengths of the arc, taking theta (angle in radians) as it's argument
# The chord is the horizontal line which separates the arc segment from the rest of the circle
# Formula for theta (radians) only, not degrees #confirmed using http://www.ambrsoft.com/Trigocalc/Sphere/Arc_.htm
# Returns floating point
def calc_chord_length(self):
return 2 * self.radius * np.sin(self.calc_radians() / 2)
# Calculates the length of arc, taking theta (angle in radians) as its argument.
# Confirmed using http://www.ambrsoft.com/Trigocalc/Sphere/Arc_.htm
# Returns floating point
def calc_arc_length(self):
return (self.calc_degrees() / 360) * self.calc_diameter() * np.pi
# Calculates the Sagitta of the arc segment. The Sagitta is the distance from the centre of the arc
# to the centre of the chord
# Confirmed correct against online calculator https://www.liutaiomottola.com/formulae/sag.htm
def calc_sagitta(self):
return self.radius - (np.sqrt((self.radius ** 2) - ((self.calc_chord_length() / 2) ** 2)))
# Calculate the distance between the chord of the segment and the centre of the circle
# Returns floating point
def calc_arc_apothem(self):
return self.radius - self.calc_sagitta()
# Calculate centre point of circle
# Returns floating point
def calc_circular_centre_x(self):
return self.calc_chord_length() / 2
# Calculate centre point of circle
# Returns a floating point number
def calc_circular_centre_y(self):
return self.calc_sagitta() - self.radius
# Calculate the diameter of the circle
# Returns a floating point number which is double the radius.
def calc_diameter(self):
return self.radius * 2
# Returns the starting angle of the circular arc as float
# Takes two arguments, starting y and x coordinates
def calc_start_angle(self, start_y, start_x):
centre_y = self.calc_circular_centre_y()
centre_x = self.calc_circular_centre_x()
return np.arctan2(start_y - centre_y, start_x - centre_x)
# Returns the ending angle of the circular arc as float
# Takes two arguments, ending y and x coordinates
def calc_end_angle(self, end_y, end_x):
centre_y = self.calc_circular_centre_y()
centre_x = self.calc_circular_centre_x()
return np.arctan2(end_y - centre_y, end_x - centre_x)
# Returns a numpy array of y-axis values for mapping on matplotlib graph. x values list is a list of distances
# in nautical miles. Each y-axis value represents the rising and falling of the earth to simulate 'curvature' which
# effects line of sight visibility.
def calculate_earth_surface_y_values(self, list_of_x_axis_values, list_of_angles):
earth_radius = self.radius
y_values_list = []
for j in range(len(list_of_x_axis_values)):
# Calculate the y axis value (height) for the corresponding x value (distance). Subtract the apothem
# of the circle to ensure the arc starts at coordinates 0,0 and ends at zero again on the y axis
y = earth_radius * np.sin(list_of_angles[j]) - self.calc_arc_apothem()
# Convert nautical miles to meters for elevation
y = y * 1852
y_values_list.append(y)
y_values_np = np.array(y_values_list)
return y_values_np
And my main.py is running as follows -
import time
import numpy as np
from haversine import haversine, Unit
from main import location_reference_converter, circle, data_handling
from main.graph_processing import create_graph
from main.kml_generator import create_kml_file
from main.validation_handling import validate_google_sample_number
def run_graphing_and_kml_process(input_file_one, input_file_two, output_folder, api_key, samples, first_file_type,
second_file_type):
# This will have options in future for different units of measurement
earth_radius = 3440.065 # in nautical miles
# Validate the number of samples requested
validate_google_sample_number(samples)
# Create a list of location objects from both .csv files
first_location_list = location_reference_converter.conversion_type(input_file_one, first_file_type)
second_location_list = location_reference_converter.conversion_type(input_file_two, second_file_type)
# Iterate through the locations in the first .csv file
for i in first_location_list:
# Retrieve the first Location object and store its coordinates
pos_1 = (i.latitude, i.longitude)
# Iterate through the locations in the second .csv file
for x in second_location_list:
# Retrieve the second Location object and store its coordinates
pos_2 = (x.latitude, x.longitude)
# Calculate the great circle distance between both locations using the haversine formula
# Currently this is returned in Nm but future functionality will allow other values
great_circle_distance = haversine(pos_1, pos_2, unit=Unit.NAUTICAL_MILES)
# Create a circle object to simulate curvature of the earth.
c1 = circle.Circle(earth_radius, great_circle_distance)
# Set start and end points for representation of the earths curvature
x1, y1 = 0, 0
x2, y2 = great_circle_distance, 0
angle_list = np.linspace(c1.calc_start_angle(0, 0), c1.calc_end_angle(0, great_circle_distance), samples)
x_values = np.linspace(x1, x2, samples)
# Create the y-axis values to draw the earth's curved surface
y_values = c1.calculate_earth_surface_y_values(x_values, angle_list)
# Construct api url and extract the elevation data
elevation_data = data_handling.send_and_receive_data(i.coordinates_string, x.coordinates_string, samples,
api_key, y_values)
create_graph(x_values, y_values, elevation_data, great_circle_distance, i, x, output_folder)
# Rest for a moment to prevent the api being bombarded with requests.
time.sleep(2)
create_kml_file(first_location_list, second_location_list, "First Locations", "Second Locations", output_folder)
return "Graphing and .kml processes complete."
Is my calculation, for graphing the assumed rise and fall of the simulated Earth surface (assuming a perfect circle) correct? I am aware that at longer distances, a small margin of error could really mess up the accuracy of the results.
Edit: link to GitHub page for full code
https://github.com/MHenderson1988/line-of-sight-analysis/tree/master/main