Examples with PVCoordinates org.orekit.utils.PVCoordinates used on opensource projects

Example 36 with PVCoordinates

use of org.orekit.utils.PVCoordinates in project Orekit by CS-SI.

the class L2TransformProvider method getTransform.

/**
 * {@inheritDoc}
 */
@Override
public Transform getTransform(final AbsoluteDate date) throws OrekitException {
    final PVCoordinates pv21 = secondaryBody.getPVCoordinates(date, frame);
    final Vector3D translation = getL2(pv21.getPosition()).negate();
    final Rotation rotation = new Rotation(pv21.getPosition(), pv21.getVelocity(), Vector3D.PLUS_I, Vector3D.PLUS_J);
    return new Transform(date, new Transform(date, translation), new Transform(date, rotation));
}

Example 37 with PVCoordinates

use of org.orekit.utils.PVCoordinates in project Orekit by CS-SI.

the class NumericalPropagator method tolerances.

/**
 * Estimate tolerance vectors for integrators.
 * <p>
 * The errors are estimated from partial derivatives properties of orbits,
 * starting from a scalar position error specified by the user.
 * Considering the energy conservation equation V = sqrt(mu (2/r - 1/a)),
 * we get at constant energy (i.e. on a Keplerian trajectory):
 * <pre>
 * V² r |dV| = mu |dr|
 * </pre>
 * <p> So we deduce a scalar velocity error consistent with the position error.
 * From here, we apply orbits Jacobians matrices to get consistent errors
 * on orbital parameters.
 *
 * <p>
 * The tolerances are only <em>orders of magnitude</em>, and integrator tolerances
 * are only local estimates, not global ones. So some care must be taken when using
 * these tolerances. Setting 1mm as a position error does NOT mean the tolerances
 * will guarantee a 1mm error position after several orbits integration.
 * </p>
 * @param dP user specified position error
 * @param orbit reference orbit
 * @param type propagation type for the meaning of the tolerance vectors elements
 * (it may be different from {@code orbit.getType()})
 * @return a two rows array, row 0 being the absolute tolerance error and row 1
 * being the relative tolerance error
 * @exception OrekitException if Jacobian is singular
 */
public static double[][] tolerances(final double dP, final Orbit orbit, final OrbitType type) throws OrekitException {
    // estimate the scalar velocity error
    final PVCoordinates pv = orbit.getPVCoordinates();
    final double r2 = pv.getPosition().getNormSq();
    final double v = pv.getVelocity().getNorm();
    final double dV = orbit.getMu() * dP / (v * r2);
    final double[] absTol = new double[7];
    final double[] relTol = new double[7];
    // we set the mass tolerance arbitrarily to 1.0e-6 kg, as mass evolves linearly
    // with trust, this often has no influence at all on propagation
    absTol[6] = 1.0e-6;
    if (type == OrbitType.CARTESIAN) {
        absTol[0] = dP;
        absTol[1] = dP;
        absTol[2] = dP;
        absTol[3] = dV;
        absTol[4] = dV;
        absTol[5] = dV;
    } else {
        // convert the orbit to the desired type
        final double[][] jacobian = new double[6][6];
        final Orbit converted = type.convertType(orbit);
        converted.getJacobianWrtCartesian(PositionAngle.TRUE, jacobian);
        for (int i = 0; i < 6; ++i) {
            final double[] row = jacobian[i];
            absTol[i] = FastMath.abs(row[0]) * dP + FastMath.abs(row[1]) * dP + FastMath.abs(row[2]) * dP + FastMath.abs(row[3]) * dV + FastMath.abs(row[4]) * dV + FastMath.abs(row[5]) * dV;
            if (Double.isNaN(absTol[i])) {
                throw new OrekitException(OrekitMessages.SINGULAR_JACOBIAN_FOR_ORBIT_TYPE, type);
            }
        }
    }
    Arrays.fill(relTol, dP / FastMath.sqrt(r2));
    return new double[][] { absTol, relTol };
}

Example 38 with PVCoordinates

use of org.orekit.utils.PVCoordinates in project Orekit by CS-SI.

the class TLEPropagator method computePVCoordinates.

/**
 * Retrieves the position and velocity.
 * @return the computed PVCoordinates.
 * @exception OrekitException if current orbit is out of supported range
 * (too large eccentricity, too low perigee ...)
 */
private PVCoordinates computePVCoordinates() throws OrekitException {
    // Long period periodics
    final double axn = e * FastMath.cos(omega);
    double temp = 1.0 / (a * (1.0 - e * e));
    final double xlcof = 0.125 * TLEConstants.A3OVK2 * sini0 * (3.0 + 5.0 * cosi0) / (1.0 + cosi0);
    final double aycof = 0.25 * TLEConstants.A3OVK2 * sini0;
    final double xll = temp * xlcof * axn;
    final double aynl = temp * aycof;
    final double xlt = xl + xll;
    final double ayn = e * FastMath.sin(omega) + aynl;
    final double elsq = axn * axn + ayn * ayn;
    final double capu = MathUtils.normalizeAngle(xlt - xnode, FastMath.PI);
    double epw = capu;
    double ecosE = 0;
    double esinE = 0;
    double sinEPW = 0;
    double cosEPW = 0;
    // Dundee changes:  items dependent on cosio get recomputed:
    final double cosi0Sq = cosi0 * cosi0;
    final double x3thm1 = 3.0 * cosi0Sq - 1.0;
    final double x1mth2 = 1.0 - cosi0Sq;
    final double x7thm1 = 7.0 * cosi0Sq - 1.0;
    if (e > (1 - 1e-6)) {
        throw new OrekitException(OrekitMessages.TOO_LARGE_ECCENTRICITY_FOR_PROPAGATION_MODEL, e);
    }
    // Solve Kepler's' Equation.
    final double newtonRaphsonEpsilon = 1e-12;
    for (int j = 0; j < 10; j++) {
        boolean doSecondOrderNewtonRaphson = true;
        sinEPW = FastMath.sin(epw);
        cosEPW = FastMath.cos(epw);
        ecosE = axn * cosEPW + ayn * sinEPW;
        esinE = axn * sinEPW - ayn * cosEPW;
        final double f = capu - epw + esinE;
        if (FastMath.abs(f) < newtonRaphsonEpsilon) {
            break;
        }
        final double fdot = 1.0 - ecosE;
        double delta_epw = f / fdot;
        if (j == 0) {
            final double maxNewtonRaphson = 1.25 * FastMath.abs(e);
            doSecondOrderNewtonRaphson = false;
            if (delta_epw > maxNewtonRaphson) {
                delta_epw = maxNewtonRaphson;
            } else if (delta_epw < -maxNewtonRaphson) {
                delta_epw = -maxNewtonRaphson;
            } else {
                doSecondOrderNewtonRaphson = true;
            }
        }
        if (doSecondOrderNewtonRaphson) {
            delta_epw = f / (fdot + 0.5 * esinE * delta_epw);
        }
        epw += delta_epw;
    }
    // Short period preliminary quantities
    temp = 1.0 - elsq;
    final double pl = a * temp;
    final double r = a * (1.0 - ecosE);
    double temp2 = a / r;
    final double betal = FastMath.sqrt(temp);
    temp = esinE / (1.0 + betal);
    final double cosu = temp2 * (cosEPW - axn + ayn * temp);
    final double sinu = temp2 * (sinEPW - ayn - axn * temp);
    final double u = FastMath.atan2(sinu, cosu);
    final double sin2u = 2.0 * sinu * cosu;
    final double cos2u = 2.0 * cosu * cosu - 1.0;
    final double temp1 = TLEConstants.CK2 / pl;
    temp2 = temp1 / pl;
    // Update for short periodics
    final double rk = r * (1.0 - 1.5 * temp2 * betal * x3thm1) + 0.5 * temp1 * x1mth2 * cos2u;
    final double uk = u - 0.25 * temp2 * x7thm1 * sin2u;
    final double xnodek = xnode + 1.5 * temp2 * cosi0 * sin2u;
    final double xinck = i + 1.5 * temp2 * cosi0 * sini0 * cos2u;
    // Orientation vectors
    final double sinuk = FastMath.sin(uk);
    final double cosuk = FastMath.cos(uk);
    final double sinik = FastMath.sin(xinck);
    final double cosik = FastMath.cos(xinck);
    final double sinnok = FastMath.sin(xnodek);
    final double cosnok = FastMath.cos(xnodek);
    final double xmx = -sinnok * cosik;
    final double xmy = cosnok * cosik;
    final double ux = xmx * sinuk + cosnok * cosuk;
    final double uy = xmy * sinuk + sinnok * cosuk;
    final double uz = sinik * sinuk;
    // Position and velocity
    final double cr = 1000 * rk * TLEConstants.EARTH_RADIUS;
    final Vector3D pos = new Vector3D(cr * ux, cr * uy, cr * uz);
    final double rdot = TLEConstants.XKE * FastMath.sqrt(a) * esinE / r;
    final double rfdot = TLEConstants.XKE * FastMath.sqrt(pl) / r;
    final double xn = TLEConstants.XKE / (a * FastMath.sqrt(a));
    final double rdotk = rdot - xn * temp1 * x1mth2 * sin2u;
    final double rfdotk = rfdot + xn * temp1 * (x1mth2 * cos2u + 1.5 * x3thm1);
    final double vx = xmx * cosuk - cosnok * sinuk;
    final double vy = xmy * cosuk - sinnok * sinuk;
    final double vz = sinik * cosuk;
    final double cv = 1000.0 * TLEConstants.EARTH_RADIUS / 60.0;
    final Vector3D vel = new Vector3D(cv * (rdotk * ux + rfdotk * vx), cv * (rdotk * uy + rfdotk * vy), cv * (rdotk * uz + rfdotk * vz));
    return new PVCoordinates(pos, vel);
}

Example 39 with PVCoordinates

use of org.orekit.utils.PVCoordinates in project Orekit by CS-SI.

the class AbstractPropagatorConverter method setSample.

/**
 * Set the states sample.
 * @param states spacecraft states sample
 * @exception OrekitException if position/velocity cannot be extracted from sample
 */
private void setSample(final List<SpacecraftState> states) throws OrekitException {
    this.sample = states;
    if (onlyPosition) {
        target = new double[states.size() * 3];
        weight = new double[states.size() * 3];
    } else {
        target = new double[states.size() * 6];
        weight = new double[states.size() * 6];
    }
    int k = 0;
    for (int i = 0; i < states.size(); i++) {
        final PVCoordinates pv = states.get(i).getPVCoordinates(frame);
        // position
        target[k] = pv.getPosition().getX();
        weight[k++] = 1;
        target[k] = pv.getPosition().getY();
        weight[k++] = 1;
        target[k] = pv.getPosition().getZ();
        weight[k++] = 1;
        // velocity
        if (!onlyPosition) {
            // velocity weight relative to position
            final double r2 = pv.getPosition().getNormSq();
            final double v = pv.getVelocity().getNorm();
            final double vWeight = v * r2 / states.get(i).getMu();
            target[k] = pv.getVelocity().getX();
            weight[k++] = vWeight;
            target[k] = pv.getVelocity().getY();
            weight[k++] = vWeight;
            target[k] = pv.getVelocity().getZ();
            weight[k++] = vWeight;
        }
    }
}

Example 40 with PVCoordinates

use of org.orekit.utils.PVCoordinates in project Orekit by CS-SI.

the class JacobianPropagatorConverter method fillRows.

/**
 * Fill up a few value rows (either 6 or 3 depending on velocities used or not).
 * @param value values array
 * @param row first row index
 * @param state spacecraft state
 * @exception OrekitException if Jacobians matrices cannot be retrieved
 */
private void fillRows(final double[] value, final int row, final SpacecraftState state) throws OrekitException {
    final PVCoordinates pv = state.getPVCoordinates(getFrame());
    value[row] = pv.getPosition().getX();
    value[row + 1] = pv.getPosition().getY();
    value[row + 2] = pv.getPosition().getZ();
    if (!isOnlyPosition()) {
        value[row + 3] = pv.getVelocity().getX();
        value[row + 4] = pv.getVelocity().getY();
        value[row + 5] = pv.getVelocity().getZ();
    }
}