Read a position from a curve

1. Find the interval containing the instant

First choose the body and quantity, such as Moon longitude. The evaluator finds the first table entry with that pair. If none exists, it returns MissingSeries. Longitude and declination need separate entries, so finding one does not prove the other is present.

Consider invented metadata with start time t₀ = 1000 s, segment length D = 100 s, and three segments. They span 1000–1100, 1100–1200 and 1200–1300 s. At t = 1175 s, subtraction gives 175 s from the start. Division by 100 gives 1.75; rounding downward gives index 1. Indices start at zero, so this is the second segment.

relative = (t − t₀) / D; index = floor(relative)

At a shared boundary such as 1200 s, index 2 selects the segment on the right. At the final endpoint 1300 s, the initial index equals the segment count, 3. The evaluator subtracts one in this special case and uses index 2. Dates before the start or beyond the final endpoint are rejected; it does not extrapolate a position.

The selected segment begins at 1100 s. Our time is 75 s into a 100-second interval, or three quarters of the way across. Mapping the start to −1 and the end to +1 turns that fraction into x = 2 × 0.75 − 1 = 0.5. This dimensionless local coordinate is the input to the stored polynomial.

segment_start = t₀ + index × D; x = 2 × (t − segment_start) / D − 1

If this series begins at byte offset 208 and stores four eight-byte coefficients per segment, index 1 starts at 208 + 1 × 4 × 8 = 240. The evaluator reads those four little-endian floating-point values. “Little-endian” means the least significant byte is stored first.

const localCoordinate = (t: number, start: number, duration: number): number => 2 * (t - start) / duration - 1;

2. Combine the stored polynomial weights

A polynomial can be built from familiar pieces. For a Chebyshev series, the first pieces are T₀(x) = 1, T₁(x) = x and T₂(x) = 2x² − 1, where x² means x multiplied by itself. With invented degree-valued coefficients [10, 2, 3] and x = 0.5, the pieces are 1, 0.5 and −0.5. Multiplying each by its weight and adding gives 10 + 1 − 1.5 = 9.5°.

p(x) = c₀T₀(x) + c₁T₁(x) + c₂T₂(x) + …

For a longer series we need not build every piece separately. The engine uses a backward calculation named Clenshaw summation. Begin with two temporary values b₁ = 0 and b₂ = 0, and visit coefficients from the last down to c₁. For each coefficient compute a new value b₀, then replace the old pair with b₁ = b₀ and b₂ = the previous b₁. The subscripts here identify temporary slots, not coefficient indices.

b₀ = 2xb₁ − b₂ + c; after the loop, p(x) = c₀ + xb₁ − b₂

In our example the c₂ = 3 step produces b₀ = 3, leaving (b₁, b₂) = (3, 0). The c₁ = 2 step produces 2 × 0.5 × 3 − 0 + 2 = 5, leaving (5, 3). The final expression is 10 + 0.5 × 5 − 3 = 9.5°. This is the same polynomial as before, evaluated by a reusable sequence of multiplication and subtraction. A one-coefficient series simply returns c₀.

The engine stores c₀ with its full weight. Some references use a half-weight convention for the constant coefficient; their final expression must be interpreted with their coefficient definition. Halving this file’s c₀ would change the answer.

const clenshaw = (coefficients: readonly number[], x: number): number => {
  const [b1, b2] = coefficients.slice(1).reduceRight<[number, number]>(([a, b], c) => [2 * x * a - b + c, a], [0, 0]);
  return coefficients[0]! + x * b1 - b2;
};

3. Give the two angular quantities their meaning

The fitter kept longitude continuous across complete turns so that a curve could pass smoothly through the zero-degree seam. Suppose evaluation returns 361.25°. On a chart circle, subtracting one full turn gives 1.25°, the same direction. Likewise −2° becomes 358°. The longitude accessor takes a nonnegative remainder after division by 360.

longitude = nonnegative remainder of p(x) divided by 360°

Declination has a different meaning: it is the angle north or south of the celestial equator, the plane through Earth’s equator extended into the sky. Its sign carries that distinction. A fitted value of −12.5° remains −12.5°; turning it into 347.5° would discard its meaning. The declination accessor returns the evaluated value unchanged. Both quantities are degrees, but only longitude is wrapped around a full circle.

These accessors recover quantities already prepared by the builder. Their use of remainder and a signed angle is a numerical convention in the current runtime; there is no separately sourced origin story for these wrapper functions. Mathematical wrap formulas describe the ideal [0, 360) interval. At extremely small negative floating-point inputs, rounding can produce exactly 360, as the source-level limitation recorded in earlier lessons explains.

const longitude = (angle: number): number => { const remainder = angle % 360; return remainder < 0 ? remainder + 360 : remainder; };
const declination = (angle: number): number => angle;