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.models.earth.troposphere;
18  
19  import org.hipparchus.util.FastMath;
20  import org.hipparchus.util.Precision;
21  import org.junit.jupiter.api.Assertions;
22  import org.junit.jupiter.api.BeforeAll;
23  import org.junit.jupiter.api.BeforeEach;
24  import org.junit.jupiter.api.Test;
25  import org.orekit.Utils;
26  import org.orekit.bodies.GeodeticPoint;
27  import org.orekit.errors.OrekitException;
28  import org.orekit.time.AbsoluteDate;
29  import org.orekit.time.TimeScalesFactory;
30  import org.orekit.utils.TrackingCoordinates;
31  
32  public class ViennaThreeTest {
33  
34      private static double epsilon = 1e-6;
35  
36      @BeforeAll
37      public static void setUpGlobal() {
38          Utils.setDataRoot("atmosphere");
39      }
40  
41      @BeforeEach
42      public void setUp() throws OrekitException {
43          Utils.setDataRoot("regular-data:potential/shm-format");
44      }
45  
46      @Test
47      public void testMappingFactors() {
48  
49          // Site:     latitude:  37.5°
50          //           longitude: 277.5°
51          //           height:    824 m
52          //
53          // Date:     25 November 2018 at 0h UT
54          //
55          // Values: ah  = 0.00123462
56          //         aw  = 0.00047101
57          //         zhd = 2.1993 m
58          //         zwd = 0.0690 m
59          //
60          // Values taken from: http://vmf.geo.tuwien.ac.at/trop_products/GRID/5x5/VMF3/VMF3_OP/2018/VMF3_20181125.H00
61          //
62          // Expected mapping factors : hydrostatic -> 1.621024
63          //                                    wet -> 1.623023
64          //
65          // Expected outputs are obtained by performing the Matlab script vmf3.m provided by TU WIEN:
66          // http://vmf.geo.tuwien.ac.at/codes/
67          //
68  
69          final AbsoluteDate date = new AbsoluteDate(2018, 11, 25, TimeScalesFactory.getUTC());
70  
71          final double latitude     = FastMath.toRadians(37.5);
72          final double longitude    = FastMath.toRadians(277.5);
73          final double height       = 824.0;
74          final GeodeticPoint point = new GeodeticPoint(latitude, longitude, height);
75  
76          final TrackingCoordinates trackingCoordinates = new TrackingCoordinates(0.0, FastMath.toRadians(38.0), 0.0);
77          final double expectedHydro = 1.621024;
78          final double expectedWet   = 1.623023;
79  
80          final ViennaThree model = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
81                                                    new ConstantAzimuthalGradientProvider(null),
82                                                    new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
83                                                    TimeScalesFactory.getUTC());
84  
85          final double[] computedMapping = model.mappingFactors(trackingCoordinates, point,
86                                                                TroposphericModelUtils.STANDARD_ATMOSPHERE,
87                                                                date);
88  
89          Assertions.assertEquals(expectedHydro, computedMapping[0], epsilon);
90          Assertions.assertEquals(expectedWet,   computedMapping[1], epsilon);
91      }
92  
93      @Test
94      public void testLowElevation() {
95  
96          // Site:     latitude:  37.5°
97          //           longitude: 277.5°
98          //           height:    824 m
99          //
100         // Date:     25 November 2018 at 0h UT
101         //
102         // Values: ah  = 0.00123462
103         //         aw  = 0.00047101
104         //         zhd = 2.1993 m
105         //         zwd = 0.0690 m
106         //
107         // Values taken from: http://vmf.geo.tuwien.ac.at/trop_products/GRID/5x5/VMF3/VMF3_OP/2018/VMF3_20181125.H00
108         //
109         // Expected mapping factors : hydrostatic -> 10.132802
110         //                                    wet -> 10.879154
111         //
112         // Expected outputs are obtained by performing the Matlab script vmf3.m provided by TU WIEN:
113         // http://vmf.geo.tuwien.ac.at/codes/
114         //
115 
116         final AbsoluteDate date = new AbsoluteDate(2018, 11, 25, TimeScalesFactory.getUTC());
117 
118         final double latitude     = FastMath.toRadians(37.5);
119         final double longitude    = FastMath.toRadians(277.5);
120         final double height       = 824.0;
121         final GeodeticPoint point = new GeodeticPoint(latitude, longitude, height);
122 
123         final TrackingCoordinates trackingCoordinates = new TrackingCoordinates(0.0, FastMath.toRadians(5.0), 0.0);
124         final double expectedHydro = 10.132802;
125         final double expectedWet   = 10.879154;
126 
127         final ViennaThree model = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
128                                                   new ConstantAzimuthalGradientProvider(null),
129                                                   new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
130                                                   TimeScalesFactory.getUTC());
131 
132         final double[] computedMapping = model.mappingFactors(trackingCoordinates, point,
133                                                               TroposphericModelUtils.STANDARD_ATMOSPHERE,
134                                                               date);
135 
136         Assertions.assertEquals(expectedHydro, computedMapping[0], epsilon);
137         Assertions.assertEquals(expectedWet,   computedMapping[1], epsilon);
138     }
139 
140     @Test
141     public void testHightElevation() {
142 
143         // Site:     latitude:  37.5°
144         //           longitude: 277.5°
145         //           height:    824 m
146         //
147         // Date:     25 November 2018 at 0h UT
148         //
149         // Values: ah  = 0.00123462
150         //         aw  = 0.00047101
151         //         zhd = 2.1993 m
152         //         zwd = 0.0690 m
153         //
154         // Values taken from: http://vmf.geo.tuwien.ac.at/trop_products/GRID/5x5/VMF3/VMF3_OP/2018/VMF3_20181125.H00
155         //
156         // Expected mapping factors : hydrostatic -> 1.003810
157         //                                    wet -> 1.003816
158         //
159         // Expected outputs are obtained by performing the Matlab script vmf3.m provided by TU WIEN:
160         // http://vmf.geo.tuwien.ac.at/codes/
161         //
162 
163         final AbsoluteDate date = new AbsoluteDate(2018, 11, 25, TimeScalesFactory.getUTC());
164 
165         final double latitude     = FastMath.toRadians(37.5);
166         final double longitude    = FastMath.toRadians(277.5);
167         final double height       = 824.0;
168         final GeodeticPoint point = new GeodeticPoint(latitude, longitude, height);
169 
170         final TrackingCoordinates trackingCoordinates = new TrackingCoordinates(0.0, FastMath.toRadians(85.0), 0.0);
171         final double expectedHydro = 1.003810;
172         final double expectedWet   = 1.003816;
173 
174         final ViennaThree model = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
175                                                   new ConstantAzimuthalGradientProvider(null),
176                                                   new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
177                                                   TimeScalesFactory.getUTC());
178 
179         final double[] computedMapping = model.mappingFactors(trackingCoordinates, point,
180                                                               TroposphericModelUtils.STANDARD_ATMOSPHERE,
181                                                               date);
182 
183         Assertions.assertEquals(expectedHydro, computedMapping[0], epsilon);
184         Assertions.assertEquals(expectedWet,   computedMapping[1], epsilon);
185     }
186 
187     @Test
188     public void testDelay() {
189         final double        azimuth   = 30.0;
190         final double        elevation = 10.0;
191         final double        height    = 100.0;
192         final AbsoluteDate  date      = new AbsoluteDate();
193         final GeodeticPoint point     = new GeodeticPoint(FastMath.toRadians(37.5), FastMath.toRadians(277.5), height);
194         final ViennaThree   model     = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
195                                                         new ConstantAzimuthalGradientProvider(null),
196                                                         new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
197                                                         TimeScalesFactory.getUTC());
198         final TroposphericDelay delay = model.pathDelay(new TrackingCoordinates(FastMath.toRadians(azimuth), FastMath.toRadians(elevation), 0.0),
199                                                         point,
200                                                         TroposphericModelUtils.STANDARD_ATMOSPHERE,
201                                                         model.getParameters(date), date);
202         Assertions.assertEquals( 2.1993, delay.getZh(),    1.0e-4);
203         Assertions.assertEquals( 0.069,  delay.getZw(),    1.0e-4);
204         Assertions.assertEquals(12.2124, delay.getSh(),    1.0e-4);
205         Assertions.assertEquals( 0.3916, delay.getSw(),    1.0e-4);
206         Assertions.assertEquals(12.6041, delay.getDelay(), 1.0e-4);
207     }
208 
209     @Test
210     public void testDelayWithAzimuthalAsymmetry() {
211         final double        azimuth   = 30.0;
212         final double        elevation = 10.0;
213         final double        height    = 100.0;
214         final AbsoluteDate  date      = new AbsoluteDate();
215         final GeodeticPoint point     = new GeodeticPoint(FastMath.toRadians(37.5), FastMath.toRadians(277.5), height);
216         final ViennaThree   model     = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
217                                                         new ConstantAzimuthalGradientProvider(new AzimuthalGradientCoefficients(12.0, 4.5,
218                                                                                                                                 0.8, 1.25)),
219                                                         new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
220                                                         TimeScalesFactory.getUTC());
221         final TroposphericDelay delay = model.pathDelay(new TrackingCoordinates(FastMath.toRadians(azimuth), FastMath.toRadians(elevation), 0.0),
222                                                         point,
223                                                         TroposphericModelUtils.STANDARD_ATMOSPHERE,
224                                                         model.getParameters(date), date);
225         Assertions.assertEquals( 2.1993,                      delay.getZh(),    1.0e-4);
226         Assertions.assertEquals( 0.069,                       delay.getZw(),    1.0e-4);
227         Assertions.assertEquals(12.2124 + 373.8241,           delay.getSh(),    1.0e-4); // second term is due to azimuthal gradient
228         Assertions.assertEquals( 0.3916 +  38.9670,           delay.getSw(),    1.0e-4); // second term is due to azimuthal gradient
229         Assertions.assertEquals(12.6041 + 373.8241 + 38.9670, delay.getDelay(), 1.0e-4);
230     }
231 
232     @Test
233     public void testFixedHeight() {
234         final AbsoluteDate date = new AbsoluteDate();
235         final GeodeticPoint point = new GeodeticPoint(FastMath.toRadians(37.5), FastMath.toRadians(277.5), 350.0);
236         ViennaThree model = new ViennaThree(new ConstantViennaAProvider(new ViennaACoefficients(0.00123462, 0.00047101)),
237                                             new ConstantAzimuthalGradientProvider(null),
238                                             new ConstantTroposphericModel(new TroposphericDelay(2.1993, 0.0690, 0, 0)),
239                                             TimeScalesFactory.getUTC());
240         double lastDelay = Double.MAX_VALUE;
241         // delay shall decline with increasing elevation angle
242         for (double elev = 10d; elev < 90d; elev += 8d) {
243             final double delay = model.pathDelay(new TrackingCoordinates(0.0, FastMath.toRadians(elev), 0.0),
244                                                  point,
245                                                  TroposphericModelUtils.STANDARD_ATMOSPHERE,
246                                                  model.getParameters(date), date).getDelay();
247             Assertions.assertTrue(Precision.compareTo(delay, lastDelay, epsilon) < 0);
248             lastDelay = delay;
249         }
250     }
251 
252 }