FundamentalNutationArguments.java

  1. /* Copyright 2002-2013 CS Systèmes d'Information
  2.  * Licensed to CS Systèmes d'Information (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;

  19. import java.io.BufferedReader;
  20. import java.io.IOException;
  21. import java.io.InputStream;
  22. import java.io.InputStreamReader;
  23. import java.io.Serializable;
  24. import java.util.ArrayList;
  25. import java.util.Arrays;
  26. import java.util.HashMap;
  27. import java.util.List;
  28. import java.util.Map;
  29. import java.util.regex.Matcher;
  30. import java.util.regex.Pattern;

  32. import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
  33. import org.apache.commons.math3.exception.util.DummyLocalizable;
  34. import org.apache.commons.math3.util.FastMath;
  35. import org.orekit.errors.OrekitException;
  36. import org.orekit.errors.OrekitMessages;
  37. import org.orekit.time.AbsoluteDate;
  38. import org.orekit.time.TimeFunction;
  39. import org.orekit.time.TimeScale;
  40. import org.orekit.utils.IERSConventions;

  42. /**
  43.  * Class computing the fundamental arguments for nutation and tides.
  44.  * <p>
  45.  * The fundamental arguments are split in two sets:
  46.  * </p>
  47.  * <ul>
  48.  *   <li>the Delaunay arguments for Moon and Sun effects</li>
  49.  *   <li>the planetary arguments for other planets</li>
  50.  * </ul>
  51.  *
  52.  * @author Luc Maisonobe
  53.  * @see SeriesTerm
  54.  * @see PoissonSeries
  55.  * @see BodiesElements
  56.  */
  57. public class FundamentalNutationArguments implements Serializable {

  59.     /** Serializable UID. */
  60.     private static final long serialVersionUID = 20131209L;

  62.     /** IERS conventions to use. */
  63.     private final IERSConventions conventions;

  65.     /** Time scale for GMST computation. */
  66.     private final TimeScale timeScale;

  68.     /** Function computing Greenwich Mean Sidereal Time. */
  69.     private final transient TimeFunction<DerivativeStructure> gmstFunction;

  71.     // luni-solar Delaunay arguments

  73.     /** Coefficients for mean anomaly of the Moon. */
  74.     private final double[] lCoefficients;

  76.     /** Coefficients for mean anomaly of the Sun. */
  77.     private final double[] lPrimeCoefficients;

  79.     /** Coefficients for L - &Omega; where L is the mean longitude of the Moon. */
  80.     private final double[] fCoefficients;

  82.     /** Coefficients for mean elongation of the Moon from the Sun. */
  83.     private final double[] dCoefficients;

  85.     /** Coefficients for mean longitude of the ascending node of the Moon. */
  86.     private final double[] omegaCoefficients;

  88.     // planetary nutation arguments

  90.     /** Coefficients for mean Mercury longitude. */
  91.     private final double[] lMeCoefficients;

  93.     /** Coefficients for mean Venus longitude. */
  94.     private final double[] lVeCoefficients;

  96.     /** Coefficients for mean Earth longitude. */
  97.     private final double[] lECoefficients;

  99.     /** Coefficients for mean Mars longitude. */
  100.     private final double[] lMaCoefficients;

  102.     /** Coefficients for mean Jupiter longitude. */
  103.     private final double[] lJCoefficients;

  105.     /** Coefficients for mean Saturn longitude. */
  106.     private final double[] lSaCoefficients;

  108.     /** Coefficients for mean Uranus longitude. */
  109.     private final double[] lUCoefficients;

  111.     /** Coefficients for mean Neptune longitude. */
  112.     private final double[] lNeCoefficients;

  114.     /** Coefficients for general accumulated precession. */
  115.     private final double[] paCoefficients;

  117.     /** Build a model of fundamental arguments from an IERS table file.
  118.      * @param conventions IERS conventions to use
  119.      * @param timeScale time scale for GMST computation
  120.      * (may be null if tide parameter γ = GMST + π is not needed)
  121.      * @param stream stream containing the IERS table
  122.      * @param name name of the resource file (for error messages only)
  123.      * @exception OrekitException if stream is null or the table cannot be parsed
  124.      */
  125.     public FundamentalNutationArguments(final IERSConventions conventions,
  126.                                         final TimeScale timeScale,
  127.                                         final InputStream stream, final String name)
  128.         throws OrekitException {
  129.         this(conventions, timeScale, parseCoefficients(stream, name));
  130.     }

  132.     /** Build a model of fundamental arguments from an IERS table file.
  133.      * @param conventions IERS conventions to use
  134.      * @param timeScale time scale for GMST computation
  135.      * (may be null if tide parameter γ = GMST + π is not needed)
  136.      * @param coefficients list of coefficients arrays (all 14 arrays must be provided,
  137.      * the 5 Delaunay first and the 9 planetary afterwards)
  138.      * @exception OrekitException if GMST function cannot be retrieved
  139.      * @since 6.1
  140.      */
  141.     public FundamentalNutationArguments(final IERSConventions conventions, final TimeScale timeScale,
  142.                                         final List<double[]> coefficients)
  143.         throws OrekitException {
  144.         this.conventions        = conventions;
  145.         this.timeScale          = timeScale;
  146.         this.gmstFunction       = (timeScale == null) ? null : conventions.getGMSTFunction(timeScale);
  147.         this.lCoefficients      = coefficients.get( 0);
  148.         this.lPrimeCoefficients = coefficients.get( 1);
  149.         this.fCoefficients      = coefficients.get( 2);
  150.         this.dCoefficients      = coefficients.get( 3);
  151.         this.omegaCoefficients  = coefficients.get( 4);
  152.         this.lMeCoefficients    = coefficients.get( 5);
  153.         this.lVeCoefficients    = coefficients.get( 6);
  154.         this.lECoefficients     = coefficients.get( 7);
  155.         this.lMaCoefficients    = coefficients.get( 8);
  156.         this.lJCoefficients     = coefficients.get( 9);
  157.         this.lSaCoefficients    = coefficients.get(10);
  158.         this.lUCoefficients     = coefficients.get(11);
  159.         this.lNeCoefficients    = coefficients.get(12);
  160.         this.paCoefficients     = coefficients.get(13);
  161.     }

  163.     /** Parse coefficients.
  164.      * @param stream stream containing the IERS table
  165.      * @param name name of the resource file (for error messages only)
  166.      * @return list of coefficients arrays
  167.      * @exception OrekitException if stream is null or the table cannot be parsed
  168.      */
  169.     private static List<double[]> parseCoefficients(final InputStream stream, final String name)
  170.         throws OrekitException {

  172.         if (stream == null) {
  173.             throw new OrekitException(OrekitMessages.UNABLE_TO_FIND_FILE, name);
  174.         }

  176.         try {

  178.             final DefinitionParser definitionParser = new DefinitionParser();

  180.             // setup the reader
  181.             final BufferedReader reader = new BufferedReader(new InputStreamReader(stream, "UTF-8"));
  182.             int lineNumber = 0;

  184.             // look for the reference date and the 14 polynomials
  185.             final int n = FundamentalName.values().length;
  186.             final Map<FundamentalName, double[]> polynomials = new HashMap<FundamentalName, double[]>(n);
  187.             for (String line = reader.readLine(); line != null; line = reader.readLine()) {
  188.                 lineNumber++;
  189.                 if (definitionParser.parseDefinition(line, lineNumber, name)) {
  190.                     polynomials.put(definitionParser.getParsedName(),
  191.                                     definitionParser.getParsedPolynomial());
  192.                 }
  193.             }

  195.             final List<double[]> coefficients = new ArrayList<double[]>(n);
  196.             coefficients.add(getCoefficients(FundamentalName.L,       polynomials, name));
  197.             coefficients.add(getCoefficients(FundamentalName.L_PRIME, polynomials, name));
  198.             coefficients.add(getCoefficients(FundamentalName.F,       polynomials, name));
  199.             coefficients.add(getCoefficients(FundamentalName.D,       polynomials, name));
  200.             coefficients.add(getCoefficients(FundamentalName.OMEGA,   polynomials, name));
  201.             if (polynomials.containsKey(FundamentalName.L_ME)) {
  202.                 // IERS conventions 2003 and later provide planetary nutation arguments
  203.                 coefficients.add(getCoefficients(FundamentalName.L_ME,    polynomials, name));
  204.                 coefficients.add(getCoefficients(FundamentalName.L_VE,    polynomials, name));
  205.                 coefficients.add(getCoefficients(FundamentalName.L_E,     polynomials, name));
  206.                 coefficients.add(getCoefficients(FundamentalName.L_MA,    polynomials, name));
  207.                 coefficients.add(getCoefficients(FundamentalName.L_J,     polynomials, name));
  208.                 coefficients.add(getCoefficients(FundamentalName.L_SA,    polynomials, name));
  209.                 coefficients.add(getCoefficients(FundamentalName.L_U,     polynomials, name));
  210.                 coefficients.add(getCoefficients(FundamentalName.L_NE,    polynomials, name));
  211.                 coefficients.add(getCoefficients(FundamentalName.PA,      polynomials, name));
  212.             } else {
  213.                 // IERS conventions 1996 and earlier don't provide planetary nutation arguments
  214.                 final double[] zero = new double[] {
  215.                     0.0
  216.                 };
  217.                 while (coefficients.size() < n) {
  218.                     coefficients.add(zero);
  219.                 }
  220.             }

  222.             return coefficients;

  224.         } catch (IOException ioe) {
  225.             throw new OrekitException(ioe, new DummyLocalizable(ioe.getMessage()));
  226.         }

  228.     }

  230.     /** Get the coefficients for a fundamental argument.
  231.      * @param argument fundamental argument
  232.      * @param polynomials map of the polynomials
  233.      * @param fileName name of the file from which the coefficients have been read
  234.      * @return polynomials coefficients (ordered from high degrees to low degrees)
  235.      * @exception OrekitException if the argument is not found
  236.      */
  237.     private static double[] getCoefficients(final FundamentalName argument,
  238.                                             final Map<FundamentalName, double[]> polynomials,
  239.                                             final String fileName)
  240.         throws OrekitException {
  241.         if (!polynomials.containsKey(argument)) {
  242.             throw new OrekitException(OrekitMessages.NOT_A_SUPPORTED_IERS_DATA_FILE, fileName);
  243.         }
  244.         return polynomials.get(argument);
  245.     }

  247.     /** Evaluate a polynomial.
  248.      * @param tc offset in Julian centuries
  249.      * @param coefficients polynomial coefficients (ordered from low degrees to high degrees)
  250.      * @return value of the polynomial
  251.      */
  252.     private double value(final double tc, final double[] coefficients) {
  253.         double value = 0;
  254.         for (int i = coefficients.length - 1; i >= 0; --i) {
  255.             value = coefficients[i] + tc * value;
  256.         }
  257.         return value;
  258.     }

  260.     /** Evaluate all fundamental arguments for the current date (Delaunay plus planetary).
  261.      * @param date current date
  262.      * @return all fundamental arguments for the current date (Delaunay plus planetary)
  263.      */
  264.     public BodiesElements evaluateAll(final AbsoluteDate date) {

  266.         final double tc = conventions.evaluateTC(date);
  267.         final double gamma = gmstFunction == null ?
  268.                              Double.NaN : gmstFunction.value(date).getValue() + FastMath.PI;

  270.         return new BodiesElements(date, tc, gamma,
  271.                                   value(tc, lCoefficients),      // mean anomaly of the Moon
  272.                                   value(tc, lPrimeCoefficients), // mean anomaly of the Sun
  273.                                   value(tc, fCoefficients),      // L - &Omega; where L is the mean longitude of the Moon
  274.                                   value(tc, dCoefficients),      // mean elongation of the Moon from the Sun
  275.                                   value(tc, omegaCoefficients),  // mean longitude of the ascending node of the Moon
  276.                                   value(tc, lMeCoefficients),    // mean Mercury longitude
  277.                                   value(tc, lVeCoefficients),    // mean Venus longitude
  278.                                   value(tc, lECoefficients),     // mean Earth longitude
  279.                                   value(tc, lMaCoefficients),    // mean Mars longitude
  280.                                   value(tc, lJCoefficients),     // mean Jupiter longitude
  281.                                   value(tc, lSaCoefficients),    // mean Saturn longitude
  282.                                   value(tc, lUCoefficients),     // mean Uranus longitude
  283.                                   value(tc, lNeCoefficients),    // mean Neptune longitude
  284.                                   value(tc, paCoefficients));    // general accumulated precession in longitude

  286.     }

  288.     /** Evaluate a polynomial.
  289.      * @param tc offset in Julian centuries
  290.      * @param coefficients polynomial coefficients (ordered from low degrees to high degrees)
  291.      * @return value of the polynomial
  292.      */
  293.     private DerivativeStructure value(final DerivativeStructure tc, final double[] coefficients) {
  294.         DerivativeStructure value = tc.getField().getZero();
  295.         for (int i = coefficients.length - 1; i >= 0; --i) {
  296.             value = value.multiply(tc).add(coefficients[i]);
  297.         }
  298.         return value;
  299.     }

  301.     /** Evaluate all fundamental arguments for the current date (Delaunay plus planetary),
  302.      * including the first time derivative.
  303.      * @param date current date
  304.      * @return all fundamental arguments for the current date (Delaunay plus planetary),
  305.      * including the first time derivative
  306.      */
  307.     public FieldBodiesElements<DerivativeStructure> evaluateDerivative(final AbsoluteDate date) {

  309.         final DerivativeStructure tc = conventions.dsEvaluateTC(date);

  311.         return new FieldBodiesElements<DerivativeStructure>(date, tc,
  312.                                                             gmstFunction.value(date).add(FastMath.PI),
  313.                                                             value(tc, lCoefficients),      // mean anomaly of the Moon
  314.                                                             value(tc, lPrimeCoefficients), // mean anomaly of the Sun
  315.                                                             value(tc, fCoefficients),      // L - &Omega; where L is the mean longitude of the Moon
  316.                                                             value(tc, dCoefficients),      // mean elongation of the Moon from the Sun
  317.                                                             value(tc, omegaCoefficients),  // mean longitude of the ascending node of the Moon
  318.                                                             value(tc, lMeCoefficients),    // mean Mercury longitude
  319.                                                             value(tc, lVeCoefficients),    // mean Venus longitude
  320.                                                             value(tc, lECoefficients),     // mean Earth longitude
  321.                                                             value(tc, lMaCoefficients),    // mean Mars longitude
  322.                                                             value(tc, lJCoefficients),     // mean Jupiter longitude
  323.                                                             value(tc, lSaCoefficients),    // mean Saturn longitude
  324.                                                             value(tc, lUCoefficients),     // mean Uranus longitude
  325.                                                             value(tc, lNeCoefficients),    // mean Neptune longitude
  326.                                                             value(tc, paCoefficients));    // general accumulated precession in longitude

  328.     }

  330.     /** Replace the instance with a data transfer object for serialization.
  331.      * <p>
  332.      * This intermediate class serializes only the frame key.
  333.      * </p>
  334.      * @return data transfer object that will be serialized
  335.      */
  336.     private Object writeReplace() {
  337.         return new DataTransferObject(conventions, timeScale,
  338.                                       Arrays.asList(lCoefficients, lPrimeCoefficients, fCoefficients,
  339.                                                     dCoefficients, omegaCoefficients,
  340.                                                     lMeCoefficients, lVeCoefficients, lECoefficients,
  341.                                                     lMaCoefficients, lJCoefficients, lSaCoefficients,
  342.                                                     lUCoefficients, lNeCoefficients, paCoefficients));
  343.     }

  345.     /** Internal class used only for serialization. */
  346.     private static class DataTransferObject implements Serializable {

  348.         /** Serializable UID. */
  349.         private static final long serialVersionUID = 20131209L;

  351.         /** IERS conventions to use. */
  352.         private final IERSConventions conventions;

  354.         /** Time scale for GMST computation. */
  355.         private final TimeScale timeScale;

  357.         /** All coefficients. */
  358.         private final List<double[]> coefficients;

  360.         /** Simple constructor.
  361.          * @param conventions IERS conventions to use
  362.          * @param timeScale time scale for GMST computation
  363.          * @param coefficients all coefficients
  364.          */
  365.         public DataTransferObject(final IERSConventions conventions, final TimeScale timeScale,
  366.                                   final List<double[]> coefficients) {
  367.             this.conventions  = conventions;
  368.             this.timeScale    = timeScale;
  369.             this.coefficients = coefficients;
  370.         }

  372.         /** Replace the deserialized data transfer object with a {@link TIRFProvider}.
  373.          * @return replacement {@link TIRFProvider}
  374.          */
  375.         private Object readResolve() {
  376.             try {
  377.                 // retrieve a managed frame
  378.                 return new FundamentalNutationArguments(conventions, timeScale, coefficients);
  379.             } catch (OrekitException oe) {
  380.                 throw OrekitException.createInternalError(oe);
  381.             }
  382.         }

  384.     }

  386.     /** Enumerate for the fundamental names. */
  387.     private enum FundamentalName {

  389.         /** Constant for Mean anomaly of the Moon. */
  390.         L() {
  391.             /** {@inheritDoc} */
  392.             public String getArgumentName() {
  393.                 return "l";
  394.             }
  395.         },

  397.         /** Constant for Mean anomaly of the Sun. */
  398.         L_PRIME() {
  399.             /** {@inheritDoc} */
  400.             public String getArgumentName() {
  401.                 return "l'";
  402.             }
  403.         },

  405.         /** Constant for L - &Omega; where L is the mean longitude of the Moon. */
  406.         F() {
  407.             /** {@inheritDoc} */
  408.             public String getArgumentName() {
  409.                 return "F";
  410.             }
  411.         },

  413.         /** Constant for mean elongation of the Moon from the Sun. */
  414.         D() {
  415.             /** {@inheritDoc} */
  416.             public String getArgumentName() {
  417.                 return "D";
  418.             }
  419.         },

  421.         /** Constant for longitude of the ascending node of the Moon. */
  422.         OMEGA() {
  423.             /** {@inheritDoc} */
  424.             public String getArgumentName() {
  425.                 return "\u03a9";
  426.             }
  427.         },

  429.         /** Constant for mean Mercury longitude. */
  430.         L_ME() {
  431.             /** {@inheritDoc} */
  432.             public String getArgumentName() {
  433.                 return "LMe";
  434.             }
  435.         },

  437.         /** Constant for mean Venus longitude. */
  438.         L_VE() {
  439.             /** {@inheritDoc} */
  440.             public String getArgumentName() {
  441.                 return "LVe";
  442.             }
  443.         },

  445.         /** Constant for mean Earth longitude. */
  446.         L_E() {
  447.             /** {@inheritDoc} */
  448.             public String getArgumentName() {
  449.                 return "LE";
  450.             }
  451.         },

  453.         /** Constant for mean Mars longitude. */
  454.         L_MA() {
  455.             /** {@inheritDoc} */
  456.             public String getArgumentName() {
  457.                 return "LMa";
  458.             }
  459.         },

  461.         /** Constant for mean Jupiter longitude. */
  462.         L_J() {
  463.             /** {@inheritDoc} */
  464.             public String getArgumentName() {
  465.                 return "LJ";
  466.             }
  467.         },

  469.         /** Constant for mean Saturn longitude. */
  470.         L_SA() {
  471.             /** {@inheritDoc} */
  472.             public String getArgumentName() {
  473.                 return "LSa";
  474.             }
  475.         },

  477.         /** Constant for mean Uranus longitude. */
  478.         L_U() {
  479.             /** {@inheritDoc} */
  480.             public String getArgumentName() {
  481.                 return "LU";
  482.             }
  483.         },

  485.         /** Constant for mean Neptune longitude. */
  486.         L_NE() {
  487.             /** {@inheritDoc} */
  488.             public String getArgumentName() {
  489.                 return "LNe";
  490.             }
  491.         },

  493.         /** Constant for general accumulated precession in longitude. */
  494.         PA() {
  495.             /** {@inheritDoc} */
  496.             public String getArgumentName() {
  497.                 return "pA";
  498.             }
  499.         };

  501.         /** Get the fundamental name.
  502.          * @return fundamental name
  503.          */
  504.         public abstract String getArgumentName();

  506.     }

  508.     /** Local parser for argument definition lines. */
  509.     private static class DefinitionParser {

  511.         /** Regular expression pattern for definitions. */
  512.         private final Pattern pattern;

  514.         /** Parser for polynomials. */
  515.         private PolynomialParser polynomialParser;

  517.         /** Last parsed fundamental name. */
  518.         private FundamentalName parsedName;

  520.         /** Last parsed polynomial. */
  521.         private double[] parsedPolynomial;

  523.         /** Simple constructor. */
  524.         public DefinitionParser() {

  526.             // the luni-solar Delaunay arguments polynomial parts should read something like:
  527.             // F5 ≡ Ω = 125.04455501° − 6962890.5431″t + 7.4722″t² + 0.007702″t³ − 0.00005939″t⁴
  528.             // whereas the planetary arguments polynomial parts should read something like:
  529.             // F14 ≡ pA  = 0.02438175 × t + 0.00000538691 × t²
  530.             final String unicodeIdenticalTo = "\u2261";

  532.             // pattern for the global line
  533.             final StringBuilder builder = new StringBuilder();
  534.             for (final FundamentalName fn : FundamentalName.values()) {
  535.                 if (builder.length() > 0) {
  536.                     builder.append('|');
  537.                 }
  538.                 builder.append(fn.getArgumentName());
  539.             }
  540.             final String fundamentalName = "\\p{Space}*((?:" + builder.toString() + ")+)";
  541.             pattern = Pattern.compile("\\p{Space}*F\\p{Digit}+\\p{Space}*" + unicodeIdenticalTo +
  542.                                       fundamentalName + "\\p{Space}*=\\p{Space}*(.*)");

  544.             polynomialParser = new PolynomialParser('t', PolynomialParser.Unit.NO_UNITS);

  546.         }

  548.         /** Parse a definition line.
  549.          * @param line line to parse
  550.          * @param lineNumber line number
  551.          * @param fileName name of the file
  552.          * @return true if a definition has been parsed
  553.          */
  554.         public boolean parseDefinition(final String line, final int lineNumber, final String fileName) {

  556.             parsedName       = null;
  557.             parsedPolynomial = null;

  559.             final Matcher matcher = pattern.matcher(line);
  560.             if (matcher.matches()) {
  561.                 for (FundamentalName fn : FundamentalName.values()) {
  562.                     if (fn.getArgumentName().equals(matcher.group(1))) {
  563.                         parsedName = fn;
  564.                     }
  565.                 }

  567.                 // parse the polynomial
  568.                 parsedPolynomial = polynomialParser.parse(matcher.group(2));

  570.                 return true;

  572.             } else {
  573.                 return false;
  574.             }

  576.         }

  578.         /** Get the last parsed fundamental name.
  579.          * @return last parsed fundamental name
  580.          */
  581.         public FundamentalName getParsedName() {
  582.             return parsedName;
  583.         }

  585.         /** Get the last parsed polynomial.
  586.          * @return last parsed polynomial
  587.          */
  588.         public double[] getParsedPolynomial() {
  589.             return parsedPolynomial.clone();
  590.         }

  592.     }

  594. }
Created with JaCoCo 0.6.3.201306030806