Distance Calculations

from geopy.distance import Distance

class Distance:
    """
    Base class for distance measurements with unit conversions.
    Supports arithmetic operations and destination calculations.
    """
    
    def __init__(self, *args, **kwargs):
        """Initialize distance measurement"""
    
    @property
    def kilometers(self) -> float:
        """Distance in kilometers"""
    
    @property
    def km(self) -> float:
        """Alias for kilometers"""
    
    @property
    def meters(self) -> float:
        """Distance in meters"""
    
    @property
    def m(self) -> float:
        """Alias for meters"""
    
    @property
    def miles(self) -> float:
        """Distance in miles"""
    
    @property
    def mi(self) -> float:
        """Alias for miles"""
    
    @property
    def feet(self) -> float:
        """Distance in feet"""
    
    @property
    def ft(self) -> float:
        """Alias for feet"""
    
    @property
    def nautical(self) -> float:
        """Distance in nautical miles"""
    
    @property
    def nm(self) -> float:
        """Alias for nautical miles"""
    
    def destination(self, point, bearing, distance=None):
        """
        Calculate destination point from bearing and distance.
        
        Parameters:
        - point (Point): Starting point
        - bearing (float): Bearing in degrees (0-360)
        - distance (Distance): Distance to travel (uses self if None)
        
        Returns:
        Point: Destination point
        """
    
    def __add__(self, other):
        """Addition with another distance"""
    
    def __sub__(self, other):
        """Subtraction with another distance"""
    
    def __mul__(self, scalar):
        """Multiplication by scalar"""
    
    def __truediv__(self, scalar):
        """Division by scalar"""
    
    def __lt__(self, other):
        """Less than comparison"""
    
    def __le__(self, other):
        """Less than or equal comparison"""
    
    def __eq__(self, other):
        """Equality comparison"""
    
    def __ne__(self, other):
        """Inequality comparison"""
    
    def __gt__(self, other):
        """Greater than comparison"""
    
    def __ge__(self, other):
        """Greater than or equal comparison"""

from geopy.distance import geodesic

class geodesic(Distance):
    """
    Calculate geodesic distance using ellipsoidal earth model.
    Most accurate distance calculation method.
    """
    
    def __init__(self, point1=None, point2=None, ellipsoid='WGS-84'):
        """
        Initialize geodesic distance calculation.
        
        Parameters:
        - point1 (Point): First point
        - point2 (Point): Second point  
        - ellipsoid (str): Reference ellipsoid name
        
        Supported ellipsoids:
        - 'WGS-84' (default) - World Geodetic System 1984
        - 'GRS-80' - Geodetic Reference System 1980
        - 'Airy (1830)' - Airy ellipsoid 1830
        - 'Intl 1924' - International ellipsoid 1924
        - 'Clarke (1880)' - Clarke ellipsoid 1880
        - 'GRS-67' - Geodetic Reference System 1967
        """
    
    def measure(self, point1, point2):
        """
        Measure distance between two points.
        
        Parameters:
        - point1 (Point): First point
        - point2 (Point): Second point
        
        Returns:
        geodesic: Distance measurement
        """
    
    def destination(self, point, bearing, distance=None):
        """
        Calculate destination point using geodesic calculation.
        
        Parameters:
        - point (Point): Starting point
        - bearing (float): Initial bearing in degrees (0-360)
        - distance (Distance): Distance to travel
        
        Returns:
        Point: Destination point with precise coordinates
        """

from geopy.distance import great_circle

class great_circle(Distance):
    """
    Calculate great-circle distance using spherical earth model.
    Faster but less accurate than geodesic calculations.
    """
    
    def __init__(self, point1=None, point2=None, radius=EARTH_RADIUS):
        """
        Initialize great-circle distance calculation.
        
        Parameters:
        - point1 (Point): First point
        - point2 (Point): Second point
        - radius (float): Earth radius in kilometers (default: 6371.009)
        """
    
    def measure(self, point1, point2):
        """
        Measure great-circle distance between two points.
        
        Parameters:
        - point1 (Point): First point
        - point2 (Point): Second point
        
        Returns:
        great_circle: Distance measurement
        """
    
    def destination(self, point, bearing, distance=None):
        """
        Calculate destination point using great-circle calculation.
        
        Parameters:
        - point (Point): Starting point
        - bearing (float): Bearing in degrees (0-360)
        - distance (Distance): Distance to travel
        
        Returns:
        Point: Destination point
        """

from geopy.distance import EARTH_RADIUS, ELLIPSOIDS, distance

# Constants
EARTH_RADIUS = 6371.009  # Mean earth radius in kilometers

ELLIPSOIDS = {
    'WGS-84': (6378137.0, 6356752.314245, 298.257223563),
    'GRS-80': (6378137.0, 6356752.314140, 298.257222101),
    'Airy (1830)': (6377563.396, 6356256.909, 299.3249646),
    'Intl 1924': (6378388.0, 6356911.946, 297.0),
    'Clarke (1880)': (6378249.145, 6356514.870, 293.465),
    'GRS-67': (6378160.0, 6356774.719, 298.25)
}

# Default distance algorithm (alias for geodesic)
distance = geodesic

from geopy.distance import lonlat

def lonlat(x, y, z=0):
    """
    Convert (longitude, latitude) to Point(latitude, longitude).
    
    Utility function to handle coordinate order differences.
    Many APIs return coordinates as [longitude, latitude] but
    Point expects (latitude, longitude).
    
    Parameters:
    - x (float): Longitude
    - y (float): Latitude
    - z (float): Altitude (optional)
    
    Returns:
    Point: Point(latitude, longitude, altitude)
    """

from geopy.distance import geodesic, great_circle
from geopy.point import Point

# Define some locations
new_york = Point(40.7128, -74.0060)
los_angeles = Point(34.0522, -118.2437)
london = Point(51.5074, -0.1278)

# Calculate distances using different algorithms
geo_distance = geodesic(new_york, los_angeles)
gc_distance = great_circle(new_york, los_angeles)

print(f"Geodesic distance: {geo_distance.miles:.2f} miles")
print(f"Great circle distance: {gc_distance.miles:.2f} miles")
print(f"Difference: {abs(geo_distance.miles - gc_distance.miles):.2f} miles")

# Access different units
print(f"Distance in kilometers: {geo_distance.km:.2f}")
print(f"Distance in meters: {geo_distance.m:.2f}")
print(f"Distance in nautical miles: {geo_distance.nautical:.2f}")
print(f"Distance in feet: {geo_distance.feet:.2f}")

from geopy.distance import geodesic

# Calculate multiple distances
ny_to_la = geodesic(new_york, los_angeles)
la_to_london = geodesic(los_angeles, london)

# Arithmetic operations
total_distance = ny_to_la + la_to_london
half_distance = ny_to_la / 2
double_distance = ny_to_la * 2

print(f"NY to LA: {ny_to_la.miles:.0f} miles")
print(f"LA to London: {la_to_london.miles:.0f} miles")
print(f"Total journey: {total_distance.miles:.0f} miles")
print(f"Half distance: {half_distance.miles:.0f} miles")

# Distance comparisons
if ny_to_la > la_to_london:
    print("NY to LA is farther than LA to London")
else:
    print("LA to London is farther than NY to LA")

# Sort locations by distance from reference point
reference = Point(41.8781, -87.6298)  # Chicago
cities = [
    ("New York", new_york),
    ("Los Angeles", los_angeles), 
    ("London", london)
]

cities_with_distance = [
    (name, point, geodesic(reference, point))
    for name, point in cities
]

# Sort by distance
cities_sorted = sorted(cities_with_distance, key=lambda x: x[2].kilometers)

print("\nCities sorted by distance from Chicago:")
for name, point, dist in cities_sorted:
    print(f"{name}: {dist.miles:.0f} miles")

from geopy.distance import geodesic
from geopy.point import Point

# Starting point
start = Point(40.7128, -74.0060)  # New York

# Calculate destinations at different bearings and distances
distance_100km = geodesic(kilometers=100)
distance_200mi = geodesic(miles=200)

# Cardinal directions (0° = North, 90° = East, 180° = South, 270° = West)
destinations = {}
bearings = [
    (0, "North"),
    (90, "East"), 
    (180, "South"),
    (270, "West"),
    (45, "Northeast"),
    (135, "Southeast"),
    (225, "Southwest"),
    (315, "Northwest")
]

print("Destinations 100km from New York:")
for bearing, direction in bearings:
    dest = distance_100km.destination(start, bearing)
    destinations[direction] = dest
    print(f"{direction} ({bearing}°): {dest.format_decimal()}")

# Calculate round trip
outbound = geodesic(start, destinations["North"])
# Bearing back is opposite direction (180° difference)
home = outbound.destination(destinations["North"], 180)

print(f"\nRound trip verification:")
print(f"Original: {start.format_decimal()}")
print(f"Returned: {home.format_decimal()}")
print(f"Difference: {geodesic(start, home).meters:.2f} meters")

from geopy.distance import geodesic
from geopy.point import Point
import math

def calculate_area_of_triangle(p1, p2, p3):
    """Calculate area of triangle using Heron's formula with geodesic distances"""
    # Calculate side lengths
    a = geodesic(p2, p3).meters
    b = geodesic(p1, p3).meters  
    c = geodesic(p1, p2).meters
    
    # Semi-perimeter
    s = (a + b + c) / 2
    
    # Heron's formula
    area = math.sqrt(s * (s - a) * (s - b) * (s - c))
    return area  # in square meters

def find_center_point(points):
    """Find geographic center (centroid) of multiple points"""
    if not points:
        return None
    
    # Convert to cartesian coordinates for averaging
    x_total = y_total = z_total = 0
    
    for point in points:
        lat_rad = math.radians(point.latitude)
        lon_rad = math.radians(point.longitude)
        
        x = math.cos(lat_rad) * math.cos(lon_rad)
        y = math.cos(lat_rad) * math.sin(lon_rad)
        z = math.sin(lat_rad)
        
        x_total += x
        y_total += y
        z_total += z
    
    # Average
    x_avg = x_total / len(points)
    y_avg = y_total / len(points)
    z_avg = z_total / len(points)
    
    # Convert back to lat/lon
    hyp = math.sqrt(x_avg**2 + y_avg**2)
    lat = math.atan2(z_avg, hyp)
    lon = math.atan2(y_avg, x_avg)
    
    return Point(math.degrees(lat), math.degrees(lon))

def create_circle_of_points(center, radius_km, num_points=8):
    """Create circle of points around a center at specified radius"""
    points = []
    distance_obj = geodesic(kilometers=radius_km)
    
    for i in range(num_points):
        bearing = (360 / num_points) * i
        point = distance_obj.destination(center, bearing)
        points.append(point)
    
    return points

# Example usage
triangle_points = [
    Point(40.7128, -74.0060),  # New York
    Point(41.8781, -87.6298),  # Chicago  
    Point(39.7392, -104.9903)  # Denver
]

area = calculate_area_of_triangle(*triangle_points)
center = find_center_point(triangle_points)
circle = create_circle_of_points(center, 500, 12)  # 500km radius, 12 points

print(f"Triangle area: {area/1e6:.2f} square kilometers")
print(f"Geographic center: {center.format_decimal()}")
print(f"Circle points around center:")
for i, point in enumerate(circle):
    bearing = (360 / len(circle)) * i
    print(f"  {bearing:3.0f}°: {point.format_decimal()}")

from geopy.distance import geodesic, great_circle
from geopy.point import Point
import time

def benchmark_distance_calculations():
    """Compare performance of geodesic vs great-circle calculations"""
    
    # Test data
    points = [
        Point(40.7128, -74.0060),  # New York
        Point(34.0522, -118.2437), # Los Angeles
        Point(51.5074, -0.1278),   # London
        Point(35.6762, 139.6503),  # Tokyo
        Point(-33.8688, 151.2093), # Sydney
    ]
    
    num_calculations = 1000
    
    # Benchmark geodesic
    start_time = time.time()
    for _ in range(num_calculations):
        for i in range(len(points)):
            for j in range(i+1, len(points)):
                geodesic(points[i], points[j])
    geodesic_time = time.time() - start_time
    
    # Benchmark great circle
    start_time = time.time()
    for _ in range(num_calculations):
        for i in range(len(points)):
            for j in range(i+1, len(points)):
                great_circle(points[i], points[j])
    great_circle_time = time.time() - start_time
    
    print(f"Geodesic calculations: {geodesic_time:.3f} seconds")
    print(f"Great circle calculations: {great_circle_time:.3f} seconds")
    print(f"Speed ratio: {geodesic_time/great_circle_time:.1f}x")
    
    # Accuracy comparison
    p1, p2 = points[0], points[1]
    geo_dist = geodesic(p1, p2)
    gc_dist = great_circle(p1, p2)
    
    print(f"\nAccuracy comparison (NY to LA):")
    print(f"Geodesic: {geo_dist.kilometers:.6f} km")
    print(f"Great circle: {gc_dist.kilometers:.6f} km")
    print(f"Difference: {abs(geo_dist.kilometers - gc_dist.kilometers):.6f} km")
    print(f"Relative error: {abs(geo_dist.kilometers - gc_dist.kilometers)/geo_dist.kilometers*100:.4f}%")

# Run benchmark
benchmark_distance_calculations()

from geopy.distance import geodesic, ELLIPSOIDS

# View available ellipsoids
print("Available ellipsoids:")
for name, (a, b, f) in ELLIPSOIDS.items():
    print(f"  {name}: a={a}m, b={b}m, f=1/{f:.1f}")

# Compare distances using different ellipsoids
p1 = Point(40.7128, -74.0060)  # New York
p2 = Point(51.5074, -0.1278)   # London

ellipsoid_results = {}
for ellipsoid_name in ELLIPSOIDS.keys():
    distance = geodesic(p1, p2, ellipsoid=ellipsoid_name)
    ellipsoid_results[ellipsoid_name] = distance.kilometers

print(f"\nDistance NY to London using different ellipsoids:")
for name, dist in ellipsoid_results.items():
    print(f"  {name}: {dist:.6f} km")

# Find the range
min_dist = min(ellipsoid_results.values())
max_dist = max(ellipsoid_results.values())
print(f"\nRange: {max_dist - min_dist:.6f} km")
print(f"Relative variation: {(max_dist - min_dist)/min_dist*100:.4f}%")