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.data; 18 19 import org.hipparchus.CalculusFieldElement; 20 21 /** Class for tide terms. 22 * <p> 23 * BEWARE! For consistency with all the other Poisson series terms, 24 * the elements in γ, l, l', F, D and Ω are ADDED together to compute 25 * the argument of the term. In classical tides series, the computed 26 * argument is cGamma * γ - (cL * l + cLPrime * l' + cF * F + cD * D 27 * + cOmega * Ω). So at parsing time, the signs of cL, cLPrime, cF, 28 * cD and cOmega must already have been reversed so the addition 29 * performed here will work. This is done automatically when the 30 * parser has been configured with a call to {@link 31 * PoissonSeriesParser#withDoodson(int, int)} as the relationship 32 * between the Doodson arguments and the traditional Delaunay 33 * arguments ensures the proper sign is known. 34 * </p> 35 * @author Luc Maisonobe 36 */ 37 class TideTerm extends SeriesTerm { 38 39 /** Coefficient for γ = GMST + π tide parameter. */ 40 private final int cGamma; 41 42 /** Coefficient for mean anomaly of the Moon. */ 43 private final int cL; 44 45 /** Coefficient for mean anomaly of the Sun. */ 46 private final int cLPrime; 47 48 /** Coefficient for L - Ω where L is the mean longitude of the Moon. */ 49 private final int cF; 50 51 /** Coefficient for mean elongation of the Moon from the Sun. */ 52 private final int cD; 53 54 /** Coefficient for mean longitude of the ascending node of the Moon. */ 55 private final int cOmega; 56 57 /** Build a tide term for nutation series. 58 * @param cGamma coefficient for γ = GMST + π tide parameter 59 * @param cL coefficient for mean anomaly of the Moon 60 * @param cLPrime coefficient for mean anomaly of the Sun 61 * @param cF coefficient for L - Ω where L is the mean longitude of the Moon 62 * @param cD coefficient for mean elongation of the Moon from the Sun 63 * @param cOmega coefficient for mean longitude of the ascending node of the Moon 64 */ 65 TideTerm(final int cGamma, 66 final int cL, final int cLPrime, final int cF, final int cD, final int cOmega) { 67 this.cGamma = cGamma; 68 this.cL = cL; 69 this.cLPrime = cLPrime; 70 this.cF = cF; 71 this.cD = cD; 72 this.cOmega = cOmega; 73 } 74 75 /** {@inheritDoc} */ 76 protected double argument(final BodiesElements elements) { 77 return cGamma * elements.getGamma() + 78 cL * elements.getL() + cLPrime * elements.getLPrime() + cF * elements.getF() + 79 cD * elements.getD() + cOmega * elements.getOmega(); 80 } 81 82 /** {@inheritDoc} */ 83 protected double argumentDerivative(final BodiesElements elements) { 84 return cGamma * elements.getGammaDot() + 85 cL * elements.getLDot() + cLPrime * elements.getLPrimeDot() + cF * elements.getFDot() + 86 cD * elements.getDDot() + cOmega * elements.getOmegaDot(); 87 } 88 89 /** {@inheritDoc} */ 90 protected <T extends CalculusFieldElement<T>> T argument(final FieldBodiesElements<T> elements) { 91 return elements.getGamma().multiply(cGamma). 92 add(elements.getL().multiply(cL)). 93 add(elements.getLPrime().multiply(cLPrime)). 94 add(elements.getF().multiply(cF)). 95 add(elements.getD().multiply(cD)). 96 add(elements.getOmega().multiply(cOmega)); 97 } 98 99 /** {@inheritDoc} */ 100 protected <T extends CalculusFieldElement<T>> T argumentDerivative(final FieldBodiesElements<T> elements) { 101 return elements.getGammaDot().multiply(cGamma). 102 add(elements.getLDot().multiply(cL)). 103 add(elements.getLPrimeDot().multiply(cLPrime)). 104 add(elements.getFDot().multiply(cF)). 105 add(elements.getDDot().multiply(cD)). 106 add(elements.getOmegaDot().multiply(cOmega)); 107 } 108 109 }