Search in sources :

Example 46 with RandomGenerator

use of org.hipparchus.random.RandomGenerator in project Orekit by CS-SI.

the class FieldTransformTest method doTestLinear.

private <T extends RealFieldElement<T>> void doTestLinear(final Field<T> field) {
    RandomGenerator random = new Well19937a(0x14f6411217b148d8l);
    for (int n = 0; n < 100; ++n) {
        FieldTransform<T> t = randomTransform(field, random);
        // build an equivalent linear transform by extracting raw translation/rotation
        FieldMatrix<T> linearA = MatrixUtils.createFieldMatrix(field, 3, 4);
        linearA.setSubMatrix(t.getRotation().getMatrix(), 0, 0);
        FieldVector3D<T> rt = t.getRotation().applyTo(t.getTranslation());
        linearA.setEntry(0, 3, rt.getX());
        linearA.setEntry(1, 3, rt.getY());
        linearA.setEntry(2, 3, rt.getZ());
        // build an equivalent linear transform by observing transformed points
        FieldMatrix<T> linearB = MatrixUtils.createFieldMatrix(field, 3, 4);
        FieldVector3D<T> p0 = t.transformPosition(FieldVector3D.getZero(field));
        FieldVector3D<T> pI = t.transformPosition(FieldVector3D.getPlusI(field)).subtract(p0);
        FieldVector3D<T> pJ = t.transformPosition(FieldVector3D.getPlusJ(field)).subtract(p0);
        FieldVector3D<T> pK = t.transformPosition(FieldVector3D.getPlusK(field)).subtract(p0);
        linearB.setEntry(0, 0, pI.getX());
        linearB.setEntry(1, 0, pI.getY());
        linearB.setEntry(2, 0, pI.getZ());
        linearB.setEntry(0, 1, pJ.getX());
        linearB.setEntry(1, 1, pJ.getY());
        linearB.setEntry(2, 1, pJ.getZ());
        linearB.setEntry(0, 2, pK.getX());
        linearB.setEntry(1, 2, pK.getY());
        linearB.setEntry(2, 2, pK.getZ());
        linearB.setEntry(0, 3, p0.getX());
        linearB.setEntry(1, 3, p0.getY());
        linearB.setEntry(2, 3, p0.getZ());
        // both linear transforms should be equal
        FieldMatrix<T> sub = linearB.subtract(linearA);
        double refMax = 0;
        double diffMax = 0;
        for (int i = 0; i < linearA.getRowDimension(); ++i) {
            for (int j = 0; j < linearA.getColumnDimension(); ++j) {
                refMax = FastMath.max(linearA.getEntry(i, j).getReal(), refMax);
                diffMax = FastMath.max(sub.getEntry(i, j).getReal(), diffMax);
            }
        }
        Assert.assertEquals(0.0, diffMax, 2.0e-12 * refMax);
        for (int i = 0; i < 100; ++i) {
            FieldVector3D<T> p = randomVector(field, 1.0e3, random);
            FieldVector3D<T> q = t.transformPosition(p);
            T[] pField = MathArrays.buildArray(field, 4);
            pField[0] = p.getX();
            pField[1] = p.getY();
            pField[2] = p.getZ();
            pField[3] = field.getOne();
            T[] qA = linearA.operate(pField);
            Assert.assertEquals(q.getX().getReal(), qA[0].getReal(), 1.0e-13 * p.getNorm().getReal());
            Assert.assertEquals(q.getY().getReal(), qA[1].getReal(), 1.0e-13 * p.getNorm().getReal());
            Assert.assertEquals(q.getZ().getReal(), qA[2].getReal(), 1.0e-13 * p.getNorm().getReal());
            T[] qB = linearB.operate(pField);
            Assert.assertEquals(q.getX().getReal(), qB[0].getReal(), 1.0e-10 * p.getNorm().getReal());
            Assert.assertEquals(q.getY().getReal(), qB[1].getReal(), 1.0e-10 * p.getNorm().getReal());
            Assert.assertEquals(q.getZ().getReal(), qB[2].getReal(), 1.0e-10 * p.getNorm().getReal());
        }
    }
}
Also used : Well19937a(org.hipparchus.random.Well19937a) RandomGenerator(org.hipparchus.random.RandomGenerator)

Example 47 with RandomGenerator

use of org.hipparchus.random.RandomGenerator in project Orekit by CS-SI.

the class FieldTransformTest method doTestTransPV.

private <T extends RealFieldElement<T>> void doTestTransPV(Field<T> field) {
    RandomGenerator rnd = new Well19937a(0x73d5554d99427af0l);
    for (int i = 0; i < 10; ++i) {
        // random position, velocity and acceleration
        FieldVector3D<T> pos = randomVector(field, 1.0e3, rnd);
        FieldVector3D<T> vel = randomVector(field, 1.0, rnd);
        FieldVector3D<T> acc = randomVector(field, 1.0e-3, rnd);
        FieldPVCoordinates<T> pvOne = new FieldPVCoordinates<>(pos, vel, acc);
        // random transform
        FieldVector3D<T> transPos = randomVector(field, 1.0e3, rnd);
        FieldVector3D<T> transVel = randomVector(field, 1.0, rnd);
        FieldVector3D<T> transAcc = randomVector(field, 1.0e-3, rnd);
        FieldTransform<T> tr = new FieldTransform<>(FieldAbsoluteDate.getJ2000Epoch(field), transPos, transVel, transAcc);
        double dt = 1;
        // we should obtain
        FieldVector3D<T> good = tr.transformPosition(pos.add(new FieldVector3D<>(dt, vel))).add(new FieldVector3D<>(dt, transVel));
        // we have
        FieldPVCoordinates<T> pvTwo = tr.transformPVCoordinates(pvOne);
        FieldVector3D<T> result = pvTwo.getPosition().add(new FieldVector3D<>(dt, pvTwo.getVelocity()));
        checkVector(good, result, 1.0e-15);
        // test inverse
        FieldVector3D<T> resultvel = tr.getInverse().transformPVCoordinates(pvTwo).getVelocity();
        checkVector(resultvel, vel, 1.0e-15);
    }
}
Also used : FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates) Well19937a(org.hipparchus.random.Well19937a) RandomGenerator(org.hipparchus.random.RandomGenerator)

Example 48 with RandomGenerator

use of org.hipparchus.random.RandomGenerator in project Orekit by CS-SI.

the class FieldTransformTest method doTestJacobianP.

private <T extends RealFieldElement<T>> void doTestJacobianP(Field<T> field) {
    // base directions for finite differences
    @SuppressWarnings("unchecked") FieldPVCoordinates<T>[] directions = (FieldPVCoordinates<T>[]) Array.newInstance(FieldPVCoordinates.class, 3);
    directions[0] = new FieldPVCoordinates<>(FieldVector3D.getPlusI(field), FieldVector3D.getZero(field), FieldVector3D.getZero(field));
    directions[1] = new FieldPVCoordinates<>(FieldVector3D.getPlusJ(field), FieldVector3D.getZero(field), FieldVector3D.getZero(field));
    directions[2] = new FieldPVCoordinates<>(FieldVector3D.getPlusK(field), FieldVector3D.getZero(field), FieldVector3D.getZero(field));
    double h = 0.01;
    RandomGenerator random = new Well19937a(0x47fd0d6809f4b173l);
    for (int i = 0; i < 20; ++i) {
        // generate a random transform
        FieldTransform<T> combined = randomTransform(field, random);
        // compute Jacobian
        T[][] jacobian = MathArrays.buildArray(field, 9, 9);
        for (int l = 0; l < jacobian.length; ++l) {
            for (int c = 0; c < jacobian[l].length; ++c) {
                jacobian[l][c] = field.getZero().add(l + 0.1 * c);
            }
        }
        combined.getJacobian(CartesianDerivativesFilter.USE_P, jacobian);
        for (int j = 0; j < 100; ++j) {
            FieldPVCoordinates<T> pv0 = new FieldPVCoordinates<>(randomVector(field, 1e3, random), randomVector(field, 1.0, random), randomVector(field, 1.0e-3, random));
            double epsilonP = 2.0e-12 * pv0.getPosition().getNorm().getReal();
            for (int c = 0; c < directions.length; ++c) {
                // eight points finite differences estimation of a Jacobian column
                FieldPVCoordinates<T> pvm4h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, -4 * h, directions[c]));
                FieldPVCoordinates<T> pvm3h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, -3 * h, directions[c]));
                FieldPVCoordinates<T> pvm2h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, -2 * h, directions[c]));
                FieldPVCoordinates<T> pvm1h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, -1 * h, directions[c]));
                FieldPVCoordinates<T> pvp1h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, +1 * h, directions[c]));
                FieldPVCoordinates<T> pvp2h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, +2 * h, directions[c]));
                FieldPVCoordinates<T> pvp3h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, +3 * h, directions[c]));
                FieldPVCoordinates<T> pvp4h = combined.transformPVCoordinates(new FieldPVCoordinates<>(1.0, pv0, +4 * h, directions[c]));
                FieldPVCoordinates<T> d4 = new FieldPVCoordinates<>(pvm4h, pvp4h);
                FieldPVCoordinates<T> d3 = new FieldPVCoordinates<>(pvm3h, pvp3h);
                FieldPVCoordinates<T> d2 = new FieldPVCoordinates<>(pvm2h, pvp2h);
                FieldPVCoordinates<T> d1 = new FieldPVCoordinates<>(pvm1h, pvp1h);
                double d = 1.0 / (840 * h);
                FieldPVCoordinates<T> estimatedColumn = new FieldPVCoordinates<>(-3 * d, d4, 32 * d, d3, -168 * d, d2, 672 * d, d1);
                // check analytical Jacobian against finite difference reference
                Assert.assertEquals(estimatedColumn.getPosition().getX().getReal(), jacobian[0][c].getReal(), epsilonP);
                Assert.assertEquals(estimatedColumn.getPosition().getY().getReal(), jacobian[1][c].getReal(), epsilonP);
                Assert.assertEquals(estimatedColumn.getPosition().getZ().getReal(), jacobian[2][c].getReal(), epsilonP);
                // check the rest of the matrix remains untouched
                for (int l = 3; l < jacobian.length; ++l) {
                    Assert.assertEquals(l + 0.1 * c, jacobian[l][c].getReal(), 1.0e-15);
                }
            }
            // check the rest of the matrix remains untouched
            for (int c = directions.length; c < jacobian[0].length; ++c) {
                for (int l = 0; l < jacobian.length; ++l) {
                    Assert.assertEquals(l + 0.1 * c, jacobian[l][c].getReal(), 1.0e-15);
                }
            }
        }
    }
}
Also used : FieldPVCoordinates(org.orekit.utils.FieldPVCoordinates) TimeStampedFieldPVCoordinates(org.orekit.utils.TimeStampedFieldPVCoordinates) Well19937a(org.hipparchus.random.Well19937a) RandomGenerator(org.hipparchus.random.RandomGenerator)

Example 49 with RandomGenerator

use of org.hipparchus.random.RandomGenerator in project Orekit by CS-SI.

the class FieldTransformTest method doTestIdentityLine.

private <T extends RealFieldElement<T>> void doTestIdentityLine(Field<T> field) {
    RandomGenerator random = new Well19937a(0x98603025df70db7cl);
    FieldVector3D<T> p1 = randomVector(field, 100.0, random);
    FieldVector3D<T> p2 = randomVector(field, 100.0, random);
    FieldLine<T> line = new FieldLine<>(p1, p2, 1.0e-6);
    FieldLine<T> transformed = FieldTransform.getIdentity(field).transformLine(line);
    Assert.assertSame(line, transformed);
}
Also used : FieldLine(org.hipparchus.geometry.euclidean.threed.FieldLine) Well19937a(org.hipparchus.random.Well19937a) RandomGenerator(org.hipparchus.random.RandomGenerator)

Example 50 with RandomGenerator

use of org.hipparchus.random.RandomGenerator in project Orekit by CS-SI.

the class EllipsoidTest method testLimb.

@Test
public void testLimb() throws OrekitException {
    final Ellipsoid ellipsoid = new Ellipsoid(FramesFactory.getEME2000(), 1, 2, 3);
    RandomGenerator random = new Well1024a(0xa69c430a67475af7l);
    for (int i = 0; i < 5000; ++i) {
        Vector3D observer = new Vector3D((random.nextDouble() - 0.5) * 5, (random.nextDouble() - 0.5) * 5, (random.nextDouble() - 0.5) * 5);
        Vector3D outside = new Vector3D((random.nextDouble() - 0.5) * 5, (random.nextDouble() - 0.5) * 5, (random.nextDouble() - 0.5) * 5);
        if (ellipsoid.isInside(observer)) {
            try {
                ellipsoid.pointOnLimb(observer, outside);
                Assert.fail("an exception should have been thrown");
            } catch (OrekitException oe) {
                Assert.assertEquals(OrekitMessages.POINT_INSIDE_ELLIPSOID, oe.getSpecifier());
            }
        } else {
            final Vector3D onLimb = ellipsoid.pointOnLimb(observer, outside);
            Assert.assertEquals(0, FastMath.sin(Vector3D.angle(Vector3D.crossProduct(observer, outside), Vector3D.crossProduct(observer, onLimb))), 2e-14);
            final double scaledX = onLimb.getX() / ellipsoid.getA();
            final double scaledY = onLimb.getY() / ellipsoid.getB();
            final double scaledZ = onLimb.getZ() / ellipsoid.getC();
            Assert.assertEquals(1.0, scaledX * scaledX + scaledY * scaledY + scaledZ * scaledZ, 9e-11);
            final Vector3D normal = new Vector3D(scaledX / ellipsoid.getA(), scaledY / ellipsoid.getB(), scaledZ / ellipsoid.getC()).normalize();
            final Vector3D lineOfSight = onLimb.subtract(observer).normalize();
            Assert.assertEquals(0.0, Vector3D.dotProduct(normal, lineOfSight), 5e-10);
        }
    }
}
Also used : Vector3D(org.hipparchus.geometry.euclidean.threed.Vector3D) OrekitException(org.orekit.errors.OrekitException) RandomGenerator(org.hipparchus.random.RandomGenerator) Well1024a(org.hipparchus.random.Well1024a) Test(org.junit.Test)

Aggregations

RandomGenerator (org.hipparchus.random.RandomGenerator)100 Test (org.junit.Test)78 Well19937a (org.hipparchus.random.Well19937a)73 Vector3D (org.hipparchus.geometry.euclidean.threed.Vector3D)33 Well1024a (org.hipparchus.random.Well1024a)27 DerivativeStructure (org.hipparchus.analysis.differentiation.DerivativeStructure)22 FieldPVCoordinates (org.orekit.utils.FieldPVCoordinates)22 GeodeticPoint (org.orekit.bodies.GeodeticPoint)20 Rotation (org.hipparchus.geometry.euclidean.threed.Rotation)19 PVCoordinates (org.orekit.utils.PVCoordinates)15 TimeStampedFieldPVCoordinates (org.orekit.utils.TimeStampedFieldPVCoordinates)15 DSFactory (org.hipparchus.analysis.differentiation.DSFactory)14 FieldVector3D (org.hipparchus.geometry.euclidean.threed.FieldVector3D)14 FieldAbsoluteDate (org.orekit.time.FieldAbsoluteDate)14 Frame (org.orekit.frames.Frame)10 GaussianRandomGenerator (org.hipparchus.random.GaussianRandomGenerator)8 UncorrelatedRandomVectorGenerator (org.hipparchus.random.UncorrelatedRandomVectorGenerator)8 FieldKeplerianOrbit (org.orekit.orbits.FieldKeplerianOrbit)8 OrbitType (org.orekit.orbits.OrbitType)8 FieldSpacecraftState (org.orekit.propagation.FieldSpacecraftState)8