src/tutorials/java/fr/cs/examples/refiningPleiades/generators/Noise.java

+/* Copyright 2013-2022 CS GROUP
+ * Licensed to CS GROUP (CS) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * CS licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *   http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package fr.cs.examples.refiningPleiades.generators;
+
+
+/** Class for adding a noise to measurements.
+ * @author Lucie Labat-Allee
+ * @author Guylaine Prat
+ * @since 2.0
+ */
+public class Noise {
+
+    /** Type of distribution. */
+    private static final int GAUSSIAN = 0;
+
+    /** Dimension. */
+    private int dimension;
+
+    /** Mean. */
+    private double[] mean;
+
+    /** Standard deviation. */
+    private double[] standardDeviation;
+
+    /** Distribution. */
+    private int distribution = GAUSSIAN;
+
+    /** Build a new instance.
+     * @param distribution noise type
+     * @param dimension noise dimension
+     */
+    public Noise(final int distribution, final int dimension) {
+
+        this.mean = new double[dimension];
+        this.standardDeviation = new double[dimension];
+        this.dimension = dimension;
+        this.distribution = distribution;
+    }
+
+    /** Get the mean.
+     * @return the mean
+     */
+    public double[] getMean() {
+        return mean.clone();
+    }
+
+    /** Set the mean.
+     * @param meanValue the mean to set
+     */
+    public void setMean(final double[] meanValue) {
+        this.mean = meanValue.clone();
+    }
+
+    /** Get the standard deviation.
+     * @return the standard deviation
+     */
+    public double[] getStandardDeviation() {
+        return standardDeviation.clone();
+    }
+
+    /** Set the standard deviation.
+     * @param standardDeviationValue the standard deviation to set
+     */
+    public void setStandardDeviation(final double[] standardDeviationValue) {
+        this.standardDeviation = standardDeviationValue.clone();
+    }
+
+    /** Get the distribution.
+     * @return the distribution
+     */
+    public int getDistribution() {
+        return distribution;
+    }
+
+    /** Get the dimension.
+     * @return the dimension
+     */
+    public int getDimension() {
+        return dimension;
+    }
+}

src/tutorials/java/fr/cs/examples/refiningPleiades/metrics/DistanceTools.java

+/* Copyright 2013-2022 CS GROUP
+ * Licensed to CS GROUP (CS) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * CS licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *   http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package fr.cs.examples.refiningPleiades.metrics;
+
+import org.hipparchus.geometry.euclidean.threed.Vector3D;
+import org.hipparchus.util.FastMath;
+import org.orekit.bodies.GeodeticPoint;
+import org.orekit.utils.Constants;
+
+/**
+ * Class for computing geodetic and anguler distance between two geodetic points.
+ * @author Jonathan Guinet
+ * @author Guylaine Prat
+ * @since 2.0
+ */
+public class DistanceTools {
+
+    /** Private constructor for utility class.
+     */
+    private DistanceTools() {
+    }
+
+    /** Choice the method for computing geodetic distance between two points.
+     * @param geoPoint1 first geodetic point (rad)
+     * @param geoPoint2 second geodetic point (rad)
+     * @param computeAngular if true, distance will be angular, otherwise will be in meters
+     * @return distance in meters or radians if flag computeAngular is true
+     */
+    public static double computeDistance(final GeodeticPoint geoPoint1, final GeodeticPoint geoPoint2,
+                                         final boolean computeAngular) {
+
+        if (computeAngular == true) {
+            return computeDistanceAngular(geoPoint1, geoPoint2);
+        } else {
+            return computeDistanceInMeter(geoPoint1, geoPoint2);
+        }
+    }
+
+    /** Compute a geodetic distance in meters between two geodetic points..
+     * @param geoPoint1 first geodetic point (rad)
+     * @param geoPoint2 second geodetic point (rad)
+     * @return distance in meters
+     */
+    public static double computeDistanceInMeter(final GeodeticPoint geoPoint1, final GeodeticPoint geoPoint2) {
+
+        // get vectors on unit sphere from angular coordinates
+        final Vector3D p1 = new Vector3D(geoPoint1.getLatitude(), geoPoint1.getLongitude()); //
+        final Vector3D p2 = new Vector3D(geoPoint2.getLatitude(), geoPoint2.getLongitude());
+        final double distance =  Constants.WGS84_EARTH_EQUATORIAL_RADIUS * Vector3D.angle(p1, p2);
+        return distance;
+    }
+
+    /** Compute an angular distance between two geodetic points.
+     * @param geoPoint1 first geodetic point (rad)
+     * @param geoPoint2 second geodetic point (rad)
+     * @return angular distance in radians
+     */
+    public static double computeDistanceAngular(final GeodeticPoint geoPoint1, final GeodeticPoint geoPoint2) {
+
+        final double lonDiff = geoPoint1.getLongitude() - geoPoint2.getLongitude();
+        final double latDiff = geoPoint1.getLatitude() - geoPoint2.getLatitude();
+        return FastMath.hypot(lonDiff, latDiff);
+    }
+}

src/tutorials/java/fr/cs/examples/refiningPleiades/metrics/LocalisationMetrics.java

+/* Copyright 2013-2022 CS GROUP
+ * Licensed to CS GROUP (CS) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * CS licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *   http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package fr.cs.examples.refiningPleiades.metrics;
+
+import java.util.Iterator;
+import java.util.Map;
+import java.util.Set;
+
+import org.hipparchus.geometry.euclidean.threed.Vector3D;
+import org.hipparchus.util.FastMath;
+import org.orekit.bodies.GeodeticPoint;
+import org.orekit.rugged.adjustment.measurements.SensorToGroundMapping;
+import org.orekit.rugged.adjustment.measurements.SensorToSensorMapping;
+import org.orekit.rugged.api.Rugged;
+import org.orekit.rugged.linesensor.LineSensor;
+import org.orekit.rugged.linesensor.SensorPixel;
+import org.orekit.rugged.utils.SpacecraftToObservedBody;
+import org.orekit.time.AbsoluteDate;
+
+/**
+ * Class for testing geometric performances in absolute location.
+ * Metrics are computed for two scenarios: ground control points and tie points.
+ * @see SensorToSensorMapping
+ * @see SensorToGroundMapping
+ * @author Jonathan Guinet
+ * @author Lucie Labat-Allee
+ * @author Guylaine Prat
+ * @since 2.0
+ */
+public class LocalisationMetrics {
+
+    /** Maximum residual distance. */
+    private double resMax;
+
+    /** Mean residual distance. */
+    private double resMean;
+
+    /** LOS distance max. */
+    private double losDistanceMax;
+
+    /** LOS distance mean. */
+    private double losDistanceMean;
+
+    /** Earth distance max. */
+    private double earthDistanceMax;
+
+    /** Earth distance mean.*/
+    private double earthDistanceMean;
+
+
+    /** Compute metrics corresponding to the Ground Control Points (GCP) study.
+     * @param measMapping Mapping of observations/measurements = the ground truth
+     * @param rugged Rugged instance
+     * @param computeAngular flag to know if distance is computed in meters (false) or with angular (true)
+     */
+    public LocalisationMetrics(final SensorToGroundMapping measMapping, final Rugged rugged, final boolean computeAngular) {
+    	
+        // Initialization
+        this.resMax = 0.0;
+        this.resMean = 0.0;
+
+        // Compute metrics - Case of Sensor to Ground mapping (Ground Control Points GCP study)
+        computeGCPmetrics(measMapping, rugged, computeAngular);
+    }
+
+    /** Compute metrics corresponding to the tie points study.
+     * @param measMapping Mapping of observations/measurements = the ground truth
+     * @param ruggedA Rugged instance corresponding to viewing model A
+     * @param ruggedB Rugged instance corresponding to viewing model B
+     * @param computeAngular flag to know if distance is computed in meters (false) or with angular (true)
+     */
+    public LocalisationMetrics(final SensorToSensorMapping measMapping, final Rugged ruggedA, final Rugged ruggedB,
+                               final boolean computeAngular) {
+
+        // Initialization
+        this.resMax = 0.0;
+        this.resMean = 0.0;
+        this.losDistanceMax = 0.0;
+        this.losDistanceMean = 0.0;
+        this.earthDistanceMax = 0.0;
+        this.earthDistanceMean = 0.0;
+
+        // Compute metrics - Case of Sensor to Sensor mapping (tie points study)
+        computeTiePointsMetrics(measMapping, ruggedA, ruggedB, computeAngular);
+    }
+
+    /**
+     * Compute metrics: case of ground control points (GCP).
+     * @param measMapping Mapping of observations/measurements = the ground truth
+     * @param rugged Rugged instance
+     * @param computeAngular flag to know if distance is computed in meters (false) or with angular (true)
+     */
+    public void computeGCPmetrics(final SensorToGroundMapping measMapping, final Rugged rugged, final boolean computeAngular) {
+
+        // Mapping of observations/measurements = the ground truth
+        final Set<Map.Entry<SensorPixel, GeodeticPoint>> measurementsMapping;
+        final LineSensor lineSensor;
+
+        // counter for compute mean distance
+        double count = 0;
+
+        // Initialization
+        measurementsMapping = measMapping.getMapping();
+        lineSensor = rugged.getLineSensor(measMapping.getSensorName());
+
+        // number of measurements
+        int nbMeas = measurementsMapping.size();
+
+        final Iterator<Map.Entry<SensorPixel, GeodeticPoint>> gtIt = measurementsMapping.iterator();
+
+        // Browse map of measurements
+        while (gtIt.hasNext()) {
+
+            final Map.Entry<SensorPixel, GeodeticPoint> gtMeas = gtIt.next();
+
+            final SensorPixel gtSP = gtMeas.getKey();
+            final GeodeticPoint gtGP = gtMeas.getValue();
+
+            final AbsoluteDate date = lineSensor.getDate(gtSP.getLineNumber());
+
+            final GeodeticPoint esGP = rugged.directLocation(date, lineSensor.getPosition(),
+                                                             lineSensor.getLOS(date,
+                                                             (int) gtSP.getPixelNumber()));
+            // Compute distance
+            double distance = DistanceTools.computeDistance(esGP, gtGP, computeAngular);
+            count += distance;
+            if (distance > resMax) {
+                resMax = distance;
+            }
+        }
+        // Mean of residuals
+        resMean = count / nbMeas;
+    }
+
+    /**
+     * Compute metrics: case of tie points.
+     * @param measMapping Mapping of observations/measurements = the ground truth
+     * @param ruggedA Rugged instance corresponding to viewing model A
+     * @param ruggedB Rugged instance corresponding to viewing model B
+     * @param computeAngular Flag to know if distance is computed in meters (false) or with angular (true)
+     */
+    public void computeTiePointsMetrics(final SensorToSensorMapping measMapping, final Rugged ruggedA, final Rugged ruggedB,
+                                      final boolean computeAngular) {
+
+        // Mapping of observations/measurements = the ground truth
+        final Set<Map.Entry<SensorPixel, SensorPixel>> measurementsMapping;
+
+        final LineSensor lineSensorA;   // line sensor corresponding to rugged A
+        final LineSensor lineSensorB;   // line sensor corresponding to rugged B
+        double count = 0;               // counter for computing remaining mean distance
+        double losDistanceCount = 0;    // counter for computing LOS distance mean
+        double earthDistanceCount = 0;  // counter for computing Earth distance mean
+        int i = 0;                      // increment of measurements
+
+        // Initialization
+        measurementsMapping = measMapping.getMapping();
+        lineSensorA = ruggedA.getLineSensor(measMapping.getSensorNameA());
+        lineSensorB = ruggedB.getLineSensor(measMapping.getSensorNameB());
+        int nbMeas = measurementsMapping.size(); // number of measurements
+
+        // Browse map of measurements
+        for (Iterator<Map.Entry<SensorPixel, SensorPixel>> gtIt = measurementsMapping.iterator();
+                        gtIt.hasNext();
+                        i++) {
+
+            if (i == measurementsMapping.size()) {
+                break;
+            }
+
+            final Map.Entry<SensorPixel, SensorPixel> gtMapping = gtIt.next();
+
+            final SensorPixel spA = gtMapping.getKey();
+            final SensorPixel spB = gtMapping.getValue();
+
+            final AbsoluteDate dateA = lineSensorA.getDate(spA.getLineNumber());
+            final AbsoluteDate dateB = lineSensorB.getDate(spB.getLineNumber());
+
+            final double pixelA = spA.getPixelNumber();
+            final double pixelB = spB.getPixelNumber();
+
+            final Vector3D losA = lineSensorA.getLOS(dateA, pixelA);
+            final Vector3D losB = lineSensorB.getLOS(dateB, pixelB);
+
+            final GeodeticPoint gpA = ruggedA.directLocation(dateA, lineSensorA.getPosition(), losA);
+            final GeodeticPoint gpB = ruggedB.directLocation(dateB, lineSensorB.getPosition(), losB);
+
+            final SpacecraftToObservedBody scToBodyA = ruggedA.getScToBody();
+
+            // Estimated distances (LOS / Earth)
+            final double[] distances = ruggedB.distanceBetweenLOS(lineSensorA, dateA, pixelA, scToBodyA, lineSensorB, dateB, pixelB);
+
+            // LOS distance control
+            final double estLosDistance =  distances[0]; // LOS distance estimation
+            if (estLosDistance > losDistanceMax) {
+                losDistanceMax = estLosDistance;
+            }
+            losDistanceCount += estLosDistance;
+
+            // Earth distance control
+            final double estEarthDistance = distances[1]; // Earth distance estimation
+            final double measEarthDistance = measMapping.getBodyDistance(i).doubleValue(); // Earth measurement distance
+            final double earthDistance =  FastMath.abs(estEarthDistance - measEarthDistance);
+
+            if (earthDistance > earthDistanceMax) {
+                earthDistanceMax = earthDistance;
+
+            }
+            earthDistanceCount += earthDistance;
+
+            // Compute remaining distance
+            double distance = DistanceTools.computeDistance(gpB, gpA, computeAngular);
+            count += distance;
+            if (distance > resMax) {
+                resMax = distance;
+            }
+        }
+
+        resMean = count / nbMeas;
+        losDistanceMean = losDistanceCount / nbMeas;
+        earthDistanceMean = earthDistanceCount / nbMeas;
+    }
+
+    /** Get the max residual.
+     * @return maximum of residuals
+     */
+    public double getMaxResidual() {
+        return resMax;
+    }
+
+    /** Get the mean of residuals.
+     * @return mean of residuals
+     */
+    public double getMeanResidual() {
+        return resMean;
+    }
+
+    /** Get the LOS maximum distance.
+     * @return LOS maximum distance
+     */
+    public double getLosMaxDistance() {
+        return losDistanceMax;
+    }
+
+    /** Get the LOS mean distance.
+     * @return LOS mean distance
+     */
+    public double getLosMeanDistance() {
+        return losDistanceMean;
+    }
+
+    /** Get the Earth distance maximum residual.
+     * @return Earth distance max residual
+     */
+    public double getEarthMaxDistance() {
+        return earthDistanceMax;
+    }
+
+    /** Get the Earth distance mean residual.
+     * @return Earth distance mean residual
+     */
+    public double getEarthMeanDistance() {
+        return earthDistanceMean;
+    }
+}

src/tutorials/java/fr/cs/examples/refiningPleiades/models/OrbitModel.java

+/* Copyright 2013-2022 CS GROUP
+ * Licensed to CS GROUP (CS) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * CS licenses this file to You under the Apache License, Version 2.0
+ * (the "License"); you may not use this file except in compliance with
+ * the License.  You may obtain a copy of the License at
+ *
+ *   http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package fr.cs.examples.refiningPleiades.models;
+
+import java.util.ArrayList;
+import java.util.List;
+
+import org.hipparchus.geometry.euclidean.threed.Rotation;
+import org.hipparchus.geometry.euclidean.threed.RotationConvention;
+import org.hipparchus.geometry.euclidean.threed.RotationOrder;
+import org.hipparchus.geometry.euclidean.threed.Vector3D;
+import org.hipparchus.util.FastMath;
+import org.orekit.attitudes.AttitudeProvider;
+import org.orekit.attitudes.LofOffset;
+import org.orekit.attitudes.NadirPointing;
+import org.orekit.attitudes.TabulatedLofOffset;
+import org.orekit.attitudes.YawCompensation;
+import org.orekit.bodies.BodyShape;
+import org.orekit.bodies.OneAxisEllipsoid;
+import org.orekit.forces.gravity.potential.GravityFieldFactory;
+import org.orekit.forces.gravity.potential.NormalizedSphericalHarmonicsProvider;
+import org.orekit.frames.Frame;
+import org.orekit.frames.FramesFactory;
+import org.orekit.frames.LOFType;
+import org.orekit.frames.Transform;
+import org.orekit.orbits.CircularOrbit;
+import org.orekit.orbits.Orbit;
+import org.orekit.orbits.PositionAngle;
+import org.orekit.propagation.Propagator;
+import org.orekit.propagation.analytical.KeplerianPropagator;
+import org.orekit.time.AbsoluteDate;
+import org.orekit.time.TimeScale;
+import org.orekit.utils.AngularDerivativesFilter;
+import org.orekit.utils.Constants;
+import org.orekit.utils.IERSConventions;
+import org.orekit.utils.PVCoordinates;
+import org.orekit.utils.TimeStampedAngularCoordinates;
+import org.orekit.utils.TimeStampedPVCoordinates;
+
+/**
+ * Orbit Model class to generate positions-velocities and attitude quaternions.
+ * <p>
+ * the aim of this class is to simulate the orbit model of an LEO satellite
+ * </p>
+ * @author Jonathan Guinet
+ * @author Guylaine Prat
+ * @since 2.0 
+ */
+
+public class OrbitModel {
+
+    /** Flag to change Attitude Law (by default Nadir Pointing Yaw Compensation). */
+    private boolean userDefinedLOFTransform;
+
+    /** User defined Roll law. */
+    private double[] lofTransformRollPoly;
+
+    /** User defined Pitch law. */
+    private double[] lofTransformPitchPoly;
+
+    /** User defined Yaw law. */
+    private double[] lofTransformYawPoly;
+
+    /** Reference date. */
+    private AbsoluteDate refDate;
+
+    /** Default constructor.
+     */
+    public OrbitModel() {
+        userDefinedLOFTransform = false;
+    }
+
+    /** Generate satellite ephemeris.
+     */
+    public static void addSatellitePV(final TimeScale gps, final Frame eme2000, final Frame itrf,
+                                      final List<TimeStampedPVCoordinates> satellitePVList,
+                                      final String absDate,
+                                      final double px, final double py, final double pz,
+                                      final double vx, final double vy, final double vz) {
+    	
+        final AbsoluteDate ephemerisDate = new AbsoluteDate(absDate, gps);
+        final Vector3D position = new Vector3D(px, py, pz);
+        final Vector3D velocity = new Vector3D(vx, vy, vz);
+        final PVCoordinates pvITRF = new PVCoordinates(position, velocity);
+        final Transform transform = itrf.getTransformTo(eme2000, ephemerisDate);
+        final Vector3D pEME2000 = transform.transformPosition(pvITRF.getPosition());
+        final Vector3D vEME2000 = transform.transformVector(pvITRF.getVelocity());
+        satellitePVList.add(new TimeStampedPVCoordinates(ephemerisDate, pEME2000, vEME2000, Vector3D.ZERO));
+    }
+
+    /** Generate satellite attitude.
+     */
+    public void addSatelliteQ(final TimeScale gps,
+                              final List<TimeStampedAngularCoordinates> satelliteQList,
+                              final String absDate,
+                              final double q0, final double q1, final double q2, final double q3) {
+    	
+        final AbsoluteDate attitudeDate = new AbsoluteDate(absDate, gps);
+        final Rotation rotation = new Rotation(q0, q1, q2, q3, true);
+        final TimeStampedAngularCoordinates pair = 
+                new TimeStampedAngularCoordinates(attitudeDate, rotation, Vector3D.ZERO, Vector3D.ZERO);
+        satelliteQList.add(pair);
+    }
+
+    /** Create an Earth.
+     * @return the Earth as the WGS84 ellipsoid
+     */
+    public BodyShape createEarth() {
+    	
+        return new OneAxisEllipsoid(Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
+                                    Constants.WGS84_EARTH_FLATTENING,
+                                    FramesFactory.getITRF(IERSConventions.IERS_2010, true));
+    }
+
+    /** Created a gravity field.
+     * @return normalized spherical harmonics coefficients
+     */
+    public NormalizedSphericalHarmonicsProvider createGravityField() {
+    	
+        return GravityFieldFactory.getNormalizedProvider(12, 12);
+    }
+
+    /** Create an orbit at a chosen date.
+     * @param mu Earth gravitational constant
+     * @return the orbit
+     */
+    public Orbit createOrbit(final double mu, final AbsoluteDate date) {
+    	
+        // the following orbital parameters have been computed using
+        // Orekit tutorial about phasing, using the following configuration:
+        //
+        // orbit.date = 2012-01-01T00:00:00.000
+        // phasing.orbits.number = 143
+        // phasing.days.number = 10
+        // sun.synchronous.reference.latitude = 0
+        // sun.synchronous.reference.ascending = false
+        // sun.synchronous.mean.solar.time = 10:30:00
+        // gravity.field.degree = 12
+        // gravity.field.order = 12
+
+        final Frame eme2000 = FramesFactory.getEME2000();
+        return new CircularOrbit(694000.0 + Constants.WGS84_EARTH_EQUATORIAL_RADIUS,
+                                 -4.029194321683225E-4,
+                                 0.0013530362644647786,
+                                 FastMath.toRadians(98.2), // Pleiades inclination 98.2 deg
+                                 FastMath.toRadians(-86.47 + 180),
+                                 FastMath.toRadians(135.9 + 0.3),
+                                 PositionAngle.TRUE,
+                                 eme2000,
+                                 date,
+                                 mu);
+    }
+
+    /** Set the Local Orbital Frame characteristics.
+     */
+    public void setLOFTransform(final double[] rollPoly, final double[] pitchPoly,
+                                final double[] yawPoly, final AbsoluteDate date) {
+
+        this.userDefinedLOFTransform = true;
+        lofTransformRollPoly = rollPoly.clone();
+        lofTransformPitchPoly = pitchPoly.clone();
+        lofTransformYawPoly = yawPoly.clone();
+        this.refDate = date;
+    }
+
+    /** Recompute the polynom coefficients with shift.
+     */
+    private double getPoly(final double[] poly, final double shift) {
+    	
+        double val = 0.0;
+        for (int coef = 0; coef < poly.length; coef++) {
+            val = val + poly[coef] * FastMath.pow(shift, coef);
+        }
+        return val;
+    }
+
+    /** Get the offset.
+     */
+   private Rotation getOffset(final BodyShape earth, final Orbit orbit, final double shift){
+    	
+        final LOFType type = LOFType.VVLH;
+        final double roll = getPoly(lofTransformRollPoly, shift);
+        final double pitch = getPoly(lofTransformPitchPoly, shift);
+        final double yaw = getPoly(lofTransformYawPoly, shift);
+
+        final LofOffset law = new LofOffset(orbit.getFrame(), type, RotationOrder.XYZ,
+                                      roll, pitch, yaw);
+        final Rotation offsetAtt = law.
+                                   getAttitude(orbit, orbit.getDate().shiftedBy(shift), orbit.getFrame()).
+                                   getRotation();
+        final NadirPointing nadirPointing = new NadirPointing(orbit.getFrame(), earth);
+        final Rotation nadirAtt = nadirPointing.
+                                  getAttitude(orbit, orbit.getDate().getDate().shiftedBy(shift), orbit.getFrame()).
+                                  getRotation();
+        final Rotation offsetProper = offsetAtt.compose(nadirAtt.revert(), RotationConvention.VECTOR_OPERATOR);
+
+        return offsetProper;
+    }
+
+   /** Create the attitude provider.
+    */
+    public AttitudeProvider createAttitudeProvider(final BodyShape earth, final Orbit orbit) {
+
+        if (userDefinedLOFTransform) {
+            final LOFType type = LOFType.VVLH;
+
+            final List<TimeStampedAngularCoordinates> list = new ArrayList<TimeStampedAngularCoordinates>();
+
+            for (double shift = -10.0; shift < +10.0; shift += 1e-2) {
+                list.add(new TimeStampedAngularCoordinates(refDate
+                    .shiftedBy(shift), getOffset(earth, orbit, shift),
+                                                           Vector3D.ZERO,
+                                                           Vector3D.ZERO));
+            }
+
+            final TabulatedLofOffset tabulated = new TabulatedLofOffset(orbit.getFrame(), type, list,
+                                                                        2, AngularDerivativesFilter.USE_R);
+
+            return tabulated;
+        } else {
+            return new YawCompensation(orbit.getFrame(), new NadirPointing(orbit.getFrame(), earth));
+        }
+    }
+
+   /** Generate the orbit.
+    */
+   public List<TimeStampedPVCoordinates> orbitToPV(final Orbit orbit, final BodyShape earth,
+    		                                        final AbsoluteDate minDate, final AbsoluteDate maxDate,
+    		                                        final double step) {
+    	
+        final Propagator propagator = new KeplerianPropagator(orbit);
+
+        propagator.setAttitudeProvider(createAttitudeProvider(earth, orbit));
+
+        propagator.propagate(minDate);
+        final List<TimeStampedPVCoordinates> list = new ArrayList<TimeStampedPVCoordinates>();
+        propagator.getMultiplexer().add(step,
+        		                        currentState ->
+                                        list.add(new TimeStampedPVCoordinates(currentState.getDate(),
+                                                                              currentState.getPVCoordinates().getPosition(),
+                                                                              currentState.getPVCoordinates().getVelocity(),
+                                                                              Vector3D.ZERO)));
+        propagator.propagate(maxDate);
+
+        return list;
+    }
+
+   /** Generate the attitude.
+    */
+   public List<TimeStampedAngularCoordinates> orbitToQ(final Orbit orbit, final BodyShape earth,
+    		                                            final AbsoluteDate minDate, final AbsoluteDate maxDate,
+    		                                            final double step) {
+    	
+        final Propagator propagator = new KeplerianPropagator(orbit);
+        propagator.setAttitudeProvider(createAttitudeProvider(earth, orbit));
+        propagator.propagate(minDate);
+        final List<TimeStampedAngularCoordinates> list = new ArrayList<>();
+        propagator.getMultiplexer().add(step,
+                                        currentState ->
+                                        list.add(new TimeStampedAngularCoordinates(currentState.getDate(),
+                                                                                   currentState.getAttitude().getRotation(),
+                                                                                   Vector3D.ZERO,
+                                                                                   Vector3D.ZERO)));
+        propagator.propagate(maxDate);
+
+        return list;
+    }
+}

src/tutorials/java/AffinagePleiades/PleiadesViewingModel.java → src/tutorials/java/fr/cs/examples/refiningPleiades/models/PleiadesViewingModel.java

-/* Copyright 2013-2016 CS Systèmes d'Information
- * Licensed to CS Systèmes d'Information (CS) under one or more
+/* Copyright 2013-2022 CS GROUP
+ * Licensed to CS GROUP (CS) under one or more
  * contributor license agreements.  See the NOTICE file distributed with
  * this work for additional information regarding copyright ownership.
  * CS licenses this file to You under the Apache License, Version 2.0
@@ -14,19 +14,18 @@
  * See the License for the specific language governing permissions and
  * limitations under the License.
  */
-package AffinagePleiades;
-
+package fr.cs.examples.refiningPleiades.models;
 
+import java.util.ArrayList;
+import java.util.List;
 
 import org.hipparchus.geometry.euclidean.threed.Rotation;
 import org.hipparchus.geometry.euclidean.threed.RotationConvention;
 import org.hipparchus.geometry.euclidean.threed.Vector3D;
-
 import org.hipparchus.util.FastMath;
-import java.util.ArrayList;
-import java.util.List;
-
-
+import org.orekit.rugged.linesensor.LineDatation;
+import org.orekit.rugged.linesensor.LineSensor;
+import org.orekit.rugged.linesensor.LinearLineDatation;
 import org.orekit.rugged.los.FixedRotation;
 import org.orekit.rugged.los.FixedZHomothety;
 import org.orekit.rugged.los.LOSBuilder;
@@ -34,157 +33,146 @@ import org.orekit.rugged.los.TimeDependentLOS;
 import org.orekit.time.AbsoluteDate;
 import org.orekit.time.TimeScale;
 import org.orekit.time.TimeScalesFactory;
-import org.orekit.rugged.linesensor.LinearLineDatation;
-import org.orekit.rugged.linesensor.LineDatation;
-import org.orekit.rugged.linesensor.LineSensor;
-
-import org.orekit.rugged.errors.RuggedException;
-import org.orekit.errors.OrekitException;
 
+import fr.cs.examples.refiningPleiades.Refining;
 
 /**
- * PleiadesViewingModel class definition 
- * @author Jonathan Guinet, Lucie LabatAllee
- *
+ * Pleiades viewing model class definition.
+ * The aim of this class is to simulate PHR sensor.
+ * @author Jonathan Guinet
+ * @author Lucie Labat-Allee
+ * @author Guylaine Prat
+ * @since 2.0
  */
 
-
 public class PleiadesViewingModel {
-	
-    
-	/** intrinsic Pleiades parameters */
-	public  double fov = 1.65; // 20km - alt 694km
-    public  int dimension = 40000;
-    private final double linePeriod =  1e-4;
-    
-	public  double angle;
-	public  LineSensor lineSensor;	
-	public  String date;
-
-	
-	
-	private String sensorName;
-	
+
+    /** Pleiades parameters. */
+    private static final double FOV = 1.65; // 20km - alt 694km
+    private static final int DIMENSION = 40000;
+    private static final double LINE_PERIOD =  1.e-4; 
+
+    private double angle;
+    private LineSensor lineSensor;
+    private String date;
+
+    private String sensorName;
+
+
     /** Simple constructor.
      * <p>
-     *  initialize PleiadesViewingModel with 
+     *  initialize PleiadesViewingModel with
      *  sensorName="line", incidenceAngle = 0.0, date = "2016-01-01T12:00:00.0"
      * </p>
      */
-    public PleiadesViewingModel() throws RuggedException, OrekitException {	
-        this("line",0.0,"2016-01-01T12:00:00.0");	
+    public PleiadesViewingModel(final String sensorName) {
+    	
+        this(sensorName, 0.0, "2016-01-01T12:00:00.0");
     }
-    
-    /** PleiadesViewingModel constructor
+
+    /** PleiadesViewingModel constructor.
      * @param sensorName sensor name
      * @param incidenceAngle incidence angle
-     * @param referenceDate reference date 
-     * @throws RuggedException
-     * @throws OrekitException
+     * @param referenceDate reference date
      */
-    public PleiadesViewingModel(final String sensorName,final double incidenceAngle,final String referenceDate) throws RuggedException, OrekitException { 
+    public PleiadesViewingModel(final String sensorName, final double incidenceAngle, final String referenceDate) {
+    	
         this.sensorName = sensorName;
         this.date = referenceDate;
-        this.angle = incidenceAngle;       
-        this.createLineSensor();        
-    }    
-    
-	
-	public   LOSBuilder rawLOS(Vector3D center, Vector3D normal, double halfAperture, int n)
-	{
-	
-        List<Vector3D> list = new ArrayList<Vector3D>(n);
+        this.angle = incidenceAngle;
+        this.createLineSensor();
+    }
+
+    /** Create raw fixed Line Of sight
+     */
+    public LOSBuilder rawLOS(final Vector3D center, final Vector3D normal, final double halfAperture, final int n) {
+
+        final List<Vector3D> list = new ArrayList<Vector3D>(n);
         for (int i = 0; i < n; ++i) {
-            double alpha = (halfAperture * (2 * i + 1 - n)) / (n - 1);
+            final double alpha = (halfAperture * (2 * i + 1 - n)) / (n - 1);
             list.add(new Rotation(normal, alpha, RotationConvention.VECTOR_OPERATOR).applyTo(center));
         }
-        
-        
+
         return new LOSBuilder(list);
-    }	
-	
-
-	
-	public  TimeDependentLOS buildLOS()
-	{	
-	LOSBuilder losBuilder = rawLOS(new Rotation(Vector3D.PLUS_I,
-            FastMath.toRadians(this.angle),
-            RotationConvention.VECTOR_OPERATOR).applyTo(Vector3D.PLUS_K), Vector3D.PLUS_I, FastMath.toRadians(fov/2), dimension); 	
-
-	losBuilder.addTransform(new FixedRotation("roll",  Vector3D.MINUS_I, 0.00));
-    losBuilder.addTransform(new FixedRotation("pitch", Vector3D.MINUS_J, 0.00));
-    losBuilder.addTransform(new FixedZHomothety("factor", 1.0));
-      
-    return  losBuilder.build();
-	}
-    
-	
-	public  AbsoluteDate getDatationReference() throws OrekitException
-	{
-	    // We use Orekit for handling time and dates, and Rugged for defining the datation model:
-	    TimeScale utc = TimeScalesFactory.getUTC();
-	    return new AbsoluteDate(date, utc);
-	}
-	
-
-	
-
-	
-	public  AbsoluteDate getMinDate() throws RuggedException
-	{
-		return lineSensor.getDate(0);
-	}
-	
-	public  AbsoluteDate  getMaxDate() throws RuggedException
-	{
-		return lineSensor.getDate(dimension);
-	}
-
-	public  LineSensor  getLineSensor() {
-			return lineSensor;
-	}
-
-	public  String getSensorName() {
-			return sensorName;
-	}
-	
-	
-	private  void  createLineSensor()
-			throws RuggedException, OrekitException {
-	
-		//System.out.println("add line sensor");
-		
-		
-	      
-	    // TBN: refining data are read only for level L1B or L1C
-		//System.out.println("refining info");
-	      
+    }
+
+    /** Build a LOS provider
+     */
+    public TimeDependentLOS buildLOS() {
+    	
+        final LOSBuilder losBuilder = rawLOS(new Rotation(Vector3D.PLUS_I,
+                                                          FastMath.toRadians(this.angle),
+                                                          RotationConvention.VECTOR_OPERATOR).applyTo(Vector3D.PLUS_K),
+                                             Vector3D.PLUS_I, FastMath.toRadians(FOV / 2), DIMENSION);
+
+        losBuilder.addTransform(new FixedRotation(sensorName + Refining.getRollsuffix(),  Vector3D.MINUS_I, 0.00));
+        losBuilder.addTransform(new FixedRotation(sensorName + Refining.getPitchsuffix(), Vector3D.MINUS_J, 0.00));
+
+        // factor is a common parameters shared between all Pleiades models
+        losBuilder.addTransform(new FixedZHomothety(Refining.getFactorname(), 1.0));
+
+        return losBuilder.build();
+    }
+
+
+    /** Get the reference date.
+     */
+    public AbsoluteDate getDatationReference() {
+
+    	// We use Orekit for handling time and dates, and Rugged for defining the datation model:
+    	final TimeScale utc = TimeScalesFactory.getUTC();
+
+    	return new AbsoluteDate(date, utc);
+    }
+
+    /** Get the min date.
+     */
+   public  AbsoluteDate getMinDate() {
+        return lineSensor.getDate(0);
+    }
+
+   /** Get the max date.
+    */
+   public  AbsoluteDate  getMaxDate() {
+        return lineSensor.getDate(DIMENSION);
+    }
+
+   /** Get the line sensor.
+    */
+   public  LineSensor  getLineSensor() {
+        return lineSensor;
+    }
+
+   /** Get the sensor name.
+    */
+   public  String getSensorName() {
+        return sensorName;
+    }
+
+   /** Get the number of lines of the sensor
+    * @return the number of lines of the sensor
+    */
+public int getDimension() {
+        return DIMENSION;
+    }
+
+   /** Create the line sensor
+    */
+   private void createLineSensor() {
+
         // Offset of the MSI from center of mass of satellite
-		//System.out.println("MSI offset from center of mass of satellite");
         // one line sensor
-        // position: 1.5m in front (+X) and 20 cm above (-Z) of the S/C center of mass
         // los: swath in the (YZ) plane, looking at 50° roll, 2.6" per pixel
-        //Vector3D msiOffset = new Vector3D(1.5, 0, -0.2);
-		Vector3D msiOffset = new Vector3D(0, 0, 0);
-        
-
-        
-        // to do build complex los 
-        TimeDependentLOS lineOfSight = buildLOS();
-        
-        
-        final double rate =  1 / linePeriod;
-        // linear datation model: at reference time we get the middle line, and the rate is one line every 1.5ms
+        final Vector3D msiOffset = new Vector3D(0, 0, 0);
 
-        
-        
-        LineDatation lineDatation = new LinearLineDatation(getDatationReference(), dimension / 2, rate);
-        //LineDatation lineDatation = new LinearLineDatation(absDate, 1d, 20);
-        lineSensor = new LineSensor(sensorName, lineDatation, msiOffset, lineOfSight);
-        
+        final TimeDependentLOS lineOfSight = buildLOS();
 
-	}
+        final double rate =  1. / LINE_PERIOD;
+        // linear datation model: at reference time we get the middle line, and the rate is one line every 1.5ms
 
+        final LineDatation lineDatation = new LinearLineDatation(getDatationReference(), DIMENSION / 2, rate);
 
+        lineSensor = new LineSensor(sensorName, lineDatation, msiOffset, lineOfSight);
+    }
 }