1 /* Copyright 2002-2021 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.models.earth.troposphere;
18
19 import org.hipparchus.Field;
20 import org.hipparchus.CalculusFieldElement;
21 import org.hipparchus.util.FastMath;
22
23 /**
24 * Utility class for tropospheric models.
25 * @author Bryan Cazabonne
26 * @since 11.0
27 */
28 public class TroposphericModelUtils {
29
30 /**
31 * Private constructor as class is a utility.
32 */
33 private TroposphericModelUtils() {
34 // Nothing to do
35 }
36
37 /** Compute the mapping function related to the coefficient values and the elevation.
38 * @param a a coefficient
39 * @param b b coefficient
40 * @param c c coefficient
41 * @param elevation the elevation of the satellite, in radians.
42 * @return the value of the function at a given elevation
43 */
44 public static double mappingFunction(final double a, final double b, final double c, final double elevation) {
45 final double sinE = FastMath.sin(elevation);
46 // Numerator
47 final double numMP = 1 + a / (1 + b / (1 + c));
48 // Denominator
49 final double denMP = sinE + a / (sinE + b / (sinE + c));
50
51 final double fElevation = numMP / denMP;
52
53 return fElevation;
54 }
55
56 /** Compute the mapping function related to the coefficient values and the elevation.
57 * @param <T> type of the elements
58 * @param a a coefficient
59 * @param b b coefficient
60 * @param c c coefficient
61 * @param elevation the elevation of the satellite, in radians.
62 * @return the value of the function at a given elevation
63 */
64 public static <T extends CalculusFieldElement<T>> T mappingFunction(final T a, final T b, final T c, final T elevation) {
65 final T sinE = FastMath.sin(elevation);
66 // Numerator
67 final T numMP = a.divide(b.divide(c.add(1.0)).add(1.0)).add(1.0);
68 // Denominator
69 final T denMP = a.divide(b.divide(c.add(sinE)).add(sinE)).add(sinE);
70
71 final T fElevation = numMP.divide(denMP);
72
73 return fElevation;
74 }
75
76 /** This method computes the height correction for the hydrostatic
77 * component of the mapping function.
78 * The formulas are given by Neill's paper, 1996:
79 *<p>
80 * Niell A. E. (1996)
81 * "Global mapping functions for the atmosphere delay of radio wavelengths,”
82 * J. Geophys. Res., 101(B2), pp. 3227–3246, doi: 10.1029/95JB03048.
83 *</p>
84 * @param elevation the elevation of the satellite, in radians.
85 * @param height the height of the station in m above sea level.
86 * @return the height correction, in m
87 */
88 public static double computeHeightCorrection(final double elevation, final double height) {
89 final double fixedHeight = FastMath.max(0.0, height);
90 final double sinE = FastMath.sin(elevation);
91 // Ref: Eq. 4
92 final double function = TroposphericModelUtils.mappingFunction(2.53e-5, 5.49e-3, 1.14e-3, elevation);
93 // Ref: Eq. 6
94 final double dmdh = (1 / sinE) - function;
95 // Ref: Eq. 7
96 final double correction = dmdh * (fixedHeight / 1000.0);
97 return correction;
98 }
99
100 /** This method computes the height correction for the hydrostatic
101 * component of the mapping function.
102 * The formulas are given by Neill's paper, 1996:
103 *<p>
104 * Niell A. E. (1996)
105 * "Global mapping functions for the atmosphere delay of radio wavelengths,”
106 * J. Geophys. Res., 101(B2), pp. 3227–3246, doi: 10.1029/95JB03048.
107 *</p>
108 * @param <T> type of the elements
109 * @param elevation the elevation of the satellite, in radians.
110 * @param height the height of the station in m above sea level.
111 * @param field field to which the elements belong
112 * @return the height correction, in m
113 */
114 public static <T extends CalculusFieldElement<T>> T computeHeightCorrection(final T elevation, final T height, final Field<T> field) {
115 final T zero = field.getZero();
116 final T fixedHeight = FastMath.max(zero, height);
117 final T sinE = FastMath.sin(elevation);
118 // Ref: Eq. 4
119 final T function = TroposphericModelUtils.mappingFunction(zero.add(2.53e-5), zero.add(5.49e-3), zero.add(1.14e-3), elevation);
120 // Ref: Eq. 6
121 final T dmdh = sinE.reciprocal().subtract(function);
122 // Ref: Eq. 7
123 final T correction = dmdh.multiply(fixedHeight.divide(1000.0));
124 return correction;
125 }
126
127 }