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.utils;
18  
19  import org.hipparchus.util.CombinatoricsUtils;
20  import org.hipparchus.util.FastMath;
21  
22  /**
23   * Computes the P<sub>nm</sub>(t) coefficients.
24   * <p>
25   * The computation of the Legendre polynomials is performed following:
26   * Heiskanen and Moritz, Physical Geodesy, 1967, eq. 1-62
27   * </p>
28   * @since 11.0
29   * @author Bryan Cazabonne
30   */
31  public class LegendrePolynomials {
32  
33      /** Array for the Legendre polynomials. */
34      private double[][] pCoef;
35  
36      /** Create Legendre polynomials for the given degree and order.
37       * @param degree degree of the spherical harmonics
38       * @param order order of the spherical harmonics
39       * @param t argument for polynomials calculation
40       */
41      public LegendrePolynomials(final int degree, final int order,
42                                 final double t) {
43  
44          // Initialize array
45          this.pCoef = new double[degree + 1][order + 1];
46  
47          final double t2 = t * t;
48  
49          for (int n = 0; n <= degree; n++) {
50  
51              // m shall be <= n (Heiskanen and Moritz, 1967, pp 21)
52              for (int m = 0; m <= FastMath.min(n, order); m++) {
53  
54                  // r = int((n - m) / 2)
55                  final int r = (int) (n - m) / 2;
56                  double sum = 0.;
57                  for (int k = 0; k <= r; k++) {
58                      final double term = FastMath.pow(-1.0, k) * CombinatoricsUtils.factorialDouble(2 * n - 2 * k) /
59                                      (CombinatoricsUtils.factorialDouble(k) * CombinatoricsUtils.factorialDouble(n - k) *
60                                       CombinatoricsUtils.factorialDouble(n - m - 2 * k)) *
61                                       FastMath.pow(t, n - m - 2 * k);
62                      sum = sum + term;
63                  }
64  
65                  pCoef[n][m] = FastMath.pow(2, -n) * FastMath.pow(1.0 - t2, 0.5 * m) * sum;
66  
67              }
68  
69          }
70  
71      }
72  
73      /** Return the coefficient P<sub>nm</sub>.
74       * @param n index
75       * @param m index
76       * @return The coefficient P<sub>nm</sub>
77       */
78      public double getPnm(final int n, final int m) {
79          return pCoef[n][m];
80      }
81  
82  }