Tools & Software

# Arch Linux
sudo pacman -S proj

# Ubuntu/Debian
sudo apt install proj-bin proj-data

# macOS
brew install proj

# Verify installation
proj -V
cs2cs -V
geod -V

+proj=PROJECTION  +datum=DATUM  +zone=ZONE  +ellps=ELLIPSOID
+towgs84=dx,dy,dz,rx,ry,rz,s  (7-parameter transformation)

# Geographic to UTM Zone 11N (Los Angeles)
echo "-118.4081 33.9425" | cs2cs \
    +proj=longlat +datum=WGS84 \
    +to +proj=utm +zone=11 +datum=WGS84
# Output: 368916.29   3757142.35 0.00

# UTM to Geographic (output in DMS)
echo "368916 3757142" | cs2cs \
    +proj=utm +zone=11 +datum=WGS84 \
    +to +proj=longlat +datum=WGS84 -f "%.6f"
# Output: -118.408102   33.942503 0.000000

# UTM Zone 38N (Iraq/Afghanistan operations area)
echo "44.3661 33.3128" | cs2cs \
    +proj=longlat +datum=WGS84 \
    +to +proj=utm +zone=38 +datum=WGS84
# Output: 500000.xx   3686797.xx

# NAD27 to WGS84 (CONUS)
echo "-98.23 26.20" | cs2cs \
    +proj=longlat +datum=NAD27 \
    +to +proj=longlat +datum=WGS84 -f "%.8f"

# ED50 to WGS84 (Europe)
echo "13.4049 52.5200" | cs2cs \
    +proj=longlat +datum=ED50 \
    +to +proj=longlat +datum=WGS84 -f "%.8f"

# Tokyo datum to WGS84 (Japan)
echo "139.6917 35.6895" | cs2cs \
    +proj=longlat +ellps=bessel +towgs84=-148,507,685,0,0,0,0 \
    +to +proj=longlat +datum=WGS84 -f "%.8f"

# Geographic to Texas South Central (EPSG:2278, feet)
echo "-98.2300 26.2034" | cs2cs \
    +proj=longlat +datum=WGS84 \
    +to +init=epsg:2278 -f "%.2f"

# List available EPSG codes
projinfo EPSG:2278

# California Zone 5 (meters)
echo "-118.25 34.05" | cs2cs \
    +proj=longlat +datum=WGS84 \
    +to +init=epsg:2229 -f "%.3f"

# Given: Start point, bearing, distance
# Find: End point, back-azimuth

# Start at Los Angeles, bearing 066°, distance 3970 km
echo "33.9425 -118.4081 66 3970" | geod +ellps=WGS84 +units=km -f "%.6f"
# Output: lat lon back-azimuth

# Military convoy movement: 10 km bearing 045°
echo "33.3128 44.3661 45 10" | geod +ellps=WGS84 +units=km -f "%.6f"

# Given: Two points
# Find: Distance, forward azimuth, back azimuth

# Los Angeles to JFK
geod +ellps=WGS84 -I +units=km <<< "33.9425 -118.4081 40.6413 -73.7781"
# Output: 66°xx'xx.xx" 250°xx'xx.xx" 3970.xxx

# Baghdad to Kabul
geod +ellps=WGS84 -I +units=km <<< "33.3128 44.3661 34.5553 69.2075"
# Output: bearing1 bearing2 distance

# Multiple points (batch processing)
cat << 'EOF' | geod +ellps=WGS84 -I +units=km
33.9425 -118.4081 40.6413 -73.7781
34.0522 -118.2437 51.5074 -0.1278
40.7128 -74.0060 48.8566 2.3522
EOF

# Generate waypoints along great circle route
# LA to JFK, 10 intermediate points
echo "33.9425 -118.4081 40.6413 -73.7781" | \
    geod +ellps=WGS84 -I +units=km -G 10 -f "%.6f"

# Create KML-compatible waypoint list
geod +ellps=WGS84 -I +units=km -G 20 <<< "33.9425 -118.4081 40.6413 -73.7781" | \
    awk 'NR>1 {printf "%.6f,%.6f,0\n", $2, $1}'

# Mercator projection
echo "-98.23 26.20" | proj +proj=merc +ellps=WGS84

# Transverse Mercator (UTM without zone)
echo "-98.23 26.20" | proj +proj=tmerc +lon_0=-99 +lat_0=0 +k=0.9996 +ellps=WGS84

# Lambert Conformal Conic (aviation charts)
echo "-98.23 26.20" | proj +proj=lcc +lat_1=33 +lat_2=45 +lat_0=39 +lon_0=-96 +ellps=WGS84

# Stereographic (polar regions)
echo "0 90" | proj +proj=stere +lat_0=90 +lon_0=0 +ellps=WGS84

# Convert file of coordinates
cs2cs +proj=longlat +datum=WGS84 +to +proj=utm +zone=11 +datum=WGS84 \
    < input_latlon.txt > output_utm.txt

# AWK integration for formatted output
cat coordinates.txt | while read lon lat; do
    echo "$lon $lat" | cs2cs +proj=longlat +datum=WGS84 \
        +to +proj=utm +zone=11 +datum=WGS84 | \
        awk -v lon="$lon" -v lat="$lat" '{printf "%s %s -> %.2f E %.2f N\n", lat, lon, $1, $2}'
done

# Parallel processing for large datasets
parallel --pipe -N1000 "cs2cs +proj=longlat +datum=WGS84 +to +proj=utm +zone=11" \
    < massive_coords.txt > output.txt

# List available timezones
timedatectl list-timezones | grep America

# Show current settings
timedatectl status

# Set timezone
sudo timedatectl set-timezone America/Chicago

# Current time in different zones
TZ=UTC date
TZ=America/Los_Angeles date
TZ=America/Chicago date
TZ=Asia/Tokyo date

# Convert between zones
TZ=America/New_York date -d "TZ=\"UTC\" 2026-03-14 18:00"

pip install geopy pyproj mgrs utm geomag

from geopy.distance import geodesic, great_circle
from geopy.geocoders import Nominatim

# Distance calculation
la = (33.9425, -118.4081)
jfk = (40.6413, -73.7781)

print(f"Geodesic: {geodesic(la, jfk).km:.2f} km")
print(f"Great circle: {great_circle(la, jfk).km:.2f} km")

# Geocoding
geolocator = Nominatim(user_agent="nav-study")
location = geolocator.geocode("McAllen, Texas")
print(f"McAllen: {location.latitude}, {location.longitude}")

from pyproj import CRS, Transformer

# Define coordinate systems
wgs84 = CRS.from_epsg(4326)      # WGS84 Geographic
utm11n = CRS.from_epsg(32611)    # UTM Zone 11N
nad27 = CRS.from_epsg(4267)      # NAD27

# Transform coordinates
to_utm = Transformer.from_crs(wgs84, utm11n, always_xy=True)
lon, lat = -118.4081, 33.9425
x, y = to_utm.transform(lon, lat)
print(f"UTM 11N: {x:.2f} E, {y:.2f} N")

# Datum transformation
nad27_to_wgs84 = Transformer.from_crs(nad27, wgs84, always_xy=True)

import mgrs

m = mgrs.MGRS()

# Geographic to MGRS
lat, lon = 33.9425, -118.4081
mgrs_10 = m.toMGRS(lat, lon, MGRSPrecision=5)  # 10-digit
mgrs_8 = m.toMGRS(lat, lon, MGRSPrecision=4)   # 8-digit
mgrs_6 = m.toMGRS(lat, lon, MGRSPrecision=3)   # 6-digit

print(f"MGRS 10-digit: {mgrs_10}")
print(f"MGRS 8-digit:  {mgrs_8}")
print(f"MGRS 6-digit:  {mgrs_6}")

# MGRS to Geographic
lat, lon = m.toLatLon("11SLA6891657142")

import utm

# Geographic to UTM
lat, lon = 33.9425, -118.4081
easting, northing, zone_number, zone_letter = utm.from_latlon(lat, lon)
print(f"{zone_number}{zone_letter} {easting:.0f}mE {northing:.0f}mN")

# UTM to Geographic
lat, lon = utm.to_latlon(368916, 3757142, 11, 'S')

import geomag
from datetime import date

# Get declination for a location
lat, lon = 26.2034, -98.2300  # McAllen, TX
declination = geomag.declination(lat, lon)

print(f"Declination at McAllen: {declination:.2f}°")
print(f"Direction: {'East' if declination > 0 else 'West'}")

# Get declination for specific date
dec_2026 = geomag.declination(lat, lon, date(2026, 3, 14))

"""
Vincenty's Formulae for Geodesic Calculations
Accuracy: ~0.5mm on the WGS84 ellipsoid
Reference: Vincenty, T. (1975). Direct and Inverse Solutions of Geodesics
           on the Ellipsoid with Application of Nested Equations.
           Survey Review, 23(176), 88-93.
"""

import math
from typing import Tuple, Optional
from dataclasses import dataclass

@dataclass
class WGS84:
    """WGS84 Ellipsoid Parameters"""
    a: float = 6378137.0              # Semi-major axis (m)
    b: float = 6356752.314245         # Semi-minor axis (m)
    f: float = 1 / 298.257223563      # Flattening

def vincenty_inverse(
    lat1: float, lon1: float,
    lat2: float, lon2: float,
    max_iterations: int = 200,
    tolerance: float = 1e-12
) -> Optional[Tuple[float, float, float]]:
    """
    Vincenty Inverse Problem: Given two points, find distance and azimuths.

    Args:
        lat1, lon1: First point (degrees)
        lat2, lon2: Second point (degrees)
        max_iterations: Maximum iterations for convergence
        tolerance: Convergence tolerance (radians)

    Returns:
        Tuple of (distance_m, azimuth1_deg, azimuth2_deg)
        or None if failed to converge (antipodal points)
    """
    e = WGS84()

    # Convert to radians
    φ1 = math.radians(lat1)
    φ2 = math.radians(lat2)
    L = math.radians(lon2 - lon1)

    # Reduced latitudes (parametric latitude)
    U1 = math.atan((1 - e.f) * math.tan(φ1))
    U2 = math.atan((1 - e.f) * math.tan(φ2))

    sin_U1 = math.sin(U1)
    cos_U1 = math.cos(U1)
    sin_U2 = math.sin(U2)
    cos_U2 = math.cos(U2)

    # Initial approximation
    λ = L

    for _ in range(max_iterations):
        sin_λ = math.sin(λ)
        cos_λ = math.cos(λ)

        sin_σ = math.sqrt(
            (cos_U2 * sin_λ) ** 2 +
            (cos_U1 * sin_U2 - sin_U1 * cos_U2 * cos_λ) ** 2
        )

        if sin_σ == 0:
            return (0.0, 0.0, 0.0)  # Coincident points

        cos_σ = sin_U1 * sin_U2 + cos_U1 * cos_U2 * cos_λ
        σ = math.atan2(sin_σ, cos_σ)

        sin_α = cos_U1 * cos_U2 * sin_λ / sin_σ
        cos²_α = 1 - sin_α ** 2

        if cos²_α == 0:
            cos_2σm = 0  # Equatorial line
        else:
            cos_2σm = cos_σ - 2 * sin_U1 * sin_U2 / cos²_α

        C = e.f / 16 * cos²_α * (4 + e.f * (4 - 3 * cos²_α))

        λ_prev = λ
        λ = L + (1 - C) * e.f * sin_α * (
            σ + C * sin_σ * (cos_2σm + C * cos_σ * (-1 + 2 * cos_2σm ** 2))
        )

        if abs(λ - λ_prev) < tolerance:
            break
    else:
        return None  # Failed to converge

    # Calculate distance
    u² = cos²_α * (e.a ** 2 - e.b ** 2) / e.b ** 2
    A = 1 + u² / 16384 * (4096 + u² * (-768 + u² * (320 - 175 * u²)))
    B = u² / 1024 * (256 + u² * (-128 + u² * (74 - 47 * u²)))

    Δσ = B * sin_σ * (
        cos_2σm + B / 4 * (
            cos_σ * (-1 + 2 * cos_2σm ** 2) -
            B / 6 * cos_2σm * (-3 + 4 * sin_σ ** 2) * (-3 + 4 * cos_2σm ** 2)
        )
    )

    s = e.b * A * (σ - Δσ)  # Distance in meters

    # Forward azimuth
    α1 = math.atan2(
        cos_U2 * sin_λ,
        cos_U1 * sin_U2 - sin_U1 * cos_U2 * cos_λ
    )

    # Back azimuth
    α2 = math.atan2(
        cos_U1 * sin_λ,
        -sin_U1 * cos_U2 + cos_U1 * sin_U2 * cos_λ
    )

    return (s, math.degrees(α1) % 360, (math.degrees(α2) + 180) % 360)


def vincenty_direct(
    lat1: float, lon1: float,
    azimuth: float, distance: float,
    max_iterations: int = 200,
    tolerance: float = 1e-12
) -> Tuple[float, float, float]:
    """
    Vincenty Direct Problem: Given start point, azimuth, and distance,
    find the destination point.

    Args:
        lat1, lon1: Start point (degrees)
        azimuth: Forward azimuth (degrees)
        distance: Distance to travel (meters)
        max_iterations: Maximum iterations
        tolerance: Convergence tolerance

    Returns:
        Tuple of (lat2, lon2, back_azimuth)
    """
    e = WGS84()

    φ1 = math.radians(lat1)
    α1 = math.radians(azimuth)
    s = distance

    # Reduced latitude
    U1 = math.atan((1 - e.f) * math.tan(φ1))

    sin_U1 = math.sin(U1)
    cos_U1 = math.cos(U1)
    sin_α1 = math.sin(α1)
    cos_α1 = math.cos(α1)

    # Azimuth of geodesic at equator
    σ1 = math.atan2(
        math.tan(U1),
        cos_α1
    )

    sin_α = cos_U1 * sin_α1
    cos²_α = 1 - sin_α ** 2

    u² = cos²_α * (e.a ** 2 - e.b ** 2) / e.b ** 2
    A = 1 + u² / 16384 * (4096 + u² * (-768 + u² * (320 - 175 * u²)))
    B = u² / 1024 * (256 + u² * (-128 + u² * (74 - 47 * u²)))

    # Initial approximation
    σ = s / (e.b * A)

    for _ in range(max_iterations):
        cos_2σm = math.cos(2 * σ1 + σ)
        sin_σ = math.sin(σ)
        cos_σ = math.cos(σ)

        Δσ = B * sin_σ * (
            cos_2σm + B / 4 * (
                cos_σ * (-1 + 2 * cos_2σm ** 2) -
                B / 6 * cos_2σm * (-3 + 4 * sin_σ ** 2) * (-3 + 4 * cos_2σm ** 2)
            )
        )

        σ_prev = σ
        σ = s / (e.b * A) + Δσ

        if abs(σ - σ_prev) < tolerance:
            break

    sin_σ = math.sin(σ)
    cos_σ = math.cos(σ)
    cos_2σm = math.cos(2 * σ1 + σ)

    # Final latitude
    φ2 = math.atan2(
        sin_U1 * cos_σ + cos_U1 * sin_σ * cos_α1,
        (1 - e.f) * math.sqrt(
            sin_α ** 2 + (sin_U1 * sin_σ - cos_U1 * cos_σ * cos_α1) ** 2
        )
    )

    # Longitude difference
    λ = math.atan2(
        sin_σ * sin_α1,
        cos_U1 * cos_σ - sin_U1 * sin_σ * cos_α1
    )

    C = e.f / 16 * cos²_α * (4 + e.f * (4 - 3 * cos²_α))
    L = λ - (1 - C) * e.f * sin_α * (
        σ + C * sin_σ * (cos_2σm + C * cos_σ * (-1 + 2 * cos_2σm ** 2))
    )

    lon2 = math.degrees(math.radians(lon1) + L)

    # Back azimuth
    α2 = math.atan2(
        sin_α,
        -sin_U1 * sin_σ + cos_U1 * cos_σ * cos_α1
    )

    return (math.degrees(φ2), lon2, (math.degrees(α2) + 180) % 360)


# Example usage
if __name__ == "__main__":
    # Los Angeles to JFK
    la = (33.9425, -118.4081)
    jfk = (40.6413, -73.7781)

    result = vincenty_inverse(*la, *jfk)
    if result:
        dist, az1, az2 = result
        print(f"Distance: {dist/1000:.3f} km")
        print(f"Forward azimuth: {az1:.4f}°")
        print(f"Back azimuth: {az2:.4f}°")

    # Direct: 100 km bearing 045° from LA
    lat2, lon2, back_az = vincenty_direct(33.9425, -118.4081, 45, 100000)
    print(f"\nDestination: {lat2:.6f}, {lon2:.6f}")
    print(f"Back azimuth: {back_az:.4f}°")

#!/usr/bin/env python3
"""
MGRS Batch Conversion Utility
For processing large coordinate datasets.
"""

import mgrs
import csv
import sys
from concurrent.futures import ProcessPoolExecutor, as_completed
from typing import Iterator, Tuple

m = mgrs.MGRS()

def convert_to_mgrs(lat: float, lon: float, precision: int = 5) -> str:
    """Convert lat/lon to MGRS with specified precision."""
    try:
        return m.toMGRS(lat, lon, MGRSPrecision=precision)
    except Exception as e:
        return f"ERROR: {e}"

def convert_from_mgrs(mgrs_str: str) -> Tuple[float, float]:
    """Convert MGRS to lat/lon."""
    try:
        lat, lon = m.toLatLon(mgrs_str)
        return (lat, lon)
    except Exception as e:
        return (None, None)

def batch_convert(
    input_file: str,
    output_file: str,
    direction: str = 'to_mgrs',
    precision: int = 5,
    workers: int = 4
) -> None:
    """
    Batch convert coordinates.

    Args:
        input_file: CSV with lat,lon or mgrs
        output_file: Output CSV
        direction: 'to_mgrs' or 'from_mgrs'
        precision: MGRS precision (1-5)
        workers: Parallel workers
    """
    with open(input_file, 'r') as infile, \
         open(output_file, 'w', newline='') as outfile:

        reader = csv.reader(infile)
        writer = csv.writer(outfile)

        if direction == 'to_mgrs':
            writer.writerow(['lat', 'lon', 'mgrs'])
            for row in reader:
                lat, lon = float(row[0]), float(row[1])
                result = convert_to_mgrs(lat, lon, precision)
                writer.writerow([lat, lon, result])
        else:
            writer.writerow(['mgrs', 'lat', 'lon'])
            for row in reader:
                mgrs_str = row[0].strip()
                lat, lon = convert_from_mgrs(mgrs_str)
                writer.writerow([mgrs_str, lat, lon])

if __name__ == "__main__":
    # Example: Convert list of coordinates
    coords = [
        (33.9425, -118.4081),   # Los Angeles
        (40.7128, -74.0060),    # New York
        (51.5074, -0.1278),     # London
        (33.3128, 44.3661),     # Baghdad
        (34.5553, 69.2075),     # Kabul
    ]

    print("Coordinate to MGRS Conversion:")
    print("-" * 60)
    for lat, lon in coords:
        mgrs_10 = convert_to_mgrs(lat, lon, 5)
        mgrs_6 = convert_to_mgrs(lat, lon, 3)
        print(f"{lat:9.4f}, {lon:10.4f} → {mgrs_10} ({mgrs_6})")

#!/bin/bash
# coord_validate.sh - Validate coordinates in multiple formats

validate_decimal() {
    local lat="$1" lon="$2"

    # Check latitude range [-90, 90]
    if ! awk -v lat="$lat" 'BEGIN { exit !(lat >= -90 && lat <= 90) }'; then
        echo "ERROR: Latitude $lat out of range [-90, 90]"
        return 1
    fi

    # Check longitude range [-180, 180]
    if ! awk -v lon="$lon" 'BEGIN { exit !(lon >= -180 && lon <= 180) }'; then
        echo "ERROR: Longitude $lon out of range [-180, 180]"
        return 1
    fi

    echo "VALID: $lat, $lon"
    return 0
}

validate_utm() {
    local zone="$1" easting="$2" northing="$3"

    # Zone must be 1-60
    zone_num="${zone%[A-Z]}"
    zone_letter="${zone#$zone_num}"

    if [[ ! "$zone_num" =~ ^[1-9]$|^[1-5][0-9]$|^60$ ]]; then
        echo "ERROR: Zone number $zone_num out of range [1, 60]"
        return 1
    fi

    # Easting typically 100000-900000
    if ! awk -v e="$easting" 'BEGIN { exit !(e >= 100000 && e <= 900000) }'; then
        echo "WARNING: Easting $easting outside typical range [100000, 900000]"
    fi

    # Northing 0-10000000 (varies by hemisphere)
    if ! awk -v n="$northing" 'BEGIN { exit !(n >= 0 && n <= 10000000) }'; then
        echo "ERROR: Northing $northing out of range [0, 10000000]"
        return 1
    fi

    echo "VALID: $zone $easting mE $northing mN"
    return 0
}

validate_mgrs() {
    local mgrs="$1"

    # Basic format check
    # Format: ZONE_NUM ZONE_LETTER 100KM_SQUARE EASTING NORTHING
    # Examples: 11SLA6891657142 (10-digit), 11SLA68915714 (8-digit)

    if [[ ! "$mgrs" =~ ^[0-9]{1,2}[A-Z][A-Z]{2}[0-9]+$ ]]; then
        echo "ERROR: Invalid MGRS format: $mgrs"
        return 1
    fi

    # Extract components
    local zone_num zone_letter grid_square coords
    zone_num=$(echo "$mgrs" | grep -oP '^\d{1,2}')
    zone_letter=$(echo "$mgrs" | grep -oP '^\d{1,2}\K[A-Z]')
    grid_square=$(echo "$mgrs" | grep -oP '^\d{1,2}[A-Z]\K[A-Z]{2}')
    coords=$(echo "$mgrs" | grep -oP '^\d{1,2}[A-Z][A-Z]{2}\K\d+')

    # Coordinate digits must be even (equal easting/northing)
    if (( ${#coords} % 2 != 0 )); then
        echo "ERROR: MGRS must have even number of coordinate digits"
        return 1
    fi

    # Valid zone letters (C-X, excluding I and O)
    if [[ ! "$zone_letter" =~ ^[C-HJ-NP-X]$ ]]; then
        echo "ERROR: Invalid zone letter: $zone_letter"
        return 1
    fi

    echo "VALID: $zone_num$zone_letter $grid_square ${coords:0:${#coords}/2} ${coords:${#coords}/2}"
    return 0
}

validate_dms() {
    local dms="$1"

    # Expected format: DD°MM'SS.S"N/S or DDD°MM'SS.S"E/W
    if [[ "$dms" =~ ^([0-9]+)°([0-9]+)\'([0-9.]+)\"([NSEW])$ ]]; then
        local deg="${BASH_REMATCH[1]}"
        local min="${BASH_REMATCH[2]}"
        local sec="${BASH_REMATCH[3]}"
        local dir="${BASH_REMATCH[4]}"

        # Validate degrees
        local max_deg=90
        [[ "$dir" =~ [EW] ]] && max_deg=180

        if (( deg > max_deg )); then
            echo "ERROR: Degrees $deg exceed maximum $max_deg"
            return 1
        fi

        # Validate minutes and seconds
        if (( min >= 60 )); then
            echo "ERROR: Minutes $min must be < 60"
            return 1
        fi

        if ! awk -v sec="$sec" 'BEGIN { exit !(sec >= 0 && sec < 60) }'; then
            echo "ERROR: Seconds $sec must be < 60"
            return 1
        fi

        echo "VALID: $deg°$min'$sec\"$dir"
        return 0
    else
        echo "ERROR: Invalid DMS format: $dms"
        return 1
    fi
}

# Round-trip validation using PROJ
round_trip_validate() {
    local lat="$1" lon="$2"
    local tolerance="${3:-0.0001}"

    # Convert to UTM and back
    local zone=$((((${lon%.*} + 180) / 6) + 1))

    # Forward transformation
    local utm=$(echo "$lon $lat" | cs2cs +proj=longlat +datum=WGS84 \
        +to +proj=utm +zone=$zone +datum=WGS84 -f "%.3f")
    local easting=$(echo "$utm" | awk '{print $1}')
    local northing=$(echo "$utm" | awk '{print $2}')

    # Inverse transformation
    local result=$(echo "$easting $northing" | cs2cs \
        +proj=utm +zone=$zone +datum=WGS84 \
        +to +proj=longlat +datum=WGS84 -f "%.8f")
    local lon2=$(echo "$result" | awk '{print $1}')
    local lat2=$(echo "$result" | awk '{print $2}')

    # Compare
    local diff_lat=$(awk -v a="$lat" -v b="$lat2" 'BEGIN { printf "%.8f", (a-b) < 0 ? -(a-b) : (a-b) }')
    local diff_lon=$(awk -v a="$lon" -v b="$lon2" 'BEGIN { printf "%.8f", (a-b) < 0 ? -(a-b) : (a-b) }')

    if awk -v d="$diff_lat" -v t="$tolerance" 'BEGIN { exit !(d < t) }' && \
       awk -v d="$diff_lon" -v t="$tolerance" 'BEGIN { exit !(d < t) }'; then
        echo "ROUND-TRIP OK: Δlat=$diff_lat° Δlon=$diff_lon°"
        return 0
    else
        echo "ROUND-TRIP FAIL: Δlat=$diff_lat° Δlon=$diff_lon° (tolerance: $tolerance°)"
        return 1
    fi
}

# Main
case "$1" in
    decimal)  validate_decimal "$2" "$3" ;;
    utm)      validate_utm "$2" "$3" "$4" ;;
    mgrs)     validate_mgrs "$2" ;;
    dms)      validate_dms "$2" ;;
    roundtrip) round_trip_validate "$2" "$3" "$4" ;;
    *)
        echo "Usage: $0 {decimal|utm|mgrs|dms|roundtrip} <args>"
        echo "  decimal <lat> <lon>"
        echo "  utm <zone> <easting> <northing>"
        echo "  mgrs <mgrs_string>"
        echo "  dms <dms_string>"
        echo "  roundtrip <lat> <lon> [tolerance]"
        ;;
esac

#!/usr/bin/env python3
"""
Calculate Circular Error Probable (CEP) from a set of positions.
Used for GPS/GNSS accuracy assessment.
"""

import math
from typing import List, Tuple
import statistics

def calculate_cep(
    points: List[Tuple[float, float]],
    true_pos: Tuple[float, float] = None
) -> dict:
    """
    Calculate CEP statistics from measured positions.

    Args:
        points: List of (x, y) positions in meters
        true_pos: True position if known, otherwise uses mean

    Returns:
        Dictionary with CEP statistics
    """
    if true_pos is None:
        # Use mean as reference
        ref_x = statistics.mean(p[0] for p in points)
        ref_y = statistics.mean(p[1] for p in points)
    else:
        ref_x, ref_y = true_pos

    # Calculate radial errors
    radial_errors = []
    for x, y in points:
        error = math.sqrt((x - ref_x)**2 + (y - ref_y)**2)
        radial_errors.append(error)

    radial_errors.sort()
    n = len(radial_errors)

    # CEP = 50th percentile of radial errors
    cep_idx = int(n * 0.5)
    cep = radial_errors[cep_idx]

    # 95th percentile (2DRMS approximation)
    p95_idx = int(n * 0.95)
    p95 = radial_errors[p95_idx]

    # Calculate standard deviations
    errors_x = [p[0] - ref_x for p in points]
    errors_y = [p[1] - ref_y for p in points]

    sigma_x = statistics.stdev(errors_x) if len(errors_x) > 1 else 0
    sigma_y = statistics.stdev(errors_y) if len(errors_y) > 1 else 0

    # Theoretical CEP (assuming circular distribution)
    sigma_r = math.sqrt(sigma_x**2 + sigma_y**2)
    cep_theoretical = 0.674 * (sigma_x + sigma_y) / 2

    # 2DRMS
    drms_2 = 2 * math.sqrt(sigma_x**2 + sigma_y**2)

    return {
        'n_samples': n,
        'reference': (ref_x, ref_y),
        'sigma_x': sigma_x,
        'sigma_y': sigma_y,
        'sigma_r': sigma_r,
        'cep_measured': cep,
        'cep_theoretical': cep_theoretical,
        '2drms': drms_2,
        'p95': p95,
        'mean_error': statistics.mean(radial_errors),
        'max_error': max(radial_errors),
    }

# Example usage
if __name__ == "__main__":
    import random
    random.seed(42)

    # Simulate GPS measurements with 5m sigma
    true_position = (0, 0)
    measurements = [
        (random.gauss(0, 5), random.gauss(0, 5))
        for _ in range(1000)
    ]

    results = calculate_cep(measurements, true_position)

    print("CEP Analysis Results")
    print("=" * 40)
    print(f"Samples:       {results['n_samples']}")
    print(f"σ_x:           {results['sigma_x']:.2f} m")
    print(f"σ_y:           {results['sigma_y']:.2f} m")
    print(f"CEP (measured):{results['cep_measured']:.2f} m")
    print(f"CEP (theory):  {results['cep_theoretical']:.2f} m")
    print(f"2DRMS:         {results['2drms']:.2f} m")
    print(f"95%:           {results['p95']:.2f} m")
    print(f"Mean error:    {results['mean_error']:.2f} m")
    print(f"Max error:     {results['max_error']:.2f} m")

# Run examples
python examples/navigation/coordinates.py

source examples/navigation/nav-calc.sh

# Distance
haversine 33.9425 -118.4081 40.6413 -73.7781 km

# Bearing
bearing 33.9425 -118.4081 40.6413 -73.7781

# Conversions
dms_to_decimal 33 45 30 N
decimal_to_dms 33.758333 lat

# Declination
apply_declination 45 12 E

# Demo
nav_demo

python examples/navigation/mgrs-convert.py

# Arch Linux
sudo pacman -S qgis

# Ubuntu
sudo apt install qgis qgis-plugin-grass

# Arch Linux (AUR)
yay -S google-earth-pro

# Or download from Google

# Using network.awk concepts for coordinates
# Parse coordinates from various formats

echo "33°45'30\"N 118°24'29\"W" | awk '
{
    # Parse DMS format
    gsub(/[°'\''"NSEW]/, " ")
    split($0, parts)
    lat = parts[1] + parts[2]/60 + parts[3]/3600
    lon = parts[4] + parts[5]/60 + parts[6]/3600
    printf "%.6f, %.6f\n", lat, -lon
}'

# ~/.local/share/bash/nav.sh
# Source this for navigation functions

source ~/.local/share/bash/math.sh  # For floating point

haversine() {
    # See examples/navigation/nav-calc.sh for implementation
    ...
}

# ~/.zshrc additions
alias nav='source ~/path/to/examples/navigation/nav-calc.sh && nav_demo'
alias declination='python -c "import geomag; print(geomag.declination($1, $2))"'