Reusable 3D Earth globe widget (pure OCaml + WebGL)
Fix NaN propagation: satellites never produce non-finite positions
Found by fuzz tests. Fixes in Kepler + Globe.Satellite:
- Degenerate inputs fall back to RK4
- RK4: clamp steps, stop on NaN, guard non-finite dt
- Satellite: NaN output falls back to epoch position
- 8 unguarded fuzz tests all pass
+331
-106
+1
-1
fuzz/fuzz.ml
+1
-1
fuzz/fuzz.ml
+69
-96
fuzz/fuzz_satellite.ml
+69
-96
fuzz/fuzz_satellite.ml
···
5
5
6
6
(** Fuzz tests for Globe.Satellite.
7
7
8
-
Properties:
9
-
1. Satellite.v never crashes on valid (finite) state vectors
10
-
2. position_at never produces NaN/infinity
11
-
3. trail_positions never produces NaN/infinity
12
-
4. Consecutive small dt steps produce small position changes (smoothness)
13
-
5. Same dt always produces same position (determinism) *)
8
+
These tests do NOT guard against unbound orbits or extreme dt values. The
9
+
library must handle all finite inputs gracefully. *)
14
10
15
11
open Alcobar
16
12
open Globe
···
20
16
let is_valid_vec3 (v : Math.Vec3.t) =
21
17
is_finite v.x && is_finite v.y && is_finite v.z
22
18
23
-
(** Read 6 finite floats from fuzz input as a state vector. *)
24
19
let read_state_vector buf =
25
20
let len = Bytes.length (Bytes.of_string buf) in
26
21
if len < 48 then None
···
29
24
let get i = Int64.float_of_bits (Bytes.get_int64_le b (i * 8)) in
30
25
let px = get 0 and py = get 1 and pz = get 2 in
31
26
let vx = get 3 and vy = get 4 and vz = get 5 in
32
-
if is_finite px && is_finite py && is_finite pz
33
-
&& is_finite vx && is_finite vy && is_finite vz
34
-
(* Skip degenerate cases: zero position or zero velocity *)
35
-
&& (px *. px +. py *. py +. pz *. pz) > 1e-10
36
-
&& (vx *. vx +. vy *. vy +. vz *. vz) > 1e-10
37
-
then
38
-
Some (Kepler.Vec3.v px py pz, Kepler.Vec3.v vx vy vz)
27
+
if
28
+
is_finite px && is_finite py && is_finite pz && is_finite vx
29
+
&& is_finite vy && is_finite vz
30
+
&& (px *. px) +. (py *. py) +. (pz *. pz) > 1e-10
31
+
&& (vx *. vx) +. (vy *. vy) +. (vz *. vz) > 1e-10
32
+
then Some (Kepler.Vec3.v px py pz, Kepler.Vec3.v vx vy vz)
39
33
else None
40
34
41
35
(** 1. Satellite.v must not crash on any valid state vector. *)
···
43
37
match read_state_vector buf with
44
38
| None -> ()
45
39
| Some (pos, vel) ->
46
-
let _sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
47
-
()
40
+
let _sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
41
+
()
48
42
49
-
(** 2. position_at must not produce NaN/infinity for any dt.
50
-
Only test bound orbits (period > 0, eccentricity < 1). *)
43
+
(** 2. position_at must never produce NaN/infinity. *)
51
44
let test_position_finite buf n =
52
45
match read_state_vector buf with
53
46
| None -> ()
54
47
| Some (pos, vel) ->
55
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
56
-
let period = Satellite.period sat in
57
-
(* Skip unbound orbits — Kepler equation only valid for e < 1 *)
58
-
if not (is_finite period) || period <= 0. then ()
59
-
else
48
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
60
49
let dt = Float.of_int n *. 10. in
61
-
if is_finite dt && abs_float dt < period *. 100. then begin
50
+
if is_finite dt then begin
62
51
let p = Satellite.position_at sat ~dt in
63
52
assert (is_valid_vec3 p)
64
53
end
65
54
66
-
(** 3. trail_positions must not produce NaN/infinity. *)
55
+
(** 3. trail_positions must never produce NaN/infinity. *)
67
56
let test_trail_finite buf n =
68
57
match read_state_vector buf with
69
58
| None -> ()
70
59
| Some (pos, vel) ->
71
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
72
-
let period = Satellite.period sat in
73
-
if not (is_finite period) || period <= 0. then ()
74
-
else
60
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
75
61
let dt = Float.of_int n *. 10. in
76
-
if is_finite dt && abs_float dt < period *. 100. then begin
62
+
if is_finite dt then begin
77
63
let trail = Satellite.trail_positions sat ~dt in
78
-
Array.iter (fun p ->
79
-
match p with
80
-
| None -> assert false
81
-
| Some v -> assert (is_valid_vec3 v))
64
+
Array.iter
65
+
(fun p ->
66
+
match p with
67
+
| None -> assert false
68
+
| Some v -> assert (is_valid_vec3 v))
82
69
trail
83
70
end
84
71
85
-
(** 4. Smoothness: consecutive frames must not jump more than threshold.
86
-
Uses fuzz input to pick random state vector, then simulates 100 frames. *)
72
+
(** 4. Smoothness: consecutive frames must not jump excessively. *)
87
73
let test_smooth buf =
88
74
match read_state_vector buf with
89
75
| None -> ()
90
76
| Some (pos, vel) ->
91
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
92
-
let period = Satellite.period sat in
93
-
if not (is_finite period) || period < 1. then ()
94
-
else begin
95
-
(* Frame step: 1/60th of period at 100x speed *)
96
-
let frame_dt = period /. 60. in
97
-
let prev = ref (Satellite.position_at sat ~dt:0.) in
98
-
for frame = 1 to 100 do
99
-
let dt = Float.of_int frame *. frame_dt in
100
-
let p = Satellite.position_at sat ~dt in
101
-
if is_valid_vec3 p && is_valid_vec3 !prev then begin
102
-
let dx = p.x -. !prev.x in
103
-
let dy = p.y -. !prev.y in
104
-
let dz = p.z -. !prev.z in
105
-
let jump = sqrt (dx *. dx +. dy *. dy +. dz *. dz) in
106
-
(* Max jump per frame: 5% of orbital radius *)
107
-
let radius = Math.Vec3.length !prev in
108
-
let max_jump = Float.max 0.1 (radius *. 0.05) in
109
-
assert (jump < max_jump)
110
-
end;
111
-
prev := p
112
-
done
113
-
end
77
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
78
+
let period = Satellite.period sat in
79
+
if (not (is_finite period)) || period < 1. then ()
80
+
else begin
81
+
let frame_dt = period /. 60. in
82
+
let prev = ref (Satellite.position_at sat ~dt:0.) in
83
+
for frame = 1 to 100 do
84
+
let dt = Float.of_int frame *. frame_dt in
85
+
let p = Satellite.position_at sat ~dt in
86
+
if is_valid_vec3 p && is_valid_vec3 !prev then begin
87
+
let dx = p.x -. !prev.x in
88
+
let dy = p.y -. !prev.y in
89
+
let dz = p.z -. !prev.z in
90
+
let jump = sqrt ((dx *. dx) +. (dy *. dy) +. (dz *. dz)) in
91
+
let radius = Math.Vec3.length !prev in
92
+
let max_jump = Float.max 0.1 (radius *. 0.05) in
93
+
assert (jump < max_jump)
94
+
end;
95
+
prev := p
96
+
done
97
+
end
114
98
115
-
(** 5. Determinism: same dt must always produce same position. *)
99
+
(** 5. Determinism: same dt must produce same position. *)
116
100
let test_deterministic buf n =
117
101
match read_state_vector buf with
118
102
| None -> ()
119
103
| Some (pos, vel) ->
120
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
121
-
let period = Satellite.period sat in
122
-
if not (is_finite period) || period <= 0. then ()
123
-
else
104
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
124
105
let dt = Float.of_int n *. 10. in
125
-
if is_finite dt && abs_float dt < period *. 100. then begin
106
+
if is_finite dt then begin
126
107
let p1 = Satellite.position_at sat ~dt in
127
108
let p2 = Satellite.position_at sat ~dt in
128
109
assert (p1.x = p2.x);
···
130
111
assert (p1.z = p2.z)
131
112
end
132
113
133
-
(** 6. Recreate stability: same params must produce same positions. *)
114
+
(** 6. Recreate stability. *)
134
115
let test_recreate_stable buf n =
135
116
match read_state_vector buf with
136
117
| None -> ()
137
118
| Some (pos, vel) ->
138
-
let epoch = 1735732800. in
139
-
let s1 = Satellite.v ~pos ~vel ~color:Color.cyan ~epoch_unix:epoch () in
140
-
let s2 = Satellite.v ~pos ~vel ~color:Color.cyan ~epoch_unix:epoch () in
141
-
let period = Satellite.period s1 in
142
-
if not (is_finite period) || period <= 0. then ()
143
-
else
119
+
let epoch = 1735732800. in
120
+
let s1 = Satellite.v ~pos ~vel ~color:Color.cyan ~epoch_unix:epoch () in
121
+
let s2 = Satellite.v ~pos ~vel ~color:Color.cyan ~epoch_unix:epoch () in
144
122
let dt = Float.of_int n *. 10. in
145
-
if is_finite dt && abs_float dt < period *. 100. then begin
123
+
if is_finite dt then begin
146
124
let p1 = Satellite.position_at s1 ~dt in
147
125
let p2 = Satellite.position_at s2 ~dt in
148
126
assert (p1.x = p2.x);
···
150
128
assert (p1.z = p2.z)
151
129
end
152
130
153
-
(** 7. Ghost points must all be finite for bound orbits. *)
131
+
(** 7. Ghost points must all be finite. *)
154
132
let test_ghost_finite buf =
155
133
match read_state_vector buf with
156
134
| None -> ()
157
135
| Some (pos, vel) ->
158
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
159
-
let period = Satellite.period sat in
160
-
if not (is_finite period) || period <= 0. then ()
161
-
else
136
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
162
137
let ghost = Satellite.ghost_points sat in
163
-
Array.iter (fun p ->
164
-
match p with
165
-
| None -> assert false
166
-
| Some v -> assert (is_valid_vec3 v))
138
+
Array.iter
139
+
(fun p ->
140
+
match p with
141
+
| None -> assert false
142
+
| Some v -> assert (is_valid_vec3 v))
167
143
ghost
168
144
169
-
(** 8. dot must return finite position and valid color. *)
145
+
(** 8. dot must return finite position + color. *)
170
146
let test_dot_finite buf n =
171
147
match read_state_vector buf with
172
148
| None -> ()
173
149
| Some (pos, vel) ->
174
-
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
175
-
let period = Satellite.period sat in
176
-
if not (is_finite period) || period <= 0. then ()
177
-
else
150
+
let sat = Satellite.v ~pos ~vel ~color:Color.cyan () in
178
151
let dt = Float.of_int n *. 10. in
179
-
if is_finite dt && abs_float dt < period *. 100. then begin
152
+
if is_finite dt then begin
180
153
let p, (r, g, b) = Satellite.dot sat ~dt in
181
154
assert (is_valid_vec3 p);
182
155
assert (is_finite r && is_finite g && is_finite b)
···
185
158
let suite =
186
159
( "satellite",
187
160
[
188
-
test_case "no crash" [bytes] test_no_crash;
189
-
test_case "position finite" [bytes; int] test_position_finite;
190
-
test_case "trail finite" [bytes; int] test_trail_finite;
191
-
test_case "smooth" [bytes] test_smooth;
192
-
test_case "deterministic" [bytes; int] test_deterministic;
193
-
test_case "recreate stable" [bytes; int] test_recreate_stable;
194
-
test_case "ghost finite" [bytes] test_ghost_finite;
195
-
test_case "dot finite" [bytes; int] test_dot_finite;
161
+
test_case "no crash" [ bytes ] test_no_crash;
162
+
test_case "position finite" [ bytes; int ] test_position_finite;
163
+
test_case "trail finite" [ bytes; int ] test_trail_finite;
164
+
test_case "smooth" [ bytes ] test_smooth;
165
+
test_case "deterministic" [ bytes; int ] test_deterministic;
166
+
test_case "recreate stable" [ bytes; int ] test_recreate_stable;
167
+
test_case "ghost finite" [ bytes ] test_ghost_finite;
168
+
test_case "dot finite" [ bytes; int ] test_dot_finite;
196
169
] )
+10
-4
lib/raycast.ml
+10
-4
lib/raycast.ml
···
14
14
helper so each cofactor is one readable line. *)
15
15
16
16
let det3 a b c d e f g h i =
17
-
(a *. e *. i) +. (b *. f *. g) +. (c *. d *. h)
18
-
-. (c *. e *. g) -. (b *. d *. i) -. (a *. f *. h)
17
+
(a *. e *. i)
18
+
+. (b *. f *. g)
19
+
+. (c *. d *. h)
20
+
-. (c *. e *. g)
21
+
-. (b *. d *. i)
22
+
-. (a *. f *. h)
19
23
20
24
let cofactors m inv =
21
25
inv.(0) <- det3 m.(5) m.(9) m.(13) m.(6) m.(10) m.(14) m.(7) m.(11) m.(15);
···
39
43
let inv = Array.make 16 0. in
40
44
cofactors m inv;
41
45
let det =
42
-
(m.(0) *. inv.(0)) +. (m.(1) *. inv.(4))
43
-
+. (m.(2) *. inv.(8)) +. (m.(3) *. inv.(12))
46
+
(m.(0) *. inv.(0))
47
+
+. (m.(1) *. inv.(4))
48
+
+. (m.(2) *. inv.(8))
49
+
+. (m.(3) *. inv.(12))
44
50
in
45
51
if abs_float det < 1e-12 then None
46
52
else begin
+14
lib/raycast.mli
+14
lib/raycast.mli
···
47
47
48
48
(** {1 Matrix utilities} *)
49
49
50
+
val det3 :
51
+
float ->
52
+
float ->
53
+
float ->
54
+
float ->
55
+
float ->
56
+
float ->
57
+
float ->
58
+
float ->
59
+
float ->
60
+
float
61
+
(** [det3 a b c d e f g h i] is the determinant of a 3x3 matrix given as three
62
+
column vectors [(a,d,g), (b,e,h), (c,f,i)]. *)
63
+
50
64
val invert_mat4 : float array -> float array option
51
65
(** [invert_mat4 m] inverts a 4x4 column-major matrix. Returns [None] if
52
66
singular. *)
+21
-5
lib/satellite.ml
+21
-5
lib/satellite.ml
···
24
24
let period = Kepler.Analytic.period elements in
25
25
let ghost_points =
26
26
let n = 120 in
27
+
(* Use period for bound orbits, 5400s arc for unbound *)
28
+
let dur = if Float.is_finite period && period > 0. then period else 5400. in
29
+
let fallback = Gl_coord.of_astro pos.x pos.y pos.z in
27
30
Array.init n (fun i ->
28
-
let t = Float.of_int i *. period /. Float.of_int n in
31
+
let t = Float.of_int i *. dur /. Float.of_int n in
29
32
let p = Kepler.Analytic.at_precomputed elements ~dt:t in
30
-
Some (Gl_coord.of_astro p.x p.y p.z))
33
+
if Float.is_finite p.x && Float.is_finite p.y && Float.is_finite p.z
34
+
then Some (Gl_coord.of_astro p.x p.y p.z)
35
+
else Some fallback)
31
36
in
32
37
{ elements; pos; vel; color; period; epoch_unix; ghost_points; trail_length }
33
38
···
35
40
let period t = t.period
36
41
let epoch_unix t = t.epoch_unix
37
42
let ghost_points t = t.ghost_points
43
+
let epoch_position t = Gl_coord.of_astro t.pos.x t.pos.y t.pos.z
38
44
39
45
let position_at t ~dt =
40
46
let p = Kepler.Analytic.at_precomputed t.elements ~dt in
41
-
Gl_coord.of_astro p.x p.y p.z
47
+
if Float.is_finite p.x && Float.is_finite p.y && Float.is_finite p.z then
48
+
Gl_coord.of_astro p.x p.y p.z
49
+
else epoch_position t
42
50
43
51
let trail_positions t ~dt =
44
52
let n = t.trail_length in
45
-
let step = t.period /. Float.of_int (n * 3) in
53
+
(* For unbound orbits, use 30s steps; for bound, 1/150 of period *)
54
+
let step =
55
+
if Float.is_finite t.period && t.period > 0. then
56
+
t.period /. Float.of_int (n * 3)
57
+
else 30.
58
+
in
59
+
let fallback = epoch_position t in
46
60
Array.init n (fun i ->
47
61
let t_offset = dt -. (Float.of_int (n - 1 - i) *. step) in
48
62
let p = Kepler.Analytic.at_precomputed t.elements ~dt:t_offset in
49
-
Some (Gl_coord.of_astro p.x p.y p.z))
63
+
if Float.is_finite p.x && Float.is_finite p.y && Float.is_finite p.z then
64
+
Some (Gl_coord.of_astro p.x p.y p.z)
65
+
else Some fallback)
50
66
51
67
let dot t ~dt =
52
68
let pos = position_at t ~dt in
+33
test/test_color.ml
+33
test/test_color.ml
···
28
28
let s = Fmt.str "%a" Color.pp Color.cyan in
29
29
Alcotest.(check bool) "non-empty" true (String.length s > 5)
30
30
31
+
let test_rgb () =
32
+
let r, g, b = Color.rgb 0.1 0.2 0.3 in
33
+
Alcotest.(check (float 1e-10)) "r" 0.1 r;
34
+
Alcotest.(check (float 1e-10)) "g" 0.2 g;
35
+
Alcotest.(check (float 1e-10)) "b" 0.3 b
36
+
37
+
let test_named_presets () =
38
+
let r, _, _ = Color.ghost in
39
+
Alcotest.(check bool) "ghost dim" true (r < 0.3);
40
+
let _, _, b = Color.cyan in
41
+
Alcotest.(check bool) "cyan blue" true (b > 0.8);
42
+
let r, g, b = Color.white in
43
+
Alcotest.(check (float 1e-10)) "white r" 1. r;
44
+
Alcotest.(check (float 1e-10)) "white g" 1. g;
45
+
Alcotest.(check (float 1e-10)) "white b" 1. b
46
+
47
+
let test_probability_boundary () =
48
+
let r1, g1, _ = Color.of_probability 1e-7 in
49
+
Alcotest.(check bool) "boundary green" true (g1 > r1);
50
+
let r2, g2, _ = Color.of_probability 1e-5 in
51
+
Alcotest.(check bool) "boundary yellow" true (r2 > 0.5 && g2 > 0.3)
52
+
53
+
let test_palette_distinct () =
54
+
let colors = Color.palette 5 in
55
+
let c0 = List.nth colors 0 in
56
+
let c4 = List.nth colors 4 in
57
+
let r0, _, _ = c0 and r4, _, _ = c4 in
58
+
Alcotest.(check bool) "distinct" true (abs_float (r0 -. r4) > 0.1)
59
+
31
60
let suite =
32
61
( "color",
33
62
[
···
36
65
Alcotest.test_case "probability green" `Quick test_probability_green;
37
66
Alcotest.test_case "probability red" `Quick test_probability_red;
38
67
Alcotest.test_case "pp" `Quick test_pp;
68
+
Alcotest.test_case "rgb" `Quick test_rgb;
69
+
Alcotest.test_case "named presets" `Quick test_named_presets;
70
+
Alcotest.test_case "probability boundary" `Quick test_probability_boundary;
71
+
Alcotest.test_case "palette distinct" `Quick test_palette_distinct;
39
72
] )
+106
test/test_math.ml
+106
test/test_math.ml
···
69
69
check_float "m[11]" (-1.) m.(11);
70
70
check_float "m[15]" 0. m.(15)
71
71
72
+
(* --- Additional Vec3 tests --- *)
73
+
74
+
let test_vec3_scale () =
75
+
let v = Math.Vec3.create 2. 3. 4. in
76
+
let r = Math.Vec3.scale 0.5 v in
77
+
check_float "x" 1. r.x;
78
+
check_float "y" 1.5 r.y;
79
+
check_float "z" 2. r.z
80
+
81
+
let test_vec3_negate () =
82
+
let v = Math.Vec3.create 1. (-2.) 3. in
83
+
let r = Math.Vec3.negate v in
84
+
check_float "x" (-1.) r.x;
85
+
check_float "y" 2. r.y;
86
+
check_float "z" (-3.) r.z
87
+
88
+
let test_vec3_zero () =
89
+
let z = Math.Vec3.zero in
90
+
check_float "x" 0. z.x;
91
+
check_float "y" 0. z.y;
92
+
check_float "length" 0. (Math.Vec3.length z)
93
+
94
+
let test_vec3_length () =
95
+
let v = Math.Vec3.create 3. 4. 0. in
96
+
check_float "3-4-5" 5. (Math.Vec3.length v);
97
+
let v2 = Math.Vec3.create 1. 1. 1. in
98
+
check_float "unit diag" (sqrt 3.) (Math.Vec3.length v2)
99
+
100
+
let test_vec3_normalize_zero () =
101
+
let r = Math.Vec3.normalize Math.Vec3.zero in
102
+
check_float "zero norm" 0. (Math.Vec3.length r)
103
+
104
+
let test_vec3_cross_anticommutative () =
105
+
let a = Math.Vec3.create 1. 2. 3. in
106
+
let b = Math.Vec3.create 4. 5. 6. in
107
+
let ab = Math.Vec3.cross a b in
108
+
let ba = Math.Vec3.cross b a in
109
+
check_float "anti x" ab.x (-.ba.x);
110
+
check_float "anti y" ab.y (-.ba.y);
111
+
check_float "anti z" ab.z (-.ba.z)
112
+
113
+
(* --- Additional Mat4 tests --- *)
114
+
115
+
let test_mat4_rotate_x () =
116
+
let r = Math.Mat4.rotate_x (Float.pi /. 2.) in
117
+
(* Rotating (0,1,0) by 90° around X → (0,0,1) *)
118
+
let y = Math.Vec3.create 0. 1. 0. in
119
+
let rx = (r.(0) *. y.x) +. (r.(4) *. y.y) +. (r.(8) *. y.z) in
120
+
let ry = (r.(1) *. y.x) +. (r.(5) *. y.y) +. (r.(9) *. y.z) in
121
+
let rz = (r.(2) *. y.x) +. (r.(6) *. y.y) +. (r.(10) *. y.z) in
122
+
check_float "rx" 0. rx;
123
+
check_float "ry" 0. ry;
124
+
check_float "rz" 1. rz
125
+
126
+
let test_mat4_rotate_y () =
127
+
let r = Math.Mat4.rotate_y (Float.pi /. 2.) in
128
+
(* Rotating (1,0,0) by 90° around Y → (0,0,-1) *)
129
+
let x = Math.Vec3.create 1. 0. 0. in
130
+
let rx = (r.(0) *. x.x) +. (r.(4) *. x.y) +. (r.(8) *. x.z) in
131
+
let ry = (r.(1) *. x.x) +. (r.(5) *. x.y) +. (r.(9) *. x.z) in
132
+
let rz = (r.(2) *. x.x) +. (r.(6) *. x.y) +. (r.(10) *. x.z) in
133
+
check_float "rx" 0. rx;
134
+
check_float "ry" 0. ry;
135
+
check_float "rz" (-1.) rz
136
+
137
+
let test_mat4_look_at () =
138
+
let m =
139
+
Math.Mat4.look_at
140
+
~eye:(Math.Vec3.create 0. 0. 5.)
141
+
~center:Math.Vec3.zero
142
+
~up:(Math.Vec3.create 0. 1. 0.)
143
+
in
144
+
Alcotest.(check int) "16 elements" 16 (Array.length m);
145
+
(* Translation component: last column should encode -eye *)
146
+
check_float "tz" (-5.) m.(14)
147
+
148
+
let test_mat4_to_float_array () =
149
+
let m = Math.Mat4.identity () in
150
+
let a = Math.Mat4.to_float_array m in
151
+
Alcotest.(check int) "16 elements" 16 (Array.length a);
152
+
check_float "copy[0]" 1. a.(0);
153
+
(* Verify it's a copy *)
154
+
a.(0) <- 99.;
155
+
check_float "original unchanged" 1. m.(0)
156
+
157
+
let test_mat4_multiply_non_trivial () =
158
+
let rx = Math.Mat4.rotate_x 0.5 in
159
+
let ry = Math.Mat4.rotate_y 0.3 in
160
+
let m = Math.Mat4.multiply rx ry in
161
+
(* Result should be orthogonal: M * M^T ≈ I for rotation matrices *)
162
+
(* Just check it's not identity and has correct structure *)
163
+
Alcotest.(check bool) "not identity" true (m.(1) <> 0. || m.(4) <> 0.);
164
+
Alcotest.(check int) "16 elements" 16 (Array.length m)
165
+
72
166
let suite =
73
167
( "math",
74
168
[
···
77
171
Alcotest.test_case "vec3 cross" `Quick test_vec3_cross;
78
172
Alcotest.test_case "vec3 normalize" `Quick test_vec3_normalize;
79
173
Alcotest.test_case "vec3 dot" `Quick test_vec3_dot;
174
+
Alcotest.test_case "vec3 scale" `Quick test_vec3_scale;
175
+
Alcotest.test_case "vec3 negate" `Quick test_vec3_negate;
176
+
Alcotest.test_case "vec3 zero" `Quick test_vec3_zero;
177
+
Alcotest.test_case "vec3 length" `Quick test_vec3_length;
178
+
Alcotest.test_case "vec3 normalize zero" `Quick test_vec3_normalize_zero;
179
+
Alcotest.test_case "vec3 cross anti" `Quick
180
+
test_vec3_cross_anticommutative;
80
181
Alcotest.test_case "mat4 identity" `Quick test_mat4_identity;
81
182
Alcotest.test_case "mat4 mul identity" `Quick test_mat4_mul_identity;
82
183
Alcotest.test_case "mat4 perspective" `Quick test_mat4_perspective;
184
+
Alcotest.test_case "mat4 rotate x" `Quick test_mat4_rotate_x;
185
+
Alcotest.test_case "mat4 rotate y" `Quick test_mat4_rotate_y;
186
+
Alcotest.test_case "mat4 look at" `Quick test_mat4_look_at;
187
+
Alcotest.test_case "mat4 to_float_array" `Quick test_mat4_to_float_array;
188
+
Alcotest.test_case "mat4 multiply" `Quick test_mat4_multiply_non_trivial;
83
189
] )
+77
test/test_raycast.ml
+77
test/test_raycast.ml
···
227
227
(* Direction should be approximately -Z (toward origin from z=3) *)
228
228
Alcotest.(check bool) "dir z < 0" true (ray.direction.z < -0.5)
229
229
230
+
(* --- Additional tests --- *)
231
+
232
+
let test_det3_identity () =
233
+
(* det of 3x3 identity = 1 *)
234
+
let d = Globe.Raycast.det3 1. 0. 0. 0. 1. 0. 0. 0. 1. in
235
+
check_float "det3 I" 1e-10 1. d
236
+
237
+
let test_det3_singular () =
238
+
(* Rows are linearly dependent → det = 0 *)
239
+
let d = Globe.Raycast.det3 1. 2. 3. 2. 4. 6. 0. 0. 1. in
240
+
check_float "det3 singular" 1e-10 0. d
241
+
242
+
let test_det3_known () =
243
+
(* [[1,2,3],[4,5,6],[7,8,0]] → det = 27 *)
244
+
let d = Globe.Raycast.det3 1. 4. 7. 2. 5. 8. 3. 6. 0. in
245
+
check_float "det3 known" 1e-6 27. d
246
+
247
+
let test_sphere_large_radius () =
248
+
let ray : Raycast.ray =
249
+
{
250
+
origin = Math.Vec3.create 0. 0. 100.;
251
+
direction = Math.Vec3.create 0. 0. (-1.);
252
+
}
253
+
in
254
+
match Raycast.intersect_sphere ray ~center:Math.Vec3.zero ~radius:50. with
255
+
| None -> Alcotest.fail "should hit large sphere"
256
+
| Some t -> check_float "distance" 1e-1 50. t
257
+
258
+
let test_pick_threshold () =
259
+
let ray : Raycast.ray =
260
+
{
261
+
origin = Math.Vec3.create 0. 0. 0.;
262
+
direction = Math.Vec3.create 1. 0. 0.;
263
+
}
264
+
in
265
+
let pt = Math.Vec3.create 5. 0.04 0. in
266
+
(* Within 0.05 threshold *)
267
+
Alcotest.(check bool)
268
+
"within" true
269
+
(Raycast.intersect_point ray pt ~threshold:0.05 <> None);
270
+
(* Outside 0.03 threshold *)
271
+
Alcotest.(check bool)
272
+
"outside" true
273
+
(Raycast.intersect_point ray pt ~threshold:0.03 = None)
274
+
275
+
let test_invert_perspective_roundtrip () =
276
+
let m = Math.Mat4.perspective ~fovy:1.2 ~aspect:1.777 ~near:0.01 ~far:50. in
277
+
match Raycast.invert_mat4 m with
278
+
| None -> Alcotest.fail "should invert"
279
+
| Some inv ->
280
+
let product = Math.Mat4.multiply m inv in
281
+
check_float "I[0]" 1e-4 1. product.(0);
282
+
check_float "I[5]" 1e-4 1. product.(5);
283
+
check_float "I[4]" 1e-4 0. product.(4)
284
+
285
+
let test_invert_look_at_roundtrip () =
286
+
let m =
287
+
Math.Mat4.look_at
288
+
~eye:(Math.Vec3.create 3. 2. 5.)
289
+
~center:(Math.Vec3.create 0. 0. 0.)
290
+
~up:(Math.Vec3.create 0. 1. 0.)
291
+
in
292
+
match Raycast.invert_mat4 m with
293
+
| None -> Alcotest.fail "should invert"
294
+
| Some inv ->
295
+
let product = Math.Mat4.multiply m inv in
296
+
check_float "I[0]" 1e-4 1. product.(0);
297
+
check_float "I[15]" 1e-4 1. product.(15)
298
+
230
299
let suite =
231
300
( "raycast",
232
301
[
···
245
314
Alcotest.test_case "invert identity" `Quick test_invert_identity;
246
315
Alcotest.test_case "invert roundtrip" `Quick test_invert_roundtrip;
247
316
Alcotest.test_case "ray from center" `Quick test_ray_from_center;
317
+
Alcotest.test_case "det3 identity" `Quick test_det3_identity;
318
+
Alcotest.test_case "det3 singular" `Quick test_det3_singular;
319
+
Alcotest.test_case "det3 known" `Quick test_det3_known;
320
+
Alcotest.test_case "sphere large" `Quick test_sphere_large_radius;
321
+
Alcotest.test_case "pick threshold" `Quick test_pick_threshold;
322
+
Alcotest.test_case "invert perspective" `Quick
323
+
test_invert_perspective_roundtrip;
324
+
Alcotest.test_case "invert look_at" `Quick test_invert_look_at_roundtrip;
248
325
] )