PotentialCoefficientsReader.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.forces.gravity.potential;

  19. import java.io.IOException;
  20. import java.io.InputStream;
  21. import java.text.ParseException;
  22. import java.util.ArrayList;
  23. import java.util.Arrays;
  24. import java.util.List;

  26. import org.apache.commons.math3.util.FastMath;
  27. import org.apache.commons.math3.util.Precision;
  28. import org.orekit.data.DataLoader;
  29. import org.orekit.errors.OrekitException;
  30. import org.orekit.errors.OrekitMessages;

  32. /**This abstract class represents a Gravitational Potential Coefficients file reader.
  33.  *
  34.  * <p> As it exits many different coefficients models and containers this
  35.  *  interface represents all the methods that should be implemented by a reader.
  36.  *  The proper way to use this interface is to call the {@link GravityFieldFactory}
  37.  *  which will determine which reader to use with the selected potential
  38.  *  coefficients file.<p>
  39.  *
  40.  * @see GravityFieldFactory
  41.  * @author Fabien Maussion
  42.  */
  43. public abstract class PotentialCoefficientsReader implements DataLoader {

  45.     /** Maximal degree to parse. */
  46.     private int maxParseDegree;

  48.     /** Maximal order to parse. */
  49.     private int maxParseOrder;

  51.     /** Regular expression for supported files names. */
  52.     private final String supportedNames;

  54.     /** Allow missing coefficients in the input data. */
  55.     private final boolean missingCoefficientsAllowed;

  57.     /** Indicator for complete read. */
  58.     private boolean readComplete;

  60.     /** Central body reference radius. */
  61.     private double ae;

  63.     /** Central body attraction coefficient. */
  64.     private double mu;

  66.     /** Raw tesseral-sectorial coefficients matrix. */
  67.     private double[][] rawC;

  69.     /** Raw tesseral-sectorial coefficients matrix. */
  70.     private double[][] rawS;

  72.     /** Indicator for normalized raw coefficients. */
  73.     private boolean normalized;

  75.     /** Tide system. */
  76.     private TideSystem tideSystem;

  78.     /** Simple constructor.
  79.      * <p>Build an uninitialized reader.</p>
  80.      * @param supportedNames regular expression for supported files names
  81.      * @param missingCoefficientsAllowed allow missing coefficients in the input data
  82.      */
  83.     protected PotentialCoefficientsReader(final String supportedNames,
  84.                                           final boolean missingCoefficientsAllowed) {
  85.         this.supportedNames             = supportedNames;
  86.         this.missingCoefficientsAllowed = missingCoefficientsAllowed;
  87.         this.maxParseDegree             = Integer.MAX_VALUE;
  88.         this.maxParseOrder              = Integer.MAX_VALUE;
  89.         this.readComplete               = false;
  90.         this.ae                         = Double.NaN;
  91.         this.mu                         = Double.NaN;
  92.         this.rawC                       = null;
  93.         this.rawS                       = null;
  94.         this.normalized                 = false;
  95.         this.tideSystem                 = TideSystem.UNKNOWN;
  96.     }

  98.     /** Get the regular expression for supported files names.
  99.      * @return regular expression for supported files names
  100.      */
  101.     public String getSupportedNames() {
  102.         return supportedNames;
  103.     }

  105.     /** Check if missing coefficients are allowed in the input data.
  106.      * @return true if missing coefficients are allowed in the input data
  107.      */
  108.     public boolean missingCoefficientsAllowed() {
  109.         return missingCoefficientsAllowed;
  110.     }

  112.     /** Set the degree limit for the next file parsing.
  113.      * @param maxParseDegree maximal degree to parse (may be safely
  114.      * set to {@link Integer#MAX_VALUE} to parse all available coefficients)
  115.      * @since 6.0
  116.      */
  117.     public void setMaxParseDegree(final int maxParseDegree) {
  118.         this.maxParseDegree = maxParseDegree;
  119.     }

  121.     /** Get the degree limit for the next file parsing.
  122.      * @return degree limit for the next file parsing
  123.      * @since 6.0
  124.      */
  125.     public int getMaxParseDegree() {
  126.         return maxParseDegree;
  127.     }

  129.     /** Set the order limit for the next file parsing.
  130.      * @param maxParseOrder maximal order to parse (may be safely
  131.      * set to {@link Integer#MAX_VALUE} to parse all available coefficients)
  132.      * @since 6.0
  133.      */
  134.     public void setMaxParseOrder(final int maxParseOrder) {
  135.         this.maxParseOrder = maxParseOrder;
  136.     }

  138.     /** Get the order limit for the next file parsing.
  139.      * @return order limit for the next file parsing
  140.      * @since 6.0
  141.      */
  142.     public int getMaxParseOrder() {
  143.         return maxParseOrder;
  144.     }

  146.     /** {@inheritDoc} */
  147.     public boolean stillAcceptsData() {
  148.         return !(readComplete &&
  149.                  getMaxAvailableDegree() >= getMaxParseDegree() &&
  150.                  getMaxAvailableOrder()  >= getMaxParseOrder());
  151.     }

  153.     /** Set the indicator for completed read.
  154.      * @param readComplete if true, a gravity field has been completely read
  155.      */
  156.     protected void setReadComplete(final boolean readComplete) {
  157.         this.readComplete = readComplete;
  158.     }

  160.     /** Set the central body reference radius.
  161.      * @param ae central body reference radius
  162.      */
  163.     protected void setAe(final double ae) {
  164.         this.ae = ae;
  165.     }

  167.     /** Get the central body reference radius.
  168.      * @return central body reference radius
  169.      */
  170.     protected double getAe() {
  171.         return ae;
  172.     }

  174.     /** Set the central body attraction coefficient.
  175.      * @param mu central body attraction coefficient
  176.      */
  177.     protected void setMu(final double mu) {
  178.         this.mu = mu;
  179.     }

  181.     /** Get the central body attraction coefficient.
  182.      * @return central body attraction coefficient
  183.      */
  184.     protected double getMu() {
  185.         return mu;
  186.     }

  188.     /** Set the {@link TideSystem} used in the gravity field.
  189.      * @param tideSystem tide system used in the gravity field
  190.      */
  191.     protected void setTideSystem(final TideSystem tideSystem) {
  192.         this.tideSystem = tideSystem;
  193.     }

  195.     /** Get the {@link TideSystem} used in the gravity field.
  196.      * @return tide system used in the gravity field
  197.      */
  198.     protected TideSystem getTideSystem() {
  199.         return tideSystem;
  200.     }

  202.     /** Set the tesseral-sectorial coefficients matrix.
  203.      * @param rawNormalized if true, raw coefficients are normalized
  204.      * @param c raw tesseral-sectorial coefficients matrix
  205.      * (a reference to the array will be stored)
  206.      * @param s raw tesseral-sectorial coefficients matrix
  207.      * (a reference to the array will be stored)
  208.      * @param name name of the file (or zip entry)
  209.      * @exception OrekitException if a coefficient is missing
  210.      */
  211.     protected void setRawCoefficients(final boolean rawNormalized,
  212.                                       final double[][] c, final double[][] s,
  213.                                       final String name)
  214.         throws OrekitException {

  216.         // normalization indicator
  217.         normalized = rawNormalized;

  219.         // cosine part
  220.         for (int i = 0; i < c.length; ++i) {
  221.             for (int j = 0; j < c[i].length; ++j) {
  222.                 if (Double.isNaN(c[i][j])) {
  223.                     throw new OrekitException(OrekitMessages.MISSING_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
  224.                                               'C', i, j, name);
  225.                 }
  226.             }
  227.         }
  228.         rawC = c;

  230.         // sine part
  231.         for (int i = 0; i < s.length; ++i) {
  232.             for (int j = 0; j < s[i].length; ++j) {
  233.                 if (Double.isNaN(s[i][j])) {
  234.                     throw new OrekitException(OrekitMessages.MISSING_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
  235.                                               'S', i, j, name);
  236.                 }
  237.             }
  238.         }
  239.         rawS = s;

  241.     }

  243.     /** Set the normalized tesseral-sectorial coefficients matrix.
  244.      * @param normalizedS tesseral-sectorial coefficients matrix
  245.      * (a reference to the array will be stored, <em>and</em>
  246.      * the elements will be un-normalized in-place)
  247.      * @param name name of the file (or zip entry)
  248.      * @exception OrekitException if a coefficient is missing
  249.      */
  250.     protected void setNormalizedS(final double[][] normalizedS, final String name)
  251.         throws OrekitException {
  252.         final int degree = normalizedS.length - 1;
  253.         final int order  = normalizedS[degree].length - 1;
  254.         final double[][] factors = GravityFieldFactory.getUnnormalizationFactors(degree, order);
  255.         for (int i = 0; i < normalizedS.length; ++i) {
  256.             for (int j = 0; j < normalizedS[i].length; ++j) {
  257.                 normalizedS[i][j] *= factors[i][j];
  258.             }
  259.         }
  260.         setUnNormalizedS(normalizedS, name);
  261.     }

  263.     /** Set the un-normalized tesseral-sectorial coefficients matrix.
  264.      * @param unNormalizedS un-normalized tesseral-sectorial coefficients matrix
  265.      * (a reference to the array will be stored)
  266.      * @param name name of the file (or zip entry)
  267.      * @exception OrekitException if a coefficient is missing
  268.      */
  269.     protected void setUnNormalizedS(final double[][] unNormalizedS, final String name)
  270.         throws OrekitException {

  272.         for (int i = 0; i < unNormalizedS.length; ++i) {
  273.             for (int j = 0; j < unNormalizedS[i].length; ++j) {
  274.                 if (Double.isNaN(unNormalizedS[i][j])) {
  275.                     throw new OrekitException(OrekitMessages.MISSING_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
  276.                                               'S', i, j, name);
  277.                 }
  278.             }
  279.         }

  281.         this.rawS = unNormalizedS;

  283.     }

  285.     /** Get the maximal degree available in the last file parsed.
  286.      * @return maximal degree available in the last file parsed
  287.      * @since 6.0
  288.      */
  289.     public int getMaxAvailableDegree() {
  290.         return rawC.length - 1;
  291.     }

  293.     /** Get the maximal order available in the last file parsed.
  294.      * @return maximal order available in the last file parsed
  295.      * @since 6.0
  296.      */
  297.     public int getMaxAvailableOrder() {
  298.         return rawC[rawC.length - 1].length - 1;
  299.     }

  301.     /** {@inheritDoc} */
  302.     public abstract void loadData(InputStream input, String name)
  303.         throws IOException, ParseException, OrekitException;

  305.     /** Get a provider for read spherical harmonics coefficients.
  306.      * @param wantNormalized if true, the provider will provide normalized coefficients,
  307.      * otherwise it will provide un-normalized coefficients
  308.      * @param degree maximal degree
  309.      * @param order maximal order
  310.      * @return a new provider
  311.      * @exception OrekitException if the requested maximal degree or order exceeds the
  312.      * available degree or order or if no gravity field has read yet
  313.      * @see #getConstantProvider(boolean, int, int)
  314.      * @since 6.0
  315.      */
  316.     public abstract RawSphericalHarmonicsProvider getProvider(boolean wantNormalized, int degree, int order)
  317.         throws OrekitException;

  319.     /** Get a time-independent provider for read spherical harmonics coefficients.
  320.      * @param wantNormalized if true, the raw provider must provide normalized coefficients,
  321.      * otherwise it will provide un-normalized coefficients
  322.      * @param degree maximal degree
  323.      * @param order maximal order
  324.      * @return a new provider, with no time-dependent parts
  325.      * @exception OrekitException if the requested maximal degree or order exceeds the
  326.      * available degree or order or if no gravity field has read yet
  327.      * @see #getProvider(boolean, int, int)
  328.      * @since 6.0
  329.      */
  330.     protected ConstantSphericalHarmonics getConstantProvider(final boolean wantNormalized,
  331.                                                              final int degree, final int order)
  332.         throws OrekitException {

  334.         if (!readComplete) {
  335.             throw new OrekitException(OrekitMessages.NO_GRAVITY_FIELD_DATA_LOADED);
  336.         }

  338.         if (degree >= rawC.length) {
  339.             throw new OrekitException(OrekitMessages.TOO_LARGE_DEGREE_FOR_GRAVITY_FIELD,
  340.                                       degree, rawC.length - 1);
  341.         }

  343.         if (order >= rawC[rawC.length - 1].length) {
  344.             throw new OrekitException(OrekitMessages.TOO_LARGE_ORDER_FOR_GRAVITY_FIELD,
  345.                                       order, rawC[rawC.length - 1].length);
  346.         }

  348.         // fix normalization
  349.         final double[][] truncatedC = buildTriangularArray(degree, order, 0.0);
  350.         final double[][] truncatedS = buildTriangularArray(degree, order, 0.0);
  351.         rescale(1.0, normalized, rawC, rawS, wantNormalized, truncatedC, truncatedS);

  353.         return new ConstantSphericalHarmonics(ae, mu, tideSystem, truncatedC, truncatedS);

  355.     }

  357.     /** Build a coefficients triangular array.
  358.      * @param degree array degree
  359.      * @param order array order
  360.      * @param value initial value to put in array elements
  361.      * @return built array
  362.      */
  363.     protected static double[][] buildTriangularArray(final int degree, final int order, final double value) {
  364.         final double[][] array = new double[degree + 1][];
  365.         for (int k = 0; k < array.length; ++k) {
  366.             array[k] = buildRow(k, order, value);
  367.         }
  368.         return array;
  369.     }

  371.     /** Build a coefficients row.
  372.      * @param degree row degree
  373.      * @param order row order
  374.      * @param value initial value to put in array elements
  375.      * @return built row
  376.      */
  377.     protected static double[] buildRow(final int degree, final int order, final double value) {
  378.         final double[] row = new double[FastMath.min(order, degree) + 1];
  379.         Arrays.fill(row, value);
  380.         return row;
  381.     }

  383.     /** Extend a list of lists of coefficients if needed.
  384.      * @param list list of lists of coefficients
  385.      * @param degree degree required to be present
  386.      * @param order order required to be present
  387.      * @param value initial value to put in list elements
  388.      */
  389.     protected void extendListOfLists(final List<List<Double>> list, final int degree, final int order,
  390.                                      final double value) {
  391.         for (int i = list.size(); i <= degree; ++i) {
  392.             list.add(new ArrayList<Double>());
  393.         }
  394.         final List<Double> listN = list.get(degree);
  395.         final Double v = Double.valueOf(value);
  396.         for (int j = listN.size(); j <= order; ++j) {
  397.             listN.add(v);
  398.         }
  399.     }

  401.     /** Convert a list of list into an array.
  402.      * @param list list of lists of coefficients
  403.      * @return a new array
  404.      */
  405.     protected double[][] toArray(final List<List<Double>> list) {
  406.         final double[][] array = new double[list.size()][];
  407.         for (int i = 0; i < array.length; ++i) {
  408.             array[i] = new double[list.get(i).size()];
  409.             for (int j = 0; j < array[i].length; ++j) {
  410.                 array[i][j] = list.get(i).get(j);
  411.             }
  412.         }
  413.         return array;
  414.     }

  416.     /** Parse a coefficient.
  417.      * @param field text field to parse
  418.      * @param list list where to put the coefficient
  419.      * @param i first index in the list
  420.      * @param j second index in the list
  421.      * @param cName name of the coefficient
  422.      * @param name name of the file
  423.      * @exception OrekitException if the coefficient is already set
  424.      */
  425.     protected void parseCoefficient(final String field, final List<List<Double>> list,
  426.                                     final int i, final int j,
  427.                                     final String cName, final String name)
  428.         throws OrekitException {
  429.         final double value    = Double.parseDouble(field.replace('D', 'E'));
  430.         final double oldValue = list.get(i).get(j);
  431.         if (Double.isNaN(oldValue) || Precision.equals(oldValue, 0.0, 1)) {
  432.             // the coefficient was not already initialized
  433.             list.get(i).set(j, value);
  434.         } else {
  435.             throw new OrekitException(OrekitMessages.DUPLICATED_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
  436.                                       name, i, j, name);
  437.         }
  438.     }

  440.     /** Parse a coefficient.
  441.      * @param field text field to parse
  442.      * @param array array where to put the coefficient
  443.      * @param i first index in the list
  444.      * @param j second index in the list
  445.      * @param cName name of the coefficient
  446.      * @param name name of the file
  447.      * @exception OrekitException if the coefficient is already set
  448.      */
  449.     protected void parseCoefficient(final String field, final double[][] array,
  450.                                     final int i, final int j,
  451.                                     final String cName, final String name)
  452.         throws OrekitException {
  453.         final double value    = Double.parseDouble(field.replace('D', 'E'));
  454.         final double oldValue = array[i][j];
  455.         if (Double.isNaN(oldValue) || Precision.equals(oldValue, 0.0, 1)) {
  456.             // the coefficient was not already initialized
  457.             array[i][j] = value;
  458.         } else {
  459.             throw new OrekitException(OrekitMessages.DUPLICATED_GRAVITY_FIELD_COEFFICIENT_IN_FILE,
  460.                                       name, i, j, name);
  461.         }
  462.     }

  464.     /** Rescale coefficients arrays.
  465.      * @param scale general scaling factor to apply to all elements
  466.      * @param normalizedOrigin if true, the origin coefficients are normalized
  467.      * @param originC cosine part of the origina coefficients
  468.      * @param originS sine part of the origin coefficients
  469.      * @param wantNormalized if true, the rescaled coefficients must be normalized
  470.      * @param rescaledC cosine part of the rescaled coefficients to fill in (may be the originC array)
  471.      * @param rescaledS sine part of the rescaled coefficients to fill in (may be the originS array)
  472.      * @exception OrekitException if normalization/unnormalization fails because of an underflow
  473.      * due to too high degree/order
  474.      */
  475.     protected static void rescale(final double scale,
  476.                                   final boolean normalizedOrigin, final double[][] originC,
  477.                                   final double[][] originS, final boolean wantNormalized,
  478.                                   final double[][] rescaledC, final double[][] rescaledS)
  479.         throws OrekitException {

  481.         if (wantNormalized == normalizedOrigin) {
  482.             // apply only the general scaling factor
  483.             for (int i = 0; i < rescaledC.length; ++i) {
  484.                 final double[] rCi = rescaledC[i];
  485.                 final double[] rSi = rescaledS[i];
  486.                 final double[] oCi = originC[i];
  487.                 final double[] oSi = originS[i];
  488.                 for (int j = 0; j < rCi.length; ++j) {
  489.                     rCi[j] = oCi[j] * scale;
  490.                     rSi[j] = oSi[j] * scale;
  491.                 }
  492.             }
  493.         } else {

  495.             // we have to re-scale the coefficients
  496.             // (we use rescaledC.length - 1 for the order instead of rescaledC[rescaledC.length - 1].length - 1
  497.             //  because typically trend or pulsation arrays are irregular, some test cases have
  498.             //  order 2 elements at degree 2, but only order 1 elements for higher degrees for example)
  499.             final double[][] factors = GravityFieldFactory.getUnnormalizationFactors(rescaledC.length - 1,
  500.                                                                                      rescaledC.length - 1);

  502.             if (wantNormalized) {
  503.                 // normalize the coefficients
  504.                 for (int i = 0; i < rescaledC.length; ++i) {
  505.                     final double[] rCi = rescaledC[i];
  506.                     final double[] rSi = rescaledS[i];
  507.                     final double[] oCi = originC[i];
  508.                     final double[] oSi = originS[i];
  509.                     final double[] fi  = factors[i];
  510.                     for (int j = 0; j < rCi.length; ++j) {
  511.                         final double factor = scale / fi[j];
  512.                         rCi[j] = oCi[j] * factor;
  513.                         rSi[j] = oSi[j] * factor;
  514.                     }
  515.                 }
  516.             } else {
  517.                 // un-normalize the coefficients
  518.                 for (int i = 0; i < rescaledC.length; ++i) {
  519.                     final double[] rCi = rescaledC[i];
  520.                     final double[] rSi = rescaledS[i];
  521.                     final double[] oCi = originC[i];
  522.                     final double[] oSi = originS[i];
  523.                     final double[] fi  = factors[i];
  524.                     for (int j = 0; j < rCi.length; ++j) {
  525.                         final double factor = scale * fi[j];
  526.                         rCi[j] = oCi[j] * factor;
  527.                         rSi[j] = oSi[j] * factor;
  528.                     }
  529.                 }
  530.             }

  532.         }
  533.     }

  535. }
Created with JaCoCo 0.6.3.201306030806