feat(globe): port to CeliumJS

tests/unit/boundingBox.spec.ts

@@ -0,0 +1,200 @@
+import { test, expect } from '@playwright/test';
+import {
+	computeBoundingBox,
+	boundingBoxRing,
+	type BoundingBox,
+} from '../../src/lib/domain/boundingBox';
+import type { LatLngTuple } from '../../src/lib/domain/geo';
+
+/**
+ * The restricted area filed with the regulator.
+ *
+ * It is a rectangle in kilometres, axis-aligned to east/north at its own centre,
+ * with four lat/lon corners joined by great circles. Not a rectangle in degrees:
+ * that form cannot work near a pole, because every meridian passes through the
+ * pole, so any lat/lon rectangle containing one spans all 360 deg of longitude.
+ * The measured cost of that was a 44 200 km^2 cap standing in for a 1 695 km^2
+ * corridor.
+ *
+ * These tests run in Node, not a browser — the module under test is pure
+ * geometry and imports only a type.
+ */
+
+const R_KM = 6371;
+const rad = (d: number) => (d * Math.PI) / 180;
+const deg = (r: number) => (r * 180) / Math.PI;
+
+/**
+ * Standard spherical destination-point formula, written out here rather than
+ * imported so the tests do not check the implementation against itself.
+ * Undefined starting exactly at a pole, which is why the polar cases use 89.99.
+ */
+function destination(lat: number, lng: number, bearingDeg: number, distKm: number): LatLngTuple {
+	const d = distKm / R_KM;
+	const br = rad(bearingDeg);
+	const p1 = rad(lat);
+	const l1 = rad(lng);
+	const p2 = Math.asin(Math.sin(p1) * Math.cos(d) + Math.cos(p1) * Math.sin(d) * Math.cos(br));
+	const l2 =
+		l1 +
+		Math.atan2(Math.sin(br) * Math.sin(d) * Math.cos(p1), Math.cos(d) - Math.sin(p1) * Math.sin(p2));
+	return [deg(p2), deg(l2)];
+}
+
+/** A straight meridional track: identical in kilometres at any latitude. */
+function meridionalTrack(lat: number, lng: number, lengthKm: number, n = 40): LatLngTuple[] {
+	return Array.from({ length: n }, (_, i) => destination(lat, lng, 180, (lengthKm * i) / (n - 1)));
+}
+
+function toVec([lat, lng]: LatLngTuple): [number, number, number] {
+	const p = rad(lat);
+	const l = rad(lng);
+	return [Math.cos(p) * Math.cos(l), Math.cos(p) * Math.sin(l), Math.sin(p)];
+}
+
+const dot = (a: number[], b: number[]) => a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
+const cross = (a: number[], b: number[]): [number, number, number] => [
+	a[1] * b[2] - a[2] * b[1],
+	a[2] * b[0] - a[0] * b[2],
+	a[0] * b[1] - a[1] * b[0],
+];
+
+/**
+ * Signed clearance in km from a point to the great circle through two corners,
+ * positive on the side the box interior is on.
+ *
+ * This is measured against the edge the regulator would draw — the great circle
+ * between the filed corners — not against the projected rectangle. The two are
+ * not the same: a great-circle edge bows toward the centre of the box relative
+ * to its chord in the projection, by roughly (half-edge)^2 / 2R. On a 500 km
+ * edge that is 4.9 km, so it silently eats a 5 km margin whole.
+ */
+function clearanceKm(point: LatLngTuple, from: LatLngTuple, to: LatLngTuple, inside: LatLngTuple) {
+	const n = cross(toVec(from), toVec(to));
+	const len = Math.hypot(...n) || 1;
+	const unit = n.map((c) => c / len);
+	const sign = Math.sign(dot(unit, toVec(inside))) || 1;
+	return sign * Math.asin(Math.max(-1, Math.min(1, dot(unit, toVec(point))))) * R_KM;
+}
+
+/** Smallest clearance from any track point to any of the four filed edges. */
+function worstClearanceKm(box: BoundingBox, path: LatLngTuple[]): number {
+	const c = box.corners;
+	const edges: [LatLngTuple, LatLngTuple][] = [
+		[c[0], c[1]],
+		[c[1], c[2]],
+		[c[2], c[3]],
+		[c[3], c[0]],
+	];
+	let worst = Infinity;
+	for (const p of path) {
+		for (const [a, b] of edges) {
+			const d = clearanceKm(p, a, b, box.centre);
+			if (d < worst) worst = d;
+		}
+	}
+	return worst;
+}
+
+/** Spherical excess area of the filed quad, km^2. */
+function areaKm2(box: BoundingBox): number {
+	const v = box.corners.map(toVec);
+	let sum = 0;
+	for (let i = 0; i < 4; i++) {
+		// Interior angle at vertex i, between the planes of its two edges.
+		const prev = v[(i + 3) % 4];
+		const next = v[(i + 1) % 4];
+		const n1 = cross(v[i], prev);
+		const n2 = cross(v[i], next);
+		const l1 = Math.hypot(...n1) || 1;
+		const l2 = Math.hypot(...n2) || 1;
+		sum += Math.acos(Math.max(-1, Math.min(1, -dot(n1, n2) / (l1 * l2))));
+	}
+	return (sum - 2 * Math.PI) * R_KM * R_KM;
+}
+
+/**
+ * The real trajectory from a launch at 89.99 N, 0 E, decimated. Under the old
+ * lat/lon-rectangle form this produced south 88.93, north 90, west -180,
+ * east 180: the entire cap, 44 200 km^2. In 1 degree of latitude it sweeps
+ * 79 degrees of longitude, because near a pole a short displacement crosses
+ * many meridians.
+ */
+const POLAR_TRACK: LatLngTuple[] = [
+	[89.99, 0],
+	[89.9564, 60.5627],
+	[89.8954, 70.1192],
+	[89.8095, 73.9615],
+	[89.6924, 76.1248],
+	[89.563, 77.1293],
+	[89.494, 77.2284],
+	[89.4553, 77.0498],
+	[89.4281, 76.8977],
+	[89.4109, 76.7301],
+	[89.3981, 76.6866],
+	[89.388, 76.9395],
+	[89.3823, 77.663],
+	[89.3771, 78.6642],
+	[89.3748, 79.0682],
+	[89.3736, 79.2029],
+	[89.3589, 79.2173],
+	[89.3177, 78.9415],
+	[89.2028, 78.9972],
+	[89.0974, 79.125],
+	[89.0269, 79.0071],
+	[88.9812, 78.8571],
+	[88.9782, 78.8563],
+];
+
+const MARGIN = 5;
+
+test('every trajectory point clears the filed edges by the full margin', () => {
+	// 500 km is where the great-circle bow matters: the east and west edges span
+	// +-255 km, and a chord-to-arc sag of 255^2/2R = 5.1 km would consume the
+	// entire 5 km margin and put the edge inside the trajectory.
+	const path = meridionalTrack(52.2, 0.1, 500);
+	const box = computeBoundingBox(path, MARGIN);
+
+	expect(box).not.toBeNull();
+	expect(worstClearanceKm(box!, path)).toBeGreaterThanOrEqual(MARGIN - 0.01);
+});
+
+test('every trajectory point clears the filed edges by the full margin near the pole', () => {
+	const box = computeBoundingBox(POLAR_TRACK, MARGIN);
+
+	expect(box).not.toBeNull();
+	expect(worstClearanceKm(box!, POLAR_TRACK)).toBeGreaterThanOrEqual(MARGIN - 0.01);
+});
+
+test('a box at the pole covers the corridor, not the whole polar cap', () => {
+	const box = computeBoundingBox(POLAR_TRACK, MARGIN);
+
+	// The lat/lon rectangle gave 44 200 km^2 for this track. The corridor it
+	// actually flies is 13.7 x 123.4 km.
+	expect(areaKm2(box!)).toBeLessThan(3000);
+});
+
+test('the same track in kilometres gives the same box at 52 N and at 89.99 N', () => {
+	// A meridional track is the one shape whose kilometre extent is independent
+	// of latitude, so any difference here is the code treating a pole specially.
+	const mid = computeBoundingBox(meridionalTrack(52.2, 0.1, 120), MARGIN)!;
+	const polar = computeBoundingBox(meridionalTrack(89.99, 0, 120), MARGIN)!;
+
+	expect(polar.heightKm).toBeCloseTo(mid.heightKm, 1);
+	expect(polar.widthKm).toBeCloseTo(mid.widthKm, 1);
+});
+
+test('the drawn ring closes and never spans the antimeridian in one segment', () => {
+	// Cesium cuts geometry at the IDL; a single segment straddling it with
+	// degenerate endpoints is what stopped the render loop before. Dense samples
+	// along each edge keep every segment short and the drawn shape equal to the
+	// filed one.
+	const ring = boundingBoxRing(computeBoundingBox(POLAR_TRACK, MARGIN)!);
+
+	expect(ring.length).toBeGreaterThan(16);
+	expect(ring[0]).toEqual(ring[ring.length - 1]);
+	for (let i = 1; i < ring.length; i++) {
+		const step = Math.abs(ring[i][1] - ring[i - 1][1]);
+		expect(Math.min(step, 360 - step)).toBeLessThan(90);
+	}
+});