1   /* Copyright 2002-2025 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  
19  import org.hipparchus.CalculusFieldElement;
20  import org.hipparchus.util.FastMath;
21  
22  /**
23   * Compute the probability of collision assuming the worst case described in : "Kyle Alfriend, Maruthi Akella, Joseph
24   * Frisbee, James Foster, Deok-Jin Lee, and Matthew Wilkins. Probability of ProbabilityOfCollision Error Analysis. Space
25   * Debris, 1(1):21–35, 1999.".
26   * <p>It assumes:
27   *     <ul>
28   *         <li>Short encounter leading to a linear relative motion.</li>
29   *         <li>Spherical collision object.</li>
30   *         <li>Uncorrelated positional covariance.</li>
31   *         <li>Gaussian distribution of the position uncertainties.</li>
32   *         <li>Deterministic velocity i.e. no velocity uncertainties.</li>
33   *         <li>Both objects are in circular orbits (eq 14).</li>
34   *         <li>Probability density function is constant over the collision circle (eq 18).</li>
35   *         <li>Covariance multiplied by a coefficient KSquared = MahalanobisDistanceSquared / 2 (eq 19-20).</li>
36   *     </ul>
37   * <p>
38   * By assuming a constant probability density function over the collision circle this method will,
39   * <b>most of the time</b>, give much higher probability of collision than other regular methods.
40   * That is why it is qualified as a maximum probability of collision computing method.
41   *
42   * @author Vincent Cucchietti
43   * @see <a href="https://en.wikipedia.org/wiki/Mahalanobis_distance">Mahalanobis distance.</a>
44   * @since 12.0
45   */
46  public class Alfriend1999Max extends AbstractAlfriend1999 {
47  
48      /** Empty constructor. */
49      public Alfriend1999Max() {
50          super(ShortTermEncounter2DPOCMethodType.ALFRIEND_1999_MAX.name());
51      }
52  
53      /**
54       * Compute the value of the probability of collision.
55       *
56       * @param radius sum of primary and secondary collision object equivalent sphere radii (m)
57       * @param squaredMahalanobisDistance squared Mahalanobis distance
58       * @param covarianceMatrixDeterminant covariance matrix determinant
59       *
60       * @return value of the probability of collision
61       */
62      @Override
63      public double computeValue(final double radius, final double squaredMahalanobisDistance,
64                                 final double covarianceMatrixDeterminant) {
65          return radius * radius / (squaredMahalanobisDistance * FastMath.sqrt(covarianceMatrixDeterminant) * FastMath.E);
66      }
67  
68      /**
69       * Compute the value of the probability of collision.
70       *
71       * @param radius sum of primary and secondary collision object equivalent sphere radii (m)
72       * @param squaredMahalanobisDistance squared Mahalanobis distance
73       * @param covarianceMatrixDeterminant covariance matrix determinant
74       * @param <T> type of the field elements
75       *
76       * @return value of the probability of collision
77       */
78      @Override
79      public <T extends CalculusFieldElement<T>> T computeValue(final T radius, final T squaredMahalanobisDistance,
80                                                                final T covarianceMatrixDeterminant) {
81          return radius.multiply(radius).divide(squaredMahalanobisDistance.multiply(covarianceMatrixDeterminant.sqrt())
82                                                                          .multiply(FastMath.E));
83      }
84  
85      /** {@inheritDoc} */
86      @Override
87      public boolean isAMaximumProbabilityOfCollisionMethod() {
88          return true;
89      }
90  
91      /** {@inheritDoc} */
92      @Override
93      public ShortTermEncounter2DPOCMethodType getType() {
94          return ShortTermEncounter2DPOCMethodType.ALFRIEND_1999_MAX;
95      }
96  }