API Reference
Types
PositionVelocityTime.SatelliteState — Type
SatelliteState{CP<:Real,D<:GNSSDecoderState,S<:AbstractGNSSSignal}Combines the GNSS decoder state with code and carrier phase measurements for a single satellite.
Fields
decoder::GNSSDecoderState: GNSS decoder state containing decoded navigation datasystem::AbstractGNSSSignal: GNSS system (e.g.,GPSL1CA(),GalileoE1B())code_phase::CP: Code phase measurement in chipscarrier_doppler: Carrier Doppler frequency in Hzcarrier_phase::CP: Carrier phase measurement in radians, matchingTracking.get_carrier_phase(default:0.0)
Constructors
SatelliteState(; decoder, system, code_phase, carrier_doppler, carrier_phase=0.0)
SatelliteState(decoder, system, sat_state)The second constructor extracts code phase, carrier Doppler, and carrier phase from a Tracking satellite state (Tracking.TrackedSat). Tracking.jl defines it (from Tracking 9.1 on), so it is available once Tracking is loaded; this package does not depend on Tracking.
PositionVelocityTime.PVTSolution — Type
PVTSolutionComplete Position, Velocity, and Time solution from GNSS measurements.
Fields
position::ECEF{Float64}: User position in ECEF coordinates (meters)velocity::ECEF{Float64}: User velocity in ECEF coordinates (m/s)course_over_ground::typeof(1.0°): Horizontal direction of travel (degrees), the azimuth of the velocity vector in the local East-North-Up frame atposition— measured clockwise from true North and wrapped to[0, 360)°, following the GNSS course-over-ground (COG) convention. Derived from the velocity alone, so this is the direction of motion, not vehicle heading (which a single-antenna receiver cannot observe). The vertical (Up) velocity component is ignored.0°when the horizontal velocity is zero (stationary or purely vertical), where course is undefined.time_correction::typeof(1.0m): Estimated receiver clock bias of the reference GNSS (meters). For a multi-GNSS solution this is the bias ofreference_system; other systems' biases aretime_correction + inter_system_biases[system].time::Union{TAITime, Nothing}: Estimated time as a TAI epoch (seeTAITime; with AstroTime loaded,TAIEpoch(time)converts it)relative_clock_drift::Float64: Relative receiver clock drift (dimensionless)dop::Union{DOP, Nothing}: Dilution of precision valuessats::Dictionary{Tuple{Symbol, Int}, SatInfo}: Maps(signal, PRN)to satellite info — position, transmit time, and the post-fit pseudorange and range-rate residuals (seeSatInfo). The signal tag (e.g.:GPSL1CA,:GalileoE1B; seeget_signal_id) keeps the same PRN apart both across constellations (GPS PRN 5 vs Galileo E05) and across signals of one constellation (a satellite tracked on GPS L1 C/A and L5 yields two entries sharing a PRN). Receiver-clock grouping is by time system, not signal — seereference_system.reference_system::Union{SupportedTimeSystem, Nothing}: GNSS time system (e.g.GNSSSignals.GPST(),GST()) thattimeandtime_correctionare referenced to.inter_system_biases::Dict{SupportedTimeSystem, typeof(1.0m)}: For each GNSS time system other thanreference_system, the offset of that system's time scale relative to the reference system's (meters) — the inter-system bias. This is a system / space-segment effect, not a receiver one (this receiver has no inter-system hardware bias): it is the GNSS system-time offset, so the Galileo entry equals−c · Δt_systems,Δt_systems = GST − GPST(the GGTO), and the BeiDou entry likewise carries the BDT steering residual (the BGTO's subject) — never the defined 14 s count offset, which is a convention already removed from the measurements (seecalc_time_scale_offsets). It is estimated directly from the geometry whenever observable (no broadcast error); the broadcast GGTO/BGTO is used to derive it only as a fallback (the offset-aided collapse), as that broadcast value may be erroneous. Reference-independent (the difference of two entries is the offset between those two systems); empty for a single-system solution.inter_frequency_biases::Dict{Symbol, InterFrequencyBias}: For each frequency band other than the reference band, the receiver inter-frequency bias relative to it — the differential hardware delay of that band's RF chain, estimated as an extra unknown when satellites are processed on more than one band. The key is the band (e.g.:L1,:L5; seeget_band_id) and the value is shared across all constellations on that band. EachInterFrequencyBiascarries both the bias (in metres) and the reference band it is measured against — the anchor whose delay is folded into the clock biases — so each reported bias's anchor is explicit rather than implicit. The reference is chosen per coverage component (seeband_ifb_layout), so a disconnected constellation can in principle yield several references. Empty for a single-band solution. A solution then needsn ≥ 3 + M + Bsatellites forMtime systems andBextra bands.
Immutable: a solution's fields cannot be reassigned once it is built. The allocation-free calc_pvt! returns a new solution too; what it reuses are the sats, inter_system_biases and inter_frequency_biases containers of the solution it is explicitly handed as its output, which the returned solution takes over. calc_pvt builds new containers every time and never modifies the prev_pvt it reads.
PositionVelocityTime.TAITimes.TAITime — Type
TAITime(second::Integer, fraction::Real)
TAITime(datetime::Dates.DateTime)An epoch on International Atomic Time (TAI): second whole seconds after J2000 (2000-01-01T12:00:00 TAI) plus fraction ∈ [0, 1) of a second — the layout of AstroTime's TAIEpoch, so the two convert without loss. A fraction outside [0, 1) is normalised into it.
TAITime(datetime) reads a DateTime as a TAI calendar label (millisecond resolution); Dates.DateTime(t) gives the TAI calendar label back, rounded to the millisecond. The difference of two TAITimes is in seconds, t + Δ shifts t by Δ seconds, and ≈ compares like AstroTime's epochs do.
Fields
second::Int64: whole TAI seconds since J2000.fraction::Float64: the fraction of the second, in[0, 1).
PositionVelocityTime.SatInfo — Type
SatInfoPer-satellite information attached to a PVTSolution (one entry per satellite used in the fix).
Fields
position::ECEF{Float64}: Satellite ECEF position at transmit time (metres).time::Float64: Satellite transmit time (system time of week, seconds).residual::typeof(1.0m): Post-fit least-squares pseudorange residual (metres) — the (atmosphere-corrected) measured minus the modeled pseudorange. A per-satellite fit-quality / outlier indicator.rate_residual::typeof(1.0m/s): Post-fit least-squares range-rate residual (metres per second) — the measured minus the modeled range rate of the carrier-Doppler velocity and clock-drift solve. The rate-domain counterpart ofresidual: it flags a satellite whose Doppler disagrees with the velocity fix (cycle slips, dynamics) independently of its pseudorange.
Both are measured minus modeled ("observed minus computed"), the orientation GNSS software reports observation residuals in — RTKLIB's rescode / resdop, and GNSS-SDR and PocketSDR through it. Note that rate_residual follows resdop in the geometric range-rate sense, positive while the satellite recedes; a receiver forming the same residual from its tracking loops' λ · carrier_doppler works in the opposite sign (see calc_user_velocity_and_clock_drift and calc_pseudorange_rate, which set that sign).
PositionVelocityTime.InterFrequencyBias — Type
InterFrequencyBiasA receiver inter-frequency bias attached to a PVTSolution: the differential hardware delay of one frequency band's RF chain, together with the band it is measured against.
Fields
value::typeof(1.0m): The bias (metres), relative toreference— how much longer this band's receiver chain delay is than the reference band's.reference::Symbol: The reference band (e.g.:L1; seeget_band_id) whose delay is folded into the clock biases and against whichvalueis measured. Chosen per coverage component (seeband_ifb_layout), so different biases in a disconnected solution can in principle carry different references.
With AstroTime loaded, a package extension converts between TAITime and AstroTime's TAIEpoch both ways — TAIEpoch(pvt.time) and TAITime(epoch) — exactly, as the two share their layout.
Signal Groups
One epoch's measurements reach calc_pvt as signal groups: a NamedTuple of PositionVelocityTime.SignalGroups, one per ranging signal. Deliberately the same names Tracking.jl uses, and deliberately unexported, so using Tracking, PositionVelocityTime does not clash and both types print fully qualified.
Grouping is what makes the solve type-stable: within a group every satellite shares one concrete SatelliteState type, so the per-satellite work of PositionVelocityTime.collect_measurements is statically dispatched. Behind that collection pass the solver sees one flat, parameter-free PositionVelocityTime.SatelliteMeasurement row per satellite and therefore compiles once for every constellation mix.
Groups are built where the satellites are tracked, not derived from a pooled vector at solve time — PositionVelocityTime.signal_groups builds a whole epoch's groups from a Tracking.TrackState. There is deliberately no conversion from a flat Vector{SatelliteState}; calc_pvt refuses one with an error saying what to build.
PositionVelocityTime.SignalGroup — Type
PositionVelocityTime.SignalGroup(signal, satellites)A per-epoch measurement snapshot of the satellites tracked on one ranging signal.
Deliberately the same name as Tracking.SignalGroup, and deliberately not exported: using Tracking, PositionVelocityTime does not clash, call sites read PositionVelocityTime.SignalGroup(...), and Julia prints the unexported type fully qualified so an error message distinguishes the two by itself. The correspondence is exact — Tracking's group holds the loop state of a set of satellites, this one holds one epoch's measurements of the same satellites, with the same group keys and, when built through Tracking.jl, the same dictionary keys.
Groups are the unit of type stability: within a group every satellite shares one concrete SatelliteState type, so the per-satellite calls of collect_measurements are statically dispatched. Pooling several groups into one vector — what calc_pvt used to take — is what made them dynamic.
Fields
signal::AbstractGNSSSignal: the ranging signal of every satellite in the group. Carried on the group itself rather than read off a satellite, so the group's signal is known even when it holds none.satellites: the satellites, as aDictionary{Int,<:SatelliteState}keyed by PRN (mirroring Tracking, and structurally preventing a duplicate PRN within the group) or as anAbstractVector{<:SatelliteState}(so a consumer can reuse one buffer per group). Iteration order is significant — seecalc_pvt.
Every satellite's own state.system must be signal. This is not checked per epoch; a group built by Tracking.jl (signal_groups) satisfies it by construction.
PositionVelocityTime.SignalGroups — Type
PositionVelocityTime.SignalGroups{N}Type alias for the storage shape calc_pvt takes: a NamedTuple of N SignalGroups. Mirrors Tracking.SignalGroups, so a receiver can carry one set of group keys through tracking, decoding and the fix.
PositionVelocityTime.signal_groups — Function
PositionVelocityTime.signal_groups(track_state, decoders) -> NamedTuple of SignalGroupBuild one epoch's SignalGroups from a receiver's own per-signal state. The method that does this is defined by Tracking.jl (from Tracking 9.1 on) and takes a Tracking.TrackState together with that receiver's decoders; see it for the details.
There is deliberately no method taking a flat Vector{SatelliteState}. Groups are not a wrapper around the old input — they are the shape a receiver already has, and converting a pooled vector back into them at every epoch would reintroduce, in the conversion, exactly the dynamic dispatch the grouping exists to remove. Build the groups where the satellites are tracked, and hand the same groups to calc_pvt.
PositionVelocityTime.SatelliteMeasurement — Type
SatelliteMeasurementOne satellite's contribution to one epoch, reduced to a single concrete type with no type parameters — the row the whole PVT solver works on.
Everything that depends on the navigation-message family, the decoder or the ranging signal has already been evaluated by collect_measurements: the transmit time, the propagated orbit, the clock drift, the group-delay-corrected observables, the classification keys and the broadcast time offsets. What remains is plain data — numbers, Symbols and the time_system singleton — stored inline in a vector, so a mixed-constellation epoch is one Vector{SatelliteMeasurement} and the solver compiles once for every mix.
Fields
prn::Int: the satellite's PRN, unique only within its constellation.signal_id::Symbol:get_signal_idof the ranging signal (:GPSL1CA, …), the per-signal half of the(signal, PRN)key ofPVTSolution.sats.band_id::Symbol:get_band_idof the ranging signal (:L1,:L5, …), the grouping key of the receiver inter-frequency biases.time_system::SupportedTimeSystem:get_time_systemof the ranging signal (GPST(),GST(),BDT()), the grouping key of the receiver clock biases. A smallUnionof singletons (seeSupportedTimeSystem), so the row stays concrete and every call on the field splits statically.count_offset_to_gpst::Float64:time_scale_offset_to_gpstoftime_system, precomputed, so the solver reads a number rather than branching on the time system per satellite.time::Float64: the corrected transmit time (calc_corrected_time) as a seconds-of-week count on this satellite's own system's scale.position::SVector{3,Float64}: satellite ECEF position attime(m).velocity::SVector{3,Float64}: satellite ECEF velocity attime(m/s).carrier_doppler::Float64: the measured carrier Doppler (Hz, unitless).center_frequency::Float64: the ranging signal's carrier (Hz, unitless) — the Doppler wavelength, and the1/f²of every ionospheric model.clock_drift::Float64:calc_satellite_clock_driftattime(s/s).week::Int: the absolute week number oftime_system, with the GPS L1 C/A 1024-week rollover already resolved againstapproximate_year.system_start_epoch::TAITime: the TAI epoch of that system's week 0, which dates the reported fix.system_start_time::DateTime: the same origin as a calendar date, forday_of_year.time_offsets::NTuple{3,BroadcastTimeOffset}: this satellite's broadcast offsets towardCANDIDATE_HUB_SYSTEMS, in that order.
PositionVelocityTime.BroadcastTimeOffset — Type
BroadcastTimeOffsetOne constellation's broadcast offset toward one target GNSS time system, flattened out of the navigation message into an isbits row so that evaluating it costs no dispatch and no decoder field access. One per entry of CANDIDATE_HUB_SYSTEMS is carried by every SatelliteMeasurement.
The polynomial coefficients are stored rather than an evaluated offset because calc_hub_range_offsets evaluates one system-wide offset at each satellite's own transmit time; an evaluated scalar per satellite would change that (negligibly in metres, but visibly against the tolerances the steering tests pin).
Fields
available::Bool: whether this satellite's message carries a usable offset toward the target at all — the flattenedtime_offset_available. Every other field is meaningless whenfalse, and the coefficients are thenNaN, so thatcalc_steering_offsetof an unavailable offset isNaN.A_0,A_1,A_2::Float64: the broadcast polynomial's constant, rate and acceleration terms.t_0::Float64: the offset's reference time of week, orNaNwhere the message carries none (BeiDou D1/D2, whose two-term offset has no reference epoch — there the polynomial argument is the time of week itself).week_term::Float64:SECONDS_PER_WEEK * (own week − WN_0), precomputed, so the polynomial argument needs no week lookup.0.0wheret_0isNaN.count_anchor::Float64:time_scale_offset_to_gpst(own) − time_scale_offset_to_gpst(target), the defined whole-second partGNSSDecoder.GNSSTimeOffset.A_0folds in andcalc_steering_offsettakes back out (see there).
PositionVelocityTime.collect_measurements — Function
collect_measurements(groups; approximate_year = year(now(UTC)))
-> (Vector{SatelliteMeasurement}, ionospheric_correction)Flatten the SignalGroups of one epoch into the solver's input: the per-satellite SatelliteMeasurement rows, and the single constellation-wide ionospheric model select_ionospheric_correction would pick for them (folded into the same pass — both need exactly the satellites that survive the health gate).
A satellite is kept only if its full navigation-data set is decoded and it reports healthy; completeness is checked first, which is what guarantees the health bit has been decoded at all. Everything else about it — orbit, clock, observables, classification keys, broadcast time offsets — is evaluated here, so nothing downstream touches a decoder again.
This is the function barrier: it specialises on the group shape (and stays small), while the solver it feeds sees one concrete element type for every constellation mix. The flat row order is group order × within-group order, and it is significant — see calc_pvt.
approximate_year resolves the GPS L1 C/A 1024-week rollover; see calc_pvt.
PositionVelocityTime.collect_measurements! — Function
collect_measurements!(measurements::Vector{SatelliteMeasurement}, groups;
approximate_year = year(now(UTC)))
-> (measurements, ionospheric_correction)collect_measurements into an existing vector: measurements is emptied and refilled with this epoch's rows. What it held is overwritten.
The returned tuple pairs the (non-isbits) vector with a Union of coefficient sets, which Julia returns boxed: one small allocation per call. A consumer that must not allocate calls collect_measurement_rows!, which returns the ionospheric correction alone.
PositionVelocityTime.collect_measurement_rows! — Function
collect_measurement_rows!(measurements::Vector{SatelliteMeasurement}, groups;
approximate_year = year(now(UTC)))
-> ionospheric_correctioncollect_measurements! without the vector in the result: measurements is emptied and refilled with this epoch's rows, and only the ionospheric correction is returned. That is a Union of isbits coefficient sets, which Julia returns unboxed, so once measurements has held as many rows before, this allocates nothing. It is the collection pass calc_pvt! runs, for a consumer that runs its own estimator over the same rows in a loop that must not allocate.
PositionVelocityTime.CANDIDATE_HUB_SYSTEMS — Constant
PositionVelocityTime.CANDIDATE_HUB_SYSTEMSThe GNSS time systems a clock collapse can be anchored on, in the fixed order decide_bias_layout tries them: (GPST(), GST(), BDT()). Each SatelliteMeasurement's time_offsets tuple is indexed by this order, so a hub's broadcast offset is a tuple index rather than a dictionary lookup.
PositionVelocityTime.SupportedTimeSystem — Type
PositionVelocityTime.SupportedTimeSystemUnion{GPST,GST,BDT}: the GNSS time systems a satellite can be on, which are exactly CANDIDATE_HUB_SYSTEMS — every system GNSSSignals defines is a candidate hub. The type of SatelliteMeasurement's time_system field, and of the keys and the reference system of a PVTSolution.
A closed Union rather than the abstract GNSSSignals.TimeSystem so that every call on a time system — hashing it into inter_system_biases, time_scale_offset_to_gpst, a hub's offset index — is split into statically dispatched branches instead of being a dynamic dispatch, which a juliac --trim build rejects. The row it types stays one concrete type, so the solver still compiles once for every constellation mix. A time system GNSSSignals adds later has to be added to CANDIDATE_HUB_SYSTEMS anyway (the time_offsets tuple is indexed by it), and this follows.
PVT Computation
PositionVelocityTime.calc_pvt — Function
calc_pvt(groups, prev_pvt::PVTSolution = PVTSolution();
approximate_year::Integer = year(now(UTC)),
enable_ionospheric_correction::Bool = true,
enable_tropospheric_correction::Bool = true) -> PVTSolutionCalculate Position, Velocity, and Time (PVT) from GNSS satellite measurements.
Satellites from different GNSS and frequency bands may be combined. Because each constellation references its broadcasts to its own system time, one receiver clock bias is estimated per GNSS time system; and because the receiver's RF chain delays each band differently, one receiver inter-frequency bias is estimated per frequency band beyond a reference band (shared across constellations on that band; see get_band_id). The state vector is therefore [x, y, z, tc₁, …, tc_M, ifb₁, …, ifb_B] for M distinct time systems and B extra bands. Position and time are found by least squares; velocity and clock drift are solved from carrier Doppler.
A solution requires n ≥ 3 + M + B healthy satellite measurements (each system needs at least one satellite, and a system contributing a single satellite spends it entirely on that system's clock bias), of which 3 + M must come from distinct satellites: a satellite tracked on several bands supplies one line of sight, and its extra measurements constrain the inter-frequency biases rather than the geometry. When either condition fails, constellations whose messages carry a broadcast offset to another tracked system — Galileo's GGTO (Galileo–GPS Time Offset), BeiDou's BGTO (BDT–GNSS Time Offset, toward GPS or Galileo), or the GGTO GPS itself broadcasts on CNAV/CNAV-2 — have their clock bias collapsed onto that hub system using the broadcast offset, which makes a 4-satellite mixed fix possible (see decide_bias_layout for the hub order). Estimating an independent bias is preferred whenever the geometry allows it, since it avoids the broadcast offset's own error.
Unless disabled via enable_ionospheric_correction, the ionospheric delay is corrected automatically using only the coefficients decoded from the navigation messages. A single model is chosen for the whole solve and applied to every satellite, preferring the more accurate global TEC models: NTCM-G if Galileo Effective Ionisation Level coefficients have been decoded, else BDGIM if a BDS-3 B-CNAV (B1C/B2a/B2b) coefficient set has, else the GPS Klobuchar model, else BeiDou's own Klobuchar variant from its legacy B1I/B3I message, else no correction. See select_ionospheric_correction and ionospheric_delay.
Unless disabled via enable_tropospheric_correction, the tropospheric delay is corrected with a blind model (no broadcast coefficients needed): Saastamoinen zenith delays mapped to the line of sight by the Niell mapping functions. See tropospheric_delay.
Arguments
groups: this epoch's measurements, asSignalGroups— aNamedTupleofSignalGroups, one per ranging signal, each holding that signal'sSatelliteStates. A bareTupleof groups is accepted (and numberedgroup1,group2, …), as is a singleSignalGroup(which becomes the:defaultgroup), so a one-constellation solve stays a one-liner:using PositionVelocityTime: SignalGroup calc_pvt(SignalGroup(GPSL1CA(), gps_states)) calc_pvt((gps = SignalGroup(GPSL1CA(), gps_states), galileo = SignalGroup(GalileoE1B(), galileo_states)))Grouping is what makes the solve type-stable: within a group every satellite shares one concrete state type, so
collect_measurementsdispatches statically, and the solver behind it compiles once for every constellation mix. Build the groups where the satellites are tracked — with Tracking.jl 9.1 or later loaded,signal_groups(track_state, decoders)builds them from a wholeTrackState. A pooledVector{SatelliteState}, which is what this function took before 6.0, is refused with an error saying what to build instead.Each
(signal, PRN)pair must appear at most once — a receiver produces one measurement per signal per satellite, and a duplicate would enter the least-squares solve twice. ADictionary-backed group prevents this within itself; two groups that share a ranging signal can still collide, and are not checked for.Order is significant, as vector order was before it: the flat measurement order is group order × within-group order, and it fixes the primary system's tie-break ("most satellites, ties by first appearance"),
band_ifb_layout's reference-band tie-break, and the insertion order ofPVTSolution.sats.prev_pvt: Previous PVT solution used as initial guess (default: origin). A previous solution whose own DOP is implausible (GDOP < 0, orPDOPaboveMAX_TRUSTED_WARM_START_PDOP= 50) is not used as a seed and the epoch is solved from cold, so one spurious fix cannot re-seed itself epoch after epoch.
Keyword Arguments
approximate_year: Calendar year of the observation, used to resolve the GPS L1 C/A 1024-week rollover ambiguity (legacy LNAV broadcasts only a 10-bit week number, so the receiver needs external information to determine which 1024-week cycle the recording is in). Anything within ±9 years of the actual observation date works. Defaults to the current UTC year, which is correct for live signals; for processing archived recordings, pass the rough year of the recording.enable_ionospheric_correction: whentrue(default), apply the automatic ionospheric correction described above. Set tofalseto skip it entirely and solve from the raw pseudoranges (e.g. for diagnostics or when an external correction is applied elsewhere).enable_tropospheric_correction: whentrue(default), apply the Saastamoinen tropospheric correction. Set tofalseto skip it.
Returns
A new PVTSolution containing position, velocity, time, DOP values, and satellite information; prev_pvt is never modified. calc_pvt! is the same solve reusing the containers of a solution passed to it, without allocating. Returns prev_pvt if the epoch cannot be solved: too few healthy satellites to solve the constellation (including the GGTO fallback and the distinct-satellite condition — a satellite tracked on several bands supplies one line of sight, so measurements alone are not enough), a geometry whose solved design matrix is rank deficient (reported as a negative GDOP). None of these throw, so a receiver can pass whatever it currently tracks each epoch and carry prev_pvt forward. Distrust of prev_pvt (see above) affects only the seed: an unsolvable epoch still returns prev_pvt exactly as passed.
Allocation-free solving
calc_pvt! is the same solve, returning a new (immutable) PVTSolution that reuses the containers of the solution it is handed, with its scratch in a reusable PVTWorkspace. Once both have held an epoch of a given size, a solve allocates nothing. The overwrite is explicit: only the containers of the solution argument are reused, and passing the same solution as prev_pvt — the receiver loop, pvt = calc_pvt!(pvt, workspace, groups, pvt) — is allowed.
PositionVelocityTime.calc_pvt! — Function
calc_pvt!(solution::PVTSolution, workspace::PVTWorkspace, groups, prev_pvt::PVTSolution;
approximate_year::Integer = year(now(UTC)),
enable_ionospheric_correction::Bool = true,
enable_tropospheric_correction::Bool = true) -> PVTSolutioncalc_pvt, reusing the containers of solution and the buffers of workspace, so that a steady stream of epochs is solved without allocating. The solve is the same — same arguments, same keywords, same numbers — and so is the result it returns.
PVTSolution is immutable, so the fix is a new solution, returned; keep the return value:
pvt = PVTSolution()
workspace = PVTWorkspace()
for groups in epochs
pvt = calc_pvt!(pvt, workspace, groups, pvt)
endsolution is the output argument, and the only thing overwritten: the returned solution takes over its sats, inter_system_biases and inter_frequency_biases containers, which are emptied and refilled with this epoch's values. solution itself is therefore spent — its own fields still describe the epoch it was computed for, but its containers now hold the new one — so use the returned solution from here on, and keep a copy of solution if you need it after the call (a Dictionary made by map over solution.sats shares its keys, and sees the change too). prev_pvt is only read. Passing the same solution as both is allowed and is the usual receiver loop above: every value the seed needs is read before anything is written.
Where calc_pvt returns prev_pvt — an epoch that cannot be solved — so does this: prev_pvt itself when it is solution (or shares its containers), and otherwise a copy of it in solution's containers, so that prev_pvt stays untouched.
It allocates nothing once workspace and solution have held an epoch with at least as many satellites and as many estimated biases; a larger epoch grows them, once. That covers everything from the collection pass to the returned solution, but not the groups themselves, which the caller builds (and can reuse too — a SignalGroup's satellites may be any vector).
The workspace holds no state between epochs that affects the result — only buffers — but it must not be used by two solves at once; give each task its own.
Base.copy — Method
copy(pvt::PVTSolution) -> PVTSolutionA copy of pvt that shares no container with it: sats, inter_system_biases and inter_frequency_biases are copied, every other field is an immutable value. Their entries are immutable too, so the copy is independent of pvt — a later calc_pvt! into pvt's containers leaves it untouched. This is how to keep a solution that a calc_pvt! loop will reuse.
PositionVelocityTime.PVTWorkspace — Type
PVTWorkspace()The reusable scratch storage of calc_pvt!: the flat measurement rows, the bias layout, the pseudoranges and their corrections, the least-squares buffers, the design and normal-equations matrices, the velocity solve's design rows and the residuals. Every buffer grows to the largest epoch it has seen and is reused from then on; its contents after a solve are not part of any interface.
Satellite positions are kept as a Vector{SVector{3,Float64}} — a fixed-size quantity per satellite, in a vector that is resized in place. The least-squares state, whose length 3 + M + B depends on the epoch's time systems and bands, is a plain vector instead of a static one, so the solver compiles once rather than once per layout; the matrices sized by it are kept at their largest and used through views.
PositionVelocityTime.empty_keeping_capacity! — Function
empty_keeping_capacity!(dict::Dictionary) -> dict
empty_keeping_capacity!(solution::PVTSolution) -> solutionEmpty a Dictionary without giving up its storage, so refilling it with as many entries as before allocates nothing. empty!(::Dictionary) hands the dictionary a fresh, minimal hash table instead, which a refill then reallocates and regrows.
The PVTSolution form empties all three of a solution's containers that way — sats, inter_system_biases and inter_frequency_biases (the last two Base Dicts, whose own empty! already keeps their capacity) — for a consumer that builds its own solutions into the containers of one it keeps, as calc_pvt! does.
Position and Velocity
PositionVelocityTime.get_LLA — Function
get_LLA(pvt::PVTSolution) -> LLAConvert the ECEF position in the PVT solution to geodetic coordinates (latitude, longitude, altitude) using the WGS84 ellipsoid.
PositionVelocityTime.calc_satellite_position — Function
calc_satellite_position(decoder::GNSSDecoder.GNSSDecoderState, t)
calc_satellite_position(state::SatelliteState)Calculate the satellite ECEF position from orbital parameters at time t.
The first method takes a decoder state and explicit time. The second method computes the corrected transmission time from the SatelliteState automatically.
Arguments
decoder: GNSS decoder state containing ephemeris datat: Transmission time in system time (seconds)state: ASatelliteStatecombining decoder, system, and phase measurements
Returns
An SVector{3, Float64} with the satellite position in ECEF coordinates (meters).
PositionVelocityTime.calc_satellite_position_and_velocity — Function
calc_satellite_position_and_velocity(decoder::GNSSDecoder.GNSSDecoderState, t)
calc_satellite_position_and_velocity(state::SatelliteState)Calculate the satellite ECEF position and velocity from orbital parameters at time t.
Uses Keplerian orbital mechanics with perturbation corrections (harmonic corrections for argument of latitude, radius, and inclination) to propagate the satellite ephemeris.
Arguments
decoder: GNSS decoder state containing ephemeris datat: Transmission time in system time (seconds)state: ASatelliteStatecombining decoder, system, and phase measurements
Returns
A named tuple (position, velocity) where each is an SVector{3, Float64} in ECEF coordinates (meters and m/s respectively).
PositionVelocityTime.get_sat_enu — Function
get_sat_enu(user_pos_ecef::ECEF, sat_pos_ecef::ECEF) -> Spherical
get_sat_enu(enu_from_ecef::ENUfromECEF, sat_pos_ecef::ECEF) -> SphericalConvert satellite position to East-North-Up (ENU) spherical coordinates (azimuth and elevation) relative to the user position.
Arguments
user_pos_ecef::ECEF: User position in ECEF coordinatesenu_from_ecef::ENUfromECEF: A precomputedENUfromECEF(user_pos_ecef, wgs84)transform. Pass this form to reuse one transform across all satellites that share a user position (one geodetic conversion per epoch instead of per satellite).sat_pos_ecef::ECEF: Satellite position in ECEF coordinates
Returns
Spherical coordinates containing azimuth and elevation of the satellite as seen from the user position.
Atmospheric Corrections
These corrections are applied automatically by calc_pvt; they are documented here for reference and for diagnostic use.
PositionVelocityTime.select_ionospheric_correction — Function
select_ionospheric_correction(groups)
-> Union{KlobucharParams,BeiDouKlobucharParams,NTCMGParams,BDGIMParams,Nothing}Scan the satellite decoders of one epoch's SignalGroups and pick the single ionospheric correction to apply to the whole solve, in a fixed order of preference:
- NTCM-G, if Galileo Effective Ionisation Level coefficients were decoded.
- BDGIM, if a BDS-3 B-CNAV (B1C/B2a/B2b) coefficient set was decoded.
- Klobuchar from GPS, if the GPS α/β were decoded.
- BeiDou's own Klobuchar variant from its legacy D1/D2 message (B1I/B3I).
nothing— no coefficients at all, so no correction.
The two global TEC models come first because they are the better ones: BDGIM is a spherical-harmonic TEC expansion, whereas Klobuchar is the eight-coefficient half-cosine approximation IS-GPS-200 only undertakes to remove ~50 % of the delay with. The rungs among equals are tie-breaks rather than claims — NTCM-G stays ahead of BDGIM so a Galileo-bearing epoch behaves exactly as it did before, and GPS stays ahead of BeiDou among the Klobuchar sets as it already did — but every rung is fixed rather than data-dependent, so the chosen model does not flip with the order the receiver happens to hand satellites over in. The coefficients are global to a constellation, so the first decoder that carries each set is used — "first" in group order × within-group order, the same order everything else about the epoch is ordered in.
Every satellite of groups is scanned, healthy or not — the caller decides which satellites to offer, exactly as it decides which to solve with. calc_pvt does not call this at all: collect_measurements folds the same scan into its own pass over the satellites, which is the one place they are already being visited, and there the scan sees only the satellites that clear the health gate. This is the standalone form, for a consumer running its own estimator over the same measurement model.
PositionVelocityTime.ionospheric_delay — Function
ionospheric_delay(correction, center_frequency, elevation, azimuth, lla, time_of_week)
-> Float64Slant ionospheric group delay in metres for one satellite, using the constellation-wide correction returned by select_ionospheric_correction:
::Nothing→0.0(no coefficients were decoded).KlobucharParams→ Klobuchar model (IS-GPS-200, Fig. 20-4).BeiDouKlobucharParams→ BeiDou's Klobuchar variant (BDS-SIS-ICD-B1I-3.0 §5.2.4.7).NTCMGParams→ NTCM-G model.BDGIMParams→ BDGIM (BDS-SIS-ICD-B1C-1.0 §7.8.2).
The line of sight is given by the satellite elevation/azimuth (radians) and the user geodetic position lla (a Geodesy.LLA); time_of_week is the measurement's system time of week in seconds. Which constellation's scale that is does not matter: GPST and GST coincide to within the sub-microsecond GGTO, and BDT trails them by 14 s, which moves Klobuchar's diurnal phase by 0.016 % of a day, NTCM-G's universal time by the same, and BDGIM's solar longitude by 0.06 milliradians — orders of magnitude below any of the three models' own error. The geometry is taken precomputed — and shared across satellites and with tropospheric_delay — so a whole-epoch correction does the user geodetic conversion only once. Derive the geometry from ECEF with LLAfromECEF(wgs84)(user) and _elevation_azimuth(ENUfromECEF(user, wgs84), sat).
center_frequency is the ranging signal's carrier in Hz as a plain Float64 (ustrip(Hz, get_center_frequency(system)), precomputed on every SatelliteMeasurement) — every model needs the carrier and nothing else about the signal, which is why the row carries the frequency rather than the signal.
PositionVelocityTime.tropospheric_delay — Function
tropospheric_delay(elevation, lla, doy; humidity = 0.7) -> Float64Slant tropospheric group delay in metres for a satellite at elevation (radians) seen from the user geodetic position lla (a Geodesy.LLA): the Saastamoinen zenith delays (see saastamoinen_zenith_delays) mapped to the line of sight by the Niell mapping functions (see niell_mapping_functions), whose seasonal term is driven by the day of year doy. The delay is non-dispersive (frequency independent). humidity is the relative humidity (0…1) used for the wet component. The delay is finite and monotone in elevation all the way down to the horizon (≈ 90 m for a grazing ray at sea level); elevations below the horizon reuse the horizon value. Returns 0.0 when the user height is outside the model's valid range (−100 m … 10 km).
The geometry is taken precomputed so a whole-epoch correction shares the user geodetic conversion and the elevation with ionospheric_delay.
PositionVelocityTime.saastamoinen_zenith_delays — Function
saastamoinen_zenith_delays(latitude, height, humidity) -> (hydrostatic, wet)Saastamoinen zenith hydrostatic and wet tropospheric delays in metres, from a standard atmosphere (15 °C at sea level) evaluated at the geodetic latitude (radians) and height (metres, clamped below at 0). humidity is the relative humidity (0…1) feeding the water-vapour partial pressure of the wet term. The formulas are those of RTKLIB's tropmodel at zenith.
PositionVelocityTime.niell_mapping_functions — Function
niell_mapping_functions(latitude, height, elevation, doy) -> (hydrostatic, wet)The Niell (1996) hydrostatic and wet mapping functions: dimensionless factors that map the zenith tropospheric delays to a satellite at elevation (radians). latitude is in radians, height is the site's ellipsoidal height in metres (it enters the hydrostatic height correction; RTKLIB likewise substitutes the ellipsoidal height for the orthometric one), and doy is the day of year driving the hydrostatic seasonal term (phase-shifted half a year in the southern hemisphere). Elevations below the horizon are evaluated at the horizon, where both factors are finite (≈ 37 hydrostatic, ≈ 57 wet); both decrease strictly to exactly 1 at zenith.
PositionVelocityTime.KlobucharParams — Type
KlobucharParams(α_0, α_1, α_2, α_3, β_0, β_1, β_2, β_3)The eight Klobuchar ionospheric coefficients decoded from a GPS L1 navigation message, in IS-GPS-200 SI units (seconds and seconds·semicircle⁻ⁿ). The field names mirror GNSSDecoder.GPSL1CAData (α_0…α_3, β_0…β_3).
PositionVelocityTime.klobuchar_group_delay — Function
klobuchar_group_delay(φ_u, λ_u, E, A, gps_time, α, β)Klobuchar single-frequency ionospheric group delay for GPS L1 (IS-GPS-200N, Fig. 20-4), returned in seconds. All angles are in semicircles: φ_u/λ_u are the user geodetic latitude/longitude, E/A the satellite elevation/azimuth. gps_time is GPS system time in seconds. α/β are the 4-element Klobuchar coefficient tuples (SI units).
PositionVelocityTime.BeiDouKlobucharParams — Type
BeiDouKlobucharParams(α_0, α_1, α_2, α_3, β_0, β_1, β_2, β_3)The eight coefficients of BeiDou's own Klobuchar variant, decoded from the legacy D1/D2 message (B1I/B3I), in BDS-SIS-ICD-B1I-3.0 §5.2.4.7 SI units (seconds and seconds·π⁻ⁿ — per semicircle, exactly like the GPS set). A type of its own rather than a KlobucharParams because the ICD prescribes a different algorithm, not merely different numbers — see beidou_klobuchar_group_delay — and defines the resulting delay at the B1I carrier (1561.098 MHz) rather than L1.
PositionVelocityTime.beidou_klobuchar_group_delay — Function
beidou_klobuchar_group_delay(φ_u, λ_u, E, A, t_E, α, β)BeiDou's Klobuchar-variant ionospheric group delay along the B1I propagation path (BDS-SIS-ICD-B1I-3.0 §5.2.4.7), returned in seconds. All angles are in radians, the ICD's working unit: φ_u/λ_u are the user geodetic latitude/longitude, E/A the satellite elevation/azimuth. t_E is the BDT second of week. α/β are the 4-element coefficient tuples (SI units; their s·π⁻ⁿ scale factors make the polynomial argument φ_M/π — semicircles, exactly like the GPS set).
The same eight-coefficient half-cosine idea as klobuchar_group_delay, but not the same algorithm — the ICD replaces IS-GPS-200's fitted approximations with exact spherical geometry and drops the geomagnetic frame:
- The Earth-central angle to the ionospheric pierce point is the exact
π/2 − E − arcsin(R/(R+h)·cos E)at R = 6378 km, h = 375 km, and the obliquity the matching exact1/√(1 − (R/(R+h)·cos E)²), instead of the GPS model's0.0137/(E+0.11) − 0.022and1 + 16(0.53−E)³fits. - The latitude feeding the A₂/A₄ polynomials is the pierce point's geographic latitude from the spherical triangle — no
+0.064·cos(λ−1.617)geomagnetic conversion, no ±0.416-semicircle clamp — and its longitude usesarcsinrather than a flat-Earth division. - The period A₄ is clamped from above at 172800 s as well as from below at 72000 s, and the cosine is evaluated exactly rather than by the two-term Taylor expansion inside |x| < 1.57.
The 5·10⁻⁹ s night floor and the |t − 50400| < A₄/4 day-time window are those of the GPS model.
PositionVelocityTime.NTCMGParams — Type
NTCMGParams(a_i0, a_i1, a_i2, week_number::Integer)The broadcast Galileo Effective Ionisation Level coefficients a_i0/a_i1/a_i2 (sfu, sfu/deg, sfu/deg²) decoded from an E1B navigation message, together with the GST week number needed to derive the day of year and universal time for NTCM-G.
PositionVelocityTime.BDGIMParams — Type
BDGIMParams(α_1, …, α_9, week_number::Integer)The nine broadcast BDGIM coefficients decoded from a BDS-3 B-CNAV message (B1C/B2a/B2b), in the ICD's units of TECu, together with the BDT week number.
The names are the ICD's own α_1 … α_9 (BDS-SIS-ICD-B1C-1.0 Table 7-10); GNSSDecoder spells them α_bdgim_1 … α_bdgim_9 on its BeiDou containers to keep them apart from the legacy D1/D2 Klobuchar α_0 … α_3, and already applies each one's scale factor and sign — including α_5's negative scale factor of −2⁻³ — so they arrive here as plain signed TECu.
week_number is needed for the same reason NTCMGParams carries one: the model's time argument is a date (a Modified Julian Date), which a time of week alone cannot supply. Unlike the Klobuchar set there is no reference frequency, because BDGIM yields TEC rather than a delay at a particular carrier.
PositionVelocityTime.klobuchar_params — Function
klobuchar_params(decoder) -> Union{KlobucharParams,BeiDouKlobucharParams,Nothing}Klobuchar α/β decoded from a GPS navigation message (LNAV GPSL1CAData, CNAV GPSCNAVData on L5/L2C, or CNAV-2 GPSL1C_DData) or from the BeiDou legacy D1/D2 message (BeiDouDNAVData on B1I/B3I), or nothing if they have not been broadcast yet or the decoder carries no Klobuchar set. The same single-frequency Klobuchar model is broadcast on all GPS civil signals; BeiDou broadcasts its own set, returned as a BeiDouKlobucharParams because its ICD prescribes its own variant of the algorithm, not just its own numbers.
PositionVelocityTime.ntcm_g_params — Function
ntcm_g_params(decoder) -> Union{NTCMGParams,Nothing}NTCM-G Effective Ionisation Level coefficients decoded from a Galileo navigation message (I/NAV on E1-B or E5b, or F/NAV on E5a — all broadcast a_i0…a_i2), or nothing if they (or the week number) have not been decoded yet or the decoder carries no such set. Galileo E6-B is in the latter group: C/NAV broadcasts no ionospheric coefficients, so it falls to the generic method rather than being dispatched here on a field it does not have.
PositionVelocityTime.bdgim_params — Function
bdgim_params(decoder) -> Union{BDGIMParams,Nothing}The nine BDGIM coefficients decoded from a BDS-3 B-CNAV message — B-CNAV1 on B1C, B-CNAV2 on B2a, B-CNAV3 on B2b, all three broadcasting the same set — or nothing if they, or the week number the model needs to date its epoch, have not been decoded yet, or the decoder carries no such set. The legacy D1/D2 message (BeiDouDNAVData, on B1I/B3I) is in the latter group: it broadcasts Klobuchar instead, so it falls to this generic method rather than being dispatched on fields it does not have.
PositionVelocityTime._elevation_azimuth — Function
_elevation_azimuth(enu_from_ecef::ENUfromECEF, sat_position) -> (elevation, azimuth)Elevation and azimuth (radians) of sat_position (ECEF) in the local East-North-Up frame defined by enu_from_ecef = ENUfromECEF(user_position, wgs84). Azimuth is measured clockwise from North. The transform is taken precomputed so it can be built once per user position and reused across satellites.
Dilution of Precision
The DOP values are read from the dop field of a PVTSolution, e.g. pvt.dop.GDOP.
PositionVelocityTime.DOP — Type
DOPDilution of Precision (DOP) values describing the geometric quality of the satellite constellation used for a PVT solution.
Fields
GDOP::Float64: Geometric DOP (overall quality)PDOP::Float64: Position DOP (3D position quality)VDOP::Float64: Vertical DOPHDOP::Float64: Horizontal DOPTDOP::Float64: Time DOP
Utilities
PositionVelocityTime.get_sat_info — Function
get_sat_info(pvt_solution::PVTSolution, signal::Symbol, prn::Integer) -> Union{SatInfo,Nothing}Return the SatInfo (position, transmit time and the post-fit pseudorange/range-rate residuals) of the satellite with the given prn on GNSS signal (e.g. :GPSL1CA, :GalileoE1B; see get_signal_id), or nothing if that satellite was not used in the fix. The signal tag is required because the same PRN can belong to different constellations or be tracked on several signals; see the sats field of PVTSolution.
The Measurement-Model Surface
The pieces calc_pvt is assembled from, documented so a consumer that runs its own estimator over the same measurement model — a navigation filter closing tracking loops, for instance — reuses this package's model as a stated contract. The names are deliberately not exported (they are solver internals, not part of every user's vocabulary); bind them explicitly with using PositionVelocityTime: collect_measurements, ….
The surface is specified on the flat measurement row. A consumer's first call is PositionVelocityTime.collect_measurements, which turns its signal groups into a Vector{PositionVelocityTime.SatelliteMeasurement}; every function below that takes per-satellite input takes those rows, in that order, and returns one value per row in the same order. A row is a single concrete type with no type parameters, so a consumer's own estimator compiles once over it as well, however many constellations it mixes.
Only calc_corrected_time, calc_satellite_clock_drift, get_week, broadcast_time_offset and the decoder form of time_offset_available still read a decoder — they are what the collection pass is built from, and are listed so a consumer can build or amend a row itself.
Times returned or taken here are seconds-of-week counts on each satellite's own system scale unless a function says otherwise.
PositionVelocityTime.SPEED_OF_LIGHT — Constant
SPEED_OF_LIGHTThe speed of light (m/s), shared by every range/time conversion in this package and its consumers.
PositionVelocityTime.EARTH_ROTATION_RATE — Constant
EARTH_ROTATION_RATEThe Earth's rotation rate ω_e (rad/s), the WGS-84 value GPS and Galileo (GTRF) both use. It sets the Earth's turn over a signal's travel time, which carries the satellite from the ECEF frame of transmission into the one of reception.
PositionVelocityTime.calc_corrected_time — Function
calc_corrected_time(state::SatelliteState) -> Float64The satellite's transmit time as a seconds-of-week count on its own system's scale: the decoder's time of week plus the bit, code-phase and carrier-phase fractions (calc_uncorrected_time), corrected by the SV clock polynomial, the relativistic term and the ranging signal's group delay (correct_clock). This is the t every ephemeris and clock evaluation takes; only pseudorange differencing converts it to another scale (see time_scale_offset_to_gpst), never the value itself.
PositionVelocityTime.calc_satellite_clock_drift — Function
calc_satellite_clock_drift(decoder::GNSSDecoderState, t) -> Float64Rate of the satellite clock correction (s/s): the time derivative of the clock polynomial correct_clock applies, so it must be evaluated about the same clock reference epoch t_0c. Using t instead would offset the a_f2 term by a spurious 2·a_f2·t_0c — zero whenever a_f2 = 0 (the usual GPS broadcast) but up to ~4e-9 s/s (≈1.3 m/s of range rate) at the edge of the broadcast range late in the week, and per-satellite, so it would leak into the estimated velocity rather than the common clock drift.
The rate of the relativistic periodic term Δt_rel = F·e·√A·sin(E) that correct_clock includes is not modelled here; it is at most ~1 mm/s.
PositionVelocityTime.get_sat_position — Function
get_sat_position(sat_pv) -> SVector{3,Float64}Satellite ECEF position (m) of one calc_satellite_position_and_velocity result.
PositionVelocityTime.get_sat_velocity — Function
get_sat_velocity(sat_pv) -> SVector{3,Float64}Satellite ECEF velocity (m/s) of one calc_satellite_position_and_velocity result.
PositionVelocityTime.fold_week_crossover — Function
fold_week_crossover(t) -> Float64Fold a difference of two seconds-of-week counts modulo the week, into (-302400, 302400]: near a week rollover the two operands sit on opposite sides of the wrap and their raw difference is off by ±604800 s. Every time difference in this package — and in a consumer differencing transmit times — must pass through this.
PositionVelocityTime.BiasColumns — Type
BiasColumnsPer-satellite assignment of the estimated bias columns of the least-squares design matrix, shared by calc_ρ_hat!, calc_H! and user_position. The state vector is [x, y, z, tc₁, …, tc_num_clock_biases, ifb₁, …, ifb_num_ifb] with num_lsq_params(bias_columns) entries. The two column kinds have different physical sources: a clock column is the receiver clock for one GNSS time system (the spacing between systems is a system/space-segment effect — the broadcast GNSS time offset, GGTO or BGTO), whereas an inter-frequency-bias column is the receiver's per-band RF-chain delay. Known per-satellite corrections (atmosphere, the broadcast time-system offset) are applied to the pseudoranges in calc_pvt, not carried here.
Fields
clock_bias_indices::Vector{Int}: per satellite, the clock column (1…num_clock_biases) of its GNSS time system; the design-matrix1.0lands at3 + clock_bias_indices[j].num_clock_biases::Int: number of clock biases (also the offset of the IFB block).ifb_indices::Vector{Int}: per satellite, the inter-frequency-bias column (1…num_ifb) of its frequency band, or0for the reference band; the1.0lands at3 + num_clock_biases + ifb_indices[j].num_ifb::Int: number of inter-frequency biases (frequency bands beyond the reference).
PositionVelocityTime.decide_bias_layout — Function
decide_bias_layout(measurements) -> Union{BiasLayout,Nothing}Decide the full least-squares bias layout — one clock column per GNSS time system plus the per-band inter-frequency-bias columns from band_ifb_layout — for the flat SatelliteMeasurement rows of one epoch, and return it as a BiasLayout, or nothing when the constellation cannot be solved. Every classification the decision needs — each satellite's time_system, band_id, prn and time_offsets — is a field of the row:
bias_columns::BiasColumns: the per-satellite column assignment (seeBiasColumns).clock_systems::Vector{SupportedTimeSystem}: the time system of each clock column, soclock_systems[k]belongs to clock columnk, in order of first appearance. A collapsed system has no column of its own; its satellites share the hub's.extra_bands::Vector{Symbol}: band of each inter-frequency-bias column, soextra_bands[i]belongs to columniof the IFB block.reference_bands::Vector{Symbol}: per IFB column, the reference band of its coverage component — the anchor that column's bias is measured against (seeband_ifb_layout).hub_system::Union{Nothing,SupportedTimeSystem}: the system the collapsed clocks were merged onto, ornothingfor a layout that estimates every clock bias independently, which is the common case.hub_rows::Vector{SatelliteMeasurement}: one row per collapsed time system — the satellite whose broadcast offset tohub_systemconverts that whole system's measurements. Each row carries its owntime_system, so it is its own key. The same vocabulary astime_offset_available, which decides membership, andcalc_steering_offset, which evaluates one entry (seecalc_hub_range_offsets). Empty exactly whenhub_systemisnothing.
Only the classification fields of each row are read, never its transmit time: the counts, the coverage graph and the availability flags are all the decision needs. The offsets its hub_rows enable are evaluated afterwards, at each satellite's own transmit time.
The decision is observability-driven, not merely count-driven:
- When the coverage graph is connected and the satellites suffice for the unknowns (
n ≥ 3 + num_systems + num_ifbmeasurements, of which3 + num_systemsmust come from distinct satellites — extra bands of an already-tracked satellite add inter-frequency-bias information, not geometry), estimate everything independently: the inter-system offset and the receiver inter-frequency biases are observed directly from the geometry, so neither inherits the broadcast-GGTO error (the satellite group delays are already removed per satellite upstream, so the per-band column carries the receiver chain). - Otherwise merge every clock that can be merged onto a hub system, using a broadcast offset to that hub (see
time_offset_available). Hubs are tried in the fixed order GPST, GST, BDT, and the first whose merged layout has enough satellites wins, so the choice is deterministic and a GPS-anchored collapse behaves exactly as it always did. Toward GPS Time the offsets on the air are Galileo's GGTO and BeiDou's BGTO; toward Galileo System Time they are BeiDou's BGTO variant and the GGTO GPS itself broadcasts on CNAV/CNAV-2 — which is what lets a GPS-free Galileo + BeiDou epoch, or a GPS + Galileo epoch whose Galileo satellites have not decoded their GGTO yet, still collapse. The merge removes a clock unknown per merged system (the scarce-satellite case) and reconnects a disjoint band split (the disconnected case — where a band's IFB column would otherwise be collinear with the stranded constellation's clock), making the inter-frequency bias observable again (it then carries the broadcast-offset error, alongside the offset-based inter-system bias). A system whose satellites carry no offset to the hub keeps its own clock column, so a mixed epoch can collapse Galileo and leave BeiDou independent, or the reverse. - Failing that, fall back to the (already observability-restricted) independent layout if the satellite count allows. No IFB column is created for a band stranded on its own constellation, so its inter-frequency bias folds into that constellation's clock and the reported inter-system bias reads out as
GGTO + (IFB difference)— the two are no longer separable. Position, residuals and DOP are still finite and non-degenerate; only the bias decomposition is ambiguous. Else returnnothing.
Because band_ifb_layout never creates an unobservable IFB column, every returned layout is structurally sound — the degenerate disjoint-band case is removed by construction, not caught after the fact. The satellite conditions above are structural too, and thus necessary but not sufficient: a returned layout can still have degenerate geometry (lines of sight that span too little), which no satellite count can see. That is left to the checks calc_pvt makes on the solved geometry — the DOP's positive-definite test and the velocity solve's own — rather than pre-screened.
PositionVelocityTime.decide_bias_layout! — Function
decide_bias_layout!(workspace::BiasLayoutWorkspace, measurements) -> Booldecide_bias_layout into workspace: true and the layout stored there (read it with bias_layout), or false when the constellation cannot be solved. The same decision, allocating nothing once the workspace has grown to the epoch.
PositionVelocityTime.BiasLayoutWorkspace — Type
BiasLayoutWorkspace()The storage decide_bias_layout! decides a layout into: the layout's own vectors, which bias_layout hands out as a BiasLayout, and the scratch the decision works in. Every vector is emptied and refilled per epoch, so a workspace reused across epochs allocates nothing once it has grown to the largest one.
PositionVelocityTime.bias_layout — Function
bias_layout(workspace::BiasLayoutWorkspace) -> BiasLayoutThe layout decide_bias_layout! last decided into workspace. Its vectors are the workspace's, so they are overwritten by the next decision.
PositionVelocityTime.BiasLayout — Type
BiasLayoutThe full least-squares bias layout of one epoch, as decide_bias_layout returns it: a NamedTuple type spelled out so the decision has a concrete return type. That matters beyond tidiness — the layout is destructured straight into the solver, which is compiled once for every constellation mix, so an unresolved field type here would become a dynamic dispatch on every epoch. See decide_bias_layout for what each field means.
PositionVelocityTime.select_primary_system — Function
select_primary_system(measurements, hub_system) -> systemThe primary time system of an epoch: the one whose clock bias, reference time, week and start epoch define the reported time (PVTSolution's reference_system).
measurements: the epoch's rows.hub_system: the hub of a collapse (BiasLayout), ornothingif there is none.
A collapse is anchored on its hub, so the hub is primary whenever there is one. Otherwise the primary is the system with the most rows — the best-conditioned reported time — with ties going to the one that appears first. Throws an ArgumentError if there is neither a hub nor a row. Allocates nothing.
PositionVelocityTime.calc_ρ_hat! — Function
Computes ̂ρ, the distance between the satellite and the estimated user position
calc_ρ_hat!(ρ, sat_positions, ξ, bias_columns)
ξ: Estimated user position, one clock correction per GNSS time system, and one receiver inter-frequency bias per frequency band beyond the reference band, i.e. [x, y, z, tc₁, …, tc_num_clock_biases, ifb₁, …, ifb_num_ifb] sat_positions: Satellite positions, a 3 × N matrix or a vector of SVector{3}s bias_columns: The per-satellite BiasColumns (clock column and inter-frequency-bias column of each satellite). Known range corrections (atmosphere, GGTO offset) are already folded into the measured ρ by calc_pvt.
PositionVelocityTime.calc_H — Function
calc_H(sat_positions, ξ, bias_columns::BiasColumns) -> Matrix{Float64}Allocating form of calc_H!: the least-squares design matrix at the state ξ, one row per satellite, with num_lsq_params(bias_columns) columns.
PositionVelocityTime.calc_DOP — Function
Calculates the dilution of precision for a given geometry matrix H
H_GEO has 3 + num_clock_biases + num_ifb columns (three position partials, one per GNSS time system, and one per frequency band beyond the reference). D = (HᵀH)⁻¹ is the parameter covariance in the units/frame of H. Because calc_H! builds H from ECEF unit vectors (the negated receiver→satellite lines of sight), the position block of D is in ECEF; it is rotated into the local ENU frame at user_pos before the horizontal/vertical split, so HDOP/VDOP are taken in the user's tangent plane (PDOP/GDOP are trace-invariant and unaffected by the rotation). GDOP spans all parameters; TDOP reports the clock variance of the primary (reference) system — see PVTSolution — while the other systems' clock (inter-system-bias) and the inter-frequency-bias variances enter GDOP only.
A rank-deficient geometry makes HᵀH singular (not positive definite, or only to within rounding); the factorisation detects this and the function returns the sentinel DOP(-1, …) instead of erroring.
calc_DOP(H_GEO, user_pos)
calc_DOP(H_GEO, user_pos, primary_clock_index)
H_GEO: Geometry matrix user_pos: User ECEF position, for the ECEF→ENU rotation of the horizontal/vertical DOPs primary_clock_index: Index (1…numclockbiases) of the clock column whose variance is reported as TDOP
PositionVelocityTime.calc_DOP! — Function
calc_DOP!(normal_matrix, H_GEO, user_pos::ECEF, primary_clock_index = 1) -> DOPcalc_DOP with the normal-equations matrix HᵀH formed, factorised and inverted in normal_matrix (n × n for the n columns of H_GEO, overwritten), so that it allocates nothing. The factorisation is LAPACK's Cholesky potrf, followed by a relative rank test on its pivots, and the inverse is potri on that factor — what inv of a Cholesky computes.
PositionVelocityTime.calc_pseudorange_rate — Function
calc_pseudorange_rate(user_pos, user_vel, user_clock_drift, sat_pos, sat_vel,
sat_clock_drift) -> Float64The modelled pseudorange rate (m/s) of one satellite: the exact time derivative, per second of receiver clock, of the pseudorange calc_ρ_hat! models — the satellite rotated by the Earth's turn over the travel time into the ECEF frame of reception, plus the receiver clock bias, minus the satellite clock offset. It is the model the velocity solve (calc_user_velocity_and_clock_drift) fits, and its rate residuals are the measured rate minus this at the solved state; an estimator that predicts the rate from its own state — a vector-tracking navigation filter, say — calls it to predict exactly what the solve fits.
Positive while the range grows (the satellite recedes), the sense of the velocity solve and of RTKLIB's resdop; the Doppler-signed rate a tracking loop works in, λ · carrier_doppler, is its negative.
user_pos, user_vel, sat_pos and sat_vel are ECEF SVector{3}s in m and m/s, the satellite's taken at transmission. user_clock_drift is the receiver clock drift in m/s, but sat_clock_drift the satellite clock drift as a rate in s/s, as calc_satellite_clock_drift returns it.
Beyond the line-of-sight projection of the relative velocity it carries two terms of order v/c, each up to about a millimetre per second:
- the rate of the Earth-rotation angle,
ω_e·τ̇withτ̇ = û·(v_sat − v_user)/c(EARTH_ROTATION_RATE), and - the light-time factor
k = 1 − ρ̇/c: one second of receiver clock spans that many seconds of the satellite's motion and clock,ρ̇being the rate of the range including the receiver clock drift.
Allocation-free and @inline.
PositionVelocityTime.calc_user_velocity_and_clock_drift — Function
Computes user velocity
calc_user_velocity_and_clock_drift(measurements, user_pos)
calc_user_velocity_and_clock_drift(
measurements,
user_pos,
start
)
Calculates the user velocity and a single receiver clock drift, returned together with the per-satellite post-fit range-rate residuals as ([vx, vy, vz, ċ], rate_residuals). Unlike the position solve — which estimates one clock bias per GNSS time system (the inter-system bias / GGTO is metre-level and must be resolved) — a single clock drift is shared by all satellites regardless of GNSS: the receiver has one oscillator, and the only per-system difference is the drift of the inter-system time offset (e.g. the GGTO rate A_1G), which is ~1e-6 m/s — far below the Doppler velocity resolution. Using one common drift lets every satellite constrain the four unknowns instead of spending a column per system.
The measured range rate is fitted to calc_pseudorange_rate, the exact time derivative of the pseudorange calc_ρ_hat! models (see calc_pseudorange_rate for its two terms of order v/c), by one Gauss–Newton step from start ([vx, vy, vz, ċ], m/s; the resting receiver by default). The model is so nearly linear that one step is the least-squares fit: within ~1e-11 m/s from a start a few m/s off, such as the previous epoch's solution, and within ~1e-7 m/s from rest even for an aircraft's speed and 30 km/s of clock drift — far below any Doppler's noise either way. That model needs the satellite's travel time, its rate and its Earth-rotation-corrected position, hence the solved user_pos (ECEF, m — the position block of the converged ξ).
The residuals are measured − modeled range rate (m/s): the measured rate minus calc_pseudorange_rate at the solved state, to the same precision, the same orientation as the pseudorange residuals of user_position and in the same satellite order as measurements. They are the range-rate analogue of the post-fit pseudorange residual: a per-satellite Doppler-consistency / outlier indicator. Measured and modeled are both taken in that model's sense, in which a receding satellite reads positive — the same quantity and sign as RTKLIB's resdop residual, and hence the negative of the Doppler-signed range rate a tracking loop works in (see the note at the residual loop).
Needs a geometry whose position design at user_pos (calc_H!) has full column rank, as calc_DOP establishes. Where the velocity design itself is not of full rank to within the tolerance calc_DOP applies, the solve does not throw but returns NaN for the state and every residual.
The allocating form of calc_user_velocity_and_clock_drift!.
PositionVelocityTime.calc_user_velocity_and_clock_drift! — Function
calc_user_velocity_and_clock_drift!(rate_residuals, design, measurements, user_pos,
start = zero(SVector{4,Float64})) -> (velocity_and_drift::SVector{4}, rate_residuals)calc_user_velocity_and_clock_drift writing the range-rate residuals into rate_residuals, with design — a Vector{SVector{4,Float64}} — as scratch for the design rows; both are resized to the satellite count. The solve itself is a fixed 4 × 4 static one and allocates nothing.
PositionVelocityTime.time_scale_offset_to_gpst — Function
time_scale_offset_to_gpst(time_system) -> Float64Signed offset of a GNSS time system's count against GPS Time's for the same instant — get_tai_offset(GPST) − get_tai_offset(time_system), so negative where the system's count reads lower (BDT: 19 − 33 = −14.0), and 0.0 for a system that counts alike. calc_time_scale_offsets adds the difference of two of these to a transmit time, which is what puts the +14 s onto a BeiDou time in a GPS-primary solve.
This is structural, not a bias: it follows from the time scales' definitions, not from either system's steering. Both GPST and GST are TAI − 19 s, so they count identically and this is 0.0 for every GPS and Galileo satellite. BDT is TAI − 33 s, so a BeiDou second-of-week reads 14 s lower than the GPS time-of-week for the same instant — BDT week 0 second 0 is GPS week 1356 time of week 14, the 14 leap seconds that had accrued between the two epochs.
Derived from GNSSSignals.get_tai_offset rather than tabulated, so a constellation added later is covered without touching this.
It would be tempting to fold the 14 s into get_time_of_week and be done. That would be wrong: the satellite position is propagated from the same transmit time against the message's own t_0e, which is on the broadcasting system's scale, so shifting the reported time would move the 14 s into the ephemeris — about 55 km of along-track error at BeiDou MEO velocities. The transmit time must stay on its own scale and the correction must land on the pseudorange.
PositionVelocityTime.calc_time_scale_offsets — Function
calc_time_scale_offsets(measurements, primary_system) -> Vector{Float64}Seconds to add to each satellite's transmit time to express it in primary_system's count, so that calc_pseudoranges may difference them. 0.0 for every satellite of a system that counts alike — which is every GPS and Galileo satellite, and all of them in a single-constellation epoch.
Relative to the primary system rather than to GPS Time, because the primary system's count is also what dates the reported epoch: reference_time comes out of the same differencing, and is combined with the primary system's week and system_start_epoch. Anchoring on GPST instead leaves a BeiDou-primary solve reporting a GPS-count time of week dated from the BDT epoch — a 14 s error in PVTSolution.time traded for the one this removes from the pseudoranges.
This is applied to the times handed to the differencing, and nowhere else. The satellite position must still be propagated from the untouched transmit time against the message's own t_0e, which is on the broadcasting system's scale, so shifting the time itself would move the offset into the ephemeris — about 55 km of along-track error at BeiDou MEO velocities. SatInfo.time and calc_steering_offset likewise keep the unconverted value.
PositionVelocityTime.time_offset_available — Function
time_offset_available(measurement::SatelliteMeasurement, target) -> Bool
time_offset_available(decoder::GNSSDecoderState, target) -> BoolWhether a usable broadcast offset from this satellite's own GNSS time system to target is available. Such an offset lets that constellation's measurements be expressed on the target system's clock, which makes a fix possible when the geometry is too weak to estimate an independent clock bias for it — see decide_bias_layout, which collapses onto a hub system this way. Toward GPS Time that is Galileo's GGTO and BeiDou's BGTO; toward Galileo System Time it is BeiDou's BGTO variant and the GGTO GPS itself broadcasts on CNAV/CNAV-2.
The measurement form reads the flag collect_measurements already flattened into measurement.time_offsets (equivalently measurement.time_offsets[i].available), and is false for any target outside CANDIDATE_HUB_SYSTEMS, which no message broadcasts an offset toward.
The decoder form is the underlying predicate, and delegates to GNSSDecoder's get_time_offset, which screens every way a signal can fail to have one: it broadcasts none toward target (GPS L1 C/A and Galileo E6-B broadcast none at all), it does but has not decoded one yet, the decoded one names a different target system, the ICD's "not available" sentinel is set (a GNSS_ID of 0; Galileo's all-ones GGTO), or — on Galileo — the reference week has not arrived yet, since word type 10 can be decoded before the week number is. Always false when target is the decoder's own system: the offset from a scale to itself is not a broadcast quantity.
PositionVelocityTime.broadcast_time_offset — Function
broadcast_time_offset(decoder, target; approximate_year = year(now(UTC)))
-> BroadcastTimeOffsetThe single BroadcastTimeOffset from decoder's own GNSS time system toward target, as collect_measurements flattens it onto every SatelliteMeasurement — the standalone form, for a consumer evaluating a broadcast time offset outside a fix. Unavailable (available == false) for a target outside CANDIDATE_HUB_SYSTEMS, and for any message that carries no such offset (see time_offset_available).
PositionVelocityTime.calc_steering_offset — Function
calc_steering_offset(offset::BroadcastTimeOffset, t) -> Float64The broadcast steering offset between a constellation's time scale and the target scale offset was built against, in seconds, evaluated at the own-scale time of week t. Subtract it to convert a transmit time to the target system's time. Only defined where offset.available (see time_offset_available).
Tens of nanoseconds: it is the residual between two atomic scales, not the whole difference between their counts. GNSSDecoder's GNSSTimeOffset.A_0 folds in the defined whole-second offset as well — so that t_target = t_own − Δt holds for the seconds — and offset.count_anchor takes that part back out, because calc_time_scale_offsets has already applied it to the transmit times before they were differenced. Both sides are the same time_scale_offset_to_gpst expression, so the two halves compose exactly rather than approximately.
Recovering a ~1e-8 s residual by subtracting 14 s from a Float64 leaves roughly 1.8e-15 s of rounding — one ULP at 14. That is 0.5 µm of range, so it is irrelevant to a fix, but it is why the tests here compare the steering term with an explicit tolerance rather than the default ≈.
PositionVelocityTime.inter_system_bias — Function
inter_system_bias(row, hub_system, t = row.time) -> Float64The range offset (metres) of row's time system against hub_system at the own-scale time of week t: −c · Δt_systems, with Δt_systems the broadcast steering offset row carries toward hub_system (see calc_steering_offset). Added to a modeled range it puts a measurement of row's system on the hub's clock. A collapse reports it as that system's entry of inter_system_biases, evaluated at the transmit time of the system's first row (other times differ by the offset's drift term). Zero when row is of hub_system itself, NaN where row carries no usable offset toward it (see time_offset_available).
Pass another row's transmit time as t to convert that row with this one's offset, as calc_hub_range_offsets does.
PositionVelocityTime.get_week — Function
get_week(decoder::GNSSDecoderState{<:GPSL1CAData}; approximate_year)Return the absolute GPS week number for a GPSL1 decoder, resolving the 1024-week rollover ambiguity using approximate_year as a calendar anchor.
The legacy GPS L1 C/A LNAV message broadcasts only a 10-bit week number (0–1023) modulo 1024, so the receiver cannot determine which 1024-week cycle the recording is in from the data alone (IS-GPS-200, §20.3.3.3). Each cycle is ~19.6 years, so any anchor within ±9 years of the true observation date selects the correct cycle.
GPS week 0 is 1980-01-06; cycle boundaries fall on 1999-08-22, 2019-04-07, 2038-11-21, 2058-07-08, …
For Galileo, the broadcast WN is 12 bits and does not need this treatment in any practical operational scenario.
PositionVelocityTime.system_start_epoch — Function
system_start_epoch(system) -> TAITimeAbsolute TAI epoch of a ranging signal's GNSS time-scale origin (week 0, time of week 0), from GNSSSignals' get_tai_system_start_time — the epoch already labelled on the atomic scale (GPS 1980-01-06T00:00:19 TAI, Galileo 1999-08-22T00:00:19 TAI, BeiDou 2006-01-01T00:00:33 TAI), so no leap-second table is consulted. This used to go through get_system_start_time (a UTC label) and AstroTime's leap-aware conversion; GNSSSignals 4.1 states the TAI labels itself, which is also what makes the value safe to derive at precompile time.
PositionVelocityTime.day_of_year — Function
day_of_year(system_start_time::DateTime, week, time_of_week) -> IntDay of year (1–366) of the GNSS system time given by the absolute week and time_of_week (s) counted from system_start_time, that system's week-0 origin as a calendar date (get_system_start_time(system), precomputed on every SatelliteMeasurement) — the seasonal argument of the Niell tropospheric mapping (see tropospheric_delay). The system time scale's offset from UTC (leap seconds, ≤ ~18 s) is negligible for an argument with a one-year period.
PositionVelocityTime.predict_atmospheric_delays — Function
predict_atmospheric_delays(ξ, measurements, correction,
reference_time, doy, enable_tropospheric_correction) -> Vector{Float64}Per-satellite ionospheric + tropospheric delay (metres), to be subtracted from the pseudoranges. The user position is the first three elements of the least-squares state vector ξ = [x, y, z, tc₁, …] (ECEF, metres); the remaining clock-bias components are ignored. correction is the constellation-wide ionospheric model from select_ionospheric_correction (nothing skips the ionosphere); the troposphere uses the blind Saastamoinen zenith delays mapped by the Niell mapping functions, whose seasonal term takes the day of year doy (see tropospheric_delay), unless enable_tropospheric_correction is false.
A single corrected solve is enough: the delays depend on position only through the satellite elevation/azimuth (and, for the troposphere, the user height), and ∂delay/∂position is negligible over the metre-level position uncertainty (a 15 m shift moves the elevation by ~1e-5°), so delays predicted at a nearby position are accurate to well under a millimetre — no iterate-to-convergence needed. The user geodetic coordinates and the ENU transform depend only on ξ, so they are computed once and reused across satellites.
PositionVelocityTime.predict_atmospheric_delays! — Function
predict_atmospheric_delays!(delays, ξ, measurements, model::IonosphericModel,
reference_time, doy, enable_tropospheric_correction)
-> delayspredict_atmospheric_delays into delays (resized to the satellite count), with the ionospheric correction carried by an IonosphericModel. Allocates nothing.
PositionVelocityTime.IonosphericModel — Type
IonosphericModel(correction)The ionospheric correction of one solve (see select_ionospheric_correction), held in a single concrete type. The correction itself is a Union of four coefficient sets and nothing; the solver takes this wrapper instead so that it compiles once whichever model an epoch selected, and so that passing the model costs no allocation — a Union value handed to an unspecialised argument is boxed, whereas a field of a concrete struct is stored inline.
PositionVelocityTime.calc_course_over_ground — Function
calc_course_over_ground(position::ECEF, velocity::ECEF) -> typeof(1.0°)Course over ground (degrees): the azimuth of velocity in the local East-North-Up frame at position, measured clockwise from true North and wrapped to [0, 360)°, following the GNSS COG convention. This is the direction of motion derived from the velocity, not vehicle heading. The vertical (Up) component is ignored. Returns 0° when the horizontal velocity is zero (stationary or purely vertical), where course is undefined.
ENUfromECEF is affine with position as its origin (which it maps to zero), so applying it to position + velocity rotates the velocity vector into the ENU frame without introducing any translation.
The unexported internals the surface above links to, documented for reference:
PositionVelocityTime.calc_H! — Function
Computes Geometry Matrix H
The matrix has num_lsq_params(layout) columns: three for the position partials, one per GNSS time system, and one per frequency band beyond the reference band. Each row carries a 1.0 in its system's clock column (3 + clock_bias_indices[j]) and, unless the satellite is on the reference band, a 1.0 in its band's inter-frequency-bias column (3 + num_clock_biases + ifb_indices[j]); the remaining bias columns are zero.
calc_H!(H, sat_positions, ξ, bias_columns)
ξ: Estimated user position, per-system clock corrections, and per-band inter-frequency biases sat_positions: Satellite positions, a 3 × N matrix or a vector of SVector{3}s bias_columns: The per-satellite BiasColumns
PositionVelocityTime.user_position — Function
Computes user position
user_position(sat_positions_mat, ρ, bias_columns)
user_position(sat_positions_mat, ρ, bias_columns, prev_ξ)
sat_positions_mat: Satellite positions, as a (3, N) matrix (xyz per satellite) or a vector of SVector{3}s. ρ: Array of pseudoranges
Calculates the user position by least squares method. The algorithm is based on the common reception method.
bias_columns: The per-satellite BiasColumns (clock and inter-frequency-bias columns).
Returns (ξ, residuals): the solved state vector ξ = [x, y, z, tc₁, …, ifb₁, …] and the per-satellite post-fit residual vector (measured minus modeled pseudorange, metres), in the same satellite order as ρ.
curve_fit (LsqFit's, mirrored by LevenbergMarquardt) reports its own residual as model - data, so the returned vector negates it. Measured − modeled ("observed minus computed") is how GNSS software reports observation residuals — RTKLIB's rescode, and GNSS-SDR and PocketSDR through it — and the negation is the whole of the difference: it is applied to the converged fit, so the solve itself is untouched.
The allocating form of user_position!, on fresh buffers.
PositionVelocityTime.user_position! — Function
user_position!(workspace::LMWorkspace, residuals, sat_positions, ρ, bias_columns, prev_ξ)
-> (ξ, residuals)user_position on the buffers of workspace, writing the post-fit residuals into residuals (resized to the satellite count). The returned ξ is workspace.x, so it is overwritten by the next solve on the same workspace; prev_ξ may be that same vector, to restart a solve from the previous one's solution.
PositionVelocityTime.num_lsq_params — Function
num_lsq_params(bias_columns::BiasColumns) -> IntLength of the least-squares state vector for bias_columns: 3 + num_clock_biases + num_ifb.
PositionVelocityTime.calc_hub_range_offsets — Function
calc_hub_range_offsets(measurements, hub_rows, hub_system) -> Vector{Float64}Per-satellite range offsets (metres) that carry a clock collapse into the measurements, as decided by decide_bias_layout: all-zero for the satellites of a system that keeps its own clock unknown — all of them when hub_rows is empty — and otherwise −c · Δt_systems for each satellite of a collapsed system, evaluated at its own transmit time.
The broadcast offset is Δt_systems = (that system's time) − (the hub system's) (see calc_steering_offset: the GGTO for Galileo toward a GPS hub, the BGTO for BeiDou toward either), so a transmit time becomes hub time by SUBTRACTING it; the modeled range therefore carries −c·Δt_systems, and the solve yields inter_system_biases[sys] = −c·Δt_systems. Which satellite of a system reported the offset does not matter — it is one constellation-wide value — so decide_bias_layout picks the first decoded copy per system (hub_rows) and it converts all of that system's measurements.
PositionVelocityTime.unique_time_systems — Function
unique_time_systems(systems) -> Vector{SupportedTimeSystem}The distinct GNSS time systems of systems (an iterable of TimeSystems), in order of first appearance — the order that fixes the clock columns of BiasColumns.
PositionVelocityTime.time_system_index — Function
time_system_index(systems, system) -> Union{Int,Nothing}Position of system in systems, or nothing if absent.
Multi-GNSS Classification
When calc_pvt combines constellations and bands it classifies each satellite by three keys, all provided by GNSSSignals (3.3+) and read from a satellite's ranging signal: GNSSSignals.get_time_system (a GNSSSignals.TimeSystem, i.e. GPST()/GST()/BDT()) drives receiver-clock grouping (one bias per time system); GNSSSignals.get_band_id (e.g. :L1, :L5) drives inter-frequency-bias grouping (one bias per band); and GNSSSignals.get_signal_id (e.g. :GPSL1CA) is the per-signal identity used in the sats key of PVTSolution, and the key a PositionVelocityTime.SignalGroup is one signal's worth of. All three are read off the ranging signal once, by the collection pass, and carried on every PositionVelocityTime.SatelliteMeasurement as time_system, band_id and signal_id.
The receiver inter-frequency biases and their reference bands (reported per InterFrequencyBias) are laid out from the constellation × band coverage graph:
PositionVelocityTime.band_ifb_layout — Function
band_ifb_layout(system_keys, bands)
-> (ifb_indices, extra_bands, reference_bands, num_components)Lay out the receiver inter-frequency biases from the (constellation × band) coverage graph. system_keys[j] and bands[j] are satellite j's constellation and frequency band; the constellation enters only through equality and first appearance, so any per-satellite key that separates the constellations will do — decide_bias_layout passes the clock column of each satellite, which is exactly that and, unlike the TimeSystem singletons themselves, concretely typed.
Two bands share a coverage component iff some constellation is tracked on both (directly or through a chain of shared constellations); within a component the clock/IFB split has a single gauge freedom, so one reference band per component (the most-populated, ties by first appearance) is fixed and an IFB column is created only for the other bands of that component. This estimates exactly the observable IFBs — a band that is the sole band of its component (its constellation lives only there) gets none, its delay folding into that constellation's clock — so the resulting design matrix is always full rank.
ifb_indices[j] is satellite j's IFB column (1…length(extra_bands)), or 0 for a per-component reference band; extra_bands[i] is the band of IFB column i (ordered deterministically) and reference_bands[i] is the reference band of that column's coverage component (the anchor its IFB is measured against); num_components is the number of coverage components (1 ⇔ the graph is connected).
PositionVelocityTime.band_ifb_layout! — Function
band_ifb_layout!(ifb_indices, extra_bands, reference_bands, scratch, system_keys, bands)
-> num_componentsband_ifb_layout writing its three vectors into the given ones (emptied first) and working in scratch, so that it allocates nothing once they have grown to the epoch. The result is identical: the same components, the same references and tie-breaks, and the same column order.
Trimmed Executables
The solver compiles into a juliac --trim=safe executable (Julia 1.12+ with JuliaC); test/trim holds an app that does so for every navigation-data type and a mixed GPS + Galileo + BeiDou epoch, and checks its output against a regular session. Signal groups are what make that possible for any constellation mix — each group is concretely typed, and the solver behind them sees one concrete row type.
Two dependencies the package would otherwise have cannot be trimmed, so it carries trim-safe stand-ins for the part of each it uses, as self-contained submodules with the original's interface — to be dropped once the originals can be trimmed:
PositionVelocityTime.TAITimesreplaces AstroTime'sTAIEpochas the type ofPVTSolution.time(AstroTime's own__init__cannot be trimmed), andPositionVelocityTime.LevenbergMarquardtreplaces LsqFit'scurve_fit, which calls the model through abstractly typed fields.
PositionVelocityTime.TAITimes — Module
TAITimesA minimal, trim-safe TAI epoch — TAITime — standing in for AstroTime.jl's TAIEpoch as the type of PositionVelocityTime.PVTSolution's time.
Loading AstroTime runs __init__ code that a juliac --trim=safe build cannot verify (it registers Dates format tokens and, through EarthOrientation and RemoteFiles, looks for data files and download tools), so the package cannot depend on it and still be trimmed. TAITime stores an epoch the way TAIEpoch does — whole TAI seconds since J2000 (2000-01-01T12:00:00 TAI) plus a fraction — so the two convert exactly: with AstroTime loaded, TAIEpoch(t) and TAITime(epoch) are provided by a package extension.
Kept as a self-contained module so that once AstroTime can be trimmed, it can be dropped in favour of TAIEpoch again.
PositionVelocityTime.LevenbergMarquardt — Module
LevenbergMarquardtA trim-safe stand-in for the part of LsqFit.jl that PositionVelocityTime.user_position uses: LsqFit's curve_fit with an in-place model and Jacobian, as curve_fit! on a reusable LMWorkspace, with the same arguments after the workspace and the param / resid fields of its result.
LsqFit keeps the model and Jacobian in abstractly typed fields of NLSolversBase's OnceDifferentiable, so every call through them is a dynamic dispatch, and curve_fit formats caught MethodErrors and can print a trace — all of which rules the solve out of a juliac --trim=safe build. Here the functions are type parameters and there is no trace, bounds, timer or error formatting.
The algorithm is LsqFit's levenberg_marquardt (after Kanzow, Yamashita & Fukushima
- step for step: the same damped normal equations with the
DtDdiagonal
floored at 1e-6, the same geodesic-acceleration correction when avv! is given, the same step-quality ratio deciding acceptance and trust-region updates, and the same gradient and relative step-size stopping rules, with the Jacobian evaluated lazily just as LsqFit's cache does — so it lands on the same iterates, to rounding (the damped normal equations are solved by Cholesky rather than LU; see curve_fit!). Running it on a reusable workspace means a fit allocates nothing.
Kept as a self-contained module, with LsqFit's interface, so that once LsqFit can be trimmed it can be dropped in favour of using LsqFit: curve_fit (on a fresh workspace, curve_fit!(LMWorkspace(), args...) is exactly that call).
PositionVelocityTime.LevenbergMarquardt.curve_fit! — Function
curve_fit!(workspace::LMWorkspace, model!, jacobian!, xdata, ydata, p0; kwargs...)
-> LMResultFit model!(out, xdata, p) to ydata by least squares, starting from p0, with jacobian!(J, xdata, p) its Jacobian in p — LsqFit's in-place curve_fit (inplace = true is the only form provided), on the buffers of workspace, allocating nothing once they have grown to the problem's size. avv!(dir_deriv, p, v) enables geodesic acceleration. Like LsqFit, throws an ArgumentError if ydata holds non-finite values.
The keyword arguments are LsqFit's, with its defaults: inplace = true, avv! = nothing, lambda = 10.0, x_tol = 1e-8, g_tol = 1e-12, maxIter = 1000, min_step_quality = 1e-3, good_step_quality = 0.75, lambda_increase = 10.0 and lambda_decrease = 0.1.
The result's param and resid are the workspace's x and f buffers: they are overwritten by the next fit on the same workspace, so copy them to keep them. p0 may itself be workspace.x (a fit restarted from the previous one's solution).
jacobian! is handed a view of the workspace's Jacobian matrix, not a Matrix.
PositionVelocityTime.LevenbergMarquardt.LMWorkspace — Type
LMWorkspace()The scratch storage of one curve_fit! — every vector and matrix of the iteration, so that repeated fits reuse them instead of allocating their own. Buffers grow to the largest problem seen and are then reused; a fit never shrinks them.
The matrices are kept at their largest size and used through views of the leading m × n block, because a Matrix cannot be resized in place. The vectors are resized to the current problem, which within the capacity they already have does not allocate.
PositionVelocityTime.LevenbergMarquardt.grown — Function
grown(matrix, rows, cols) -> Matrix{Float64}matrix if it holds at least rows × cols elements in that shape, else a new matrix that does. The caller stores the result back and works on view(matrix, 1:rows, 1:cols).
PositionVelocityTime.LevenbergMarquardt.LMResult — Type
LMResultThe part of LsqFit.LsqFitResult that is used: the fitted parameters param, the residual resid = model(xdata, param) − ydata at them, and whether the iteration converged (rather than running out of iterations).