1   /* Copyright 2002-2025 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.files.rinex.navigation;
18  
19  import java.io.BufferedReader;
20  import java.io.IOException;
21  import java.io.Reader;
22  import java.util.Arrays;
23  import java.util.Collections;
24  import java.util.InputMismatchException;
25  import java.util.function.Function;
26  import java.util.function.Predicate;
27  
28  import org.hipparchus.util.FastMath;
29  import org.orekit.annotation.DefaultDataContext;
30  import org.orekit.data.DataContext;
31  import org.orekit.data.DataSource;
32  import org.orekit.errors.OrekitException;
33  import org.orekit.errors.OrekitInternalError;
34  import org.orekit.errors.OrekitMessages;
35  import org.orekit.files.rinex.utils.parsing.RinexUtils;
36  import org.orekit.gnss.PredefinedGnssSignal;
37  import org.orekit.gnss.SatelliteSystem;
38  import org.orekit.gnss.TimeSystem;
39  import org.orekit.propagation.analytical.gnss.data.AbstractNavigationMessage;
40  import org.orekit.propagation.analytical.gnss.data.BeidouCivilianNavigationMessage;
41  import org.orekit.propagation.analytical.gnss.data.BeidouLegacyNavigationMessage;
42  import org.orekit.propagation.analytical.gnss.data.BeidouSatelliteType;
43  import org.orekit.propagation.analytical.gnss.data.CivilianNavigationMessage;
44  import org.orekit.propagation.analytical.gnss.data.GLONASSNavigationMessage;
45  import org.orekit.propagation.analytical.gnss.data.GPSCivilianNavigationMessage;
46  import org.orekit.propagation.analytical.gnss.data.GPSLegacyNavigationMessage;
47  import org.orekit.propagation.analytical.gnss.data.GalileoNavigationMessage;
48  import org.orekit.propagation.analytical.gnss.data.NavICL1NVNavigationMessage;
49  import org.orekit.propagation.analytical.gnss.data.NavICLegacyNavigationMessage;
50  import org.orekit.propagation.analytical.gnss.data.LegacyNavigationMessage;
51  import org.orekit.propagation.analytical.gnss.data.QZSSCivilianNavigationMessage;
52  import org.orekit.propagation.analytical.gnss.data.QZSSLegacyNavigationMessage;
53  import org.orekit.propagation.analytical.gnss.data.SBASNavigationMessage;
54  import org.orekit.time.AbsoluteDate;
55  import org.orekit.time.GNSSDate;
56  import org.orekit.time.TimeScale;
57  import org.orekit.time.TimeScales;
58  import org.orekit.utils.Constants;
59  import org.orekit.utils.units.Unit;
60  
61  /**
62   * Parser for RINEX navigation messages files.
63   * <p>
64   * This parser handles RINEX version from 2 to 4.02.
65   * </p>
66   * @see <a href="https://files.igs.org/pub/data/format/rinex2.txt">rinex 2.0</a>
67   * @see <a href="https://files.igs.org/pub/data/format/rinex210.txt">rinex 2.10</a>
68   * @see <a href="https://files.igs.org/pub/data/format/rinex211.pdf">rinex 2.11</a>
69   * @see <a href="https://files.igs.org/pub/data/format/rinex301.pdf"> 3.01 navigation messages file format</a>
70   * @see <a href="https://files.igs.org/pub/data/format/rinex302.pdf"> 3.02 navigation messages file format</a>
71   * @see <a href="https://files.igs.org/pub/data/format/rinex303.pdf"> 3.03 navigation messages file format</a>
72   * @see <a href="https://files.igs.org/pub/data/format/rinex304.pdf"> 3.04 navigation messages file format</a>
73   * @see <a href="https://files.igs.org/pub/data/format/rinex305.pdf"> 3.05 navigation messages file format</a>
74   * @see <a href="https://files.igs.org/pub/data/format/rinex_4.00.pdf"> 4.00 navigation messages file format</a>
75   * @see <a href="https://files.igs.org/pub/data/format/rinex_4.01.pdf"> 4.01 navigation messages file format</a>
76   * @see <a href="https://files.igs.org/pub/data/format/rinex_4.02.pdf"> 4.02 navigation messages file format</a>
77   *
78   * @author Bryan Cazabonne
79   * @since 11.0
80   *
81   */
82  public class RinexNavigationParser {
83  
84      /** Converter for positions. */
85      private static final Unit KM = Unit.KILOMETRE;
86  
87      /** Converter for velocities. */
88      private static final Unit KM_PER_S = Unit.parse("km/s");
89  
90      /** Converter for accelerations. */
91      private static final Unit KM_PER_S2 = Unit.parse("km/s²");
92  
93      /** Converter for velocities. */
94      private static final Unit M_PER_S = Unit.parse("m/s");
95  
96      /** Converter for clock drift. */
97      private static final Unit S_PER_S = Unit.parse("s/s");
98  
99      /** Converter for clock drift rate. */
100     private static final Unit S_PER_S2 = Unit.parse("s/s²");
101 
102     /** Converter for ΔUT₁ first derivative. */
103     private static final Unit S_PER_DAY = Unit.parse("s/d");
104 
105     /** Converter for ΔUT₁ second derivative. */
106     private static final Unit S_PER_DAY2 = Unit.parse("s/d²");
107 
108     /** Converter for square root of semi-major axis. */
109     private static final Unit SQRT_M = Unit.parse("√m");
110 
111     /** Converter for angular rates. */
112     private static final Unit RAD_PER_S = Unit.parse("rad/s");
113 
114     /** Converter for angular accelerations. */
115     private static final Unit RAD_PER_S2 = Unit.parse("rad/s²");
116 
117     /** Converter for rates of small angle. */
118     private static final Unit AS_PER_DAY = Unit.parse("as/d");
119 
120     /** Converter for accelerations of small angles. */
121     private static final Unit AS_PER_DAY2 = Unit.parse("as/d²");
122 
123     /** System initials. */
124     private static final String INITIALS = "GRECIJS";
125 
126     /** URA index to URA mapping (table 23 of NavIC ICD). */
127     // CHECKSTYLE: stop Indentation check
128     private static final double[] NAVIC_URA = {
129            2.40,    3.40,    4.85,   6.85,
130            9.65,   13.65,   24.00,  48.00,
131           96.00,  192.00,  384.00, 768.00,
132         1536.00, 3072.00, 6144.00, Double.NaN
133     };
134     // CHECKSTYLE: resume Indentation check
135 
136     /** Set of time scales. */
137     private final TimeScales timeScales;
138 
139     /**
140      * Constructor.
141      * <p>This constructor uses the {@link DataContext#getDefault() default data context}.</p>
142      * @see #RinexNavigationParser(TimeScales)
143      *
144      */
145     @DefaultDataContext
146     public RinexNavigationParser() {
147         this(DataContext.getDefault().getTimeScales());
148     }
149 
150     /**
151      * Constructor.
152      * @param timeScales the set of time scales used for parsing dates.
153      */
154     public RinexNavigationParser(final TimeScales timeScales) {
155         this.timeScales = timeScales;
156     }
157 
158     /**
159      * Parse RINEX navigation messages.
160      * @param source source providing the data to parse
161      * @return a parsed  RINEX navigation messages file
162      * @throws IOException if {@code reader} throws one
163      */
164     public RinexNavigation parse(final DataSource source) throws IOException {
165 
166         // initialize internal data structures
167         final ParseInfo pi = new ParseInfo(source.getName());
168 
169         Iterable<LineParser> candidateParsers = Collections.singleton(LineParser.HEADER_VERSION);
170         try (Reader reader = source.getOpener().openReaderOnce();
171              BufferedReader br = new BufferedReader(reader)) {
172             nextLine:
173                 for (String line = br.readLine(); line != null; line = br.readLine()) {
174                     ++pi.lineNumber;
175                     for (final LineParser candidate : candidateParsers) {
176                         if (candidate.canHandle.test(line)) {
177                             try {
178                                 candidate.parsingMethod.parse(line, pi);
179                                 candidateParsers = candidate.allowedNextProvider.apply(pi);
180                                 continue nextLine;
181                             } catch (StringIndexOutOfBoundsException | NumberFormatException | InputMismatchException e) {
182                                 throw new OrekitException(e,
183                                                           OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
184                                                           pi.lineNumber, source.getName(), line);
185                             }
186                         }
187                     }
188                     throw new OrekitException(OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
189                                               pi.lineNumber, source.getName(), line);
190                 }
191         }
192 
193         if (!pi.headerParsed) {
194             throw new OrekitException(OrekitMessages.UNEXPECTED_END_OF_FILE, source.getName());
195         }
196 
197         pi.closePendingMessage();
198 
199         return pi.file;
200 
201     }
202 
203     /** Transient data used for parsing a RINEX navigation messages file. */
204     private class ParseInfo {
205 
206         /** Name of the data source. */
207         private final String name;
208 
209         /** Set of time scales for parsing dates. */
210         private final TimeScales timeScales;
211 
212         /** The corresponding navigation messages file object. */
213         private final RinexNavigation file;
214 
215         /** Number of initial spaces in broadcast orbits lines. */
216         private int initialSpaces;
217 
218         /** Flag indicating header has been completely parsed. */
219         private boolean headerParsed;
220 
221         /** Flag indicating the distinction between "alpha" and "beta" ionospheric coefficients. */
222         private boolean isIonosphereAlphaInitialized;
223 
224         /** Satellite system line parser. */
225         private SatelliteSystemLineParser systemLineParser;
226 
227         /** Current global line number. */
228         private int lineNumber;
229 
230         /** Current line number within the navigation message. */
231         private int messageLineNumber;
232 
233         /** Container for GPS navigation message. */
234         private GPSLegacyNavigationMessage gpsLNav;
235 
236         /** Container for GPS navigation message. */
237         private GPSCivilianNavigationMessage gpsCNav;
238 
239         /** Container for Galileo navigation message. */
240         private GalileoNavigationMessage galileoNav;
241 
242         /** Container for Beidou navigation message. */
243         private BeidouLegacyNavigationMessage beidouLNav;
244 
245         /** Container for Beidou navigation message. */
246         private BeidouCivilianNavigationMessage beidouCNav;
247 
248         /** Container for QZSS navigation message. */
249         private QZSSLegacyNavigationMessage qzssLNav;
250 
251         /** Container for QZSS navigation message. */
252         private QZSSCivilianNavigationMessage qzssCNav;
253 
254         /** Container for NavIC navigation message. */
255         private NavICLegacyNavigationMessage navicLNav;
256 
257         /** Container for NavIC navigation message.
258          * @since 13.0
259          */
260         private NavICL1NVNavigationMessage navicL1NV;
261 
262         /** Container for GLONASS navigation message. */
263         private GLONASSNavigationMessage glonassNav;
264 
265         /** Container for SBAS navigation message. */
266         private SBASNavigationMessage sbasNav;
267 
268         /** Container for System Time Offset message. */
269         private SystemTimeOffsetMessage sto;
270 
271         /** Container for Earth Orientation Parameter message. */
272         private EarthOrientationParameterMessage eop;
273 
274         /** Container for ionosphere Klobuchar message. */
275         private IonosphereKlobucharMessage klobuchar;
276 
277         /** Container for ionosphere Nequick-G message. */
278         private IonosphereNequickGMessage nequickG;
279 
280         /** Container for ionosphere BDGIM message. */
281         private IonosphereBDGIMMessage bdgim;
282 
283         /** Constructor, build the ParseInfo object.
284          * @param name name of the data source
285          */
286         ParseInfo(final String name) {
287             // Initialize default values for fields
288             this.name                         = name;
289             this.timeScales                   = RinexNavigationParser.this.timeScales;
290             this.isIonosphereAlphaInitialized = false;
291             this.file                         = new RinexNavigation();
292 
293         }
294 
295         /** Ensure navigation message has been closed.
296          */
297         void closePendingMessage() {
298             if (systemLineParser != null) {
299                 systemLineParser.closeMessage(this);
300                 systemLineParser = null;
301             }
302 
303         }
304 
305     }
306 
307     /** Parsers for specific lines. */
308     private enum LineParser {
309 
310         /** Parser for version, file type and satellite system. */
311         HEADER_VERSION(line -> RinexUtils.matchesLabel(line, "RINEX VERSION / TYPE"),
312                        (line, pi) -> {
313                            RinexUtils.parseVersionFileTypeSatelliteSystem(line, pi.name, pi.file.getHeader(),
314                                                                           2.0, 2.01, 2.10, 2.11,
315                                                                           3.01, 3.02, 3.03, 3.04, 3.05,
316                                                                           4.00, 4.01, 4.02);
317                            pi.initialSpaces = pi.file.getHeader().getFormatVersion() < 3.0 ? 3 : 4;
318                        },
319                        LineParser::headerNext),
320 
321         /** Parser for generating program and emitting agency. */
322         HEADER_PROGRAM(line -> RinexUtils.matchesLabel(line, "PGM / RUN BY / DATE"),
323                        (line, pi) -> RinexUtils.parseProgramRunByDate(line, pi.lineNumber, pi.name, pi.timeScales, pi.file.getHeader()),
324                        LineParser::headerNext),
325 
326         /** Parser for comments. */
327         HEADER_COMMENT(line -> RinexUtils.matchesLabel(line, "COMMENT"),
328                        (line, pi) -> RinexUtils.parseComment(pi.lineNumber, line, pi.file),
329                        LineParser::headerNext),
330 
331         /** Parser for ionospheric correction parameters. */
332         HEADER_ION_ALPHA(line -> RinexUtils.matchesLabel(line, "ION ALPHA"),
333                          (line, pi) -> {
334 
335                              pi.file.getHeader().setIonosphericCorrectionType(IonosphericCorrectionType.GPS);
336 
337                              // Read coefficients
338                              final double[] parameters = new double[4];
339                              parameters[0] = RinexUtils.parseDouble(line, 2,  12);
340                              parameters[1] = RinexUtils.parseDouble(line, 14, 12);
341                              parameters[2] = RinexUtils.parseDouble(line, 26, 12);
342                              parameters[3] = RinexUtils.parseDouble(line, 38, 12);
343                              pi.file.setKlobucharAlpha(parameters);
344                              pi.isIonosphereAlphaInitialized = true;
345 
346                          },
347                          LineParser::headerNext),
348 
349         /** Parser for ionospheric correction parameters. */
350         HEADER_ION_BETA(line -> RinexUtils.matchesLabel(line, "ION BETA"),
351                         (line, pi) -> {
352 
353                             pi.file.getHeader().setIonosphericCorrectionType(IonosphericCorrectionType.GPS);
354 
355                             // Read coefficients
356                             final double[] parameters = new double[4];
357                             parameters[0] = RinexUtils.parseDouble(line, 2,  12);
358                             parameters[1] = RinexUtils.parseDouble(line, 14, 12);
359                             parameters[2] = RinexUtils.parseDouble(line, 26, 12);
360                             parameters[3] = RinexUtils.parseDouble(line, 38, 12);
361                             pi.file.setKlobucharBeta(parameters);
362 
363                         },
364                         LineParser::headerNext),
365 
366         /** Parser for ionospheric correction parameters. */
367         HEADER_IONOSPHERIC(line -> RinexUtils.matchesLabel(line, "IONOSPHERIC CORR"),
368                            (line, pi) -> {
369 
370                                // ionospheric correction type
371                                final IonosphericCorrectionType ionoType =
372                                                IonosphericCorrectionType.valueOf(RinexUtils.parseString(line, 0, 3));
373                                pi.file.getHeader().setIonosphericCorrectionType(ionoType);
374 
375                                // Read coefficients
376                                final double[] parameters = new double[4];
377                                parameters[0] = RinexUtils.parseDouble(line, 5,  12);
378                                parameters[1] = RinexUtils.parseDouble(line, 17, 12);
379                                parameters[2] = RinexUtils.parseDouble(line, 29, 12);
380                                parameters[3] = RinexUtils.parseDouble(line, 41, 12);
381 
382                                // Verify if we are parsing Galileo ionospheric parameters
383                                if (ionoType == IonosphericCorrectionType.GAL) {
384 
385                                    // We are parsing Galileo ionospheric parameters
386                                    pi.file.setNeQuickAlpha(parameters);
387 
388                                } else {
389                                    // We are parsing Klobuchar ionospheric parameters
390 
391                                    // Verify if we are parsing "alpha" or "beta" ionospheric parameters
392                                    if (pi.isIonosphereAlphaInitialized) {
393 
394                                        // Ionospheric "beta" parameters
395                                        pi.file.setKlobucharBeta(parameters);
396 
397                                    } else {
398 
399                                        // Ionospheric "alpha" parameters
400                                        pi.file.setKlobucharAlpha(parameters);
401 
402                                        // Set the flag to true
403                                        pi.isIonosphereAlphaInitialized = true;
404 
405                                    }
406 
407                                }
408 
409                            },
410                            LineParser::headerNext),
411 
412         /** Parser for corrections to transform the system time to UTC or to other time systems. */
413         HEADER_DELTA_UTC(line -> RinexUtils.matchesLabel(line, "DELTA-UTC: A0,A1,T,W"),
414                          (line, pi) -> {
415                              // Read fields
416                              final double a0      = RinexUtils.parseDouble(line, 3,  19);
417                              final double a1      = RinexUtils.parseDouble(line, 22, 19);
418                              final int    refTime = RinexUtils.parseInt(line, 41, 9);
419                              final int    refWeek = RinexUtils.parseInt(line, 50, 9);
420 
421                              // convert date
422                              final SatelliteSystem satSystem = pi.file.getHeader().getSatelliteSystem();
423                              final AbsoluteDate    date      = new GNSSDate(refWeek, refTime, satSystem, pi.timeScales).getDate();
424 
425                              // Add to the list
426                              final TimeSystemCorrection tsc = new TimeSystemCorrection("GPUT", date, a0, a1);
427                              pi.file.getHeader().addTimeSystemCorrections(tsc);
428                          },
429                          LineParser::headerNext),
430 
431         /** Parser for corrections to transform the GLONASS system time to UTC or to other time systems. */
432         HEADER_CORR_SYSTEM_TIME(line -> RinexUtils.matchesLabel(line, "CORR TO SYSTEM TIME"),
433                          (line, pi) -> {
434                              // Read fields
435                              final int year        = RinexUtils.parseInt(line,  0, 6);
436                              final int month       = RinexUtils.parseInt(line,  6, 6);
437                              final int day         = RinexUtils.parseInt(line, 12, 6);
438                              final double minusTau = RinexUtils.parseDouble(line, 21, 19);
439 
440                              // convert date
441                              final SatelliteSystem satSystem = pi.file.getHeader().getSatelliteSystem();
442                              final TimeScale       timeScale = satSystem.getObservationTimeScale().getTimeScale(pi.timeScales);
443                              final AbsoluteDate    date      = new AbsoluteDate(year, month, day, timeScale);
444 
445                              // Add to the list
446                              final TimeSystemCorrection tsc = new TimeSystemCorrection("GLUT", date, minusTau, 0.0);
447                              pi.file.getHeader().addTimeSystemCorrections(tsc);
448 
449                          },
450                          LineParser::headerNext),
451 
452         /** Parser for corrections to transform the system time to UTC or to other time systems. */
453         HEADER_TIME(line -> RinexUtils.matchesLabel(line, "TIME SYSTEM CORR"),
454                     (line, pi) -> {
455 
456                         // Read fields
457                         final String type    = RinexUtils.parseString(line, 0,  4);
458                         final double a0      = RinexUtils.parseDouble(line, 5,  17);
459                         final double a1      = RinexUtils.parseDouble(line, 22, 16);
460                         final int    refTime = RinexUtils.parseInt(line, 38, 7);
461                         final int    refWeek = RinexUtils.parseInt(line, 46, 5);
462 
463                         // convert date
464                         final SatelliteSystem satSystem = pi.file.getHeader().getSatelliteSystem();
465                         final AbsoluteDate    date;
466                         if (satSystem == SatelliteSystem.GLONASS) {
467                             date = null;
468                         } else if (satSystem == SatelliteSystem.BEIDOU) {
469                             date = new GNSSDate(refWeek, refTime, satSystem, pi.timeScales).getDate();
470                         } else {
471                             // all other systems are converted to GPS week in Rinex files!
472                             date = new GNSSDate(refWeek, refTime, SatelliteSystem.GPS, pi.timeScales).getDate();
473                         }
474 
475                         // Add to the list
476                         final TimeSystemCorrection tsc = new TimeSystemCorrection(type, date, a0, a1);
477                         pi.file.getHeader().addTimeSystemCorrections(tsc);
478 
479                     },
480                     LineParser::headerNext),
481 
482         /** Parser for leap seconds. */
483         HEADER_LEAP_SECONDS(line -> RinexUtils.matchesLabel(line, "LEAP SECONDS"),
484                             (line, pi) -> pi.file.getHeader().setNumberOfLeapSeconds(RinexUtils.parseInt(line, 0, 6)),
485                             LineParser::headerNext),
486 
487         /** Parser for DOI.
488          * @since 12.0
489          */
490         HEADER_DOI(line -> RinexUtils.matchesLabel(line, "DOI"),
491                    (line, pi) -> pi.file.getHeader().setDoi(RinexUtils.parseString(line, 0, RinexUtils.LABEL_INDEX)),
492                    LineParser::headerNext),
493 
494         /** Parser for license.
495          * @since 12.0
496          */
497         HEADER_LICENSE(line -> RinexUtils.matchesLabel(line, "LICENSE OF USE"),
498                        (line, pi) -> pi.file.getHeader().setLicense(RinexUtils.parseString(line, 0, RinexUtils.LABEL_INDEX)),
499                        LineParser::headerNext),
500 
501         /** Parser for stationInformation.
502          * @since 12.0
503          */
504         HEADER_STATION_INFORMATION(line -> RinexUtils.matchesLabel(line, "STATION INFORMATION"),
505                                    (line, pi) -> pi.file.getHeader().setStationInformation(RinexUtils.parseString(line, 0, RinexUtils.LABEL_INDEX)),
506                                    LineParser::headerNext),
507 
508         /** Parser for merged files.
509          * @since 12.0
510          */
511         HEADER_MERGED_FILE(line -> RinexUtils.matchesLabel(line, "MERGED FILE"),
512                            (line, pi) -> pi.file.getHeader().setMergedFiles(RinexUtils.parseInt(line, 0, 9)),
513                            LineParser::headerNext),
514 
515        /** Parser for the end of header. */
516         HEADER_END(line -> RinexUtils.matchesLabel(line, "END OF HEADER"),
517                    (line, pi) -> {
518                        // get rinex format version
519                        final RinexNavigationHeader header = pi.file.getHeader();
520                        final double version = header.getFormatVersion();
521 
522                        // check mandatory header fields
523                        if (header.getRunByName() == null ||
524                            version >= 4 && header.getNumberOfLeapSeconds() < 0) {
525                            throw new OrekitException(OrekitMessages.INCOMPLETE_HEADER, pi.name);
526                        }
527 
528                        pi.headerParsed = true;
529 
530                    },
531                    LineParser::navigationNext),
532 
533         /** Parser for navigation message space vehicle epoch and clock. */
534         NAVIGATION_SV_EPOCH_CLOCK_RINEX_2(line -> true,
535                                           (line, pi) -> {
536 
537                                               // Set the line number to 0
538                                               pi.messageLineNumber = 0;
539 
540                                               // Initialize parser
541                                               pi.closePendingMessage();
542                                               pi.systemLineParser = SatelliteSystemLineParser.getParser(pi.file.getHeader().getSatelliteSystem(),
543                                                                                                         null, pi, line);
544 
545                                               pi.systemLineParser.parseSvEpochSvClockLine(line, pi);
546 
547                                           },
548                                           LineParser::navigationNext),
549 
550         /** Parser for navigation message space vehicle epoch and clock. */
551         NAVIGATION_SV_EPOCH_CLOCK(line -> INITIALS.indexOf(line.charAt(0)) >= 0,
552                                   (line, pi) -> {
553 
554                                       // Set the line number to 0
555                                       pi.messageLineNumber = 0;
556 
557                                       if (pi.file.getHeader().getFormatVersion() < 4) {
558                                           // Current satellite system
559                                           final SatelliteSystem system = SatelliteSystem.parseSatelliteSystem(RinexUtils.parseString(line, 0, 1));
560 
561                                           // Initialize parser
562                                           pi.closePendingMessage();
563                                           pi.systemLineParser = SatelliteSystemLineParser.getParser(system, null, pi, line);
564                                       }
565 
566                                       // Read first line
567                                       pi.systemLineParser.parseSvEpochSvClockLine(line, pi);
568 
569                                   },
570                                   LineParser::navigationNext),
571 
572         /** Parser for navigation message type. */
573         EPH_TYPE(line -> line.startsWith("> EPH"),
574                  (line, pi) -> {
575                      final SatelliteSystem system = SatelliteSystem.parseSatelliteSystem(RinexUtils.parseString(line, 6, 1));
576                      final String          type   = RinexUtils.parseString(line, 10, 4);
577                      pi.closePendingMessage();
578                      pi.systemLineParser = SatelliteSystemLineParser.getParser(system, type, pi, line);
579                  },
580                  pi -> Collections.singleton(NAVIGATION_SV_EPOCH_CLOCK)),
581 
582         /** Parser for broadcast orbit. */
583         BROADCAST_ORBIT(line -> line.startsWith("   "),
584                         (line, pi) -> {
585                             switch (++pi.messageLineNumber) {
586                                 case 1: pi.systemLineParser.parseFirstBroadcastOrbit(line, pi);
587                                 break;
588                                 case 2: pi.systemLineParser.parseSecondBroadcastOrbit(line, pi);
589                                 break;
590                                 case 3: pi.systemLineParser.parseThirdBroadcastOrbit(line, pi);
591                                 break;
592                                 case 4: pi.systemLineParser.parseFourthBroadcastOrbit(line, pi);
593                                 break;
594                                 case 5: pi.systemLineParser.parseFifthBroadcastOrbit(line, pi);
595                                 break;
596                                 case 6: pi.systemLineParser.parseSixthBroadcastOrbit(line, pi);
597                                 break;
598                                 case 7: pi.systemLineParser.parseSeventhBroadcastOrbit(line, pi);
599                                 break;
600                                 case 8: pi.systemLineParser.parseEighthBroadcastOrbit(line, pi);
601                                 break;
602                                 case 9: pi.systemLineParser.parseNinthBroadcastOrbit(line, pi);
603                                 break;
604                                 default:
605                                     // this should never happen
606                                     throw new OrekitException(OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
607                                                               pi.lineNumber, pi.name, line);
608                             }
609 
610                         },
611                         LineParser::navigationNext),
612 
613         /** Parser for system time offset message model. */
614         STO_LINE_1(line -> true,
615                    (line, pi) -> {
616                        pi.sto.setTransmissionTime(Unit.SECOND.toSI(RinexUtils.parseDouble(line,  4, 19)));
617                        pi.sto.setA0(Unit.SECOND.toSI(RinexUtils.parseDouble(line, 23, 19)));
618                        pi.sto.setA1(S_PER_S.toSI(RinexUtils.parseDouble(line, 42, 19)));
619                        pi.sto.setA2(S_PER_S2.toSI(RinexUtils.parseDouble(line, 61, 19)));
620                        pi.file.addSystemTimeOffset(pi.sto);
621                        pi.sto = null;
622                    },
623                    LineParser::navigationNext),
624 
625         /** Parser for system time offset message space vehicle epoch and clock. */
626         STO_SV_EPOCH_CLOCK(line -> true,
627                            (line, pi) -> {
628 
629                                pi.sto.setDefinedTimeSystem(TimeSystem.parseTwoLettersCode(RinexUtils.parseString(line, 24, 2)));
630                                pi.sto.setReferenceTimeSystem(TimeSystem.parseTwoLettersCode(RinexUtils.parseString(line, 26, 2)));
631                                final String sbas = RinexUtils.parseString(line, 43, 18);
632                                pi.sto.setSbasId(!sbas.isEmpty() ? SbasId.valueOf(sbas) : null);
633                                final String utc = RinexUtils.parseString(line, 62, 18);
634                                pi.sto.setUtcId(!utc.isEmpty() ? UtcId.parseUtcId(utc) : null);
635 
636                                // TODO is the reference date relative to one or the other time scale?
637                                final int year  = RinexUtils.parseInt(line, 4, 4);
638                                final int month = RinexUtils.parseInt(line, 9, 2);
639                                final int day   = RinexUtils.parseInt(line, 12, 2);
640                                final int hours = RinexUtils.parseInt(line, 15, 2);
641                                final int min   = RinexUtils.parseInt(line, 18, 2);
642                                final int sec   = RinexUtils.parseInt(line, 21, 2);
643                                pi.sto.setReferenceEpoch(new AbsoluteDate(year, month, day, hours, min, sec,
644                                                                          pi.sto.getDefinedTimeSystem().getTimeScale(pi.timeScales)));
645 
646                            },
647                            pi -> Collections.singleton(STO_LINE_1)),
648 
649         /** Parser for system time offset message type. */
650         STO_TYPE(line -> line.startsWith("> STO"),
651                  (line, pi) -> {
652                      pi.closePendingMessage();
653                      pi.sto = new SystemTimeOffsetMessage(SatelliteSystem.parseSatelliteSystem(RinexUtils.parseString(line, 6, 1)),
654                                                           RinexUtils.parseInt(line, 7, 2),
655                                                           RinexUtils.parseString(line, 10, 4));
656                  },
657                  pi -> Collections.singleton(STO_SV_EPOCH_CLOCK)),
658 
659         /** Parser for Earth orientation parameter message model. */
660         EOP_LINE_2(line -> true,
661                    (line, pi) -> {
662                        pi.eop.setTransmissionTime(Unit.SECOND.toSI(RinexUtils.parseDouble(line,  4, 19)));
663                        pi.eop.setDut1(Unit.SECOND.toSI(RinexUtils.parseDouble(line, 23, 19)));
664                        pi.eop.setDut1Dot(S_PER_DAY.toSI(RinexUtils.parseDouble(line, 42, 19)));
665                        pi.eop.setDut1DotDot(S_PER_DAY2.toSI(RinexUtils.parseDouble(line, 61, 19)));
666                        pi.file.addEarthOrientationParameter(pi.eop);
667                        pi.eop = null;
668                    },
669                    LineParser::navigationNext),
670 
671         /** Parser for Earth orientation parameter message model. */
672         EOP_LINE_1(line -> true,
673                    (line, pi) -> {
674                        pi.eop.setYp(Unit.ARC_SECOND.toSI(RinexUtils.parseDouble(line, 23, 19)));
675                        pi.eop.setYpDot(AS_PER_DAY.toSI(RinexUtils.parseDouble(line, 42, 19)));
676                        pi.eop.setYpDotDot(AS_PER_DAY2.toSI(RinexUtils.parseDouble(line, 61, 19)));
677                    },
678                    pi -> Collections.singleton(EOP_LINE_2)),
679 
680         /** Parser for Earth orientation parameter message model. */
681         EOP_LINE_0(line -> true,
682                            (line, pi) -> {
683                                final int year  = RinexUtils.parseInt(line, 4, 4);
684                                final int month = RinexUtils.parseInt(line, 9, 2);
685                                final int day   = RinexUtils.parseInt(line, 12, 2);
686                                final int hours = RinexUtils.parseInt(line, 15, 2);
687                                final int min   = RinexUtils.parseInt(line, 18, 2);
688                                final int sec   = RinexUtils.parseInt(line, 21, 2);
689                                pi.eop.setReferenceEpoch(new AbsoluteDate(year, month, day, hours, min, sec,
690                                                                          pi.eop.getSystem().getObservationTimeScale().getTimeScale(pi.timeScales)));
691                                pi.eop.setXp(Unit.ARC_SECOND.toSI(RinexUtils.parseDouble(line, 23, 19)));
692                                pi.eop.setXpDot(AS_PER_DAY.toSI(RinexUtils.parseDouble(line, 42, 19)));
693                                pi.eop.setXpDotDot(AS_PER_DAY2.toSI(RinexUtils.parseDouble(line, 61, 19)));
694                            },
695                            pi -> Collections.singleton(EOP_LINE_1)),
696 
697         /** Parser for Earth orientation parameter message type. */
698         EOP_TYPE(line -> line.startsWith("> EOP"),
699                  (line, pi) -> {
700                      pi.closePendingMessage();
701                      pi.eop = new EarthOrientationParameterMessage(SatelliteSystem.parseSatelliteSystem(RinexUtils.parseString(line, 6, 1)),
702                                                                    RinexUtils.parseInt(line, 7, 2),
703                                                                    RinexUtils.parseString(line, 10, 4));
704                  },
705                  pi -> Collections.singleton(EOP_LINE_0)),
706 
707         /** Parser for ionosphere Klobuchar message model. */
708         KLOBUCHAR_LINE_2(line -> true,
709                          (line, pi) -> {
710                              pi.klobuchar.setBetaI(3, IonosphereKlobucharMessage.S_PER_SC_N[3].toSI(RinexUtils.parseDouble(line,  4, 19)));
711                              pi.klobuchar.setRegionCode(RinexUtils.parseDouble(line, 23, 19) < 0.5 ?
712                                                         RegionCode.WIDE_AREA : RegionCode.JAPAN);
713                              pi.file.addKlobucharMessage(pi.klobuchar);
714                              pi.klobuchar = null;
715                          },
716                          LineParser::navigationNext),
717 
718         /** Parser for ionosphere Klobuchar message model. */
719         KLOBUCHAR_LINE_1(line -> true,
720                          (line, pi) -> {
721                              pi.klobuchar.setAlphaI(3, IonosphereKlobucharMessage.S_PER_SC_N[3].toSI(RinexUtils.parseDouble(line,  4, 19)));
722                              pi.klobuchar.setBetaI(0, IonosphereKlobucharMessage.S_PER_SC_N[0].toSI(RinexUtils.parseDouble(line, 23, 19)));
723                              pi.klobuchar.setBetaI(1, IonosphereKlobucharMessage.S_PER_SC_N[1].toSI(RinexUtils.parseDouble(line, 42, 19)));
724                              pi.klobuchar.setBetaI(2, IonosphereKlobucharMessage.S_PER_SC_N[2].toSI(RinexUtils.parseDouble(line, 61, 19)));
725                          },
726                          pi -> Collections.singleton(KLOBUCHAR_LINE_2)),
727 
728         /** Parser for ionosphere Klobuchar message model. */
729         KLOBUCHAR_LINE_0(line -> true,
730                          (line, pi) -> {
731                              final int year  = RinexUtils.parseInt(line, 4, 4);
732                              final int month = RinexUtils.parseInt(line, 9, 2);
733                              final int day   = RinexUtils.parseInt(line, 12, 2);
734                              final int hours = RinexUtils.parseInt(line, 15, 2);
735                              final int min   = RinexUtils.parseInt(line, 18, 2);
736                              final int sec   = RinexUtils.parseInt(line, 21, 2);
737                              pi.klobuchar.setTransmitTime(new AbsoluteDate(year, month, day, hours, min, sec,
738                                                                            pi.klobuchar.getSystem().getObservationTimeScale().getTimeScale(pi.timeScales)));
739                              pi.klobuchar.setAlphaI(0, IonosphereKlobucharMessage.S_PER_SC_N[0].toSI(RinexUtils.parseDouble(line, 23, 19)));
740                              pi.klobuchar.setAlphaI(1, IonosphereKlobucharMessage.S_PER_SC_N[1].toSI(RinexUtils.parseDouble(line, 42, 19)));
741                              pi.klobuchar.setAlphaI(2, IonosphereKlobucharMessage.S_PER_SC_N[2].toSI(RinexUtils.parseDouble(line, 61, 19)));
742                          },
743                          pi -> Collections.singleton(KLOBUCHAR_LINE_1)),
744 
745         /** Parser for ionosphere Nequick-G message model. */
746         NEQUICK_LINE_1(line -> true,
747                        (line, pi) -> {
748                            pi.nequickG.setFlags((int) FastMath.rint(RinexUtils.parseDouble(line, 4, 19)));
749                            pi.file.addNequickGMessage(pi.nequickG);
750                            pi.nequickG = null;
751                        },
752                        LineParser::navigationNext),
753 
754         /** Parser for ionosphere Nequick-G message model. */
755         NEQUICK_LINE_0(line -> true,
756                        (line, pi) -> {
757                            final int year  = RinexUtils.parseInt(line, 4, 4);
758                            final int month = RinexUtils.parseInt(line, 9, 2);
759                            final int day   = RinexUtils.parseInt(line, 12, 2);
760                            final int hours = RinexUtils.parseInt(line, 15, 2);
761                            final int min   = RinexUtils.parseInt(line, 18, 2);
762                            final int sec   = RinexUtils.parseInt(line, 21, 2);
763                            pi.nequickG.setTransmitTime(new AbsoluteDate(year, month, day, hours, min, sec,
764                                                                         pi.nequickG.getSystem().getObservationTimeScale().getTimeScale(pi.timeScales)));
765                            pi.nequickG.setAi0(IonosphereNequickGMessage.SFU.toSI(RinexUtils.parseDouble(line, 23, 19)));
766                            pi.nequickG.setAi1(IonosphereNequickGMessage.SFU_PER_DEG.toSI(RinexUtils.parseDouble(line, 42, 19)));
767                            pi.nequickG.setAi2(IonosphereNequickGMessage.SFU_PER_DEG2.toSI(RinexUtils.parseDouble(line, 61, 19)));
768                        },
769                        pi -> Collections.singleton(NEQUICK_LINE_1)),
770 
771         /** Parser for ionosphere BDGIM message model. */
772         BDGIM_LINE_2(line -> true,
773                      (line, pi) -> {
774                          pi.bdgim.setAlphaI(7, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line,  4, 19)));
775                          pi.bdgim.setAlphaI(8, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 23, 19)));
776                          pi.file.addBDGIMMessage(pi.bdgim);
777                          pi.bdgim = null;
778                      },
779                      LineParser::navigationNext),
780 
781         /** Parser for ionosphere BDGIM message model. */
782         BDGIM_LINE_1(line -> true,
783                      (line, pi) -> {
784                          pi.bdgim.setAlphaI(3, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line,  4, 19)));
785                          pi.bdgim.setAlphaI(4, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 23, 19)));
786                          pi.bdgim.setAlphaI(5, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 42, 19)));
787                          pi.bdgim.setAlphaI(6, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 61, 19)));
788                      },
789                      pi -> Collections.singleton(BDGIM_LINE_2)),
790 
791         /** Parser for ionosphere BDGIM message model. */
792         BDGIM_LINE_0(line -> true,
793                      (line, pi) -> {
794                          final int year  = RinexUtils.parseInt(line, 4, 4);
795                          final int month = RinexUtils.parseInt(line, 9, 2);
796                          final int day   = RinexUtils.parseInt(line, 12, 2);
797                          final int hours = RinexUtils.parseInt(line, 15, 2);
798                          final int min   = RinexUtils.parseInt(line, 18, 2);
799                          final int sec   = RinexUtils.parseInt(line, 21, 2);
800                          pi.bdgim.setTransmitTime(new AbsoluteDate(year, month, day, hours, min, sec,
801                                                                    pi.bdgim.getSystem().getObservationTimeScale().getTimeScale(pi.timeScales)));
802                          pi.bdgim.setAlphaI(0, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 23, 19)));
803                          pi.bdgim.setAlphaI(1, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 42, 19)));
804                          pi.bdgim.setAlphaI(2, Unit.TOTAL_ELECTRON_CONTENT_UNIT.toSI(RinexUtils.parseDouble(line, 61, 19)));
805                      },
806                      pi -> Collections.singleton(BDGIM_LINE_1)),
807 
808         /** Parser for ionosphere message type. */
809         IONO_TYPE(line -> line.startsWith("> ION"),
810                   (line, pi) -> {
811                       pi.closePendingMessage();
812                       final SatelliteSystem system = SatelliteSystem.parseSatelliteSystem(RinexUtils.parseString(line, 6, 1));
813                       final int             prn    = RinexUtils.parseInt(line, 7, 2);
814                       final String          type   = RinexUtils.parseString(line, 10, 4);
815                       if (system == SatelliteSystem.GALILEO) {
816                           pi.nequickG = new IonosphereNequickGMessage(system, prn, type);
817                       } else {
818                           // in Rinex 4.00, tables A32 and A34 are ambiguous as both seem to apply
819                           // to Beidou CNVX messages, we consider BDGIM is the proper model in this case
820                           if (system == SatelliteSystem.BEIDOU && "CNVX".equals(type)) {
821                               pi.bdgim = new IonosphereBDGIMMessage(system, prn, type);
822                           } else {
823                               pi.klobuchar = new IonosphereKlobucharMessage(system, prn, type);
824                           }
825                       }
826                   },
827                   pi -> Collections.singleton(pi.nequickG != null ? NEQUICK_LINE_0 : (pi.bdgim != null ? BDGIM_LINE_0 : KLOBUCHAR_LINE_0)));
828 
829         /** Predicate for identifying lines that can be parsed. */
830         private final Predicate<String> canHandle;
831 
832         /** Parsing method. */
833         private final ParsingMethod parsingMethod;
834 
835         /** Provider for next line parsers. */
836         private final Function<ParseInfo, Iterable<LineParser>> allowedNextProvider;
837 
838         /** Simple constructor.
839          * @param canHandle predicate for identifying lines that can be parsed
840          * @param parsingMethod parsing method
841          * @param allowedNextProvider supplier for allowed parsers for next line
842          */
843         LineParser(final Predicate<String> canHandle, final ParsingMethod parsingMethod,
844                    final Function<ParseInfo, Iterable<LineParser>> allowedNextProvider) {
845             this.canHandle           = canHandle;
846             this.parsingMethod       = parsingMethod;
847             this.allowedNextProvider = allowedNextProvider;
848         }
849 
850         /** Get the allowed parsers for next lines while parsing Rinex header.
851          * @param parseInfo holder for transient data
852          * @return allowed parsers for next line
853          */
854         private static Iterable<LineParser> headerNext(final ParseInfo parseInfo) {
855             if (parseInfo.file.getHeader().getFormatVersion() < 3) {
856                 // Rinex 2.x header entries
857                 return Arrays.asList(HEADER_COMMENT, HEADER_PROGRAM,
858                                      HEADER_ION_ALPHA, HEADER_ION_BETA,
859                                      HEADER_DELTA_UTC, HEADER_CORR_SYSTEM_TIME,
860                                      HEADER_LEAP_SECONDS, HEADER_END);
861             } else if (parseInfo.file.getHeader().getFormatVersion() < 4) {
862                 // Rinex 3.x header entries
863                 return Arrays.asList(HEADER_COMMENT, HEADER_PROGRAM,
864                                      HEADER_IONOSPHERIC, HEADER_TIME,
865                                      HEADER_LEAP_SECONDS, HEADER_END);
866             } else {
867                 // Rinex 4.x header entries
868                 return Arrays.asList(HEADER_COMMENT, HEADER_PROGRAM,
869                                      HEADER_DOI, HEADER_LICENSE, HEADER_STATION_INFORMATION, HEADER_MERGED_FILE,
870                                      HEADER_LEAP_SECONDS, HEADER_END);
871             }
872         }
873 
874         /** Get the allowed parsers for next lines while parsing navigation date.
875          * @param parseInfo holder for transient data
876          * @return allowed parsers for next line
877          */
878         private static Iterable<LineParser> navigationNext(final ParseInfo parseInfo) {
879             if (parseInfo.gpsLNav    != null || parseInfo.gpsCNav    != null || parseInfo.galileoNav != null ||
880                 parseInfo.beidouLNav != null || parseInfo.beidouCNav != null || parseInfo.qzssLNav   != null ||
881                 parseInfo.qzssCNav   != null || parseInfo.navicLNav  != null || parseInfo.navicL1NV  != null ||
882                 parseInfo.sbasNav    != null) {
883                 return Collections.singleton(BROADCAST_ORBIT);
884             } else if (parseInfo.glonassNav != null) {
885                 if (parseInfo.messageLineNumber < 3) {
886                     return Collections.singleton(BROADCAST_ORBIT);
887                 } else {
888                     // workaround for some invalid files that should nevertheless be parsed
889                     // we have encountered in the wild merged files that claimed to be in 3.05 version
890                     // and hence needed at least 4 broadcast GLONASS orbit lines (the fourth line was
891                     // introduced in 3.05), but in fact only had 3 broadcast lines. We think they were
892                     // merged from files in 3.04 or earlier format. In order to parse these files,
893                     // we accept after the third line either another broadcast orbit line or a new message
894                     if (parseInfo.file.getHeader().getFormatVersion() < 4) {
895                         return Arrays.asList(BROADCAST_ORBIT, NAVIGATION_SV_EPOCH_CLOCK);
896                     } else {
897                         return Arrays.asList(BROADCAST_ORBIT, EPH_TYPE, STO_TYPE, EOP_TYPE, IONO_TYPE);
898                     }
899                 }
900             } else if (parseInfo.file.getHeader().getFormatVersion() < 3) {
901                 return Collections.singleton(NAVIGATION_SV_EPOCH_CLOCK_RINEX_2);
902             } else if (parseInfo.file.getHeader().getFormatVersion() < 4) {
903                 return Collections.singleton(NAVIGATION_SV_EPOCH_CLOCK);
904             } else {
905                 return Arrays.asList(EPH_TYPE, STO_TYPE, EOP_TYPE, IONO_TYPE);
906             }
907         }
908 
909     }
910 
911     /** Parsers for satellite system specific lines. */
912     private enum SatelliteSystemLineParser {
913 
914         /** GPS legacy. */
915         GPS_LNAV() {
916 
917             /** {@inheritDoc} */
918             @Override
919             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
920                 if (pi.file.getHeader().getFormatVersion() < 3.0) {
921                     parseSvEpochSvClockLineRinex2(line, pi.timeScales.getGPS(), pi.gpsLNav);
922                 } else {
923                     parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.gpsLNav);
924                 }
925             }
926 
927             /** {@inheritDoc} */
928             @Override
929             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
930                 pi.gpsLNav.setIODE(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
931                 pi.gpsLNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
932                 pi.gpsLNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
933                 pi.gpsLNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
934             }
935 
936             /** {@inheritDoc} */
937             @Override
938             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
939                 pi.gpsLNav.setCuc(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
940                 pi.gpsLNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
941                 pi.gpsLNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
942                 pi.gpsLNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
943             }
944 
945             /** {@inheritDoc} */
946             @Override
947             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
948                 pi.gpsLNav.setTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
949                 pi.gpsLNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
950                 pi.gpsLNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
951                 pi.gpsLNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
952             }
953 
954             /** {@inheritDoc} */
955             @Override
956             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
957                 pi.gpsLNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
958                 pi.gpsLNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
959                 pi.gpsLNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
960                 pi.gpsLNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
961             }
962 
963             /** {@inheritDoc} */
964             @Override
965             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
966                 // iDot
967                 pi.gpsLNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
968                 // Codes on L2 channel (ignored)
969                 // RinexUtils.parseDouble(line, 23, 19)
970                 // GPS week (to go with Toe)
971                 pi.gpsLNav.setWeek((int) RinexUtils.parseDouble(line, 42, 19));
972             }
973 
974             /** {@inheritDoc} */
975             @Override
976             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
977                 pi.gpsLNav.setSvAccuracy(parseBroadcastDouble1(line, pi.initialSpaces, Unit.METRE));
978                 pi.gpsLNav.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
979                 pi.gpsLNav.setTGD(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.SECOND));
980                 pi.gpsLNav.setIODC(parseBroadcastInt4(line,     pi.initialSpaces));
981             }
982 
983             /** {@inheritDoc} */
984             @Override
985             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
986                 pi.gpsLNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
987                 pi.gpsLNav.setFitInterval(parseBroadcastInt2(line, pi.initialSpaces));
988                 pi.closePendingMessage();
989             }
990 
991             /** {@inheritDoc} */
992             @Override
993             public void closeMessage(final ParseInfo pi) {
994                 pi.file.addGPSLegacyNavigationMessage(pi.gpsLNav);
995                 pi.gpsLNav = null;
996             }
997 
998         },
999 
1000         /** GPS civilian.
1001          * @since 12.0
1002          */
1003         GPS_CNAV() {
1004 
1005             /** {@inheritDoc} */
1006             @Override
1007             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1008                 parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.gpsCNav);
1009             }
1010 
1011             /** {@inheritDoc} */
1012             @Override
1013             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1014                 pi.gpsCNav.setADot(parseBroadcastDouble1(line, pi.initialSpaces, M_PER_S));
1015                 pi.gpsCNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
1016                 pi.gpsCNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1017                 pi.gpsCNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
1018             }
1019 
1020             /** {@inheritDoc} */
1021             @Override
1022             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1023                 pi.gpsCNav.setCuc(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.RADIAN));
1024                 pi.gpsCNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1025                 pi.gpsCNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1026                 pi.gpsCNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1027             }
1028 
1029             /** {@inheritDoc} */
1030             @Override
1031             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1032                 pi.gpsCNav.setTime(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.SECOND));
1033                 pi.gpsCNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1034                 pi.gpsCNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1035                 pi.gpsCNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1036             }
1037 
1038             /** {@inheritDoc} */
1039             @Override
1040             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1041                 pi.gpsCNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1042                 pi.gpsCNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1043                 pi.gpsCNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1044                 pi.gpsCNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1045             }
1046 
1047             /** {@inheritDoc} */
1048             @Override
1049             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1050                 pi.gpsCNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1051                 pi.gpsCNav.setDeltaN0Dot(parseBroadcastDouble2(line, pi.initialSpaces, RAD_PER_S2));
1052                 pi.gpsCNav.setUraiNed0(parseBroadcastInt3(line, pi.initialSpaces));
1053                 pi.gpsCNav.setUraiNed1(parseBroadcastInt4(line, pi.initialSpaces));
1054             }
1055 
1056             /** {@inheritDoc} */
1057             @Override
1058             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1059                 pi.gpsCNav.setUraiEd(parseBroadcastInt1(line, pi.initialSpaces));
1060                 pi.gpsCNav.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
1061                 pi.gpsCNav.setTGD(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1062                 pi.gpsCNav.setUraiNed2(parseBroadcastInt4(line, pi.initialSpaces));
1063             }
1064 
1065             /** {@inheritDoc} */
1066             @Override
1067             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1068                 pi.gpsCNav.setIscL1CA(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1069                 pi.gpsCNav.setIscL2C(parseBroadcastDouble2(line,  pi.initialSpaces, Unit.SECOND));
1070                 pi.gpsCNav.setIscL5I5(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1071                 pi.gpsCNav.setIscL5Q5(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1072             }
1073 
1074             /** {@inheritDoc} */
1075             @Override
1076             public void parseEighthBroadcastOrbit(final String line, final ParseInfo pi) {
1077                 if (pi.gpsCNav.isCnv2()) {
1078                     // in CNAV2 messages, there is an additional line for L1 CD and L1 CP inter signal delay
1079                     pi.gpsCNav.setIscL1CD(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1080                     pi.gpsCNav.setIscL1CP(parseBroadcastDouble2(line, pi.initialSpaces, Unit.SECOND));
1081                 } else {
1082                     parseTransmissionTimeLine(line, pi);
1083                 }
1084             }
1085 
1086             /** {@inheritDoc} */
1087             @Override
1088             public void parseNinthBroadcastOrbit(final String line, final ParseInfo pi) {
1089                 parseTransmissionTimeLine(line, pi);
1090             }
1091 
1092             /** Parse transmission time line.
1093              * @param line line to parse
1094              * @param pi holder for transient data
1095              */
1096             private void parseTransmissionTimeLine(final String line, final ParseInfo pi) {
1097                 pi.gpsCNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1098                 pi.gpsCNav.setWeek(parseBroadcastInt2(line, pi.initialSpaces));
1099                 pi.closePendingMessage();
1100             }
1101 
1102             /** {@inheritDoc} */
1103             @Override
1104             public void closeMessage(final ParseInfo pi) {
1105                 pi.file.addGPSCivilianNavigationMessage(pi.gpsCNav);
1106                 pi.gpsCNav = null;
1107             }
1108 
1109         },
1110 
1111         /** Galileo. */
1112         GALILEO() {
1113 
1114             /** {@inheritDoc} */
1115             @Override
1116             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1117                 parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.galileoNav);
1118             }
1119 
1120             /** {@inheritDoc} */
1121             @Override
1122             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1123                 pi.galileoNav.setIODNav(parseBroadcastInt1(line, pi.initialSpaces));
1124                 pi.galileoNav.setCrs(parseBroadcastDouble2(line,       pi.initialSpaces, Unit.METRE));
1125                 pi.galileoNav.setDeltaN0(parseBroadcastDouble3(line,    pi.initialSpaces, RAD_PER_S));
1126                 pi.galileoNav.setM0(parseBroadcastDouble4(line,        pi.initialSpaces, Unit.RADIAN));
1127             }
1128 
1129             /** {@inheritDoc} */
1130             @Override
1131             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1132                 pi.galileoNav.setCuc(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.RADIAN));
1133                 pi.galileoNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1134                 pi.galileoNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1135                 pi.galileoNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1136             }
1137 
1138             /** {@inheritDoc} */
1139             @Override
1140             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1141                 pi.galileoNav.setTime(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.SECOND));
1142                 pi.galileoNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1143                 pi.galileoNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1144                 pi.galileoNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1145             }
1146 
1147             /** {@inheritDoc} */
1148             @Override
1149             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1150                 pi.galileoNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1151                 pi.galileoNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1152                 pi.galileoNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1153                 pi.galileoNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1154             }
1155 
1156             /** {@inheritDoc} */
1157             @Override
1158             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1159                 // iDot
1160                 pi.galileoNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1161                 pi.galileoNav.setDataSource(parseBroadcastInt2(line, pi.initialSpaces));
1162                 // GAL week (to go with Toe)
1163                 pi.galileoNav.setWeek(parseBroadcastInt3(line, pi.initialSpaces));
1164             }
1165 
1166             /** {@inheritDoc} */
1167             @Override
1168             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1169                 pi.galileoNav.setSisa(parseBroadcastDouble1(line, pi.initialSpaces, Unit.METRE));
1170                 pi.galileoNav.setSvHealth(parseBroadcastDouble2(line, pi.initialSpaces, Unit.NONE));
1171                 pi.galileoNav.setBGDE1E5a(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1172                 pi.galileoNav.setBGDE5bE1(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1173             }
1174 
1175             /** {@inheritDoc} */
1176             @Override
1177             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1178                 pi.galileoNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1179                 pi.closePendingMessage();
1180             }
1181 
1182             /** {@inheritDoc} */
1183             @Override
1184             public void closeMessage(final ParseInfo pi) {
1185                 pi.file.addGalileoNavigationMessage(pi.galileoNav);
1186                 pi.galileoNav = null;
1187             }
1188 
1189         },
1190 
1191         /** Glonass. */
1192         GLONASS() {
1193 
1194             /** {@inheritDoc} */
1195             @Override
1196             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1197 
1198                 if (pi.file.getHeader().getFormatVersion() < 3.0) {
1199 
1200                     pi.glonassNav.setPRN(RinexUtils.parseInt(line, 0, 2));
1201 
1202                     // Toc
1203                     final int    year  = RinexUtils.convert2DigitsYear(RinexUtils.parseInt(line,  3, 2));
1204                     final int    month = RinexUtils.parseInt(line,  6, 2);
1205                     final int    day   = RinexUtils.parseInt(line,  9, 2);
1206                     final int    hours = RinexUtils.parseInt(line, 12, 2);
1207                     final int    min   = RinexUtils.parseInt(line, 15, 2);
1208                     final double sec   = RinexUtils.parseDouble(line, 17, 5);
1209                     pi.glonassNav.setEpochToc(new AbsoluteDate(year, month, day, hours, min, sec,
1210                                                                pi.timeScales.getUTC()));
1211 
1212                     // clock
1213                     pi.glonassNav.setTauN(-RinexUtils.parseDouble(line, 22, 19));
1214                     pi.glonassNav.setGammaN(RinexUtils.parseDouble(line, 41, 19));
1215                     pi.glonassNav.setTime(fmod(RinexUtils.parseDouble(line, 60, 19), Constants.JULIAN_DAY));
1216 
1217                     // Set the ephemeris epoch (same as time of clock epoch)
1218                     pi.glonassNav.setDate(pi.glonassNav.getEpochToc());
1219 
1220                 } else {
1221                     pi.glonassNav.setPRN(RinexUtils.parseInt(line, 1, 2));
1222 
1223                     // Toc
1224                     pi.glonassNav.setEpochToc(parsePrnSvEpochClock(line, pi.timeScales.getUTC()));
1225 
1226                     // clock
1227                     pi.glonassNav.setTauN(-RinexUtils.parseDouble(line, 23, 19));
1228                     pi.glonassNav.setGammaN(RinexUtils.parseDouble(line, 42, 19));
1229                     pi.glonassNav.setTime(fmod(RinexUtils.parseDouble(line, 61, 19), Constants.JULIAN_DAY));
1230 
1231                     // Set the ephemeris epoch (same as time of clock epoch)
1232                     pi.glonassNav.setDate(pi.glonassNav.getEpochToc());
1233                 }
1234 
1235             }
1236 
1237             /** {@inheritDoc} */
1238             @Override
1239             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1240                 pi.glonassNav.setX(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1241                 pi.glonassNav.setXDot(parseBroadcastDouble2(line,    pi.initialSpaces, KM_PER_S));
1242                 pi.glonassNav.setXDotDot(parseBroadcastDouble3(line, pi.initialSpaces, KM_PER_S2));
1243                 pi.glonassNav.setHealth(parseBroadcastDouble4(line,  pi.initialSpaces, Unit.NONE));
1244             }
1245 
1246             /** {@inheritDoc} */
1247             @Override
1248             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1249                 pi.glonassNav.setY(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1250                 pi.glonassNav.setYDot(parseBroadcastDouble2(line,            pi.initialSpaces, KM_PER_S));
1251                 pi.glonassNav.setYDotDot(parseBroadcastDouble3(line,         pi.initialSpaces, KM_PER_S2));
1252                 pi.glonassNav.setFrequencyNumber(parseBroadcastDouble4(line, pi.initialSpaces, Unit.NONE));
1253             }
1254 
1255             /** {@inheritDoc} */
1256             @Override
1257             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1258                 pi.glonassNav.setZ(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1259                 pi.glonassNav.setZDot(parseBroadcastDouble2(line,    pi.initialSpaces, KM_PER_S));
1260                 pi.glonassNav.setZDotDot(parseBroadcastDouble3(line, pi.initialSpaces, KM_PER_S2));
1261                 if (pi.file.getHeader().getFormatVersion() < 3.045) {
1262                     pi.closePendingMessage();
1263                 }
1264             }
1265 
1266             /** {@inheritDoc} */
1267             @Override
1268             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1269                 pi.glonassNav.setStatusFlags(parseBroadcastDouble1(line, pi.initialSpaces, Unit.NONE));
1270                 pi.glonassNav.setGroupDelayDifference(parseBroadcastDouble2(line, pi.initialSpaces, Unit.NONE));
1271                 pi.glonassNav.setURA(parseBroadcastDouble3(line,                  pi.initialSpaces, Unit.NONE));
1272                 pi.glonassNav.setHealthFlags(parseBroadcastDouble4(line,          pi.initialSpaces, Unit.NONE));
1273                 pi.closePendingMessage();
1274             }
1275 
1276             /** {@inheritDoc} */
1277             @Override
1278             public void closeMessage(final ParseInfo pi) {
1279                 pi.file.addGlonassNavigationMessage(pi.glonassNav);
1280                 pi.glonassNav = null;
1281             }
1282 
1283         },
1284 
1285         /** QZSS legacy. */
1286         QZSS_LNAV() {
1287 
1288             /** {@inheritDoc} */
1289             @Override
1290             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1291                 parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.qzssLNav);
1292             }
1293 
1294             /** {@inheritDoc} */
1295             @Override
1296             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1297                 pi.qzssLNav.setIODE(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1298                 pi.qzssLNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
1299                 pi.qzssLNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1300                 pi.qzssLNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
1301             }
1302 
1303             /** {@inheritDoc} */
1304             @Override
1305             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1306                 pi.qzssLNav.setCuc(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1307                 pi.qzssLNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1308                 pi.qzssLNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1309                 pi.qzssLNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1310             }
1311 
1312             /** {@inheritDoc} */
1313             @Override
1314             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1315                 pi.qzssLNav.setTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1316                 pi.qzssLNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1317                 pi.qzssLNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1318                 pi.qzssLNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1319             }
1320 
1321             /** {@inheritDoc} */
1322             @Override
1323             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1324                 pi.qzssLNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1325                 pi.qzssLNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1326                 pi.qzssLNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1327                 pi.qzssLNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1328             }
1329 
1330             /** {@inheritDoc} */
1331             @Override
1332             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1333                 // iDot
1334                 pi.qzssLNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1335                 // Codes on L2 channel (ignored)
1336                 // RinexUtils.parseDouble(line, 23, 19)
1337                 // GPS week (to go with Toe)
1338                 pi.qzssLNav.setWeek(parseBroadcastInt3(line, pi.initialSpaces));
1339             }
1340 
1341             /** {@inheritDoc} */
1342             @Override
1343             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1344                 pi.qzssLNav.setSvAccuracy(parseBroadcastDouble1(line,  pi.initialSpaces, Unit.METRE));
1345                 pi.qzssLNav.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
1346                 pi.qzssLNav.setTGD(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.SECOND));
1347                 pi.qzssLNav.setIODC(parseBroadcastInt4(line,     pi.initialSpaces));
1348             }
1349 
1350             /** {@inheritDoc} */
1351             @Override
1352             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1353                 pi.qzssLNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1354                 pi.qzssLNav.setFitInterval(parseBroadcastInt2(line, pi.initialSpaces));
1355                 pi.closePendingMessage();
1356             }
1357 
1358             /** {@inheritDoc} */
1359             @Override
1360             public void closeMessage(final ParseInfo pi) {
1361                 pi.file.addQZSSLegacyNavigationMessage(pi.qzssLNav);
1362                 pi.qzssLNav = null;
1363             }
1364 
1365         },
1366 
1367         /** QZSS civilian.
1368          * @since 12.0
1369          */
1370         QZSS_CNAV() {
1371 
1372             /** {@inheritDoc} */
1373             @Override
1374             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1375                 parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.qzssCNav);
1376             }
1377 
1378             /** {@inheritDoc} */
1379             @Override
1380             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1381                 pi.qzssCNav.setADot(parseBroadcastDouble1(line, pi.initialSpaces, M_PER_S));
1382                 pi.qzssCNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
1383                 pi.qzssCNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1384                 pi.qzssCNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
1385             }
1386 
1387             /** {@inheritDoc} */
1388             @Override
1389             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1390                 pi.qzssCNav.setCuc(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1391                 pi.qzssCNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1392                 pi.qzssCNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1393                 pi.qzssCNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1394             }
1395 
1396             /** {@inheritDoc} */
1397             @Override
1398             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1399                 pi.qzssCNav.setTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1400                 pi.qzssCNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1401                 pi.qzssCNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1402                 pi.qzssCNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1403             }
1404 
1405             /** {@inheritDoc} */
1406             @Override
1407             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1408                 pi.qzssCNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1409                 pi.qzssCNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1410                 pi.qzssCNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1411                 pi.qzssCNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1412             }
1413 
1414             /** {@inheritDoc} */
1415             @Override
1416             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1417                 pi.qzssCNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1418                 pi.qzssCNav.setDeltaN0Dot(parseBroadcastDouble2(line, pi.initialSpaces, RAD_PER_S2));
1419                 pi.qzssCNav.setUraiNed0(parseBroadcastInt3(line, pi.initialSpaces));
1420                 pi.qzssCNav.setUraiNed1(parseBroadcastInt4(line, pi.initialSpaces));
1421             }
1422 
1423             /** {@inheritDoc} */
1424             @Override
1425             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1426                 pi.qzssCNav.setUraiEd(parseBroadcastInt1(line, pi.initialSpaces));
1427                 pi.qzssCNav.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
1428                 pi.qzssCNav.setTGD(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1429                 pi.qzssCNav.setUraiNed2(parseBroadcastInt4(line, pi.initialSpaces));
1430             }
1431 
1432             /** {@inheritDoc} */
1433             @Override
1434             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1435                 pi.qzssCNav.setIscL1CA(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1436                 pi.qzssCNav.setIscL2C(parseBroadcastDouble2(line,  pi.initialSpaces, Unit.SECOND));
1437                 pi.qzssCNav.setIscL5I5(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1438                 pi.qzssCNav.setIscL5Q5(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1439             }
1440 
1441             /** {@inheritDoc} */
1442             @Override
1443             public void parseEighthBroadcastOrbit(final String line, final ParseInfo pi) {
1444                 if (pi.qzssCNav.isCnv2()) {
1445                     // in CNAV2 messages, there is an additional line for L1 CD and L1 CP inter signal delay
1446                     pi.qzssCNav.setIscL1CD(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1447                     pi.qzssCNav.setIscL1CP(parseBroadcastDouble2(line, pi.initialSpaces, Unit.SECOND));
1448                 } else {
1449                     parseTransmissionTimeLine(line, pi);
1450                 }
1451             }
1452 
1453             /** {@inheritDoc} */
1454             @Override
1455             public void parseNinthBroadcastOrbit(final String line, final ParseInfo pi) {
1456                 parseTransmissionTimeLine(line, pi);
1457             }
1458 
1459             /** Parse transmission time line.
1460              * @param line line to parse
1461              * @param pi holder for transient data
1462              */
1463             private void parseTransmissionTimeLine(final String line, final ParseInfo pi) {
1464                 pi.qzssCNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1465                 pi.qzssCNav.setWeek(parseBroadcastInt2(line, pi.initialSpaces));
1466                 pi.closePendingMessage();
1467             }
1468 
1469             /** {@inheritDoc} */
1470             @Override
1471             public void closeMessage(final ParseInfo pi) {
1472                 pi.file.addQZSSCivilianNavigationMessage(pi.qzssCNav);
1473                 pi.qzssCNav = null;
1474             }
1475 
1476         },
1477 
1478         /** Beidou legacy. */
1479         BEIDOU_D1_D2() {
1480 
1481             /** {@inheritDoc} */
1482             @Override
1483             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1484                 parseSvEpochSvClockLine(line, pi.timeScales.getBDT(), pi.beidouLNav);
1485             }
1486 
1487             /** {@inheritDoc} */
1488             @Override
1489             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1490                 pi.beidouLNav.setAODE(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1491                 pi.beidouLNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
1492                 pi.beidouLNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1493                 pi.beidouLNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
1494             }
1495 
1496             /** {@inheritDoc} */
1497             @Override
1498             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1499                 pi.beidouLNav.setCuc(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.RADIAN));
1500                 pi.beidouLNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1501                 pi.beidouLNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1502                 pi.beidouLNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1503             }
1504 
1505             /** {@inheritDoc} */
1506             @Override
1507             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1508                 pi.beidouLNav.setTime(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.SECOND));
1509                 pi.beidouLNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1510                 pi.beidouLNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1511                 pi.beidouLNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1512             }
1513 
1514             /** {@inheritDoc} */
1515             @Override
1516             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1517                 pi.beidouLNav.setI0(parseBroadcastDouble1(line,       pi.initialSpaces, Unit.RADIAN));
1518                 pi.beidouLNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1519                 pi.beidouLNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1520                 pi.beidouLNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1521             }
1522 
1523             /** {@inheritDoc} */
1524             @Override
1525             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1526                 // iDot
1527                 pi.beidouLNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1528                 // BDT week (to go with Toe)
1529                 pi.beidouLNav.setWeek(parseBroadcastInt3(line, pi.initialSpaces));
1530             }
1531 
1532             /** {@inheritDoc} */
1533             @Override
1534             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1535                 pi.beidouLNav.setSvAccuracy(parseBroadcastDouble1(line, pi.initialSpaces, Unit.METRE));
1536                 // TODO SatH1
1537                 pi.beidouLNav.setTGD1(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.SECOND));
1538                 pi.beidouLNav.setTGD2(parseBroadcastDouble4(line,       pi.initialSpaces, Unit.SECOND));
1539             }
1540 
1541             /** {@inheritDoc} */
1542             @Override
1543             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1544                 pi.beidouLNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1545                 pi.beidouLNav.setAODC(parseBroadcastDouble2(line,             pi.initialSpaces, Unit.SECOND));
1546                 pi.closePendingMessage();
1547             }
1548 
1549             /** {@inheritDoc} */
1550             @Override
1551             public void closeMessage(final ParseInfo pi) {
1552                 pi.file.addBeidouLegacyNavigationMessage(pi.beidouLNav);
1553                 pi.beidouLNav = null;
1554             }
1555 
1556         },
1557 
1558         /** Beidou-3 CNAV. */
1559         BEIDOU_CNV_123() {
1560 
1561             /** {@inheritDoc} */
1562             @Override
1563             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1564                 parseSvEpochSvClockLine(line, pi.timeScales.getBDT(), pi.beidouCNav);
1565             }
1566 
1567             /** {@inheritDoc} */
1568             @Override
1569             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1570                 pi.beidouCNav.setADot(parseBroadcastDouble1(line, pi.initialSpaces, M_PER_S));
1571                 pi.beidouCNav.setCrs(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.METRE));
1572                 pi.beidouCNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1573                 pi.beidouCNav.setM0(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.RADIAN));
1574             }
1575 
1576             /** {@inheritDoc} */
1577             @Override
1578             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1579                 pi.beidouCNav.setCuc(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1580                 pi.beidouCNav.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1581                 pi.beidouCNav.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1582                 pi.beidouCNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1583             }
1584 
1585             /** {@inheritDoc} */
1586             @Override
1587             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1588                 pi.beidouCNav.setTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1589                 pi.beidouCNav.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1590                 pi.beidouCNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1591                 pi.beidouCNav.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1592             }
1593 
1594             /** {@inheritDoc} */
1595             @Override
1596             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1597                 pi.beidouCNav.setI0(parseBroadcastDouble1(line,       pi.initialSpaces, Unit.RADIAN));
1598                 pi.beidouCNav.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1599                 pi.beidouCNav.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1600                 pi.beidouCNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1601             }
1602 
1603             /** {@inheritDoc} */
1604             @Override
1605             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1606                 pi.beidouCNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1607                 pi.beidouCNav.setDeltaN0Dot(parseBroadcastDouble2(line, pi.initialSpaces, RAD_PER_S2));
1608                 switch (parseBroadcastInt3(line, pi.initialSpaces)) {
1609                     case 0 :
1610                         pi.beidouCNav.setSatelliteType(BeidouSatelliteType.RESERVED);
1611                         break;
1612                     case 1 :
1613                         pi.beidouCNav.setSatelliteType(BeidouSatelliteType.GEO);
1614                         break;
1615                     case 2 :
1616                         pi.beidouCNav.setSatelliteType(BeidouSatelliteType.IGSO);
1617                         break;
1618                     case 3 :
1619                         pi.beidouCNav.setSatelliteType(BeidouSatelliteType.MEO);
1620                         break;
1621                     default:
1622                         throw new OrekitException(OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
1623                                                   pi.lineNumber, pi.name, line);
1624                 }
1625                 pi.beidouCNav.setTime(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.SECOND));
1626             }
1627 
1628             /** {@inheritDoc} */
1629             @Override
1630             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1631                 pi.beidouCNav.setSisaiOe(parseBroadcastInt1(line, pi.initialSpaces));
1632                 pi.beidouCNav.setSisaiOcb(parseBroadcastInt2(line, pi.initialSpaces));
1633                 pi.beidouCNav.setSisaiOc1(parseBroadcastInt3(line, pi.initialSpaces));
1634                 pi.beidouCNav.setSisaiOc2(parseBroadcastInt4(line, pi.initialSpaces));
1635             }
1636 
1637             /** {@inheritDoc} */
1638             @Override
1639             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1640                 if (pi.beidouCNav.getRadioWave().closeTo(PredefinedGnssSignal.B1C)) {
1641                     pi.beidouCNav.setIscB1CD(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1642                     // field 2 is spare
1643                     pi.beidouCNav.setTgdB1Cp(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1644                     pi.beidouCNav.setTgdB2ap(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1645                 } else if (pi.beidouCNav.getRadioWave().closeTo(PredefinedGnssSignal.B2A)) {
1646                     // field 1 is spare
1647                     pi.beidouCNav.setIscB2AD(parseBroadcastDouble2(line, pi.initialSpaces, Unit.SECOND));
1648                     pi.beidouCNav.setTgdB1Cp(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1649                     pi.beidouCNav.setTgdB2ap(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1650                 } else {
1651                     parseSismaiHealthIntegrity(line, pi);
1652                 }
1653             }
1654 
1655             /** {@inheritDoc} */
1656             @Override
1657             public void parseEighthBroadcastOrbit(final String line, final ParseInfo pi) {
1658                 if (pi.beidouCNav.getRadioWave().closeTo(PredefinedGnssSignal.B2B)) {
1659                     pi.beidouCNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1660                     pi.closePendingMessage();
1661                 } else {
1662                     parseSismaiHealthIntegrity(line, pi);
1663                 }
1664             }
1665 
1666             /** {@inheritDoc} */
1667             @Override
1668             public void parseNinthBroadcastOrbit(final String line, final ParseInfo pi) {
1669                 pi.beidouCNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1670                 // field 2 is spare
1671                 // field 3 is spare
1672                 pi.beidouCNav.setIODE(parseBroadcastInt4(line, pi.initialSpaces));
1673                 pi.closePendingMessage();
1674             }
1675 
1676             /** {@inheritDoc} */
1677             @Override
1678             public void closeMessage(final ParseInfo pi) {
1679                 pi.file.addBeidouCivilianNavigationMessage(pi.beidouCNav);
1680                 pi.beidouCNav = null;
1681             }
1682 
1683             /**
1684              * Parse the SISMAI/Health/integrity line.
1685              * @param line line to read
1686              * @param pi holder for transient data
1687              */
1688             private void parseSismaiHealthIntegrity(final String line, final ParseInfo pi) {
1689                 pi.beidouCNav.setSismai(parseBroadcastInt1(line, pi.initialSpaces));
1690                 pi.beidouCNav.setHealth(parseBroadcastInt2(line, pi.initialSpaces));
1691                 pi.beidouCNav.setIntegrityFlags(parseBroadcastInt3(line, pi.initialSpaces));
1692                 pi.beidouCNav.setIODC(parseBroadcastInt4(line, pi.initialSpaces));
1693             }
1694 
1695         },
1696 
1697         /** SBAS. */
1698         SBAS() {
1699 
1700             /** {@inheritDoc} */
1701             @Override
1702             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1703 
1704                 // parse PRN
1705                 pi.sbasNav.setPRN(RinexUtils.parseInt(line, 1, 2));
1706 
1707                 // Time scale (UTC for Rinex 3.01 and GPS for other RINEX versions)
1708                 final int       version100 = (int) FastMath.rint(pi.file.getHeader().getFormatVersion() * 100);
1709                 final TimeScale timeScale  = (version100 == 301) ? pi.timeScales.getUTC() : pi.timeScales.getGPS();
1710 
1711                 pi.sbasNav.setEpochToc(parsePrnSvEpochClock(line, timeScale));
1712                 pi.sbasNav.setAGf0(parseBroadcastDouble2(line, pi.initialSpaces, Unit.SECOND));
1713                 pi.sbasNav.setAGf1(parseBroadcastDouble3(line, pi.initialSpaces, S_PER_S));
1714                 pi.sbasNav.setTime(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1715 
1716                 // Set the ephemeris epoch (same as time of clock epoch)
1717                 pi.sbasNav.setDate(pi.sbasNav.getEpochToc());
1718 
1719             }
1720 
1721             /** {@inheritDoc} */
1722             @Override
1723             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1724                 pi.sbasNav.setX(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1725                 pi.sbasNav.setXDot(parseBroadcastDouble2(line,    pi.initialSpaces, KM_PER_S));
1726                 pi.sbasNav.setXDotDot(parseBroadcastDouble3(line, pi.initialSpaces, KM_PER_S2));
1727                 pi.sbasNav.setHealth(parseBroadcastDouble4(line,  pi.initialSpaces, Unit.NONE));
1728             }
1729 
1730             /** {@inheritDoc} */
1731             @Override
1732             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1733                 pi.sbasNav.setY(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1734                 pi.sbasNav.setYDot(parseBroadcastDouble2(line,    pi.initialSpaces, KM_PER_S));
1735                 pi.sbasNav.setYDotDot(parseBroadcastDouble3(line, pi.initialSpaces, KM_PER_S2));
1736                 pi.sbasNav.setURA(parseBroadcastDouble4(line,     pi.initialSpaces, Unit.NONE));
1737             }
1738 
1739             /** {@inheritDoc} */
1740             @Override
1741             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1742                 pi.sbasNav.setZ(parseBroadcastDouble1(line, pi.initialSpaces, KM));
1743                 pi.sbasNav.setZDot(parseBroadcastDouble2(line,    pi.initialSpaces, KM_PER_S));
1744                 pi.sbasNav.setZDotDot(parseBroadcastDouble3(line, pi.initialSpaces, KM_PER_S2));
1745                 pi.sbasNav.setIODN(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.NONE));
1746                 pi.closePendingMessage();
1747             }
1748 
1749             /** {@inheritDoc} */
1750             @Override
1751             public void closeMessage(final ParseInfo pi) {
1752                 pi.file.addSBASNavigationMessage(pi.sbasNav);
1753                 pi.sbasNav = null;
1754             }
1755 
1756         },
1757 
1758         /** NavIC. */
1759         NAVIC_LNAV() {
1760 
1761             /** {@inheritDoc} */
1762             @Override
1763             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1764                 parseSvEpochSvClockLine(line, pi.timeScales.getNavIC(), pi.navicLNav);
1765             }
1766 
1767             /** {@inheritDoc} */
1768             @Override
1769             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1770                 pi.navicLNav.setIODC(parseBroadcastInt1(line, pi.initialSpaces));
1771                 pi.navicLNav.setCrs(parseBroadcastDouble2(line, pi.initialSpaces, Unit.METRE));
1772                 pi.navicLNav.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1773                 pi.navicLNav.setM0(parseBroadcastDouble4(line, pi.initialSpaces, Unit.RADIAN));
1774             }
1775 
1776             /** {@inheritDoc} */
1777             @Override
1778             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1779                 pi.navicLNav.setCuc(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1780                 pi.navicLNav.setE(parseBroadcastDouble2(line, pi.initialSpaces, Unit.NONE));
1781                 pi.navicLNav.setCus(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1782                 pi.navicLNav.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1783             }
1784 
1785             /** {@inheritDoc} */
1786             @Override
1787             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1788                 pi.navicLNav.setTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1789                 pi.navicLNav.setCic(parseBroadcastDouble2(line, pi.initialSpaces, Unit.RADIAN));
1790                 pi.navicLNav.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1791                 pi.navicLNav.setCis(parseBroadcastDouble4(line, pi.initialSpaces, Unit.RADIAN));
1792             }
1793 
1794             /** {@inheritDoc} */
1795             @Override
1796             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1797                 pi.navicLNav.setI0(parseBroadcastDouble1(line, pi.initialSpaces, Unit.RADIAN));
1798                 pi.navicLNav.setCrc(parseBroadcastDouble2(line, pi.initialSpaces, Unit.METRE));
1799                 pi.navicLNav.setPa(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1800                 pi.navicLNav.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1801             }
1802 
1803             /** {@inheritDoc} */
1804             @Override
1805             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1806                 // iDot
1807                 pi.navicLNav.setIDot(parseBroadcastDouble1(line, pi.initialSpaces, RAD_PER_S));
1808                 // NavIC week (to go with Toe)
1809                 pi.navicLNav.setWeek(parseBroadcastInt3(line, pi.initialSpaces));
1810             }
1811 
1812             /** {@inheritDoc} */
1813             @Override
1814             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1815                 final int uraIndex = parseBroadcastInt1(line, pi.initialSpaces);
1816                 pi.navicLNav.setSvAccuracy(NAVIC_URA[FastMath.min(uraIndex, NAVIC_URA.length - 1)]);
1817                 pi.navicLNav.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
1818                 pi.navicLNav.setTGD(parseBroadcastDouble3(line, pi.initialSpaces, Unit.SECOND));
1819             }
1820 
1821             /** {@inheritDoc} */
1822             @Override
1823             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1824                 pi.navicLNav.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1825                 pi.closePendingMessage();
1826             }
1827 
1828             /** {@inheritDoc} */
1829             @Override
1830             public void closeMessage(final ParseInfo pi) {
1831                 pi.file.addNavICLegacyNavigationMessage(pi.navicLNav);
1832                 pi.navicLNav = null;
1833             }
1834 
1835         },
1836 
1837         /** NavIC L1NV.
1838          * @since 12.0
1839          */
1840         NAVIC_L1NV() {
1841 
1842             /** {@inheritDoc} */
1843             @Override
1844             public void parseSvEpochSvClockLine(final String line, final ParseInfo pi) {
1845                 parseSvEpochSvClockLine(line, pi.timeScales.getGPS(), pi.navicL1NV);
1846             }
1847 
1848             /** {@inheritDoc} */
1849             @Override
1850             public void parseFirstBroadcastOrbit(final String line, final ParseInfo pi) {
1851                 pi.navicL1NV.setADot(parseBroadcastDouble1(line,    pi.initialSpaces, M_PER_S));
1852                 pi.navicL1NV.setCrs(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.METRE));
1853                 pi.navicL1NV.setDeltaN0(parseBroadcastDouble3(line, pi.initialSpaces, RAD_PER_S));
1854                 pi.navicL1NV.setM0(parseBroadcastDouble4(line,      pi.initialSpaces, Unit.RADIAN));
1855             }
1856 
1857             /** {@inheritDoc} */
1858             @Override
1859             public void parseSecondBroadcastOrbit(final String line, final ParseInfo pi) {
1860                 pi.navicL1NV.setCuc(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.RADIAN));
1861                 pi.navicL1NV.setE(parseBroadcastDouble2(line,     pi.initialSpaces, Unit.NONE));
1862                 pi.navicL1NV.setCus(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.RADIAN));
1863                 pi.navicL1NV.setSqrtA(parseBroadcastDouble4(line, pi.initialSpaces, SQRT_M));
1864             }
1865 
1866             /** {@inheritDoc} */
1867             @Override
1868             public void parseThirdBroadcastOrbit(final String line, final ParseInfo pi) {
1869                 pi.navicL1NV.setTime(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.SECOND));
1870                 pi.navicL1NV.setCic(parseBroadcastDouble2(line,    pi.initialSpaces, Unit.RADIAN));
1871                 pi.navicL1NV.setOmega0(parseBroadcastDouble3(line, pi.initialSpaces, Unit.RADIAN));
1872                 pi.navicL1NV.setCis(parseBroadcastDouble4(line,    pi.initialSpaces, Unit.RADIAN));
1873             }
1874 
1875             /** {@inheritDoc} */
1876             @Override
1877             public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
1878                 pi.navicL1NV.setI0(parseBroadcastDouble1(line,       pi.initialSpaces, Unit.RADIAN));
1879                 pi.navicL1NV.setCrc(parseBroadcastDouble2(line,      pi.initialSpaces, Unit.METRE));
1880                 pi.navicL1NV.setPa(parseBroadcastDouble3(line,       pi.initialSpaces, Unit.RADIAN));
1881                 pi.navicL1NV.setOmegaDot(parseBroadcastDouble4(line, pi.initialSpaces, RAD_PER_S));
1882             }
1883 
1884             /** {@inheritDoc} */
1885             @Override
1886             public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
1887                 pi.navicL1NV.setIDot(parseBroadcastDouble1(line,             pi.initialSpaces, RAD_PER_S));
1888                 pi.navicL1NV.setDeltaN0Dot(parseBroadcastDouble2(line,       pi.initialSpaces, RAD_PER_S2));
1889                 pi.navicL1NV.setReferenceSignalFlag(parseBroadcastInt4(line, pi.initialSpaces));
1890             }
1891 
1892             /** {@inheritDoc} */
1893             @Override
1894             public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
1895                 final int uraIndex = parseBroadcastInt1(line, pi.initialSpaces);
1896                 pi.navicL1NV.setSvAccuracy(NAVIC_URA[FastMath.min(uraIndex, NAVIC_URA.length - 1)]);
1897                 pi.navicL1NV.setSvHealth(parseBroadcastInt2(line, pi.initialSpaces));
1898                 pi.navicL1NV.setTGD(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.SECOND));
1899                 pi.navicL1NV.setTGDSL5(parseBroadcastDouble4(line, pi.initialSpaces, Unit.SECOND));
1900             }
1901 
1902             /** {@inheritDoc} */
1903             @Override
1904             public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
1905                 pi.navicL1NV.setIscSL1P(parseBroadcastDouble1(line,   pi.initialSpaces, Unit.SECOND));
1906                 pi.navicL1NV.setIscL1DL1P(parseBroadcastDouble2(line, pi.initialSpaces, Unit.SECOND));
1907                 pi.navicL1NV.setIscL1PS(parseBroadcastDouble3(line,   pi.initialSpaces, Unit.SECOND));
1908                 pi.navicL1NV.setIscL1DS(parseBroadcastDouble4(line,   pi.initialSpaces, Unit.SECOND));
1909             }
1910 
1911             /** {@inheritDoc} */
1912             @Override
1913             public void parseEighthBroadcastOrbit(final String line, final ParseInfo pi) {
1914                 pi.navicL1NV.setTransmissionTime(parseBroadcastDouble1(line, pi.initialSpaces, Unit.SECOND));
1915                 pi.navicL1NV.setWeek(parseBroadcastInt2(line, pi.initialSpaces));
1916                 pi.closePendingMessage();
1917             }
1918 
1919             /** {@inheritDoc} */
1920             @Override
1921             public void closeMessage(final ParseInfo pi) {
1922                 pi.file.addNavICL1NVNavigationMessage(pi.navicL1NV);
1923                 pi.navicL1NV = null;
1924             }
1925 
1926         };
1927 
1928         /** Get the parse for navigation message.
1929          * @param system satellite system
1930          * @param type message type (null for Rinex 3.x)
1931          * @param parseInfo container for transient data
1932          * @param line line being parsed
1933          * @return the satellite system line parser
1934          */
1935         private static SatelliteSystemLineParser getParser(final SatelliteSystem system, final String type,
1936                                                            final ParseInfo parseInfo, final String line) {
1937             switch (system) {
1938                 case GPS :
1939                     if (type == null || type.equals(LegacyNavigationMessage.LNAV)) {
1940                         // in Rinex, week number is aligned to GPS week!
1941                         parseInfo.gpsLNav = new GPSLegacyNavigationMessage(parseInfo.timeScales,
1942                                                                            SatelliteSystem.GPS);
1943                         return GPS_LNAV;
1944                     } else if (type.equals(CivilianNavigationMessage.CNAV)) {
1945                         // in Rinex, week number is aligned to GPS week!
1946                         parseInfo.gpsCNav = new GPSCivilianNavigationMessage(false,
1947                                                                              parseInfo.timeScales,
1948                                                                              SatelliteSystem.GPS);
1949                         return GPS_CNAV;
1950                     } else if (type.equals(CivilianNavigationMessage.CNV2)) {
1951                         // in Rinex, week number is aligned to GPS week!
1952                         parseInfo.gpsCNav = new GPSCivilianNavigationMessage(true,
1953                                                                              parseInfo.timeScales,
1954                                                                              SatelliteSystem.GPS);
1955                         return GPS_CNAV;
1956                     }
1957                     break;
1958                 case GALILEO :
1959                     if (type == null || type.equals("INAV") || type.equals("FNAV")) {
1960                         // in Rinex, week number is aligned to GPS week!
1961                         parseInfo.galileoNav = new GalileoNavigationMessage(parseInfo.timeScales,
1962                                                                             SatelliteSystem.GPS);
1963                         return GALILEO;
1964                     }
1965                     break;
1966                 case GLONASS :
1967                     if (type == null || type.equals("FDMA")) {
1968                         parseInfo.glonassNav = new GLONASSNavigationMessage();
1969                         return GLONASS;
1970                     }
1971                     break;
1972                 case QZSS :
1973                     if (type == null || type.equals(LegacyNavigationMessage.LNAV)) {
1974                         // in Rinex, week number is aligned to GPS week!
1975                         parseInfo.qzssLNav = new QZSSLegacyNavigationMessage(parseInfo.timeScales,
1976                                                                              SatelliteSystem.GPS);
1977                         return QZSS_LNAV;
1978                     } else if (type.equals(CivilianNavigationMessage.CNAV)) {
1979                         // in Rinex, week number is aligned to GPS week!
1980                         parseInfo.qzssCNav = new QZSSCivilianNavigationMessage(false,
1981                                                                                parseInfo.timeScales,
1982                                                                                SatelliteSystem.GPS);
1983                         return QZSS_CNAV;
1984                     } else if (type.equals(CivilianNavigationMessage.CNV2)) {
1985                         // in Rinex, week number is aligned to GPS week!
1986                         parseInfo.qzssCNav = new QZSSCivilianNavigationMessage(true,
1987                                                                                parseInfo.timeScales,
1988                                                                                SatelliteSystem.GPS);
1989                         return QZSS_CNAV;
1990                     }
1991                     break;
1992                 case BEIDOU :
1993                     if (type == null ||
1994                         type.equals(BeidouLegacyNavigationMessage.D1) ||
1995                         type.equals(BeidouLegacyNavigationMessage.D2)) {
1996                         // in Rinex, week number for Beidou is really aligned to Beidou week!
1997                         parseInfo.beidouLNav = new BeidouLegacyNavigationMessage(parseInfo.timeScales,
1998                                                                                  SatelliteSystem.BEIDOU);
1999                         return BEIDOU_D1_D2;
2000                     } else if (type.equals(BeidouCivilianNavigationMessage.CNV1)) {
2001                         // in Rinex, week number for Beidou is really aligned to Beidou week!
2002                         parseInfo.beidouCNav = new BeidouCivilianNavigationMessage(PredefinedGnssSignal.B1C,
2003                                                                                    parseInfo.timeScales,
2004                                                                                    SatelliteSystem.BEIDOU);
2005                         return BEIDOU_CNV_123;
2006                     } else if (type.equals(BeidouCivilianNavigationMessage.CNV2)) {
2007                         // in Rinex, week number for Beidou is really aligned to Beidou week!
2008                         parseInfo.beidouCNav = new BeidouCivilianNavigationMessage(PredefinedGnssSignal.B2A,
2009                                                                                    parseInfo.timeScales,
2010                                                                                    SatelliteSystem.BEIDOU);
2011                         return BEIDOU_CNV_123;
2012                     } else if (type.equals(BeidouCivilianNavigationMessage.CNV3)) {
2013                         // in Rinex, week number for Beidou is really aligned to Beidou week!
2014                         parseInfo.beidouCNav = new BeidouCivilianNavigationMessage(PredefinedGnssSignal.B2B,
2015                                                                                    parseInfo.timeScales,
2016                                                                                    SatelliteSystem.BEIDOU);
2017                         return BEIDOU_CNV_123;
2018                     }
2019                     break;
2020                 case NAVIC:
2021                     if (type == null || type.equals("LNAV")) {
2022                         // in Rinex, week number is aligned to GPS week!
2023                         parseInfo.navicLNav = new NavICLegacyNavigationMessage(parseInfo.timeScales,
2024                                                                                SatelliteSystem.GPS);
2025                         return NAVIC_LNAV;
2026                     } else if (type.equals(CivilianNavigationMessage.L1NV)) {
2027                         // in Rinex, week number is aligned to GPS week!
2028                         parseInfo.navicL1NV = new NavICL1NVNavigationMessage(parseInfo.timeScales,
2029                                                                              SatelliteSystem.GPS);
2030                         return NAVIC_L1NV;
2031                     }
2032                     break;
2033                 case SBAS :
2034                     if (type == null || type.equals("SBAS")) {
2035                         parseInfo.sbasNav = new SBASNavigationMessage();
2036                         return SBAS;
2037                     }
2038                     break;
2039                 default:
2040                     // do nothing, handle error after the switch
2041             }
2042             throw new OrekitException(OrekitMessages.UNABLE_TO_PARSE_LINE_IN_FILE,
2043                                       parseInfo.lineNumber, parseInfo.name, line);
2044         }
2045 
2046         /**
2047          * Parse the SV/Epoch/Sv clock of the navigation message.
2048          * @param line line to read
2049          * @param timeScale time scale to use
2050          * @param message navigation message
2051          */
2052         protected void parseSvEpochSvClockLineRinex2(final String line, final TimeScale timeScale,
2053                                                      final AbstractNavigationMessage<?> message) {
2054             // PRN
2055             message.setPRN(RinexUtils.parseInt(line, 0, 2));
2056 
2057             // Toc
2058             final int    year  = RinexUtils.convert2DigitsYear(RinexUtils.parseInt(line,  2, 3));
2059             final int    month = RinexUtils.parseInt(line,  5, 3);
2060             final int    day   = RinexUtils.parseInt(line,  8, 3);
2061             final int    hours = RinexUtils.parseInt(line, 11, 3);
2062             final int    min   = RinexUtils.parseInt(line, 14, 3);
2063             final double sec   = RinexUtils.parseDouble(line, 17, 5);
2064             message.setEpochToc(new AbsoluteDate( year, month, day, hours, min, sec, timeScale));
2065 
2066             // clock
2067             message.setAf0(RinexUtils.parseDouble(line, 22, 19));
2068             message.setAf1(RinexUtils.parseDouble(line, 41, 19));
2069             message.setAf2(RinexUtils.parseDouble(line, 60, 19));
2070 
2071         }
2072 
2073         /**
2074          * Parse the SV/Epoch/Sv clock of the navigation message.
2075          * @param line line to read
2076          * @param timeScale time scale to use
2077          * @param message navigation message
2078          */
2079         protected void parseSvEpochSvClockLine(final String line, final TimeScale timeScale,
2080                                                final AbstractNavigationMessage<?> message) {
2081             // PRN
2082             message.setPRN(RinexUtils.parseInt(line, 1, 2));
2083 
2084             // Toc
2085             message.setEpochToc(parsePrnSvEpochClock(line, timeScale));
2086 
2087             // clock
2088             message.setAf0(RinexUtils.parseDouble(line, 23, 19));
2089             message.setAf1(RinexUtils.parseDouble(line, 42, 19));
2090             message.setAf2(RinexUtils.parseDouble(line, 61, 19));
2091 
2092         }
2093 
2094         /** Parse epoch field of a Sv/epoch/clock line.
2095          * @param line line to parse
2096          * @param timeScale time scale to use
2097          * @return parsed field
2098          */
2099         protected AbsoluteDate parsePrnSvEpochClock(final String line, final TimeScale timeScale) {
2100             final int year  = RinexUtils.parseInt(line, 4, 4);
2101             final int month = RinexUtils.parseInt(line, 9, 2);
2102             final int day   = RinexUtils.parseInt(line, 12, 2);
2103             final int hours = RinexUtils.parseInt(line, 15, 2);
2104             final int min   = RinexUtils.parseInt(line, 18, 2);
2105             final int sec   = RinexUtils.parseInt(line, 21, 2);
2106             return new AbsoluteDate(year, month, day, hours, min, sec, timeScale);
2107         }
2108 
2109         /** Parse double field 1 of a broadcast orbit line.
2110          * @param line line to parse
2111          * @param initialSpaces number of initial spaces in the line
2112          * @param unit unit to used for parsing the field
2113          * @return parsed field
2114          */
2115         protected double parseBroadcastDouble1(final String line, final int initialSpaces, final Unit unit) {
2116             return unit.toSI(RinexUtils.parseDouble(line, initialSpaces, 19));
2117         }
2118 
2119         /** Parse integer field 1 of a broadcast orbit line.
2120          * @param line line to parse
2121          * @param initialSpaces number of initial spaces in the line
2122          * @return parsed field
2123          */
2124         protected int parseBroadcastInt1(final String line, final int initialSpaces) {
2125             return (int) FastMath.rint(RinexUtils.parseDouble(line, initialSpaces, 19));
2126         }
2127 
2128         /** Parse double field 2 of a broadcast orbit line.
2129          * @param line line to parse
2130          * @param initialSpaces number of initial spaces in the line
2131          * @param unit unit to used for parsing the field
2132          * @return parsed field
2133          */
2134         protected double parseBroadcastDouble2(final String line, final int initialSpaces, final Unit unit) {
2135             return unit.toSI(RinexUtils.parseDouble(line, initialSpaces + 19, 19));
2136         }
2137 
2138         /** Parse integer field 2 of a broadcast orbit line.
2139          * @param line line to parse
2140          * @param initialSpaces number of initial spaces in the line
2141          * @return parsed field
2142          */
2143         protected int parseBroadcastInt2(final String line, final int initialSpaces) {
2144             return (int) FastMath.rint(RinexUtils.parseDouble(line, initialSpaces + 19, 19));
2145         }
2146 
2147         /** Parse double field 3 of a broadcast orbit line.
2148          * @param line line to parse
2149          * @param initialSpaces number of initial spaces in the line
2150          * @param unit unit to used for parsing the field
2151          * @return parsed field
2152          */
2153         protected double parseBroadcastDouble3(final String line, final int initialSpaces, final Unit unit) {
2154             return unit.toSI(RinexUtils.parseDouble(line, initialSpaces + 38, 19));
2155         }
2156 
2157         /** Parse integer field 3 of a broadcast orbit line.
2158          * @param line line to parse
2159          * @param initialSpaces number of initial spaces in the line
2160          * @return parsed field
2161          */
2162         protected int parseBroadcastInt3(final String line, final int initialSpaces) {
2163             return (int) FastMath.rint(RinexUtils.parseDouble(line, initialSpaces + 38, 19));
2164         }
2165 
2166         /** Parse double field 4 of a broadcast orbit line.
2167          * @param line line to parse
2168          * @param initialSpaces number of initial spaces in the line
2169          * @param unit unit to used for parsing the field
2170          * @return parsed field
2171          */
2172         protected double parseBroadcastDouble4(final String line, final int initialSpaces, final Unit unit) {
2173             return unit.toSI(RinexUtils.parseDouble(line, initialSpaces + 57, 19));
2174         }
2175 
2176         /** Parse integer field 4 of a broadcast orbit line.
2177          * @param line line to parse
2178          * @param initialSpaces number of initial spaces in the line
2179          * @return parsed field
2180          */
2181         protected int parseBroadcastInt4(final String line, final int initialSpaces) {
2182             return (int) FastMath.rint(RinexUtils.parseDouble(line, initialSpaces + 57, 19));
2183         }
2184 
2185         /**
2186          * Parse the SV/Epoch/Sv clock of the navigation message.
2187          * @param line line to read
2188          * @param pi holder for transient data
2189          */
2190         public abstract void parseSvEpochSvClockLine(String line, ParseInfo pi);
2191 
2192         /**
2193          * Parse the "BROADCASTORBIT - 1" line.
2194          * @param line line to read
2195          * @param pi holder for transient data
2196          */
2197         public abstract void parseFirstBroadcastOrbit(String line, ParseInfo pi);
2198 
2199         /**
2200          * Parse the "BROADCASTORBIT - 2" line.
2201          * @param line line to read
2202          * @param pi holder for transient data
2203          */
2204         public abstract void parseSecondBroadcastOrbit(String line, ParseInfo pi);
2205 
2206         /**
2207          * Parse the "BROADCASTORBIT - 3" line.
2208          * @param line line to read
2209          * @param pi holder for transient data
2210          */
2211         public abstract void parseThirdBroadcastOrbit(String line, ParseInfo pi);
2212 
2213         /**
2214          * Parse the "BROADCASTORBIT - 4" line.
2215          * @param line line to read
2216          * @param pi holder for transient data
2217          */
2218         public void parseFourthBroadcastOrbit(final String line, final ParseInfo pi) {
2219             // this should never be called (except by some tests that use reflection)
2220             throw new OrekitInternalError(null);
2221         }
2222 
2223         /**
2224          * Parse the "BROADCASTORBIT - 5" line.
2225          * @param line line to read
2226          * @param pi holder for transient data
2227          */
2228         public void parseFifthBroadcastOrbit(final String line, final ParseInfo pi) {
2229             // this should never be called (except by some tests that use reflection)
2230             throw new OrekitInternalError(null);
2231         }
2232 
2233         /**
2234          * Parse the "BROADCASTORBIT - 6" line.
2235          * @param line line to read
2236          * @param pi holder for transient data
2237          */
2238         public void parseSixthBroadcastOrbit(final String line, final ParseInfo pi) {
2239             // this should never be called (except by some tests that use reflection)
2240             throw new OrekitInternalError(null);
2241         }
2242 
2243         /**
2244          * Parse the "BROADCASTORBIT - 7" line.
2245          * @param line line to read
2246          * @param pi holder for transient data
2247          */
2248         public void parseSeventhBroadcastOrbit(final String line, final ParseInfo pi) {
2249             // this should never be called (except by some tests that use reflection)
2250             throw new OrekitInternalError(null);
2251         }
2252 
2253         /**
2254          * Parse the "BROADCASTORBIT - 8" line.
2255          * @param line line to read
2256          * @param pi holder for transient data
2257          */
2258         public void parseEighthBroadcastOrbit(final String line, final ParseInfo pi) {
2259             // this should never be called (except by some tests that use reflection)
2260             throw new OrekitInternalError(null);
2261         }
2262 
2263         /**
2264          * Parse the "BROADCASTORBIT - 9" line.
2265          * @param line line to read
2266          * @param pi holder for transient data
2267          */
2268         public void parseNinthBroadcastOrbit(final String line, final ParseInfo pi) {
2269             // this should never be called (except by some tests that use reflection)
2270             throw new OrekitInternalError(null);
2271         }
2272 
2273         /**
2274          * Close a message as last line was parsed.
2275          * @param pi holder for transient data
2276          */
2277         public abstract void closeMessage(ParseInfo pi);
2278 
2279         /**
2280          * Calculates the floating-point remainder of a / b.
2281          * <p>
2282          * fmod = a - x * b
2283          * where x = (int) a / b
2284          * </p>
2285          * @param a numerator
2286          * @param b denominator
2287          * @return the floating-point remainder of a / b
2288          */
2289         private static double fmod(final double a, final double b) {
2290             final double x = (int) (a / b);
2291             return a - x * b;
2292         }
2293 
2294     }
2295 
2296     /** Parsing method. */
2297     @FunctionalInterface
2298     private interface ParsingMethod {
2299         /** Parse a line.
2300          * @param line line to parse
2301          * @param parseInfo holder for transient data
2302          */
2303         void parse(String line, ParseInfo parseInfo);
2304     }
2305 
2306 }