QGIS » QGIS Application

#!/usr/bin/env python

# -*- coding: utf-8 -*-

from math import cos, radians, sin, sqrt

## Input

lat = 43.21009

lon = -78.120123

h = 124

## Ellipsoid Parameters: semi-major axis, reciprocal flattening

grs80 = 6378137, 298.257222100882711

wgs84 = 6378137, 298.257223563

def geodetic_to_geocentric(spheroid, latitude, longitude, height):

    """

    Return the Geocentric (Cartesian) Coordinates X, Y, Z

    given:

    - Datum spheroid given as a pair: semi-major axis, reciprocal flattening

    - Geodetic Coordinates latitude, longitude (in degrees)

    - Ellipsoid Height h (in metres)

    """

    a, rf = spheroid

    ϕ = radians(latitude)

    λ = radians(longitude)

    sin_ϕ = sin(ϕ)

    e2 = (2 - 1 / rf) / rf            # eccentricity squared

    N = a / sqrt(1 - e2 * sin_ϕ ** 2) # prime vertical radius of curvature

    r = (N + h) * cos(ϕ)              # distance from z axis

    X = r * cos(λ)

    Y = r * sin(λ)

    Z = ((1 - e2) * N + h) * sin_ϕ

    return X, Y, Z

print(geodetic_to_geocentric(wgs84, lat, lon, h))