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.Field;
21 import org.hipparchus.util.FastMath;
22 import org.hipparchus.util.MathArrays;
23 import org.orekit.ssa.metrics.FieldProbabilityOfCollision;
24 import org.orekit.ssa.metrics.ProbabilityOfCollision;
25
26 /**
27 * Compute the probability of collision using the method described in : <br> "Chan, K. “Collision Probability Analyses for
28 * Earth Orbiting Satellites.” In Space Cooperation into the 21st Century: 7th AAS/JRS/CSA Symposium, International Space
29 * Conference of Pacific-Basin Societies (ISCOPS; formerly PISSTA) (July 15-18, 1997, Nagasaki, Japan), edited by Peter M.
30 * Bainum, et al., 1033-1048. Advances in the Astronautical Sciences Series 96. San Diego, California: Univelt, 1997. (Zeroth
31 * order analytical expression).
32 * <p>
33 * This method is also described in depth in : "CHAN, F. Kenneth, et al. Spacecraft collision probability. El Segundo, CA :
34 * Aerospace Press, 2008."
35 * <p>
36 * It assumes :
37 * <ul>
38 * <li>Short encounter leading to a linear relative motion.</li>
39 * <li>Spherical collision object.</li>
40 * <li>Uncorrelated positional covariance.</li>
41 * <li>Gaussian distribution of the position uncertainties.</li>
42 * <li>Deterministic velocity i.e. no velocity uncertainties.</li>
43 * <li>Approximate ellipse by a disk</li>
44 * </ul>
45 *
46 * @author Vincent Cucchietti
47 * @since 12.0
48 */
49 public class Chan1997 extends AbstractShortTermEncounter2DPOCMethod {
50
51 /** Empty constructor. */
52 public Chan1997() {
53 super(ShortTermEncounter2DPOCMethodType.CHAN_1997.name());
54 }
55
56 /** {@inheritDoc} */
57 public ProbabilityOfCollision compute(final double xm, final double ym,
58 final double sigmaX, final double sigmaY,
59 final double radius) {
60
61 // Intermediary terms u and v
62 final double u = radius * radius / (sigmaX * sigmaY);
63 final double v = (xm * xm / (sigmaX * sigmaX)) + (ym * ym / (sigmaY * sigmaY));
64
65 // Number of terms M recommended by Chan
66 final int M;
67
68 if (u <= 0.01 || v <= 1) {
69 M = 3;
70 } else if (u > 0.01 && u <= 1 || v > 1 && v <= 9) {
71 M = 10;
72 } else if (u > 1 && u <= 25 || v > 9 && v <= 25) {
73 M = 20;
74 } else {
75 M = 60;
76 }
77
78 double t = 1.0;
79 double s = 1.0;
80 double sum = 1.0;
81
82 // first iteration
83 double value = FastMath.exp(-v * 0.5) * t - FastMath.exp(-(u + v) * 0.5) * t * sum;
84
85 // iterative expression
86 for (int i = 1; i < M; i++) {
87 t = (v * 0.5) / i * t;
88 s = (u * 0.5) / i * s;
89 sum = sum + s;
90 value = MathArrays.linearCombination(1, value, FastMath.exp(-v * 0.5), t) -
91 FastMath.exp(-(u + v) * 0.5) * t * sum;
92 }
93
94 return new ProbabilityOfCollision(value, getName(), isAMaximumProbabilityOfCollisionMethod());
95 }
96
97 /** {@inheritDoc} */
98 public <T extends CalculusFieldElement<T>> FieldProbabilityOfCollision<T> compute(final T xm, final T ym,
99 final T sigmaX, final T sigmaY,
100 final T radius) {
101
102 // Intermediary terms u and v
103 final T u = radius.pow(2).divide(sigmaX.multiply(sigmaY));
104 final T v = xm.divide(sigmaX).pow(2).add(ym.divide(sigmaY).pow(2));
105
106 // Number of terms M recommended by Chan
107 final int M;
108
109 if (u.getReal() <= 0.01 || v.getReal() <= 1) {
110 M = 3;
111 } else if (u.getReal() > 0.01 && u.getReal() <= 1 || v.getReal() > 1 && v.getReal() <= 9) {
112 M = 10;
113 } else if (u.getReal() > 1 && u.getReal() <= 25 || v.getReal() > 9 && v.getReal() <= 25) {
114 M = 20;
115 } else {
116 M = 60;
117 }
118
119 final Field<T> field = radius.getField();
120
121 T t = field.getOne();
122 T s = field.getOne();
123 T sum = field.getOne();
124
125 // first iteration
126 T value = v.multiply(-0.5).exp().multiply(t)
127 .subtract(u.add(v).multiply(-0.5).exp().multiply(t).multiply(sum));
128
129 // iterative expression
130 for (int i = 1; i < M; i++) {
131 t = v.multiply(0.5).divide(i).multiply(t);
132 s = u.multiply(0.5).divide(i).multiply(s);
133 sum = sum.add(s);
134
135 value = value.add(v.multiply(-0.5).exp().multiply(t)
136 .subtract(u.add(v).multiply(-0.5).exp().multiply(t).multiply(sum)));
137 }
138
139 return new FieldProbabilityOfCollision<>(value, getName(), isAMaximumProbabilityOfCollisionMethod());
140 }
141
142 /** {@inheritDoc} */
143 @Override
144 public ShortTermEncounter2DPOCMethodType getType() {
145 return ShortTermEncounter2DPOCMethodType.CHAN_1997;
146 }
147
148 }