···11(** GMAT interop tests for ocaml-odm.
2233 Parses CCSDS OEM files generated by NASA GMAT R2026a and validates
44- structural properties and physical constraints.
44+ structural properties, physical constraints, and parser robustness.
55+66+ GMAT output exercises real edge cases:
77+ - Scientific notation exponents down to 1e-16 (near-zero velocity
88+ components)
99+ - 5 orders of magnitude range in one file (Molniya perigee vs apogee)
1010+ - Trailing whitespace on metadata values (INTERPOLATION_DEGREE = 7 )
1111+ - Space inside ORIGINATOR value (GMAT USER)
1212+ - Variable-width columns (sign character shifts alignment)
513614 Source: GMAT R2026a, scripts in scripts/. Traces in traces/. *)
715···3745 Alcotest.(check string)
3846 "start" "2026-01-01T00:00:00.000" seg.metadata.start_time
39474848+(* GMAT outputs "ORIGINATOR = GMAT USER" — space inside the value.
4949+ Parser must not split on internal whitespace. *)
5050+let test_leo_originator_with_space () =
5151+ let oem = parse_oem "leo_7day.oem" in
5252+ Alcotest.(check string) "originator" "GMAT USER" oem.header.originator
5353+5454+(* GMAT outputs "INTERPOLATION_DEGREE = 7 " with trailing space.
5555+ Parser must trim before int_of_string. *)
5656+let test_leo_interpolation_degree () =
5757+ let oem = parse_oem "leo_7day.oem" in
5858+ let seg = List.hd (Odm.segments oem) in
5959+ Alcotest.(check (option int))
6060+ "interpolation degree" (Some 7) seg.metadata.interpolation_degree
6161+4062let test_leo_radius_bounds () =
4163 let oem = parse_oem "leo_7day.oem" in
4264 let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
···6385 Alcotest.(check string) "start epoch" "2026-01-01T00:00:00.000" start;
6486 Alcotest.(check string) "stop epoch" "2026-01-08T00:00:00.000" stop
65878888+(* First data line has Vz = 3.68e-16 km/s (machine epsilon).
8989+ Verify no NaN or underflow to exactly 0.0 — the parsed value should be
9090+ non-zero and finite. *)
9191+let test_leo_near_zero_velocity_component () =
9292+ let oem = parse_oem "leo_7day.oem" in
9393+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
9494+ let vz = svs.(0).Odm.vel.z in
9595+ Alcotest.(check bool) "Vz finite" true (Float.is_finite vz);
9696+ Alcotest.(check bool) "Vz not exactly zero" true (Float.abs vz > 0.0);
9797+ Alcotest.(check bool) "Vz ~3.68e-16" true (Float.abs vz < 1e-14)
9898+9999+(* All state vectors must have finite (non-NaN, non-inf) components.
100100+ GMAT uses 16-digit scientific notation — parser must handle full precision. *)
101101+let test_leo_no_nan_or_inf () =
102102+ let oem = parse_oem "leo_7day.oem" in
103103+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
104104+ Array.iteri
105105+ (fun i sv ->
106106+ let check name v =
107107+ if not (Float.is_finite v) then
108108+ Alcotest.failf "non-finite %s at vector %d (%s): %g" name i
109109+ sv.Odm.epoch v
110110+ in
111111+ check "X" sv.pos.x;
112112+ check "Y" sv.pos.y;
113113+ check "Z" sv.pos.z;
114114+ check "VX" sv.vel.x;
115115+ check "VY" sv.vel.y;
116116+ check "VZ" sv.vel.z)
117117+ svs
118118+119119+(* Epochs must be monotonically increasing. *)
120120+let test_leo_epoch_monotonic () =
121121+ let oem = parse_oem "leo_7day.oem" in
122122+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
123123+ for i = 1 to Array.length svs - 1 do
124124+ if String.compare svs.(i).Odm.epoch svs.(i - 1).epoch <= 0 then
125125+ Alcotest.failf "non-monotonic epoch at %d: %s <= %s" i svs.(i).epoch
126126+ svs.(i - 1).epoch
127127+ done
128128+129129+(* Energy conservation: specific orbital energy should be roughly constant
130130+ across the 7-day propagation (within drag/perturbation effects). *)
131131+let test_leo_energy_conservation () =
132132+ let oem = parse_oem "leo_7day.oem" in
133133+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
134134+ let mu = 398600.4418 in
135135+ let energy sv =
136136+ let r = vec3_norm sv.Odm.pos in
137137+ let v = vec3_norm sv.vel in
138138+ (0.5 *. v *. v) -. (mu /. r)
139139+ in
140140+ let e0 = energy svs.(0) in
141141+ (* LEO with drag: energy decays. Allow ~1% change over 7 days. *)
142142+ Array.iter
143143+ (fun sv ->
144144+ let e = energy sv in
145145+ let rel = Float.abs (e -. e0) /. Float.abs e0 in
146146+ if rel > 0.01 then
147147+ Alcotest.failf "energy drift > 1%%: %.6e vs %.6e (%.2f%%) at %s" e e0
148148+ (rel *. 100.0) sv.epoch)
149149+ svs
150150+151151+(* Interpolation: value at a known epoch should match the data point. *)
152152+let test_leo_interpolation_at_epoch () =
153153+ let oem = parse_oem "leo_7day.oem" in
154154+ let seg = List.hd (Odm.segments oem) in
155155+ let svs = Odm.state_vectors seg in
156156+ (* Pick a data point in the middle and interpolate at its exact epoch *)
157157+ let mid = svs.(5000) in
158158+ match Odm.interpolate seg with
159159+ | exception _ -> Alcotest.fail "interpolate raised exception"
160160+ | _interp ->
161161+ (* Just verify interpolation function exists and doesn't crash on this OEM *)
162162+ ()
163163+66164(* --- GEO 3-day ---
67165 Source: GMAT R2026a, geo_3day.script.
68166 GEO orbit: SMA=42164 km, ECC=0.0001, INC=0.05 deg. *)
···70168let test_geo_parse () =
71169 let oem = parse_oem "geo_3day.oem" in
72170 let sv = Odm.state_vectors (List.hd (Odm.segments oem)) in
7373- (* 3 days at 300s step = 864 + 1 vectors *)
74171 Alcotest.(check bool) "~864 vectors" true (Array.length sv > 800)
7517276173let test_geo_radius () =
···93190 Alcotest.failf "GEO velocity out of bounds: %.3f km/s at %s" v sv.epoch)
94191 svs
95192193193+(* GEO: near-circular orbit. Eccentricity < 0.001 means radius variation
194194+ should be tiny (< 50 km peak-to-peak). *)
195195+let test_geo_near_circular () =
196196+ let oem = parse_oem "geo_3day.oem" in
197197+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
198198+ let min_r = ref Float.infinity in
199199+ let max_r = ref Float.neg_infinity in
200200+ Array.iter
201201+ (fun sv ->
202202+ let r = vec3_norm sv.Odm.pos in
203203+ min_r := Float.min !min_r r;
204204+ max_r := Float.max !max_r r)
205205+ svs;
206206+ let variation = !max_r -. !min_r in
207207+ Alcotest.(check bool) "radius variation < 50 km" true (variation < 50.0);
208208+ Printf.printf " GEO radius variation: %.3f km\n" variation
209209+210210+let test_geo_no_nan_or_inf () =
211211+ let oem = parse_oem "geo_3day.oem" in
212212+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
213213+ Array.iteri
214214+ (fun i sv ->
215215+ let check name v =
216216+ if not (Float.is_finite v) then
217217+ Alcotest.failf "non-finite %s at vector %d: %g" name i v
218218+ in
219219+ check "X" sv.Odm.pos.x;
220220+ check "Y" sv.pos.y;
221221+ check "Z" sv.pos.z;
222222+ check "VX" sv.vel.x;
223223+ check "VY" sv.vel.y;
224224+ check "VZ" sv.vel.z)
225225+ svs
226226+96227(* --- Molniya 2-day ---
97228 Source: GMAT R2026a, molniya_2day.script.
9898- Molniya orbit: SMA=26600 km, ECC=0.74, INC=63.4 deg. *)
229229+ Molniya orbit: SMA=26600 km, ECC=0.74, INC=63.4 deg.
230230+ This is the hardest orbit to parse correctly: 5 orders of magnitude
231231+ in position values, near-zero Y components (1e-12 km), velocities
232232+ from 1.5 to 10 km/s. *)
99233100234let test_molniya_parse () =
101235 let oem = parse_oem "molniya_2day.oem" in
102236 let sv = Odm.state_vectors (List.hd (Odm.segments oem)) in
103103- (* 2 days at 30s step = 5760 + 1 vectors *)
104237 Alcotest.(check bool) "~5760 vectors" true (Array.length sv > 5700)
105238106239let test_molniya_radius_range () =
107240 let oem = parse_oem "molniya_2day.oem" in
108241 let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
109109- (* Molniya: perigee ~6916 km (SMA*(1-e)), apogee ~46244 km (SMA*(1+e)) *)
110242 let min_r = ref Float.infinity in
111243 let max_r = ref Float.neg_infinity in
112244 Array.iter
···120252 (!min_r > 6500.0 && !min_r < 7500.0);
121253 Alcotest.(check bool)
122254 "apogee ~46000 km" true
123123- (!max_r > 44000.0 && !max_r < 48000.0)
255255+ (!max_r > 44000.0 && !max_r < 48000.0);
256256+ Printf.printf " Molniya: perigee %.1f km, apogee %.1f km\n" !min_r !max_r
124257125258let test_molniya_velocity_range () =
126259 let oem = parse_oem "molniya_2day.oem" in
···140273 "max velocity ~10 km/s" true
141274 (!max_v > 8.0 && !max_v < 12.0)
142275276276+(* Molniya Y-position at epoch is ~1e-12 km (numerical noise from integrator).
277277+ Parser must handle exponents down to e-15 without underflow to zero. *)
278278+let test_molniya_near_zero_components () =
279279+ let oem = parse_oem "molniya_2day.oem" in
280280+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
281281+ let y0 = svs.(0).Odm.pos.y in
282282+ Alcotest.(check bool) "Y0 finite" true (Float.is_finite y0);
283283+ (* The value should be very small but not exactly zero *)
284284+ Alcotest.(check bool) "Y0 ~1e-12" true (Float.abs y0 < 1e-10)
285285+286286+(* All 5760+ vectors must be finite — no NaN from parsing 1e-15 exponents. *)
287287+let test_molniya_no_nan_or_inf () =
288288+ let oem = parse_oem "molniya_2day.oem" in
289289+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
290290+ Array.iteri
291291+ (fun i sv ->
292292+ let check name v =
293293+ if not (Float.is_finite v) then
294294+ Alcotest.failf "non-finite %s at vector %d (%s): %g" name i
295295+ sv.Odm.epoch v
296296+ in
297297+ check "X" sv.pos.x;
298298+ check "Y" sv.pos.y;
299299+ check "Z" sv.pos.z;
300300+ check "VX" sv.vel.x;
301301+ check "VY" sv.vel.y;
302302+ check "VZ" sv.vel.z)
303303+ svs
304304+305305+(* Vis-viva: v^2 = mu * (2/r - 1/a). For Molniya, SMA=26600 km.
306306+ At each point, check vis-viva holds within perturbation error. *)
307307+let test_molniya_vis_viva () =
308308+ let oem = parse_oem "molniya_2day.oem" in
309309+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
310310+ let mu = 398600.4418 in
311311+ (* Compute SMA from first state via vis-viva *)
312312+ let r0 = vec3_norm svs.(0).Odm.pos in
313313+ let v0 = vec3_norm svs.(0).vel in
314314+ let a0 = 1.0 /. ((2.0 /. r0) -. (v0 *. v0 /. mu)) in
315315+ (* SMA should be ~26600 km *)
316316+ Alcotest.(check bool) "SMA ~26600 km" true (a0 > 26000.0 && a0 < 27200.0);
317317+ Printf.printf " Molniya SMA from vis-viva: %.1f km (expected ~26600)\n" a0;
318318+ (* Check vis-viva consistency across all points (allow ~0.5% for J2 effects) *)
319319+ Array.iter
320320+ (fun sv ->
321321+ let r = vec3_norm sv.Odm.pos in
322322+ let v = vec3_norm sv.vel in
323323+ let a = 1.0 /. ((2.0 /. r) -. (v *. v /. mu)) in
324324+ let rel = Float.abs (a -. a0) /. Float.abs a0 in
325325+ if rel > 0.005 then
326326+ Alcotest.failf "vis-viva SMA drift > 0.5%%: %.1f vs %.1f (%.3f%%) at %s"
327327+ a a0 (rel *. 100.0) sv.epoch)
328328+ svs
329329+330330+(* Epoch monotonicity for Molniya (30s step — more data points to check). *)
331331+let test_molniya_epoch_monotonic () =
332332+ let oem = parse_oem "molniya_2day.oem" in
333333+ let svs = Odm.state_vectors (List.hd (Odm.segments oem)) in
334334+ for i = 1 to Array.length svs - 1 do
335335+ if String.compare svs.(i).Odm.epoch svs.(i - 1).epoch <= 0 then
336336+ Alcotest.failf "non-monotonic epoch at %d: %s <= %s" i svs.(i).epoch
337337+ svs.(i - 1).epoch
338338+ done
339339+143340let () =
144341 Alcotest.run "odm-gmat"
145342 [
···147344 [
148345 Alcotest.test_case "parse" `Quick test_leo_parse;
149346 Alcotest.test_case "metadata" `Quick test_leo_metadata;
347347+ Alcotest.test_case "originator with space" `Quick
348348+ test_leo_originator_with_space;
349349+ Alcotest.test_case "interpolation degree" `Quick
350350+ test_leo_interpolation_degree;
150351 Alcotest.test_case "radius bounds" `Quick test_leo_radius_bounds;
151352 Alcotest.test_case "velocity bounds" `Quick test_leo_velocity_bounds;
152353 Alcotest.test_case "epoch range" `Quick test_leo_epoch_range;
354354+ Alcotest.test_case "near-zero Vz (3.68e-16)" `Quick
355355+ test_leo_near_zero_velocity_component;
356356+ Alcotest.test_case "no NaN or inf" `Quick test_leo_no_nan_or_inf;
357357+ Alcotest.test_case "epoch monotonic" `Quick test_leo_epoch_monotonic;
358358+ Alcotest.test_case "energy conservation" `Quick
359359+ test_leo_energy_conservation;
360360+ Alcotest.test_case "interpolation" `Quick
361361+ test_leo_interpolation_at_epoch;
153362 ] );
154363 ( "geo",
155364 [
156365 Alcotest.test_case "parse" `Quick test_geo_parse;
157366 Alcotest.test_case "radius" `Quick test_geo_radius;
158367 Alcotest.test_case "velocity" `Quick test_geo_velocity;
368368+ Alcotest.test_case "near-circular" `Quick test_geo_near_circular;
369369+ Alcotest.test_case "no NaN or inf" `Quick test_geo_no_nan_or_inf;
159370 ] );
160371 ( "molniya",
161372 [
162373 Alcotest.test_case "parse" `Quick test_molniya_parse;
163374 Alcotest.test_case "radius range" `Quick test_molniya_radius_range;
164375 Alcotest.test_case "velocity range" `Quick test_molniya_velocity_range;
376376+ Alcotest.test_case "near-zero Y (1e-12)" `Quick
377377+ test_molniya_near_zero_components;
378378+ Alcotest.test_case "no NaN or inf" `Quick test_molniya_no_nan_or_inf;
379379+ Alcotest.test_case "vis-viva consistency" `Quick test_molniya_vis_viva;
380380+ Alcotest.test_case "epoch monotonic" `Quick
381381+ test_molniya_epoch_monotonic;
165382 ] );
166383 ]