Innovare

We used GNSS and Python-based coordinate transformations to decode a campus clue chain, validating accuracy and revealing the Big Dipper constellation.

Comment

Geodetic Coordinate Transformation System

Project Overview

This project implements a comprehensive geodetic coordinate transformation system capable of converting between multiple coordinate reference frames and performing complex spatial calculations. The system is built on WGS84/NAD83 ellipsoid parameters and supports professional surveying and geospatial analysis workflows.

What We Built

Core Conversion Engine (GeoConverter Class)

A unified coordinate transformation system supporting:

  1. Geodetic ↔ ECEF Transformations

    • Geodetic (Latitude, Longitude, Height) to Earth-Centered Earth-Fixed (X, Y, Z)
    • ECEF to Geodetic with iterative refinement algorithm
    • Based on WGS84 ellipsoid (a = 6378137.0 m, f = 1/298.257223563)

  2. ECEF ↔ ENU Transformations

    • Earth-Centered Earth-Fixed to East-North-Up (local tangent plane)
    • ENU to ECEF using rotation matrices
    • Reference point-based local coordinate systems

  3. Geodesic Calculations

    • Distance and bearing between two geographic points
    • Forward projection (point + bearing + distance → new point)
    • Powered by PyProj's Geod implementation

  4. Advanced Point Calculations

    • Point B: Bearing-distance projection with elevation changes
    • Point C: ECEF coordinate offsets
    • Point D: ENU local coordinate offsets
    • Point E: Bearing intersection solution (2D line intersection in ENU frame)

Technical Implementation

Coordinate Systems Explained

1. Geodetic Coordinates

  • Format: Latitude (φ), Longitude (λ), Height (h)
  • Use Case: Standard GPS coordinates, mapping applications
  • Reference: WGS84/NAD83 ellipsoid surface

2. ECEF (Earth-Centered Earth-Fixed)

  • Format: X, Y, Z in meters
  • Origin: Earth's center of mass
  • Use Case: Satellite positioning, global transformations
  • Axes:
    • X: Through prime meridian
    • Y: 90° east of prime meridian
    • Z: Through North Pole

3. ENU (East-North-Up)

  • Format: East, North, Up in meters
  • Origin: User-defined reference point
  • Use Case: Local surveying, navigation, robotics
  • Axes:
    • East: Tangent to ellipsoid pointing east
    • North: Tangent to ellipsoid pointing north
    • Up: Normal to ellipsoid pointing up

Point Calculation Algorithms

Point B: Geodesic Projection

Input: Starting point A (lat, lon, h), bearing, distance, height change Process:

  1. Use ellipsoidal geodesic equations to project along bearing
  2. Calculate new latitude/longitude at specified distance
  3. Apply height change to starting elevation Output: New point B in geodetic and DMS formats

Point C: ECEF Offset

Input: Point B coordinates, ECEF offsets (ΔX, ΔY, ΔZ) Process:

  1. Convert Point B to ECEF coordinates
  2. Apply Cartesian offsets in global frame
  3. Convert result back to geodetic Output: Point C in ECEF, geodetic, and DMS formats

Point D: Local ENU Offset

Input: Point C coordinates, local offsets (ΔEast, ΔNorth, ΔUp) Process:

  1. Convert Point C to ECEF, then to ENU using C as reference
  2. Apply local offsets in tangent plane
  3. Convert back through ECEF to geodetic Output: Point D in ENU, ECEF, geodetic, and DMS formats

Point E: Bearing Intersection (Most Complex)

Input:

  • Point D (ECEF coordinates)
  • Reference point for ENU frame
  • Target Point P (ENU coordinates)
  • Bearing from D to E
  • Bearing from E to P
  • Height change

Mathematical Process:

  1. Convert D to Local Frame
   D_enu = ECEF_to_ENU(D_ecef, ref_point)
  1. Create Direction Vectors
   v1 = [sin(bearing_D_to_E), cos(bearing_D_to_E)]  # D→E direction
   v2 = [sin(bearing_E_to_P), cos(bearing_E_to_P)]  # E→P direction
  1. Set Up Linear System Point E satisfies both:
   E = D + t₁ · v1  (lies on line from D in direction of bearing_D_to_E)
   E = P - t₂ · v2  (lies on line to P from direction of bearing_E_to_P)
  1. Solve 2×2 Matrix Equation
   [v1_east   v2_east ] [t₁]   [P_east  - D_east ]
   [v1_north  v2_north] [t₂] = [P_north - D_north]

Using Cramer's rule:

   det = v1_east·v2_north - v1_north·v2_east
   t₁ = [(P_east - D_east)·v2_north - (P_north - D_north)·v2_east] / det
  1. Calculate Intersection Point
   E_east = D_east + t₁ · v1_east
   E_north = D_north + t₁ · v1_north
   E_up = D_up + height_change
  1. Convert to All Coordinate Systems E_ecef = ENU_to_ECEF(E_enu, ref_point) E_geodetic = ECEF_to_Geodetic(E_ecef)

Why This Works:

  • Two bearing lines create two equations in 2D (east-north plane)
  • Intersection is the unique point satisfying both bearing constraints
  • Determinant check ensures lines aren't parallel
  • Height is handled independently in the vertical (up) component

Real-World Application:

  • Survey triangulation
  • Navigation waypoint calculation
  • Position fixing with bearing measurements

Final Results (NAD83)

Decimal Degrees

Point Latitude (°) Longitude (°) Elevation (m)
A 51.0790180556 -114.1325483333 1114.70
B 51.0779852778 -114.1317241667 1110.99
C 51.0769152778 -114.1323066667 1109.78
D 51.0757341667 -114.1320875000 1108.22
E 51.0745880556 -114.1361938889 1109.35

DMS Format

Point Latitude Longitude Elevation (m)
A 51°04'44.465"N 114°07'57.174"W 1114.70
B 51°04'40.747"N 114°07'54.207"W 1110.99
C 51°04'36.895"N 114°07'56.304"W 1109.78
D 51°04'32.643"N 114°07'55.515"W 1108.22
E 51°04'28.517"N 114°08'10.298"W 1109.35

Project Components

1. GeoConverter Class

File: convert_pts.py Features:

  • WGS84/NAD83 ellipsoid support
  • Multiple coordinate system conversions
  • Geodesic calculations via PyProj
  • DMS formatting utilities
  • Five specialized point calculation methods

2. Interactive CLI Interface

Options:

  • Basic conversions (Geodetic ↔ ECEF ↔ ENU)
  • Distance and bearing calculations
  • Point projection
  • Advanced multi-step point calculations (B, C, D, E)

3. Coordinate Conversion Script

File: coordinate_pts.py Purpose: Batch conversion of DMS to decimal degrees Output: Formatted tables with high precision (10 decimal places)

4. Origin Point Statistics Calculator

File: origin_statistics.py Purpose: Calculate averaged origin coordinates from multiple BESTPOSA measurements Features:

  • Parses BESTPOSA format logs (NovAtel GPS receivers)
  • Extracts latitude, longitude, and height from all data points
  • Computes statistical mean from all measurements
  • Calculates precision metrics (standard deviation, range)
  • Validates measurement consistency across all points Output: High-precision averaged coordinates with statistical analysis (10+ decimal places)

Point A Results from 14 PPP Measurements:

  • Average Latitude: 51.07900053607°
  • Average Longitude: -114.13253483192°
  • Average Height: 1113.9210 m
  • Horizontal Precision: ~2 mm (0.0019-0.0021 m standard deviation)
  • Vertical Precision: 2.4 mm (0.0024 m standard deviation)
  • Total Measurement Spread: ~6 mm horizontal, 8.8 mm vertical

This represents excellent PPP (Precise Point Positioning) quality with sub-centimeter accuracy.

5. Raw Data Parser

Purpose: Extract BESTPOSA logs from raw GPS console output Process:

  • Reads raw console log files
  • Splits data by delimiter (#)
  • Extracts first 208 characters (BESTPOSA message length)
  • Outputs clean CSV format for further processing Input: Raw NovAtel console text file Output: Parsed CSV with BESTPOSA messages ready for statistical analysis

Mathematical Foundations

Ellipsoid Parameters (WGS84)

  • Semi-major axis (a): 6,378,137.0 meters
  • Flattening (f): 1/298.257223563
  • Eccentricity² (e²): f(2 - f) ≈ 0.00669438

Key Formulas

Geodetic to ECEF

N = a / √(1 - e²·sin²φ)
X = (N + h)·cos(φ)·cos(λ)
Y = (N + h)·cos(φ)·sin(λ)
Z = ((1 - e²)·N + h)·sin(φ)

ECEF to ENU Rotation Matrix

     [-sin(λ)              cos(λ)               0        ]
T =  [-sin(φ)cos(λ)  -sin(φ)sin(λ)   cos(φ)    ]
     [ cos(φ)cos(λ)   cos(φ)sin(λ)   sin(φ)    ]

Usage Examples

Example 1: Calculate Point B from Point A

geo = GeoConverter()
result = geo.calculate_point_b(
    lat_a=51.079,
    lon_a=-114.132,
    h_a=1114.70,
    bearing_deg=180.0,    # South
    distance_m=100.0,     # 100 meters
    delta_h=-3.71         # Drop 3.71 meters
)
print(result['geodetic'])  # (lat, lon, height)
print(result['dms'])       # DMS formatted strings

Example 2: Find Intersection Point E

result = geo.calculate_point_e(
    point_d_x=-1634567.89,  # ECEF coordinates of D
    point_d_y=-3664321.12,
    point_d_z=4940123.45,
    ref_lat=51.0,           # Reference for ENU frame
    ref_lon=-114.0,
    ref_h=1100.0,
    point_p_east=150.0,     # Target P in ENU
    point_p_north=200.0,
    bearing_d_to_e=45.0,    # Northeast from D
    bearing_e_to_p=135.0,   # Southeast to P
    height_change=1.13      # Rise 1.13 meters
)

Dependencies

  • Python 3.x
  • pyproj: Geodesic calculations and coordinate reference systems
  • math: Standard mathematical operations

Install dependencies:

pip install pyproj

Accuracy and Precision

  • Decimal Precision: 10 decimal places (~1.1 cm at equator)
  • ECEF Conversion: Iterative algorithm with 1e-12 convergence threshold
  • Geodesic Calculations: PyProj uses Karney's algorithms (accurate to ~15 nanometers)
  • DMS Display: 4 decimal places for seconds (sub-meter precision)

Validation and Error Handling

  • Parallel bearing check: Determinant threshold of 1e-10
  • Convergence verification: Iterative ECEF→Geodetic conversion
  • Input validation: Type checking in interactive mode
  • Coordinate consistency: All transformations preserve point relationships

Applications

  1. Surveying: Precise point positioning and coordinate transformations
  2. Navigation: Waypoint calculation and route planning
  3. GIS Analysis: Multi-reference frame spatial analysis
  4. Robotics: Local navigation in ENU coordinates
  5. Satellite Systems: ECEF-based positioning
  6. Civil Engineering: Construction layout and measurements

Project Achievement Summary

What We Accomplished

  1. ✅ Built complete multi-coordinate system transformation engine
  2. ✅ Implemented 5 specialized point calculation algorithms
  3. ✅ Created interactive CLI for all operations
  4. ✅ Converted real-world survey data from DMS to decimal degrees
  5. ✅ Achieved professional-grade accuracy (sub-centimeter precision)
  6. ✅ Documented all mathematical foundations and algorithms
  7. ✅ Provided working examples for all use cases

Technical Highlights

  • Zero external API dependencies for core math (except PyProj for geodesics)
  • Fully invertible transformations (bidirectional conversions)
  • Numerically stable algorithms with convergence guarantees
  • Clean object-oriented architecture
  • Comprehensive error handling

Project Completed: 7 February 2026
Coordinate Reference System: NAD83 / WGS84 Precision Level: Survey-grade (centimeter accuracy)

Built With

Share this project:

Updates