GlobalIonosphereMapModel.java

/* Copyright 2002-2022 CS GROUP
 * Licensed to CS GROUP (CS) under one or more
 * contributor license agreements.  See the NOTICE file distributed with
 * this work for additional information regarding copyright ownership.
 * CS licenses this file to You under the Apache License, Version 2.0
 * (the "License"); you may not use this file except in compliance with
 * the License.  You may obtain a copy of the License at
 *
 *   http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
package org.orekit.models.earth.ionosphere;

import java.io.BufferedReader;
import java.io.IOException;
import java.io.InputStream;
import java.io.InputStreamReader;
import java.io.Serializable;
import java.nio.charset.StandardCharsets;
import java.text.ParseException;
import java.util.ArrayList;
import java.util.Collections;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
import java.util.regex.Pattern;

import org.hipparchus.Field;
import org.hipparchus.CalculusFieldElement;
import org.hipparchus.analysis.interpolation.BilinearInterpolatingFunction;
import org.hipparchus.geometry.euclidean.threed.FieldVector3D;
import org.hipparchus.geometry.euclidean.threed.Vector3D;
import org.hipparchus.util.FastMath;
import org.hipparchus.util.MathUtils;
import org.orekit.annotation.DefaultDataContext;
import org.orekit.bodies.GeodeticPoint;
import org.orekit.data.AbstractSelfFeedingLoader;
import org.orekit.data.DataContext;
import org.orekit.data.DataLoader;
import org.orekit.data.DataProvidersManager;
import org.orekit.errors.OrekitException;
import org.orekit.errors.OrekitInternalError;
import org.orekit.errors.OrekitMessages;
import org.orekit.frames.TopocentricFrame;
import org.orekit.propagation.FieldSpacecraftState;
import org.orekit.propagation.SpacecraftState;
import org.orekit.time.AbsoluteDate;
import org.orekit.time.DateTimeComponents;
import org.orekit.time.FieldAbsoluteDate;
import org.orekit.time.TimeComponents;
import org.orekit.time.TimeScale;
import org.orekit.utils.ParameterDriver;

/**
 * Global Ionosphere Map (GIM) model.
 * The ionospheric delay is computed according to the formulas:
 * <pre>
 *           40.3
 *    δ =  --------  *  STEC      with, STEC = VTEC * F(elevation)
 *            f²
 * </pre>
 * With:
 * <ul>
 * <li>f: The frequency of the signal in Hz.</li>
 * <li>STEC: The Slant Total Electron Content in TECUnits.</li>
 * <li>VTEC: The Vertical Total Electron Content in TECUnits.</li>
 * <li>F(elevation): A mapping function which depends on satellite elevation.</li>
 * </ul>
 * The VTEC is read from a IONEX file. A stream contains, for a given day, the values of the TEC for each hour of the day.
 * Values are given on a global 2.5° x 5.0° (latitude x longitude) grid.
 * <p>
 * A bilinear interpolation is performed the case of the user initialize the latitude and the
 * longitude with values that are not contained in the stream.
 * </p><p>
 * A temporal interpolation is also performed to compute the VTEC at the desired date.
 * </p><p>
 * IONEX files are obtained from
 * <a href="ftp://cddis.nasa.gov/gnss/products/ionex/"> The Crustal Dynamics Data Information System</a>.
 * </p><p>
 * The files have to be extracted to UTF-8 text files before being read by this loader.
 * </p><p>
 * Example of file:
 * </p>
 * <pre>
 *      1.0            IONOSPHERE MAPS     GPS                 IONEX VERSION / TYPE
 * BIMINX V5.3         AIUB                16-JAN-19 07:26     PGM / RUN BY / DATE
 * BROADCAST IONOSPHERE MODEL FOR DAY 015, 2019                COMMENT
 *   2019     1    15     0     0     0                        EPOCH OF FIRST MAP
 *   2019     1    16     0     0     0                        EPOCH OF LAST MAP
 *   3600                                                      INTERVAL
 *     25                                                      # OF MAPS IN FILE
 *   NONE                                                      MAPPING FUNCTION
 *      0.0                                                    ELEVATION CUTOFF
 *                                                             OBSERVABLES USED
 *   6371.0                                                    BASE RADIUS
 *      2                                                      MAP DIMENSION
 *    350.0 350.0   0.0                                        HGT1 / HGT2 / DHGT
 *     87.5 -87.5  -2.5                                        LAT1 / LAT2 / DLAT
 *   -180.0 180.0   5.0                                        LON1 / LON2 / DLON
 *     -1                                                      EXPONENT
 * TEC/RMS values in 0.1 TECU; 9999, if no value available     COMMENT
 *                                                             END OF HEADER
 *      1                                                      START OF TEC MAP
 *   2019     1    15     0     0     0                        EPOCH OF CURRENT MAP
 *     87.5-180.0 180.0   5.0 350.0                            LAT/LON1/LON2/DLON/H
 *    92   92   92   92   92   92   92   92   92   92   92   92   92   92   92   92
 *    92   92   92   92   92   92   92   92   92   92   92   92   92   92   92   92
 *    92   92   92   92   92   92   92   92   92   92   92   92   92   92   92   92
 *    92   92   92   92   92   92   92   92   92   92   92   92   92   92   92   92
 *    92   92   92   92   92   92   92   92   92
 *    ...
 * </pre>
 *
 * @see "Schaer, S., W. Gurtner, and J. Feltens, 1998, IONEX: The IONosphere Map EXchange
 *       Format Version 1, February 25, 1998, Proceedings of the IGS AC Workshop
 *       Darmstadt, Germany, February 9–11, 1998"
 *
 * @author Bryan Cazabonne
 *
 */
public class GlobalIonosphereMapModel extends AbstractSelfFeedingLoader
        implements IonosphericModel {

    /** Serializable UID. */
    private static final long serialVersionUID = 201928052L;

    /** Pattern for delimiting regular expressions. */
    private static final Pattern SEPARATOR = Pattern.compile("\\s+");

    /** Threshold for latitude and longitude difference. */
    private static final double THRESHOLD = 0.001;

    /** Geodetic site latitude, radians.*/
    private double latitude;

    /** Geodetic site longitude, radians.*/
    private double longitude;

    /** Mean earth radius [m]. */
    private double r0;

    /** Height of the ionospheric single layer [m]. */
    private double h;

    /** Time interval between two TEC maps [s]. */
    private double dt;

    /** Number of TEC maps as read on the header of the file. */
    private int nbMaps;

    /** Flag for mapping function computation. */
    private boolean mapping;

    /** Epoch of the first TEC map as read in the header of the IONEX file. */
    private AbsoluteDate startDate;

    /** Epoch of the last TEC map as read in the header of the IONEX file. */
    private AbsoluteDate endDate;

    /** Map of interpolated TEC at a specific date. */
    private Map<AbsoluteDate, Double> tecMap;

    /** UTC time scale. */
    private final TimeScale utc;

    /**
     * Constructor with supported names given by user. This constructor uses the {@link
     * DataContext#getDefault() default data context}.
     *
     * @param supportedNames regular expression that matches the names of the IONEX files
     *                       to be loaded. See {@link DataProvidersManager#feed(String,
     *                       DataLoader)}.
     * @see #GlobalIonosphereMapModel(String, DataProvidersManager, TimeScale)
     */
    @DefaultDataContext
    public GlobalIonosphereMapModel(final String supportedNames) {
        this(supportedNames,
                DataContext.getDefault().getDataProvidersManager(),
                DataContext.getDefault().getTimeScales().getUTC());
    }

    /**
     * Constructor that uses user defined supported names and data context.
     *
     * @param supportedNames       regular expression that matches the names of the IONEX
     *                             files to be loaded. See {@link DataProvidersManager#feed(String,
     *                             DataLoader)}.
     * @param dataProvidersManager provides access to auxiliary data files.
     * @param utc                  UTC time scale.
     * @since 10.1
     */
    public GlobalIonosphereMapModel(final String supportedNames,
                                    final DataProvidersManager dataProvidersManager,
                                    final TimeScale utc) {
        super(supportedNames, dataProvidersManager);
        this.latitude       = Double.NaN;
        this.longitude      = Double.NaN;
        this.tecMap         = new HashMap<>();
        this.utc = utc;
    }

    /**
     * Calculates the ionospheric path delay for the signal path from a ground
     * station to a satellite.
     * <p>
     * The path delay can be computed for any elevation angle.
     * </p>
     * @param date current date
     * @param geo geodetic point of receiver/station
     * @param elevation elevation of the satellite in radians
     * @param frequency frequency of the signal in Hz
     * @return the path delay due to the ionosphere in m
     */
    public double pathDelay(final AbsoluteDate date, final GeodeticPoint geo,
                            final double elevation, final double frequency) {
        // TEC in TECUnits
        final double tec = getTEC(date, geo);
        // Square of the frequency
        final double freq2 = frequency * frequency;
        // "Slant" Total Electron Content
        final double stec;
        // Check if a mapping factor is needed
        if (mapping) {
            stec = tec;
        } else {
            // Mapping factor
            final double fz = mappingFunction(elevation);
            stec = tec * fz;
        }
        // Delay computation
        final double alpha  = 40.3e16 / freq2;
        return alpha * stec;
    }

    @Override
    public double pathDelay(final SpacecraftState state, final TopocentricFrame baseFrame,
                            final double frequency, final double[] parameters) {

        // Elevation in radians
        final Vector3D position  = state.getPVCoordinates(baseFrame).getPosition();
        final double   elevation = position.getDelta();

        // Only consider measures above the horizon
        if (elevation > 0.0) {
            // Date
            final AbsoluteDate date = state.getDate();
            // Geodetic point
            final GeodeticPoint geo = baseFrame.getPoint();
            // Delay
            return pathDelay(date, geo, elevation, frequency);
        }

        return 0.0;

    }

    /**
     * Calculates the ionospheric path delay for the signal path from a ground
     * station to a satellite.
     * <p>
     * The path delay can be computed for any elevation angle.
     * </p>
     * @param <T> type of the elements
     * @param date current date
     * @param geo geodetic point of receiver/station
     * @param elevation elevation of the satellite in radians
     * @param frequency frequency of the signal in Hz
     * @return the path delay due to the ionosphere in m
     */
    public <T extends CalculusFieldElement<T>> T pathDelay(final FieldAbsoluteDate<T> date, final GeodeticPoint geo,
                                                       final T elevation, final double frequency) {
        // TEC in TECUnits
        final T tec = getTEC(date, geo);
        // Square of the frequency
        final double freq2 = frequency * frequency;
        // "Slant" Total Electron Content
        final T stec;
        // Check if a mapping factor is needed
        if (mapping) {
            stec = tec;
        } else {
            // Mapping factor
            final T fz = mappingFunction(elevation);
            stec = tec.multiply(fz);
        }
        // Delay computation
        final double alpha  = 40.3e16 / freq2;
        return stec.multiply(alpha);
    }

    @Override
    public <T extends CalculusFieldElement<T>> T pathDelay(final FieldSpacecraftState<T> state, final TopocentricFrame baseFrame,
                                                       final double frequency, final T[] parameters) {

        // Elevation in radians
        final FieldVector3D<T> position = state.getPVCoordinates(baseFrame).getPosition();
        final T elevation = position.getDelta();

        // Only consider measures above the horizon
        if (elevation.getReal() > 0.0) {
            // Date
            final FieldAbsoluteDate<T> date = state.getDate();
            // Geodetic point
            final GeodeticPoint geo = baseFrame.getPoint();
            // Delay
            return pathDelay(date, geo, elevation, frequency);
        }

        return elevation.getField().getZero();

    }

    /**
     * Computes the Total Electron Content (TEC) at a given date by performing a
     * temporal interpolation with the two closest date in the IONEX file.
     * @param date current date
     * @param recPoint geodetic point of receiver/station
     * @return the TEC after a temporal interpolation, in TECUnits
     */
    public double getTEC(final AbsoluteDate date, final GeodeticPoint recPoint) {

        // Load TEC data only if needed
        loadsIfNeeded(recPoint);

        // Check if the date is out of range
        checkDate(date);

        // Date and Time components
        final DateTimeComponents dateTime = date.getComponents(utc);
        // Find the two closest dates of the current date
        final double secInDay   = dateTime.getTime().getSecondsInLocalDay();
        final double ratio      = FastMath.floor(secInDay / dt) * dt;
        final AbsoluteDate tI   = new AbsoluteDate(dateTime.getDate(),
                                                   new TimeComponents(ratio),
                                                   utc);
        final AbsoluteDate tIp1 = tI.shiftedBy(dt);

        // Get the TEC values at the two closest dates
        final double tecI   = tecMap.get(tI);
        final double tecIp1 = tecMap.get(tIp1);

        // Perform temporal interpolation (Ref, Eq. 2)
        final double tec = (tIp1.durationFrom(date) / dt) * tecI + (date.durationFrom(tI) / dt) * tecIp1;
        return tec;
    }

    /**
     * Computes the Total Electron Content (TEC) at a given date by performing a
     * temporal interpolation with the two closest date in the IONEX file.
     * @param <T> type of the elements
     * @param date current date
     * @param recPoint geodetic point of receiver/station
     * @return the TEC after a temporal interpolation, in TECUnits
     */
    public <T extends CalculusFieldElement<T>> T getTEC(final FieldAbsoluteDate<T> date, final GeodeticPoint recPoint) {

        // Load TEC data only if needed
        loadsIfNeeded(recPoint);

        // Check if the date is out of range
        checkDate(date.toAbsoluteDate());

        // Field
        final Field<T> field = date.getField();

        // Date and Time components
        final DateTimeComponents dateTime = date.getComponents(utc);
        // Find the two closest dates of the current date
        final double secInDay           = dateTime.getTime().getSecondsInLocalDay();
        final double ratio              = FastMath.floor(secInDay / dt) * dt;
        final FieldAbsoluteDate<T> tI   = new FieldAbsoluteDate<>(field, dateTime.getDate(),
                                                                  new TimeComponents(ratio),
                                                                  utc);
        final FieldAbsoluteDate<T> tIp1 = tI.shiftedBy(dt);

        // Get the TEC values at the two closest dates
        final double tecI   = tecMap.get(tI.toAbsoluteDate());
        final double tecIp1 = tecMap.get(tIp1.toAbsoluteDate());

        // Perform temporal interpolation (Ref, Eq. 2)
        final T tec = tIp1.durationFrom(date).divide(dt).multiply(tecI).add(date.durationFrom(tI).divide(dt).multiply(tecIp1));
        return tec;
    }

    @Override
    public List<ParameterDriver> getParametersDrivers() {
        return Collections.emptyList();
    }

    /**
     * Computes the ionospheric mapping function.
     * @param elevation the elevation of the satellite in radians
     * @return the mapping function
     */
    private double mappingFunction(final double elevation) {
        // Calculate the zenith angle from the elevation
        final double z = FastMath.abs(0.5 * FastMath.PI - elevation);
        // Distance ratio
        final double ratio = r0 / (r0 + h);
        // Mapping function
        final double coef = FastMath.sin(z) * ratio;
        final double fz = 1.0 / FastMath.sqrt(1.0 - coef * coef);
        return fz;
    }

    /**
     * Computes the ionospheric mapping function.
     * @param <T> type of the elements
     * @param elevation the elevation of the satellite in radians
     * @return the mapping function
     */
    private <T extends CalculusFieldElement<T>> T mappingFunction(final T elevation) {
        // Calculate the zenith angle from the elevation
        final T z = FastMath.abs(elevation.getPi().multiply(0.5).subtract(elevation));
        // Distance ratio
        final double ratio = r0 / (r0 + h);
        // Mapping function
        final T coef = FastMath.sin(z).multiply(ratio);
        final T fz = FastMath.sqrt(coef.multiply(coef).negate().add(1.0)).reciprocal();
        return fz;
    }

    /**
     * Lazy loading of TEC data.
     * @param recPoint geodetic point of receiver/station
     */
    private void loadsIfNeeded(final GeodeticPoint recPoint) {

        // Current latitude and longitude of the geodetic point
        final double lat = recPoint.getLatitude();
        final double lon = MathUtils.normalizeAngle(recPoint.getLongitude(), 0.0);

        // Read the file only if the TEC map is empty or if the geodetic point displacement is
        // greater than 0.001 radians (in latitude or longitude)
        if (tecMap.isEmpty() || FastMath.abs(lat - latitude) > THRESHOLD ||  FastMath.abs(lon - longitude) > THRESHOLD) {
            this.latitude  = lat;
            this.longitude = lon;

            // Read file
            final Parser parser = new Parser();
            feed(parser);

            // File header
            final IONEXHeader top = parser.getIONEXHeader();
            this.startDate  = top.getFirstDate();
            this.endDate    = top.getLastDate();
            this.dt         = top.getInterval();
            this.nbMaps     = top.getTECMapsNumer();
            this.r0         = top.getEarthRadius();
            this.h          = top.getHIon();
            this.mapping    = top.isMappingFunction();

            // TEC map
            for (TECMap map : parser.getTECMaps()) {
                tecMap.put(map.getDate(), map.getTEC());
            }
        }
        checkSize();
    }

    /**
     * Check if the current date is between the startDate and
     * the endDate of the IONEX file.
     * @param date current date
     */
    private void checkDate(final AbsoluteDate date) {
        if (startDate.durationFrom(date) > 0 || date.durationFrom(endDate) > 0) {
            throw new OrekitException(OrekitMessages.NO_TEC_DATA_IN_FILE_FOR_DATE,
                    getSupportedNames(), date);
        }
    }

    /**
     * Check if the number of parsed TEC maps is consistent with the header specification.
     */
    private void checkSize() {
        if (tecMap.size() != nbMaps) {
            throw new OrekitException(OrekitMessages.INCONSISTENT_NUMBER_OF_TEC_MAPS_IN_FILE, tecMap.size(), nbMaps);
        }
    }

    /** Replace the instance with a data transfer object for serialization.
     * @return data transfer object that will be serialized
     */
    @DefaultDataContext
    private Object writeReplace() {
        return new DataTransferObject(getSupportedNames());
    }

    /** Parser for IONEX files. */
    private class Parser implements DataLoader {

        /** String for the end of a TEC map. */
        private static final String END = "END OF TEC MAP";

        /** String for the epoch of a TEC map. */
        private static final String EPOCH = "EPOCH OF CURRENT MAP";

        /** Index of label in data lines. */
        private static final int LABEL_START = 60;

        /** Kilometers to meters conversion factor. */
        private static final double KM_TO_M = 1000.0;

        /** Header of the IONEX file. */
        private IONEXHeader header;

        /** List of TEC Maps. */
        private List<TECMap> maps;

        @Override
        public boolean stillAcceptsData() {
            return true;
        }

        @Override
        public void loadData(final InputStream input, final String name)
            throws IOException,  ParseException {

            maps = new ArrayList<>();

            // Open stream and parse data
            int   lineNumber = 0;
            String line      = null;
            try (InputStreamReader isr = new InputStreamReader(input, StandardCharsets.UTF_8);
                 BufferedReader    br = new BufferedReader(isr)) {

                // Placeholders for parsed data
                int               interval    = 3600;
                int               nbOfMaps    = 1;
                int               exponent    = -1;
                double            baseRadius  = 6371.0e3;
                double            hIon        = 350e3;
                boolean           mappingF    = false;
                boolean           inTEC       = false;
                double[]          latitudes   = null;
                double[]          longitudes  = null;
                AbsoluteDate      firstEpoch  = null;
                AbsoluteDate      lastEpoch   = null;
                AbsoluteDate      epoch       = firstEpoch;
                ArrayList<Double> values      = new ArrayList<>();

                for (line = br.readLine(); line != null; line = br.readLine()) {
                    ++lineNumber;
                    if (line.length() > LABEL_START) {
                        switch(line.substring(LABEL_START).trim()) {
                            case "EPOCH OF FIRST MAP" :
                                firstEpoch = parseDate(line);
                                break;
                            case "EPOCH OF LAST MAP" :
                                lastEpoch = parseDate(line);
                                break;
                            case "INTERVAL" :
                                interval = parseInt(line, 2, 4);
                                break;
                            case "# OF MAPS IN FILE" :
                                nbOfMaps = parseInt(line, 2, 4);
                                break;
                            case "BASE RADIUS" :
                                // Value is in kilometers
                                baseRadius = parseDouble(line, 2, 6) * KM_TO_M;
                                break;
                            case "MAPPING FUNCTION" :
                                mappingF = !parseString(line, 2, 4).equals("NONE");
                                break;
                            case "EXPONENT" :
                                exponent = parseInt(line, 4, 2);
                                break;
                            case "HGT1 / HGT2 / DHGT" :
                                if (parseDouble(line, 17, 3) == 0.0) {
                                    // Value is in kilometers
                                    hIon = parseDouble(line, 3, 5) * KM_TO_M;
                                }
                                break;
                            case "LAT1 / LAT2 / DLAT" :
                                latitudes = parseCoordinate(line);
                                break;
                            case "LON1 / LON2 / DLON" :
                                longitudes = parseCoordinate(line);
                                break;
                            case "END OF HEADER" :
                                // Check that latitude and longitude bondaries were found
                                if (latitudes == null || longitudes == null) {
                                    throw new OrekitException(OrekitMessages.NO_LATITUDE_LONGITUDE_BONDARIES_IN_IONEX_HEADER, getSupportedNames());
                                }
                                // Check that first and last epochs were found
                                if (firstEpoch == null || lastEpoch == null) {
                                    throw new OrekitException(OrekitMessages.NO_EPOCH_IN_IONEX_HEADER, getSupportedNames());
                                }
                                // At the end of the header, we build the IONEXHeader object
                                header = new IONEXHeader(firstEpoch, lastEpoch, interval, nbOfMaps,
                                                         baseRadius, hIon, mappingF);
                                break;
                            case "START OF TEC MAP" :
                                inTEC = true;
                                break;
                            case END :
                                final double tec = interpolateTEC(values, exponent, latitudes, longitudes);
                                final TECMap map = new TECMap(epoch, tec);
                                maps.add(map);
                                // Reset parameters
                                inTEC  = false;
                                values = new ArrayList<>();
                                epoch  = null;
                                break;
                            default :
                                if (inTEC) {
                                    // Date
                                    if (line.endsWith(EPOCH)) {
                                        epoch = parseDate(line);
                                    }
                                    // Fill TEC values list
                                    if (!line.endsWith("LAT/LON1/LON2/DLON/H") &&
                                        !line.endsWith(END) &&
                                        !line.endsWith(EPOCH)) {
                                        line = line.trim();
                                        final String[] readLine = SEPARATOR.split(line);
                                        for (final String s : readLine) {
                                            values.add(Double.valueOf(s));
                                        }
                                    }
                                }
                                break;
                        }
                    } else {
                        if (inTEC) {
                            // Here, we are parsing the last line of TEC data for a given latitude
                            // The size of this line is lower than 60.
                            line = line.trim();
                            final String[] readLine = SEPARATOR.split(line);
                            for (final String s : readLine) {
                                values.add(Double.valueOf(s));
                            }
                        }
                    }

                }

                // Close the stream after reading
                input.close();

            } catch (NumberFormatException nfe) {
                throw new OrekitException(OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
                                          lineNumber, name, line);
            }

        }

        /**
         * Get the header of the IONEX file.
         * @return the header of the IONEX file
         */
        public IONEXHeader getIONEXHeader() {
            return header;
        }

        /**
         * Get the list of the TEC maps.
         * @return the list of TEC maps.
         */
        public List<TECMap> getTECMaps() {
            return maps;
        }

        /** Extract a string from a line.
         * @param line to parse
         * @param start start index of the string
         * @param length length of the string
         * @return parsed string
         */
        private String parseString(final String line, final int start, final int length) {
            return line.substring(start, FastMath.min(line.length(), start + length)).trim();
        }

        /** Extract an integer from a line.
         * @param line to parse
         * @param start start index of the integer
         * @param length length of the integer
         * @return parsed integer
         */
        private int parseInt(final String line, final int start, final int length) {
            return Integer.parseInt(parseString(line, start, length));
        }

        /** Extract a double from a line.
         * @param line to parse
         * @param start start index of the real
         * @param length length of the real
         * @return parsed real
         */
        private double parseDouble(final String line, final int start, final int length) {
            return Double.parseDouble(parseString(line, start, length));
        }

        /** Extract a date from a parsed line.
         * @param line to parse
         * @return an absolute date
         */
        private AbsoluteDate parseDate(final String line) {
            return new AbsoluteDate(parseInt(line, 0, 6),
                                    parseInt(line, 6, 6),
                                    parseInt(line, 12, 6),
                                    parseInt(line, 18, 6),
                                    parseInt(line, 24, 6),
                                    parseDouble(line, 30, 13),
                                    utc);
        }

        /** Build the coordinate array from a parsed line.
         * @param line to parse
         * @return an array of coordinates in radians
         */
        private double[] parseCoordinate(final String line) {
            final double a = parseDouble(line, 2, 6);
            final double b = parseDouble(line, 8, 6);
            final double c = parseDouble(line, 14, 6);
            final double[] coordinate = new double[((int) FastMath.abs((a - b) / c)) + 1];
            int i = 0;
            for (double cor = FastMath.min(a, b); cor <= FastMath.max(a, b); cor += FastMath.abs(c)) {
                coordinate[i] = FastMath.toRadians(cor);
                i++;
            }
            return coordinate;
        }

        /** Interpolate the TEC in latitude and longitude.
         * @param exponent exponent defining the unit of the values listed in the data blocks
         * @param values TEC values
         * @param latitudes array containing the different latitudes in radians
         * @param longitudes array containing the different latitudes in radians
         * @return the interpolated TEC in TECUnits
         */
        private double interpolateTEC(final ArrayList<Double> values, final double exponent,
                                      final double[] latitudes, final double[] longitudes) {
            // Array dimensions
            final int dimLat = latitudes.length;
            final int dimLon = longitudes.length;

            // Build the array of TEC data
            final double[][] fvalTEC = new double[dimLat][dimLon];
            int index = dimLon * dimLat;
            for (int x = 0; x < dimLat; x++) {
                for (int y = dimLon - 1; y >= 0; y--) {
                    index = index - 1;
                    fvalTEC[x][y] = values.get(index);
                }
            }

            // Build Bilinear Interpolation function
            final BilinearInterpolatingFunction functionTEC = new BilinearInterpolatingFunction(latitudes, longitudes, fvalTEC);
            final double tec = functionTEC.value(latitude, longitude) * FastMath.pow(10.0, exponent);
            return tec;
        }
    }

    /**
     * Container for IONEX data.
     * <p>
     * The TEC contained in the map is previously interpolated
     * according to the latitude and the longitude given by the user.
     * </p>
     */
    private static class TECMap {

        /** Date of the TEC Map. */
        private AbsoluteDate date;

        /** Interpolated TEC [TECUnits]. */
        private double tec;

        /**
         * Constructor.
         * @param date date of the TEC map
         * @param tec interpolated tec
         */
        TECMap(final AbsoluteDate date, final double tec) {
            this.date = date;
            this.tec  = tec;
        }

        /**
         * Get the date of the TEC map.
         * @return the date
         */
        public AbsoluteDate getDate() {
            return date;
        }

        /**
         * Get the value of the interpolated TEC.
         * @return the TEC in TECUnits
         */
        public double getTEC() {
            return tec;
        }

    }

    /** Container for IONEX header. */
    private static class IONEXHeader {

        /** Epoch of the first TEC map. */
        private AbsoluteDate firstDate;

        /** Epoch of the last TEC map. */
        private AbsoluteDate lastDate;

        /** Interval between two maps [s]. */
        private int interval;

        /** Number of maps contained in the IONEX file. */
        private int nbOfMaps;

        /** Mean earth radius [m]. */
        private double baseRadius;

        /** Height of the ionospheric single layer [m]. */
        private double hIon;

        /** Flag for mapping function adopted for TEC determination. */
        private boolean isMappingFunction;

        /**
         * Constructor.
         * @param firstDate epoch of the first TEC map.
         * @param lastDate epoch of the last TEC map.
         * @param nbOfMaps number of TEC maps contained in the file
         * @param interval number of seconds between two tec maps.
         * @param baseRadius mean earth radius in meters
         * @param hIon height of the ionospheric single layer in meters
         * @param mappingFunction flag for mapping function adopted for TEC determination
         */
        IONEXHeader(final AbsoluteDate firstDate, final AbsoluteDate lastDate,
                    final int interval, final int nbOfMaps,
                    final double baseRadius, final double hIon,
                    final boolean mappingFunction) {
            this.firstDate         = firstDate;
            this.lastDate          = lastDate;
            this.interval          = interval;
            this.nbOfMaps          = nbOfMaps;
            this.baseRadius        = baseRadius;
            this.hIon              = hIon;
            this.isMappingFunction = mappingFunction;
        }

        /**
         * Get the first date of the IONEX file.
         * @return the first date of the IONEX file
         */
        public AbsoluteDate getFirstDate() {
            return firstDate;
        }

        /**
         * Get the last date of the IONEX file.
         * @return the last date of the IONEX file
         */
        public AbsoluteDate getLastDate() {
            return lastDate;
        }

        /**
         * Get the time interval between two TEC maps.
         * @return the interval between two TEC maps
         */
        public int getInterval() {
            return interval;
        }

        /**
         * Get the number of TEC maps contained in the file.
         * @return the number of TEC maps
         */
        public int getTECMapsNumer() {
            return nbOfMaps;
        }

        /**
         * Get the mean earth radius in meters.
         * @return the mean earth radius
         */
        public double getEarthRadius() {
            return baseRadius;
        }

        /**
         * Get the height of the ionospheric single layer in meters.
         * @return the height of the ionospheric single layer
         */
        public double getHIon() {
            return hIon;
        }

        /**
         * Get the mapping function flag.
         * @return false if mapping function computation is needed
         */
        public boolean isMappingFunction() {
            return isMappingFunction;
        }

    }

    /** Internal class used only for serialization. */
    @DefaultDataContext
    private static class DataTransferObject implements Serializable {

        /** Serializable UID. */
        private static final long serialVersionUID = 201928052L;

        /** Regular expression that matches the names of the IONEX files. */
        private final String supportedNames;

        /** Simple constructor.
         * @param supportedNames regular expression that matches the names of the IONEX files
         */
        DataTransferObject(final String supportedNames) {
            this.supportedNames = supportedNames;
        }

        /** Replace the deserialized data transfer object with a {@link GlobalIonosphereMapModel}.
         * @return replacement {@link GlobalIonosphereMapModel}
         */
        private Object readResolve() {
            try {
                return new GlobalIonosphereMapModel(supportedNames);
            } catch (OrekitException oe) {
                throw new OrekitInternalError(oe);
            }
        }

    }

}