1use super::{ArrayVec, Point, PointDot, PointIndex, PointNorm};
4use num_traits::{Float, NumCast};
5
6const DEFAULT_LENGTH_STEPS: usize = 64;
7
8#[derive(Copy, Clone, Debug, PartialEq)]
16pub struct LineSegment<P> {
17 pub(crate) start: P,
18 pub(crate) end: P,
19}
20
21impl<P> LineSegment<P>
22where
23 P: Point,
24{
25 pub fn new(start: P, end: P) -> Self {
27 LineSegment { start, end }
28 }
29
30 pub fn start(&self) -> P {
32 let zero = <P::Scalar as NumCast>::from(0.0).unwrap();
33 self.eval(zero)
34 }
35
36 pub fn end(&self) -> P {
38 let one = <P::Scalar as NumCast>::from(1.0).unwrap();
39 self.eval(one)
40 }
41
42 pub fn reverse(&self) -> Self {
44 LineSegment {
45 start: self.end,
46 end: self.start,
47 }
48 }
49
50 pub fn eval(&self, t: P::Scalar) -> P {
52 self.start + (self.end - self.start) * t
53 }
54
55 pub fn split(&self, t: P::Scalar) -> (Self, Self) {
57 let ctrl_ab = self.start + (self.end - self.start) * t;
59
60 (
61 LineSegment {
62 start: self.start,
63 end: ctrl_ab,
64 },
65 LineSegment {
66 start: ctrl_ab,
67 end: self.end,
68 },
69 )
70 }
71
72 pub fn distance_to_point(&self, p: P) -> P::Scalar
74 where
75 P: PointNorm,
76 {
77 let l2 = (self.end - self.start).squared_norm();
78 if l2 < P::Scalar::epsilon() {
80 (self.start - p).squared_norm().sqrt()
81 } else {
82 let v1 = p - self.start;
83 let v2 = self.end - self.start;
84 let dot = PointDot::dot(&v1, &v2);
85 let t = (dot / l2).clamp(
91 <P::Scalar as NumCast>::from(0.0).unwrap(),
92 <P::Scalar as NumCast>::from(1.0).unwrap(),
93 );
94 let projection = self.start + (self.end - self.start) * t; (p - projection).squared_norm().sqrt()
97 }
98 }
99
100 pub fn derivative(&self) -> P {
102 self.end - self.start
103 }
104
105 pub fn tangent(&self, _t: P::Scalar) -> P
107 where
108 P: PointNorm,
109 {
110 let one = <P::Scalar as NumCast>::from(1.0).unwrap();
111 let dir = self.derivative();
112 let len = dir.squared_norm().sqrt();
113 if len <= P::Scalar::epsilon() {
114 dir
115 } else {
116 dir * (one / len)
117 }
118 }
119
120 pub fn curvature(&self, _t: P::Scalar) -> P::Scalar
124 where
125 P: PointNorm,
126 {
127 <P::Scalar as NumCast>::from(0.0).unwrap()
128 }
129
130 pub fn normal(&self, _t: P::Scalar) -> Option<P>
134 where
135 P: PointNorm,
136 {
137 None
138 }
139
140 pub fn t_at_length_approx(&self, s: P::Scalar, _nsteps: usize) -> P::Scalar
142 where
143 P: PointNorm,
144 {
145 let zero = <P::Scalar as NumCast>::from(0.0).unwrap();
146 let one = <P::Scalar as NumCast>::from(1.0).unwrap();
147 let length = (self.end - self.start).squared_norm().sqrt();
148 if length <= P::Scalar::epsilon() {
149 zero
150 } else {
151 (s / length).clamp(zero, one)
152 }
153 }
154
155 pub fn t_at_length(&self, s: P::Scalar) -> P::Scalar
157 where
158 P: PointNorm,
159 {
160 self.t_at_length_approx(s, DEFAULT_LENGTH_STEPS)
161 }
162
163 pub fn point_at_length_approx(&self, s: P::Scalar, nsteps: usize) -> P
165 where
166 P: PointNorm,
167 {
168 self.eval(self.t_at_length_approx(s, nsteps))
169 }
170
171 pub fn point_at_length(&self, s: P::Scalar) -> P
173 where
174 P: PointNorm,
175 {
176 self.point_at_length_approx(s, DEFAULT_LENGTH_STEPS)
177 }
178
179 pub(crate) fn root(&self, a: P::Scalar, b: P::Scalar) -> ArrayVec<[P::Scalar; 1]>
180 where
181 [P::Scalar; 1]: tinyvec::Array<Item = P::Scalar>,
182 {
183 let mut r = ArrayVec::new();
184 if a.abs() < P::Scalar::epsilon() {
185 return r;
186 }
187 r.push(-b / a);
188 r
189 }
190
191 pub fn bounding_box(&self) -> [(P::Scalar, P::Scalar); P::DIM]
193 where
194 P: PointIndex,
195 {
196 let mut bounds = [(
197 <P::Scalar as NumCast>::from(0.0).unwrap(),
198 <P::Scalar as NumCast>::from(0.0).unwrap(),
199 ); P::DIM];
200
201 for i in 0..P::DIM {
203 if self.start[i] < self.end[i] {
204 bounds[i] = (self.start[i], self.end[i]);
205 } else {
206 bounds[i] = (self.end[i], self.start[i]);
207 }
208 }
209
210 bounds
211 }
212}
213
214#[cfg(test)]
215mod tests {
216 use super::*;
217 use crate::{EPSILON, PointN};
218 #[test]
221 fn line_segment_interpolation() {
222 let line = LineSegment {
223 start: PointN::new([0f64, 1.77f64]),
224 end: PointN::new([4.3f64, 3f64]),
225 };
226
227 let mid = line.eval(0.5);
228 let diff = ((mid - line.start).squared_norm() - (mid - line.end).squared_norm()).abs();
229 assert!(diff < EPSILON * 10.0);
230 }
231
232 #[test]
234 fn line_segment_distance_to_point() {
235 let line = LineSegment {
237 start: PointN::new([0f64, 1f64, 0f64]),
238 end: PointN::new([3f64, 1f64, 0f64]),
239 };
240 let p1 = PointN::new([0f64, 5f64, 0f64]);
242 assert!((line.distance_to_point(p1) - 4.0).abs() < EPSILON);
243 assert!(((p1 - line.start).squared_norm().sqrt() - 4.0).abs() < EPSILON);
244 assert!(((p1 - line.end).squared_norm().sqrt() - 5.0).abs() < EPSILON);
245 let p2 = PointN::new([1.5f64, 2f64, 0f64]);
247 assert!(((p2 - line.eval(0.5)).squared_norm().sqrt() - 1.0).abs() < EPSILON);
248 }
249
250 #[test]
251 fn line_segment_split_midpoint() {
252 let line = LineSegment {
253 start: PointN::new([0f64, 0f64]),
254 end: PointN::new([4f64, 0f64]),
255 };
256
257 let (left, right) = line.split(0.5);
258 assert_eq!(left.start, PointN::new([0.0, 0.0]));
259 assert_eq!(left.end, PointN::new([2.0, 0.0]));
260 assert_eq!(right.start, PointN::new([2.0, 0.0]));
261 assert_eq!(right.end, PointN::new([4.0, 0.0]));
262 }
263
264 #[test]
265 fn line_segment_distance_to_point_after_end() {
266 let line = LineSegment {
267 start: PointN::new([0f64, 0f64]),
268 end: PointN::new([1f64, 0f64]),
269 };
270 let p = PointN::new([2f64, 0f64]);
271 assert!((line.distance_to_point(p) - 1.0).abs() < EPSILON);
272 }
273
274 #[test]
275 fn line_segment_api_parity() {
276 let line = LineSegment::new(PointN::new([0.0, 0.0]), PointN::new([2.0, 0.0]));
277
278 assert_eq!(line.start(), line.eval(0.0));
279 assert_eq!(line.end(), line.eval(1.0));
280
281 let reversed = line.reverse();
282 assert_eq!(reversed.start(), line.end());
283 assert_eq!(reversed.end(), line.start());
284
285 let tangent = line.tangent(0.3);
286 assert!((tangent[0] - 1.0).abs() < EPSILON);
287 assert!(tangent[1].abs() < EPSILON);
288
289 let curvature = line.curvature(0.3);
290 assert!(curvature.abs() < EPSILON);
291 assert!(line.normal(0.3).is_none());
292
293 let t = line.t_at_length(1.0);
294 assert!((t - 0.5).abs() < EPSILON);
295
296 let p = line.point_at_length(1.0);
297 assert!((p - PointN::new([1.0, 0.0])).squared_norm() < EPSILON);
298 }
299}