1   /* Copyright 2022-2025 Thales Alenia Space
2    * Licensed to CS Communication & Systèmes (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.weather.water;
18  
19  import org.hipparchus.CalculusFieldElement;
20  import org.hipparchus.util.FastMath;
21  import org.orekit.models.earth.troposphere.TroposphericModelUtils;
22  
23  /** Official model CIPM-2007 (identical to CIPM-1981/91) from Comité International des Poids et Mesures.
24   * <p>
25   * This water vapor model is the one from Giacomo and Davis as indicated in IERS TN 32, chap. 9.
26   * </p>
27   * @see <a href="https://www.nist.gov/system/files/documents/calibrations/CIPM-2007.pdf">Revised
28   * formula for the density of moist air (CIPM-2007), Metrologia 45 (2008) 149–155</a>
29   *
30   * @author Luc Maisonobe
31   * @since 12.1
32   */
33  public class CIPM2007 implements WaterVaporPressureProvider {
34  
35      /** Laurent series coefficient for degree +2. */
36      private static final double L_P2 = 1.2378847e-5;
37  
38      /** Laurent series coefficient for degree +1. */
39      private static final double L_P1 = -1.9121316e-2;
40  
41      /** Laurent series coefficient for degree 0. */
42      private static final double L_0 = 33.93711047;
43  
44      /** Laurent series coefficient for degree -1. */
45      private static final double L_M1 = -6343.1645;
46  
47      /** Celsius temperature offset. */
48      private static final double CELSIUS = 273.15;
49  
50      /** Constant enhancement factor. */
51      private static final double F_0 = 1.00062;
52  
53      /** Pressure enhancement factor. */
54      private static final double F_P = 3.14e-6;
55  
56      /** Temperature enhancement factor. */
57      private static final double F_T2 = 5.6e-7;
58  
59      /** {@inheritDoc} */
60      @Override
61      public double waterVaporPressure(final double p, final double t, final double rh) {
62  
63          // saturation water vapor, equation A1.1 (now in Pa, not hPa)
64          final double psv = FastMath.exp(t * (t * L_P2 + L_P1) + L_0 + L_M1 / t);
65  
66          // enhancement factor, equation A1.2
67          final double tC = t - CELSIUS;
68          final double fw = TroposphericModelUtils.HECTO_PASCAL.fromSI(p) * F_P + tC * tC * F_T2 + F_0;
69  
70          return rh * fw * psv;
71  
72      }
73  
74      /** {@inheritDoc} */
75      @Override
76      public <T extends CalculusFieldElement<T>> T waterVaporPressure(final T p, final T t, final T rh) {
77  
78          // saturation water vapor, equation A1.1 (now in Pa, not hPa)
79          final T psv = FastMath.exp(t.multiply(t.multiply(L_P2).add(L_P1)).add(L_0).add(t.reciprocal().multiply(L_M1)));
80  
81          // enhancement factor, equation A1.2
82          final T tC = t.subtract(CELSIUS);
83          final T fw = TroposphericModelUtils.HECTO_PASCAL.fromSI(p).multiply(F_P).add(tC.multiply(tC).multiply(F_T2)).add(F_0);
84  
85          return rh.multiply(fw).multiply(psv);
86  
87      }
88  
89  }