src/Datum.php
public function geocentric_to_geodetic_noniter(Point $p)
    {
        /*
        $Longitude;
        $Latitude;
        $Height;

        $W;        // distance from Z axis 
        $W2;       // square of distance from Z axis 
        $T0;       // initial estimate of vertical component 
        $T1;       // corrected estimate of vertical component 
        $S0;       // initial estimate of horizontal component 
        $S1;       // corrected estimate of horizontal component
        $Sin_B0;   // sin(B0), B0 is estimate of Bowring aux variable 
        $Sin3_B0;  // cube of sin(B0) 
        $Cos_B0;   // cos(B0)
        $Sin_p1;   // sin(phi1), phi1 is estimated latitude 
        $Cos_p1;   // cos(phi1) 
        $Rn;       // Earth radius at location 
        $Sum;      // numerator of cos(phi1) 
        $AtPole;  // indicates location is in polar region 
        */

        // Cast from string to float.
        // Since we are accepting the Point class only, then we can already
        // guarantee we have floats. A simple list($x, $y $Z) = $p->toArray() will
        // give us our values.

        $X = floatval($p->x);
        $Y = floatval($p->y);
        $Z = floatval($p->z ? $p->z : 0);

        $AtPole = false;

        if ($X <> 0.0) {
            $Longitude = atan2($Y, $X);
        } else {
            if ($Y > 0) {
                $Longitude = Common::HALF_PI;
            } elseif (Y < 0) {
                $Longitude = -Common::HALF_PI;
            } else {
                $AtPole = true;
                $Longitude = 0.0;
                if ($Z > 0.0) { /* north pole */
                    $Latitude = Common::HALF_PI;
                } elseif (Z < 0.0) { /* south pole */
                    $Latitude = -Common::HALF_PI;
                } else { /* center of earth */
                    $Latitude = Common::HALF_PI;
                    $Height = -$this->b;
                    return;
                }
            }
        }

        $W2 = $X * $X + $Y * $Y;
        $W = sqrt($W2);
        $T0 = $Z * Common::AD_C;
        $S0 = sqrt($T0 * $T0 + $W2);
        $Sin_B0 = $T0 / $S0;
        $Cos_B0 = $W / $S0;
        $Sin3_B0 = $Sin_B0 * $Sin_B0 * $Sin_B0;
        $T1 = $Z + $this->b * $this->ep2 * $Sin3_B0;
        $Sum = $W - $this->a * $this->es * $Cos_B0 * $Cos_B0 * $Cos_B0;
        $S1 = sqrt($T1 * $T1 + $Sum * $Sum);
        $Sin_p1 = $T1 / $S1;
        $Cos_p1 = $Sum / $S1;
        $Rn = $this->a / sqrt( 1.0 - $this->es * $Sin_p1 * $Sin_p1);

        if ($Cos_p1 >= Common::COS_67P5) {
            $Height = $W / $Cos_p1 - $Rn;
        } elseif ($Cos_p1 <= -Common::COS_67P5) {
            $Height = $W / -$Cos_p1 - $Rn;
        } else {
            $Height = $Z / $Sin_p1 + $Rn * ($this->es - 1.0);
        }

        if ($AtPole == false) {
            $Latitude = atan($Sin_p1 / $Cos_p1);
        }

        $p->x = $Longitude;
        $p->y = $Latitude;
        $p->z = $Height;

        return $p;
    }
public function geocentric_to_geodetic_noniter(Point $p)
    {
        /*
        $Longitude;
        $Latitude;
        $Height;

        $W;        // distance from Z axis 
        $W2;       // square of distance from Z axis 
        $T0;       // initial estimate of vertical component 
        $T1;       // corrected estimate of vertical component 
        $S0;       // initial estimate of horizontal component 
        $S1;       // corrected estimate of horizontal component
        $Sin_B0;   // sin(B0), B0 is estimate of Bowring aux variable 
        $Sin3_B0;  // cube of sin(B0) 
        $Cos_B0;   // cos(B0)
        $Sin_p1;   // sin(phi1), phi1 is estimated latitude 
        $Cos_p1;   // cos(phi1) 
        $Rn;       // Earth radius at location 
        $Sum;      // numerator of cos(phi1) 
        $AtPole;  // indicates location is in polar region 
        */

        // Cast from string to float.
        // Since we are accepting the Point class only, then we can already
        // guarantee we have floats. A simple list($x, $y $Z) = $p->toArray() will
        // give us our values.

        $X = floatval($p->x);
        $Y = floatval($p->y);
        $Z = floatval($p->z ? $p->z : 0);

        $AtPole = false;

        if ($X <> 0.0) {
            $Longitude = atan2($Y, $X);
        } else {
            if ($Y > 0) {
                $Longitude = Common::HALF_PI;
            } elseif (Y < 0) {
                $Longitude = -Common::HALF_PI;
            } else {
                $AtPole = true;
                $Longitude = 0.0;
                if ($Z > 0.0) { /* north pole */
                    $Latitude = Common::HALF_PI;
                } elseif (Z < 0.0) { /* south pole */
                    $Latitude = -Common::HALF_PI;
                } else { /* center of earth */
                    $Latitude = Common::HALF_PI;
                    $Height = -$this->b;
                    return;
                }
            }
        }

        $W2 = $X * $X + $Y * $Y;
        $W = sqrt($W2);
        $T0 = $Z * Common::AD_C;
        $S0 = sqrt($T0 * $T0 + $W2);
        $Sin_B0 = $T0 / $S0;
        $Cos_B0 = $W / $S0;
        $Sin3_B0 = $Sin_B0 * $Sin_B0 * $Sin_B0;
        $T1 = $Z + $this->b * $this->ep2 * $Sin3_B0;
        $Sum = $W - $this->a * $this->es * $Cos_B0 * $Cos_B0 * $Cos_B0;
        $S1 = sqrt($T1 * $T1 + $Sum * $Sum);
        $Sin_p1 = $T1 / $S1;
        $Cos_p1 = $Sum / $S1;
        $Rn = $this->a / sqrt( 1.0 - $this->es * $Sin_p1 * $Sin_p1);

        if ($Cos_p1 >= Common::COS_67P5) {
            $Height = $W / $Cos_p1 - $Rn;
        } elseif ($Cos_p1 <= -Common::COS_67P5) {
            $Height = $W / -$Cos_p1 - $Rn;
        } else {
            $Height = $Z / $Sin_p1 + $Rn * ($this->es - 1.0);
        }

        if ($AtPole == false) {
            $Latitude = atan($Sin_p1 / $Cos_p1);
        }

        $p->x = $Longitude;
        $p->y = $Latitude;
        $p->z = $Height;

        return $p;
    }
public function geocentric_to_geodetic_noniter ( Point $ p ) { /*
        $Longitude;
        $Latitude;
        $Height;

        $W;        // distance from Z axis 
        $W2;       // square of distance from Z axis 
        $T0;       // initial estimate of vertical component 
        $T1;       // corrected estimate of vertical component 
        $S0;       // initial estimate of horizontal component 
        $S1;       // corrected estimate of horizontal component
        $Sin_B0;   // sin(B0), B0 is estimate of Bowring aux variable 
        $Sin3_B0;  // cube of sin(B0) 
        $Cos_B0;   // cos(B0)
        $Sin_p1;   // sin(phi1), phi1 is estimated latitude 
        $Cos_p1;   // cos(phi1) 
        $Rn;       // Earth radius at location 
        $Sum;      // numerator of cos(phi1) 
        $AtPole;  // indicates location is in polar region 
        */ // Cast from string to float.
 // Since we are accepting the Point class only, then we can already
 // guarantee we have floats. A simple list($x, $y $Z) = $p->toArray() will
 // give us our values.
 $ X = floatval ( $ p -> x ) ; $ Y = floatval ( $ p -> y ) ; $ Z = floatval ( $ p -> z ? $ p -> z : 0 ) ; $ AtPole = false ; if ( $ X <> 0.0 ) { $ Longitude = atan2 ( $ Y , $ X ) ; } else { if ( $ Y > 0 ) { $ Longitude = Common :: HALF_PI ; } elseif ( Y < 0 ) { $ Longitude = - Common :: HALF_PI ; } else { $ AtPole = true ; $ Longitude = 0.0 ; if ( $ Z > 0.0 ) { /* north pole */ $ Latitude = Common :: HALF_PI ; } elseif ( Z < 0.0 ) { /* south pole */ $ Latitude = - Common :: HALF_PI ; } else { /* center of earth */ $ Latitude = Common :: HALF_PI ; $ Height = - $ this -> b ; return ; } } } $ W2 = $ X * $ X + $ Y * $ Y ; $ W = sqrt ( $ W2 ) ; $ T0 = $ Z * Common :: AD_C ; $ S0 = sqrt ( $ T0 * $ T0 + $ W2 ) ; $ Sin_B0 = $ T0 / $ S0 ; $ Cos_B0 = $ W / $ S0 ; $ Sin3_B0 = $ Sin_B0 * $ Sin_B0 * $ Sin_B0 ; $ T1 = $ Z + $ this -> b * $ this -> ep2 * $ Sin3_B0 ; $ Sum = $ W - $ this -> a * $ this -> es * $ Cos_B0 * $ Cos_B0 * $ Cos_B0 ; $ S1 = sqrt ( $ T1 * $ T1 + $ Sum * $ Sum ) ; $ Sin_p1 = $ T1 / $ S1 ; $ Cos_p1 = $ Sum / $ S1 ; $ Rn = $ this -> a / sqrt ( 1.0 - $ this -> es * $ Sin_p1 * $ Sin_p1 ) ; if ( $ Cos_p1 >= Common :: COS_67P5 ) { $ Height = $ W / $ Cos_p1 - $ Rn ; } elseif ( $ Cos_p1 <= - Common :: COS_67P5 ) { $ Height = $ W / - $ Cos_p1 - $ Rn ; } else { $ Height = $ Z / $ Sin_p1 + $ Rn * ( $ this -> es - 1.0 ) ; } if ( $ AtPole == false ) { $ Latitude = atan ( $ Sin_p1 / $ Cos_p1 ) ; } $ p -> x = $ Longitude ; $ p -> y = $ Latitude ; $ p -> z = $ Height ; return $ p ; }
Convert_Geocentric_To_Geodetic
The method used here is derived from 'An Improved Algorithm for
Geocentric to Geodetic Coordinate Conversion', by Ralph Toms, Feb 1996

@param object Point $p
@return object Point $p
Convert_Geocentric_To_Geodetic The method used here is derived from An Improved Algorithm for Geocentric to Geodetic Coordinate Conversion by Ralph Toms Feb 1996