| 1 | #include "erfa.h"
|
|---|
| 2 |
|
|---|
| 3 | int eraGc2gde ( double a, double f, double xyz[3],
|
|---|
| 4 | double *elong, double *phi, double *height )
|
|---|
| 5 | /*
|
|---|
| 6 | ** - - - - - - - - - -
|
|---|
| 7 | ** e r a G c 2 g d e
|
|---|
| 8 | ** - - - - - - - - - -
|
|---|
| 9 | **
|
|---|
| 10 | ** Transform geocentric coordinates to geodetic for a reference
|
|---|
| 11 | ** ellipsoid of specified form.
|
|---|
| 12 | **
|
|---|
| 13 | ** Given:
|
|---|
| 14 | ** a double equatorial radius (Notes 2,4)
|
|---|
| 15 | ** f double flattening (Note 3)
|
|---|
| 16 | ** xyz double[3] geocentric vector (Note 4)
|
|---|
| 17 | **
|
|---|
| 18 | ** Returned:
|
|---|
| 19 | ** elong double longitude (radians, east +ve)
|
|---|
| 20 | ** phi double latitude (geodetic, radians)
|
|---|
| 21 | ** height double height above ellipsoid (geodetic, Note 4)
|
|---|
| 22 | **
|
|---|
| 23 | ** Returned (function value):
|
|---|
| 24 | ** int status: 0 = OK
|
|---|
| 25 | ** -1 = illegal f
|
|---|
| 26 | ** -2 = illegal a
|
|---|
| 27 | **
|
|---|
| 28 | ** Notes:
|
|---|
| 29 | **
|
|---|
| 30 | ** 1) This function is based on the GCONV2H Fortran subroutine by
|
|---|
| 31 | ** Toshio Fukushima (see reference).
|
|---|
| 32 | **
|
|---|
| 33 | ** 2) The equatorial radius, a, can be in any units, but meters is
|
|---|
| 34 | ** the conventional choice.
|
|---|
| 35 | **
|
|---|
| 36 | ** 3) The flattening, f, is (for the Earth) a value around 0.00335,
|
|---|
| 37 | ** i.e. around 1/298.
|
|---|
| 38 | **
|
|---|
| 39 | ** 4) The equatorial radius, a, and the geocentric vector, xyz,
|
|---|
| 40 | ** must be given in the same units, and determine the units of
|
|---|
| 41 | ** the returned height, height.
|
|---|
| 42 | **
|
|---|
| 43 | ** 5) If an error occurs (status < 0), elong, phi and height are
|
|---|
| 44 | ** unchanged.
|
|---|
| 45 | **
|
|---|
| 46 | ** 6) The inverse transformation is performed in the function
|
|---|
| 47 | ** eraGd2gce.
|
|---|
| 48 | **
|
|---|
| 49 | ** 7) The transformation for a standard ellipsoid (such as ERFA_WGS84) can
|
|---|
| 50 | ** more conveniently be performed by calling eraGc2gd, which uses a
|
|---|
| 51 | ** numerical code to identify the required A and F values.
|
|---|
| 52 | **
|
|---|
| 53 | ** Reference:
|
|---|
| 54 | **
|
|---|
| 55 | ** Fukushima, T., "Transformation from Cartesian to geodetic
|
|---|
| 56 | ** coordinates accelerated by Halley's method", J.Geodesy (2006)
|
|---|
| 57 | ** 79: 689-693
|
|---|
| 58 | **
|
|---|
| 59 | ** Copyright (C) 2013-2017, NumFOCUS Foundation.
|
|---|
| 60 | ** Derived, with permission, from the SOFA library. See notes at end of file.
|
|---|
| 61 | */
|
|---|
| 62 | {
|
|---|
| 63 | double aeps2, e2, e4t, ec2, ec, b, x, y, z, p2, absz, p, s0, pn, zc,
|
|---|
| 64 | c0, c02, c03, s02, s03, a02, a0, a03, d0, f0, b0, s1,
|
|---|
| 65 | cc, s12, cc2;
|
|---|
| 66 |
|
|---|
| 67 |
|
|---|
| 68 | /* ------------- */
|
|---|
| 69 | /* Preliminaries */
|
|---|
| 70 | /* ------------- */
|
|---|
| 71 |
|
|---|
| 72 | /* Validate ellipsoid parameters. */
|
|---|
| 73 | if ( f < 0.0 || f >= 1.0 ) return -1;
|
|---|
| 74 | if ( a <= 0.0 ) return -2;
|
|---|
| 75 |
|
|---|
| 76 | /* Functions of ellipsoid parameters (with further validation of f). */
|
|---|
| 77 | aeps2 = a*a * 1e-32;
|
|---|
| 78 | e2 = (2.0 - f) * f;
|
|---|
| 79 | e4t = e2*e2 * 1.5;
|
|---|
| 80 | ec2 = 1.0 - e2;
|
|---|
| 81 | if ( ec2 <= 0.0 ) return -1;
|
|---|
| 82 | ec = sqrt(ec2);
|
|---|
| 83 | b = a * ec;
|
|---|
| 84 |
|
|---|
| 85 | /* Cartesian components. */
|
|---|
| 86 | x = xyz[0];
|
|---|
| 87 | y = xyz[1];
|
|---|
| 88 | z = xyz[2];
|
|---|
| 89 |
|
|---|
| 90 | /* Distance from polar axis squared. */
|
|---|
| 91 | p2 = x*x + y*y;
|
|---|
| 92 |
|
|---|
| 93 | /* Longitude. */
|
|---|
| 94 | *elong = p2 > 0.0 ? atan2(y, x) : 0.0;
|
|---|
| 95 |
|
|---|
| 96 | /* Unsigned z-coordinate. */
|
|---|
| 97 | absz = fabs(z);
|
|---|
| 98 |
|
|---|
| 99 | /* Proceed unless polar case. */
|
|---|
| 100 | if ( p2 > aeps2 ) {
|
|---|
| 101 |
|
|---|
| 102 | /* Distance from polar axis. */
|
|---|
| 103 | p = sqrt(p2);
|
|---|
| 104 |
|
|---|
| 105 | /* Normalization. */
|
|---|
| 106 | s0 = absz / a;
|
|---|
| 107 | pn = p / a;
|
|---|
| 108 | zc = ec * s0;
|
|---|
| 109 |
|
|---|
| 110 | /* Prepare Newton correction factors. */
|
|---|
| 111 | c0 = ec * pn;
|
|---|
| 112 | c02 = c0 * c0;
|
|---|
| 113 | c03 = c02 * c0;
|
|---|
| 114 | s02 = s0 * s0;
|
|---|
| 115 | s03 = s02 * s0;
|
|---|
| 116 | a02 = c02 + s02;
|
|---|
| 117 | a0 = sqrt(a02);
|
|---|
| 118 | a03 = a02 * a0;
|
|---|
| 119 | d0 = zc*a03 + e2*s03;
|
|---|
| 120 | f0 = pn*a03 - e2*c03;
|
|---|
| 121 |
|
|---|
| 122 | /* Prepare Halley correction factor. */
|
|---|
| 123 | b0 = e4t * s02 * c02 * pn * (a0 - ec);
|
|---|
| 124 | s1 = d0*f0 - b0*s0;
|
|---|
| 125 | cc = ec * (f0*f0 - b0*c0);
|
|---|
| 126 |
|
|---|
| 127 | /* Evaluate latitude and height. */
|
|---|
| 128 | *phi = atan(s1/cc);
|
|---|
| 129 | s12 = s1 * s1;
|
|---|
| 130 | cc2 = cc * cc;
|
|---|
| 131 | *height = (p*cc + absz*s1 - a * sqrt(ec2*s12 + cc2)) /
|
|---|
| 132 | sqrt(s12 + cc2);
|
|---|
| 133 | } else {
|
|---|
| 134 |
|
|---|
| 135 | /* Exception: pole. */
|
|---|
| 136 | *phi = ERFA_DPI / 2.0;
|
|---|
| 137 | *height = absz - b;
|
|---|
| 138 | }
|
|---|
| 139 |
|
|---|
| 140 | /* Restore sign of latitude. */
|
|---|
| 141 | if ( z < 0 ) *phi = -*phi;
|
|---|
| 142 |
|
|---|
| 143 | /* OK status. */
|
|---|
| 144 | return 0;
|
|---|
| 145 |
|
|---|
| 146 | }
|
|---|
| 147 | /*----------------------------------------------------------------------
|
|---|
| 148 | **
|
|---|
| 149 | **
|
|---|
| 150 | ** Copyright (C) 2013-2017, NumFOCUS Foundation.
|
|---|
| 151 | ** All rights reserved.
|
|---|
| 152 | **
|
|---|
| 153 | ** This library is derived, with permission, from the International
|
|---|
| 154 | ** Astronomical Union's "Standards of Fundamental Astronomy" library,
|
|---|
| 155 | ** available from http://www.iausofa.org.
|
|---|
| 156 | **
|
|---|
| 157 | ** The ERFA version is intended to retain identical functionality to
|
|---|
| 158 | ** the SOFA library, but made distinct through different function and
|
|---|
| 159 | ** file names, as set out in the SOFA license conditions. The SOFA
|
|---|
| 160 | ** original has a role as a reference standard for the IAU and IERS,
|
|---|
| 161 | ** and consequently redistribution is permitted only in its unaltered
|
|---|
| 162 | ** state. The ERFA version is not subject to this restriction and
|
|---|
| 163 | ** therefore can be included in distributions which do not support the
|
|---|
| 164 | ** concept of "read only" software.
|
|---|
| 165 | **
|
|---|
| 166 | ** Although the intent is to replicate the SOFA API (other than
|
|---|
| 167 | ** replacement of prefix names) and results (with the exception of
|
|---|
| 168 | ** bugs; any that are discovered will be fixed), SOFA is not
|
|---|
| 169 | ** responsible for any errors found in this version of the library.
|
|---|
| 170 | **
|
|---|
| 171 | ** If you wish to acknowledge the SOFA heritage, please acknowledge
|
|---|
| 172 | ** that you are using a library derived from SOFA, rather than SOFA
|
|---|
| 173 | ** itself.
|
|---|
| 174 | **
|
|---|
| 175 | **
|
|---|
| 176 | ** TERMS AND CONDITIONS
|
|---|
| 177 | **
|
|---|
| 178 | ** Redistribution and use in source and binary forms, with or without
|
|---|
| 179 | ** modification, are permitted provided that the following conditions
|
|---|
| 180 | ** are met:
|
|---|
| 181 | **
|
|---|
| 182 | ** 1 Redistributions of source code must retain the above copyright
|
|---|
| 183 | ** notice, this list of conditions and the following disclaimer.
|
|---|
| 184 | **
|
|---|
| 185 | ** 2 Redistributions in binary form must reproduce the above copyright
|
|---|
| 186 | ** notice, this list of conditions and the following disclaimer in
|
|---|
| 187 | ** the documentation and/or other materials provided with the
|
|---|
| 188 | ** distribution.
|
|---|
| 189 | **
|
|---|
| 190 | ** 3 Neither the name of the Standards Of Fundamental Astronomy Board,
|
|---|
| 191 | ** the International Astronomical Union nor the names of its
|
|---|
| 192 | ** contributors may be used to endorse or promote products derived
|
|---|
| 193 | ** from this software without specific prior written permission.
|
|---|
| 194 | **
|
|---|
| 195 | ** THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
|
|---|
| 196 | ** "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
|
|---|
| 197 | ** LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
|
|---|
| 198 | ** FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
|
|---|
| 199 | ** COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
|
|---|
| 200 | ** INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
|
|---|
| 201 | ** BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
|
|---|
| 202 | ** LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
|
|---|
| 203 | ** CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
|
|---|
| 204 | ** LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
|
|---|
| 205 | ** ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
|
|---|
| 206 | ** POSSIBILITY OF SUCH DAMAGE.
|
|---|
| 207 | **
|
|---|
| 208 | */
|
|---|