Search in sources :

Example 11 with NormalizedSphericalHarmonics

use of org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics in project Orekit by CS-SI.

the class HolmesFeatherstoneAttractionModel method gradient.

/**
 * Compute the gradient of the non-central part of the gravity field.
 * @param date current date
 * @param position position at which gravity field is desired in body frame
 * @param mu central attraction coefficient to use
 * @return gradient of the non-central part of the gravity field
 * @exception OrekitException if position cannot be converted to central body frame
 */
public double[] gradient(final AbsoluteDate date, final Vector3D position, final double mu) throws OrekitException {
    final int degree = provider.getMaxDegree();
    final int order = provider.getMaxOrder();
    final NormalizedSphericalHarmonics harmonics = provider.onDate(date);
    // allocate the columns for recursion
    double[] pnm0Plus2 = new double[degree + 1];
    double[] pnm0Plus1 = new double[degree + 1];
    double[] pnm0 = new double[degree + 1];
    final double[] pnm1 = new double[degree + 1];
    // compute polar coordinates
    final double x = position.getX();
    final double y = position.getY();
    final double z = position.getZ();
    final double x2 = x * x;
    final double y2 = y * y;
    final double z2 = z * z;
    final double r2 = x2 + y2 + z2;
    final double r = FastMath.sqrt(r2);
    final double rho2 = x2 + y2;
    final double rho = FastMath.sqrt(rho2);
    // cos(theta), where theta is the polar angle
    final double t = z / r;
    // sin(theta), where theta is the polar angle
    final double u = rho / r;
    final double tOu = z / rho;
    // compute distance powers
    final double[] aOrN = createDistancePowersArray(provider.getAe() / r);
    // compute longitude cosines/sines
    final double[][] cosSinLambda = createCosSinArrays(position.getX() / rho, position.getY() / rho);
    // outer summation over order
    int index = 0;
    double value = 0;
    final double[] gradient = new double[3];
    for (int m = degree; m >= 0; --m) {
        // compute tesseral terms with derivatives
        index = computeTesseral(m, degree, index, t, u, tOu, pnm0Plus2, pnm0Plus1, null, pnm0, pnm1, null);
        if (m <= order) {
            // compute contribution of current order to field (equation 5 of the paper)
            // inner summation over degree, for fixed order
            double sumDegreeS = 0;
            double sumDegreeC = 0;
            double dSumDegreeSdR = 0;
            double dSumDegreeCdR = 0;
            double dSumDegreeSdTheta = 0;
            double dSumDegreeCdTheta = 0;
            for (int n = FastMath.max(2, m); n <= degree; ++n) {
                final double qSnm = aOrN[n] * harmonics.getNormalizedSnm(n, m);
                final double qCnm = aOrN[n] * harmonics.getNormalizedCnm(n, m);
                final double nOr = n / r;
                final double s0 = pnm0[n] * qSnm;
                final double c0 = pnm0[n] * qCnm;
                final double s1 = pnm1[n] * qSnm;
                final double c1 = pnm1[n] * qCnm;
                sumDegreeS += s0;
                sumDegreeC += c0;
                dSumDegreeSdR -= nOr * s0;
                dSumDegreeCdR -= nOr * c0;
                dSumDegreeSdTheta += s1;
                dSumDegreeCdTheta += c1;
            }
            // contribution to outer summation over order
            // beware that we need to order gradient using the mathematical conventions
            // compliant with the SphericalCoordinates class, so our lambda is its theta
            // (and hence at index 1) and our theta is its phi (and hence at index 2)
            final double sML = cosSinLambda[1][m];
            final double cML = cosSinLambda[0][m];
            value = value * u + sML * sumDegreeS + cML * sumDegreeC;
            gradient[0] = gradient[0] * u + sML * dSumDegreeSdR + cML * dSumDegreeCdR;
            gradient[1] = gradient[1] * u + m * (cML * sumDegreeS - sML * sumDegreeC);
            gradient[2] = gradient[2] * u + sML * dSumDegreeSdTheta + cML * dSumDegreeCdTheta;
        }
        // rotate the recursion arrays
        final double[] tmp = pnm0Plus2;
        pnm0Plus2 = pnm0Plus1;
        pnm0Plus1 = pnm0;
        pnm0 = tmp;
    }
    // scale back
    value = FastMath.scalb(value, SCALING);
    gradient[0] = FastMath.scalb(gradient[0], SCALING);
    gradient[1] = FastMath.scalb(gradient[1], SCALING);
    gradient[2] = FastMath.scalb(gradient[2], SCALING);
    // apply the global mu/r factor
    final double muOr = mu / r;
    value *= muOr;
    gradient[0] = muOr * gradient[0] - value / r;
    gradient[1] *= muOr;
    gradient[2] *= muOr;
    // convert gradient from spherical to Cartesian
    return new SphericalCoordinates(position).toCartesianGradient(gradient);
}
Also used : SphericalCoordinates(org.hipparchus.geometry.euclidean.threed.SphericalCoordinates) NormalizedSphericalHarmonics(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics)

Example 12 with NormalizedSphericalHarmonics

use of org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics in project Orekit by CS-SI.

the class SHMFormatReaderTest method testRegular03cNormalized.

@Test
public void testRegular03cNormalized() throws OrekitException {
    Utils.setDataRoot("potential");
    GravityFieldFactory.addPotentialCoefficientsReader(new SHMFormatReader("eigen_cg03c_coef", false));
    NormalizedSphericalHarmonicsProvider provider = GravityFieldFactory.getNormalizedProvider(5, 5);
    Assert.assertEquals(TideSystem.TIDE_FREE, provider.getTideSystem());
    AbsoluteDate refDate = new AbsoluteDate("1997-01-01T12:00:00", TimeScalesFactory.getTT());
    Assert.assertEquals(refDate, provider.getReferenceDate());
    AbsoluteDate date = new AbsoluteDate("2011-05-01T01:02:03", TimeScalesFactory.getTT());
    Assert.assertEquals(date.durationFrom(refDate), provider.getOffset(date), Precision.SAFE_MIN);
    NormalizedSphericalHarmonics harmonics = provider.onDate(date);
    double offset = date.durationFrom(refDate);
    double offsetYear = offset / Constants.JULIAN_YEAR;
    Assert.assertEquals(0.957201462136e-06 + offsetYear * 0.490000000000e-11, harmonics.getNormalizedCnm(3, 0), 1.0e-15);
    Assert.assertEquals(0.174786174485e-06, harmonics.getNormalizedCnm(5, 5), 1.0e-15);
    Assert.assertEquals(0.0, harmonics.getNormalizedSnm(4, 0), 1.0e-15);
    Assert.assertEquals(0.308834784975e-06, harmonics.getNormalizedSnm(4, 4), 1.0e-15);
    Assert.assertEquals(0.3986004415E+15, provider.getMu(), 0);
    Assert.assertEquals(0.6378136460E+07, provider.getAe(), 0);
}
Also used : NormalizedSphericalHarmonics(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics) AbsoluteDate(org.orekit.time.AbsoluteDate) Test(org.junit.Test)

Example 13 with NormalizedSphericalHarmonics

use of org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics in project Orekit by CS-SI.

the class OceanTidesFieldTest method testDeltaCnmSnm.

@Test
public void testDeltaCnmSnm() throws OrekitException {
    // this is an arbitrarily truncated model, limited to 4x4 and with only a few waves
    List<OceanTidesWave> waves = getWaves(4, 4, 55565, 56554, 85455, 135655, 273555);
    UT1Scale ut1 = TimeScalesFactory.getUT1(IERSConventions.IERS_2010, true);
    TimeScale utc = TimeScalesFactory.getUTC();
    AbsoluteDate date = new AbsoluteDate(2003, 5, 6, 13, 43, 32.125, utc);
    OceanTidesField tidesField = new OceanTidesField(Constants.EIGEN5C_EARTH_EQUATORIAL_RADIUS, Constants.EIGEN5C_EARTH_MU, waves, IERSConventions.IERS_2010.getNutationArguments(ut1), null);
    NormalizedSphericalHarmonics harmonics = tidesField.onDate(date);
    double[][] refDeltaCnm = new double[][] { { 0.0, 0.0, 0.0, 0.0, 0.0 }, { 0.0, 0.0, 0.0, 0.0, 0.0 }, { -4.812565797928061E-11, -4.1748378190052583E-11, 7.013273986245356E-11, 0.0, 0.0 }, { -2.5341227608443308E-11, 9.76515813742254E-11, -1.21931214469994E-10, 1.3179722429471184E-10, 0.0 }, { -2.7496974839179478E-11, 8.419627031293907E-11, 6.56546217101275E-11, -3.375298928713117E-11, -7.588006744166988E-11 } };
    double[][] refDeltaSnm = new double[][] { { 0.0, 0.0, 0.0, 0.0, 0.0 }, { 0.0, 0.0, 0.0, 0.0, 0.0 }, { -1.168129177701461E-10, 5.646187590518608E-12, 1.742233297668071E-10, 0.0, 0.0 }, { -6.586546350227345E-11, -8.032186864783105E-11, -3.118910148495339E-11, 1.0566857199592183E-10, 0.0 }, { 7.665313525684617E-11, 7.37884528812169E-11, -1.3085142873419844E-10, -1.5813709543115768E-10, 1.770903634801541E-10 } };
    for (int n = 0; n < refDeltaCnm.length; ++n) {
        double threshold = 4.0e-17;
        for (int m = 0; m <= n; ++m) {
            Assert.assertEquals(refDeltaCnm[n][m], harmonics.getNormalizedCnm(n, m), threshold);
            Assert.assertEquals(refDeltaSnm[n][m], harmonics.getNormalizedSnm(n, m), threshold);
        }
    }
}
Also used : UT1Scale(org.orekit.time.UT1Scale) OceanTidesWave(org.orekit.forces.gravity.potential.OceanTidesWave) NormalizedSphericalHarmonics(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics) TimeScale(org.orekit.time.TimeScale) AbsoluteDate(org.orekit.time.AbsoluteDate) Test(org.junit.Test)

Example 14 with NormalizedSphericalHarmonics

use of org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics in project Orekit by CS-SI.

the class CachedNormalizedSphericalHarmonicsProviderTest method testReverseEntryGeneration.

@Test
public void testReverseEntryGeneration() throws OrekitException {
    // setup
    // generate points on grid with date as the origin
    cache.onDate(date);
    // sample before the current cached values
    AbsoluteDate sampleDate = date.shiftedBy(-step * 3);
    NormalizedSphericalHarmonics expected = raw.onDate(sampleDate);
    // action
    NormalizedSphericalHarmonics actual = cache.onDate(sampleDate);
    // verify
    double tol = Precision.EPSILON;
    for (int n = 0; n < raw.getMaxDegree(); n++) {
        for (int m = 0; m < n; m++) {
            Assert.assertEquals(expected.getNormalizedCnm(n, m), actual.getNormalizedCnm(n, m), tol);
            Assert.assertEquals(expected.getNormalizedSnm(n, m), actual.getNormalizedSnm(n, m), tol);
        }
    }
}
Also used : NormalizedSphericalHarmonics(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics) AbsoluteDate(org.orekit.time.AbsoluteDate) Test(org.junit.Test)

Example 15 with NormalizedSphericalHarmonics

use of org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics in project Orekit by CS-SI.

the class GravityFieldFactoryTest method testNormalizer.

@Test
public void testNormalizer() throws OrekitException {
    Utils.setDataRoot("potential/icgem-format");
    final double shift = 1.23456e8;
    NormalizedSphericalHarmonicsProvider ref = GravityFieldFactory.getNormalizedProvider(5, 5);
    NormalizedSphericalHarmonics refHarmonics = ref.onDate(ref.getReferenceDate().shiftedBy(shift));
    UnnormalizedSphericalHarmonicsProvider unnormalized = GravityFieldFactory.getUnnormalizedProvider(5, 5);
    NormalizedSphericalHarmonicsProvider normalized = GravityFieldFactory.getNormalizedProvider(unnormalized);
    NormalizedSphericalHarmonics normalizedHarmonics = normalized.onDate(normalized.getReferenceDate().shiftedBy(shift));
    Assert.assertEquals(ref.getMaxDegree(), normalized.getMaxDegree());
    Assert.assertEquals(ref.getMaxOrder(), normalized.getMaxOrder());
    Assert.assertEquals(ref.getReferenceDate(), normalized.getReferenceDate());
    Assert.assertEquals(ref.getAe(), normalized.getAe(), FastMath.ulp(ref.getAe()));
    Assert.assertEquals(ref.getMu(), normalized.getMu(), FastMath.ulp(ref.getMu()));
    Assert.assertEquals(ref.getOffset(AbsoluteDate.GPS_EPOCH), normalized.getOffset(AbsoluteDate.GPS_EPOCH), FastMath.ulp(ref.getOffset(AbsoluteDate.GPS_EPOCH)));
    for (int i = 0; i <= 5; ++i) {
        for (int j = 0; j <= i; ++j) {
            double cRef = refHarmonics.getNormalizedCnm(i, j);
            double cTest = normalizedHarmonics.getNormalizedCnm(i, j);
            Assert.assertEquals(cRef, cTest, FastMath.ulp(cRef));
            double sRef = refHarmonics.getNormalizedSnm(i, j);
            double sTest = normalizedHarmonics.getNormalizedSnm(i, j);
            Assert.assertEquals(sRef, sTest, FastMath.ulp(sRef));
        }
    }
}
Also used : NormalizedSphericalHarmonics(org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics) Test(org.junit.Test)

Aggregations

NormalizedSphericalHarmonics (org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider.NormalizedSphericalHarmonics)15 Test (org.junit.Test)10 AbsoluteDate (org.orekit.time.AbsoluteDate)8 SphericalCoordinates (org.hipparchus.geometry.euclidean.threed.SphericalCoordinates)3 NormalizedSphericalHarmonicsProvider (org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider)3 TimeScale (org.orekit.time.TimeScale)3 UT1Scale (org.orekit.time.UT1Scale)3 CachedNormalizedSphericalHarmonicsProvider (org.orekit.forces.gravity.potential.CachedNormalizedSphericalHarmonicsProvider)2 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)2 InvocationTargetException (java.lang.reflect.InvocationTargetException)1 StreamingStatistics (org.hipparchus.stat.descriptive.StreamingStatistics)1 OceanTidesWave (org.orekit.forces.gravity.potential.OceanTidesWave)1 Frame (org.orekit.frames.Frame)1 IERSConventions (org.orekit.utils.IERSConventions)1