AbstractAlfriend1999.java

  1. /* Copyright 2002-2024 CS GROUP
  2.  * Licensed to CS GROUP (CS) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * CS licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *   http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.orekit.ssa.collision.shorttermencounter.probability.twod;

  18. import org.hipparchus.CalculusFieldElement;
  19. import org.hipparchus.linear.BlockFieldMatrix;
  20. import org.hipparchus.linear.BlockRealMatrix;
  21. import org.hipparchus.linear.FieldLUDecomposition;
  22. import org.hipparchus.linear.FieldMatrix;
  23. import org.hipparchus.linear.LUDecomposition;
  24. import org.hipparchus.linear.RealMatrix;
  25. import org.hipparchus.util.MathArrays;
  26. import org.orekit.ssa.metrics.FieldProbabilityOfCollision;
  27. import org.orekit.ssa.metrics.ProbabilityOfCollision;

  28. /**
  29.  * Abstract class for Alfriend1999 normal and maximised methods as they share lots of similarities.
  30.  *
  31.  * @author Vincent Cucchietti
  32.  * @since 12.0
  33.  */
  34. public abstract class AbstractAlfriend1999 extends AbstractShortTermEncounter2DPOCMethod {

  35.     /**
  36.      * Default constructor.
  37.      *
  38.      * @param name name of the method
  39.      */
  40.     protected AbstractAlfriend1999(final String name) {
  41.         super(name);
  42.     }

  43.     /** {@inheritDoc} */
  44.     @Override
  45.     public ProbabilityOfCollision compute(final double xm, final double ym, final double sigmaX, final double sigmaY,
  46.                                           final double radius) {
  47.         // Reconstruct necessary values from inputs
  48.         final double squaredMahalanobisDistance =
  49.                 ShortTermEncounter2DDefinition.computeSquaredMahalanobisDistance(xm, ym, sigmaX, sigmaY);

  50.         // Compute covariance matrix determinant
  51.         final double covarianceMatrixDeterminant = computeCovarianceMatrixDeterminant(sigmaX, sigmaY);

  52.         // Compute probability of collision
  53.         final double value = computeValue(radius, squaredMahalanobisDistance, covarianceMatrixDeterminant);

  54.         return new ProbabilityOfCollision(value, getName(), this.isAMaximumProbabilityOfCollisionMethod());
  55.     }

  56.     /** {@inheritDoc} */
  57.     @Override
  58.     public <T extends CalculusFieldElement<T>> FieldProbabilityOfCollision<T> compute(final T xm, final T ym, final T sigmaX,
  59.                                                                                       final T sigmaY,
  60.                                                                                       final T radius) {
  61.         // Reconstruct necessary values from inputs
  62.         final T squaredMahalanobisDistance =
  63.                 FieldShortTermEncounter2DDefinition.computeSquaredMahalanobisDistance(xm, ym, sigmaX, sigmaY);

  64.         // Compute covariance matrix determinant
  65.         final T covarianceMatrixDeterminant = computeCovarianceMatrixDeterminant(sigmaX, sigmaY);

  66.         // Compute probability of collision
  67.         final T value = computeValue(radius, squaredMahalanobisDistance, covarianceMatrixDeterminant);

  68.         return new FieldProbabilityOfCollision<>(value, getName(), this.isAMaximumProbabilityOfCollisionMethod());
  69.     }

  70.     /**
  71.      * Compute the covariance matrix determinant.
  72.      *
  73.      * @param sigmaX square root of the x-axis eigen value of the diagonalized combined covariance matrix projected onto the
  74.      * collision plane (m)
  75.      * @param sigmaY square root of the y-axis eigen value of the diagonalized combined covariance matrix projected onto the
  76.      * collision plane (m)
  77.      *
  78.      * @return covariance matrix determinant
  79.      */
  80.     private double computeCovarianceMatrixDeterminant(final double sigmaX, final double sigmaY) {
  81.         // Rebuild covariance matrix
  82.         final RealMatrix covarianceMatrix = new BlockRealMatrix(new double[][] {
  83.                 { sigmaX * sigmaX, 0 },
  84.                 { 0, sigmaY * sigmaY }
  85.         });

  86.         // Compute determinant
  87.         return new LUDecomposition(covarianceMatrix).getDeterminant();
  88.     }

  89.     /**
  90.      * Compute the covariance matrix determinant.
  91.      *
  92.      * @param sigmaX square root of the x-axis eigen value of the diagonalized combined covariance matrix projected onto the
  93.      * collision plane (m)
  94.      * @param sigmaY square root of the y-axis eigen value of the diagonalized combined covariance matrix projected onto the
  95.      * collision plane (m)
  96.      * @param <T> type of the field elements
  97.      *
  98.      * @return covariance matrix determinant
  99.      */
  100.     private <T extends CalculusFieldElement<T>> T computeCovarianceMatrixDeterminant(final T sigmaX, final T sigmaY) {
  101.         // Rebuild covariance matrix
  102.         final T[][] covarianceMatrixData = MathArrays.buildArray(sigmaX.getField(), 2, 2);
  103.         covarianceMatrixData[0][0] = sigmaX.multiply(sigmaX);
  104.         covarianceMatrixData[1][1] = sigmaY.multiply(sigmaY);

  105.         final FieldMatrix<T> covarianceMatrix = new BlockFieldMatrix<>(covarianceMatrixData);

  106.         // Compute determinant
  107.         return new FieldLUDecomposition<>(covarianceMatrix).getDeterminant();
  108.     }

  109.     /**
  110.      * Compute the value of the probability of collision.
  111.      *
  112.      * @param radius sum of primary and secondary collision object equivalent sphere radii (m)
  113.      * @param squaredMahalanobisDistance squared Mahalanobis distance
  114.      * @param covarianceMatrixDeterminant covariance matrix determinant
  115.      *
  116.      * @return value of the probability of collision
  117.      */
  118.     abstract double computeValue(double radius, double squaredMahalanobisDistance, double covarianceMatrixDeterminant);

  119.     /**
  120.      * Compute the value of the probability of collision.
  121.      *
  122.      * @param radius sum of primary and secondary collision object equivalent sphere radii (m)
  123.      * @param squaredMahalanobisDistance squared Mahalanobis distance
  124.      * @param covarianceMatrixDeterminant covariance matrix determinant
  125.      * @param <T> type of the field elements
  126.      *
  127.      * @return value of the probability of collision
  128.      */
  129.     abstract <T extends CalculusFieldElement<T>> T computeValue(T radius, T squaredMahalanobisDistance,
  130.                                                                 T covarianceMatrixDeterminant);
  131. }