s2json_core/geometry/
vector_point.rs

1use crate::*;
2use core::{
3    cmp::Ordering,
4    f64::consts::PI,
5    fmt::Debug,
6    ops::{Add, AddAssign, Div, DivAssign, Mul, MulAssign, Neg, Rem, RemAssign, Sub, SubAssign},
7};
8use libm::{atan, atan2, fabs, fmod, log, sin, sinh, sqrt};
9use serde::{Deserialize, Serialize};
10
11/// A Vector Point uses a structure for 2D or 3D points
12#[derive(Serialize, Deserialize, Debug, Clone, Default)]
13#[repr(C)]
14pub struct VectorPoint<M: Clone = MValue> {
15    /// X coordinate
16    pub x: f64,
17    /// Y coordinate
18    pub y: f64,
19    /// Z coordinate or "altitude". May be None
20    #[serde(skip_serializing_if = "Option::is_none")]
21    pub z: Option<f64>,
22    /// M-Value
23    #[serde(skip_serializing_if = "Option::is_none")]
24    pub m: Option<M>,
25    /// T for tolerance. A tmp value used for simplification
26    #[serde(skip)]
27    pub t: Option<f64>,
28}
29// GET
30impl<M: Clone> GetXY for VectorPoint<M> {
31    fn x(&self) -> f64 {
32        self.x
33    }
34    fn y(&self) -> f64 {
35        self.y
36    }
37}
38impl<M: Clone> GetZ for VectorPoint<M> {
39    fn z(&self) -> Option<f64> {
40        self.z
41    }
42}
43impl<M: Clone> GetM<M> for VectorPoint<M> {
44    fn m(&self) -> Option<&M> {
45        self.m.as_ref()
46    }
47}
48// SET
49impl<M: Clone> SetXY for VectorPoint<M> {
50    fn set_xy(&mut self, x: f64, y: f64) {
51        self.x = x;
52        self.y = y;
53    }
54
55    fn set_x(&mut self, x: f64) {
56        self.x = x;
57    }
58
59    fn set_y(&mut self, y: f64) {
60        self.y = y;
61    }
62}
63impl<M: Clone> SetZ for VectorPoint<M> {
64    fn set_z(&mut self, z: f64) {
65        self.z = Some(z);
66    }
67}
68impl<M: Clone> SetM<M> for VectorPoint<M> {
69    fn set_m(&mut self, m: M) {
70        self.m = Some(m);
71    }
72}
73// NEW
74impl<M: Clone> NewXY for VectorPoint<M> {
75    fn new_xy(x: f64, y: f64) -> Self {
76        Self { x, y, z: None, m: None, t: None }
77    }
78}
79impl<M: Clone> NewXYZ for VectorPoint<M> {
80    fn new_xyz(x: f64, y: f64, z: f64) -> Self {
81        Self { x, y, z: Some(z), m: None, t: None }
82    }
83}
84impl<M: Clone> NewXYZM<M> for VectorPoint<M> {
85    fn new_xyzm(x: f64, y: f64, z: f64, m: M) -> Self {
86        Self { x, y, z: Some(z), m: Some(m), t: None }
87    }
88}
89// MValue METHODS
90impl VectorPoint<MValue> {
91    /// Helper function for tests. Create a new point quickly from an xy coordinate
92    pub fn from_xy(x: f64, y: f64) -> Self {
93        Self { x, y, z: None, m: None, t: None }
94    }
95
96    /// Helper function for tests. Create a new point quickly from an xyz coordinate
97    pub fn from_xyz(x: f64, y: f64, z: f64) -> Self {
98        Self { x, y, z: Some(z), m: None, t: None }
99    }
100}
101// M METHODS
102impl<M: Clone> VectorPoint<M> {
103    /// Create a new point
104    pub fn new(x: f64, y: f64, z: Option<f64>, m: Option<M>) -> Self {
105        Self { x, y, z, m, t: None }
106    }
107
108    /// Create a new point with xy
109    pub fn new_xy(x: f64, y: f64, m: Option<M>) -> Self {
110        Self { x, y, z: None, m, t: None }
111    }
112
113    /// Create a new point with xyz
114    pub fn new_xyz(x: f64, y: f64, z: f64, m: Option<M>) -> Self {
115        Self { x, y, z: Some(z), m, t: None }
116    }
117
118    /// Convert to an MValue VectorPoint
119    pub fn to_m_value(&self) -> VectorPoint<MValue>
120    where
121        M: Into<MValue>,
122    {
123        VectorPoint {
124            x: self.x,
125            y: self.y,
126            z: self.z,
127            m: self.m.clone().map(|m| m.into()),
128            t: None,
129        }
130    }
131
132    /// Project the point into the 0->1 coordinate system
133    pub fn project(&mut self, bbox: Option<&mut BBox3D>) {
134        let y = self.y;
135        let x = self.x;
136        let sin = sin((y * PI) / 180.);
137        let y2 = 0.5 - (0.25 * log((1. + sin) / (1. - sin))) / PI;
138        self.x = x / 360. + 0.5;
139        self.y = y2.clamp(0., 1.);
140
141        if let Some(bbox) = bbox {
142            bbox.extend_from_point(self)
143        };
144    }
145
146    /// Unproject the point from the 0->1 coordinate system back to a lon-lat coordinate
147    pub fn unproject(&mut self) {
148        let lon = (self.x - 0.5) * 360.;
149        let y2 = 0.5 - self.y;
150        let lat = atan(sinh(PI * (y2 * 2.))).to_degrees();
151
152        self.x = lon;
153        self.y = lat;
154    }
155
156    /// Returns true if the point is the zero vector.
157    pub fn is_empty(&self) -> bool {
158        let zero = f64::default();
159        self.x == zero && self.y == zero && (self.z.is_none() || self.z.unwrap() == zero)
160    }
161
162    /// Returns the S2 face assocated with this point
163    pub fn face(&self, f: u8) -> f64 {
164        if f == 0 {
165            self.x
166        } else if f == 1 {
167            self.y
168        } else {
169            self.z.unwrap_or(0.)
170        }
171    }
172
173    /// Apply modular arithmetic to x, y, and z using `modulus`
174    pub fn modulo(self, modulus: f64) -> Self {
175        let modulus = fabs(modulus); // Ensure positive modulus
176        Self {
177            x: fmod(self.x, modulus),
178            y: fmod(self.y, modulus),
179            z: self.z.map(|z| fmod(z, modulus)),
180            m: self.m,
181            t: None,
182        }
183    }
184
185    /// Returns the angle between "this" and v in radians, in the range [0, pi]. If
186    /// either vector is zero-length, or nearly zero-length, the result will be
187    /// zero, regardless of the other value.
188    pub fn angle<M2: Clone>(&self, b: &VectorPoint<M2>) -> f64 {
189        atan2(self.cross(b).norm(), self.dot(b))
190    }
191
192    /// Get the cross product of two Vector Points
193    pub fn cross<M2: Clone>(&self, b: &VectorPoint<M2>) -> Self {
194        if let (Some(z), Some(bz)) = (self.z, b.z) {
195            Self::new_xyz(
196                self.y * bz - z * b.y,
197                z * b.x - self.x * bz,
198                self.x * b.y - self.y * b.x,
199                None,
200            )
201        } else {
202            Self::new_xy(self.x * b.y - self.y * b.x, self.y * b.x - self.x * b.y, None)
203        }
204    }
205
206    /// dot returns the standard dot product of v and ov.
207    pub fn dot<M2: Clone>(&self, b: &VectorPoint<M2>) -> f64 {
208        if let (Some(z), Some(bz)) = (self.z, b.z) {
209            self.x * b.x + self.y * b.y + z * bz
210        } else {
211            self.x * b.x + self.y * b.y
212        }
213    }
214
215    /// Returns the absolute value of the point.
216    pub fn abs(&self) -> Self {
217        Self::new(fabs(self.x), fabs(self.y), self.z.map(fabs), None)
218    }
219
220    /// Returns the inverse of the point
221    pub fn invert(&self) -> Self {
222        -self
223    }
224
225    /// Returns the length of the point.
226    pub fn len(&self) -> f64 {
227        self.norm()
228    }
229
230    /// norm returns the vector's norm.
231    pub fn norm(&self) -> f64 {
232        sqrt(self.norm2())
233    }
234
235    /// norm2 returns the vector's squared norm.
236    pub fn norm2(&self) -> f64 {
237        self.dot(self)
238    }
239
240    /// Normalize this point to unit length.
241    pub fn normalize(&mut self) {
242        let len = self.len();
243        if len > 0.0 {
244            self.x /= len;
245            self.y /= len;
246            if let Some(z) = self.z {
247                self.z = Some(z / len);
248            }
249        }
250    }
251
252    /// return the distance from this point to the other point in radians
253    pub fn distance<M2: Clone>(&self, b: &VectorPoint<M2>) -> f64 {
254        let d: VectorPoint<M> = self - b;
255        d.len()
256    }
257
258    /// Returns the largest absolute component of the point.
259    pub fn largest_abs_component(&self) -> u8 {
260        let tmp = self.abs();
261        let tmp_z = tmp.z.unwrap_or(-2.);
262        if tmp.x > tmp.y {
263            if tmp.x > tmp_z { 0 } else { 2 }
264        } else if tmp.y > tmp_z {
265            1
266        } else {
267            2
268        }
269    }
270
271    /// Returns the intermediate point between this and the other point.
272    pub fn intermediate<M2: Clone>(&self, b: &VectorPoint<M2>, t: f64) -> Self {
273        if let (Some(z), Some(bz)) = (self.z, b.z) {
274            Self::new_xyz(
275                self.x + ((b.x - self.x) * (1.0 - t)),
276                self.y + ((b.y - self.y) * (1.0 - t)),
277                z + ((bz - z) * (1.0 - t)),
278                None,
279            )
280        } else {
281            Self::new_xy(
282                self.x + ((b.x - self.x) * (1.0 - t)),
283                self.y + ((b.y - self.y) * (1.0 - t)),
284                None,
285            )
286        }
287    }
288
289    /// Returns the perpendicular vector
290    pub fn perpendicular(&self) -> Self {
291        if let Some(z) = self.z {
292            let ref_point = if fabs(self.x) > fabs(z) {
293                Self::new_xyz(0., 0., 1., None)
294            } else {
295                Self::new_xyz(1., 0., 0., None)
296            };
297            let cross = self.cross(&ref_point);
298            Self::new_xyz(-self.y, self.x, cross.z.unwrap(), None)
299        } else {
300            Self::new_xy(-self.y, self.x, None)
301        }
302    }
303}
304impl<M: Clone> From<Point> for VectorPoint<M> {
305    fn from(p: Point) -> Self {
306        Self { x: p.0, y: p.1, z: None, m: None, t: None }
307    }
308}
309impl<M: Clone> From<&Point> for VectorPoint<M> {
310    fn from(p: &Point) -> Self {
311        Self { x: p.0, y: p.1, z: None, m: None, t: None }
312    }
313}
314impl<M: Clone> From<Point3D> for VectorPoint<M> {
315    fn from(p: Point3D) -> Self {
316        Self { x: p.0, y: p.1, z: Some(p.2), m: None, t: None }
317    }
318}
319impl<M: Clone> From<&Point3D> for VectorPoint<M> {
320    fn from(p: &Point3D) -> Self {
321        Self { x: p.0, y: p.1, z: Some(p.2), m: None, t: None }
322    }
323}
324impl<M1: Clone, M2: Clone> Add<&VectorPoint<M2>> for &VectorPoint<M1> {
325    type Output = VectorPoint<M1>;
326    fn add(self, other: &VectorPoint<M2>) -> Self::Output {
327        VectorPoint {
328            x: self.x + other.x,
329            y: self.y + other.y,
330            // Only add `z` if both `self.z` and `other.z` are `Some`
331            z: match (self.z, other.z) {
332                (Some(z1), Some(z2)) => Some(z1 + z2),
333                _ => None, // If either `z` is None, the result is `None`
334            },
335            m: None, // Combine m values
336            t: None, // Handle `t` as necessary
337        }
338    }
339}
340impl<M1: Clone, M2: Clone> AddAssign<&VectorPoint<M2>> for VectorPoint<M1> {
341    fn add_assign(&mut self, other: &VectorPoint<M2>) {
342        self.x += other.x;
343        self.y += other.y;
344        if let (Some(z), Some(other_z)) = (self.z, other.z) {
345            self.z = Some(z + other_z);
346        }
347    }
348}
349impl<M: Clone> Add<f64> for &VectorPoint<M> {
350    type Output = VectorPoint<M>;
351    fn add(self, other: f64) -> Self::Output {
352        VectorPoint {
353            x: self.x + other,
354            y: self.y + other,
355            z: self.z.map(|z| z + other),
356            m: None,
357            t: None,
358        }
359    }
360}
361impl<M: Clone> AddAssign<f64> for VectorPoint<M> {
362    fn add_assign(&mut self, other: f64) {
363        self.x += other;
364        self.y += other;
365        if let Some(z) = self.z {
366            self.z = Some(z + other);
367        }
368    }
369}
370// Implementing the Sub trait for VectorPoint
371impl<M1: Clone, M2: Clone> Sub<&VectorPoint<M2>> for &VectorPoint<M1> {
372    type Output = VectorPoint<M1>;
373    fn sub(self, other: &VectorPoint<M2>) -> Self::Output {
374        VectorPoint {
375            x: self.x - other.x,
376            y: self.y - other.y,
377            z: match (self.z, other.z) {
378                (Some(z1), Some(z2)) => Some(z1 - z2),
379                _ => None, // If either `z` is None, the result is `None`
380            },
381            m: None, // Combine m values
382            t: None, // Handle `t` as necessary
383        }
384    }
385}
386impl<M1: Clone, M2: Clone> SubAssign<&VectorPoint<M2>> for VectorPoint<M1> {
387    fn sub_assign(&mut self, other: &VectorPoint<M2>) {
388        self.x -= other.x;
389        self.y -= other.y;
390        if let (Some(z), Some(other_z)) = (self.z, other.z) {
391            self.z = Some(z - other_z);
392        }
393    }
394}
395impl<M: Clone> Sub<f64> for &VectorPoint<M> {
396    type Output = VectorPoint<M>;
397    fn sub(self, other: f64) -> Self::Output {
398        VectorPoint {
399            x: self.x - other,
400            y: self.y - other,
401            z: self.z.map(|z| z - other),
402            m: None,
403            t: None,
404        }
405    }
406}
407impl<M: Clone> SubAssign<f64> for VectorPoint<M> {
408    fn sub_assign(&mut self, other: f64) {
409        self.x -= other;
410        self.y -= other;
411        if let Some(z) = self.z {
412            self.z = Some(z - other);
413        }
414    }
415}
416// Implementing the Neg trait for VectorPoint
417impl<M: Clone> Neg for &VectorPoint<M> {
418    type Output = VectorPoint<M>;
419    fn neg(self) -> Self::Output {
420        VectorPoint { x: -self.x, y: -self.y, z: self.z.map(|z| -z), m: None, t: None }
421    }
422}
423// Implementing the Div trait for VectorPoint
424impl<M1: Clone, M2: Clone> Div<&VectorPoint<M2>> for &VectorPoint<M1> {
425    type Output = VectorPoint<M1>;
426    fn div(self, other: &VectorPoint<M2>) -> Self::Output {
427        VectorPoint {
428            x: self.x / other.x,
429            y: self.y / other.y,
430            z: match (self.z, other.z) {
431                (Some(z1), Some(z2)) => Some(z1 / z2),
432                _ => None, // If either `z` is None, the result is `None`
433            },
434            m: None, // Combine m values
435            t: None, // Handle `t` as necessary
436        }
437    }
438}
439impl<M1: Clone, M2: Clone> DivAssign<&VectorPoint<M2>> for VectorPoint<M1> {
440    fn div_assign(&mut self, other: &VectorPoint<M2>) {
441        self.x /= other.x;
442        self.y /= other.y;
443        if let (Some(z), Some(other_z)) = (self.z, other.z) {
444            self.z = Some(z / other_z);
445        }
446    }
447}
448impl<M: Clone> Div<f64> for &VectorPoint<M> {
449    type Output = VectorPoint<M>;
450    fn div(self, other: f64) -> Self::Output {
451        VectorPoint {
452            x: self.x / other,
453            y: self.y / other,
454            z: self.z.map(|z| z / other),
455            m: None,
456            t: None,
457        }
458    }
459}
460impl<M: Clone> DivAssign<f64> for VectorPoint<M> {
461    fn div_assign(&mut self, other: f64) {
462        self.x /= other;
463        self.y /= other;
464        if let Some(z) = self.z {
465            self.z = Some(z / other);
466        }
467    }
468}
469// Implementing the Mul trait for VectorPoint
470impl<M1: Clone, M2: Clone> Mul<&VectorPoint<M2>> for &VectorPoint<M1> {
471    type Output = VectorPoint<M1>;
472    fn mul(self, other: &VectorPoint<M2>) -> Self::Output {
473        VectorPoint {
474            x: self.x * other.x,
475            y: self.y * other.y,
476            z: match (self.z, other.z) {
477                (Some(z1), Some(z2)) => Some(z1 * z2),
478                _ => None, // If either `z` is None, the result is `None`
479            },
480            m: None, // Combine m values
481            t: None, // Handle `t` as necessary
482        }
483    }
484}
485impl<M1: Clone, M2: Clone> MulAssign<&VectorPoint<M2>> for VectorPoint<M1> {
486    fn mul_assign(&mut self, other: &VectorPoint<M2>) {
487        self.x *= other.x;
488        self.y *= other.y;
489        if let (Some(z), Some(other_z)) = (self.z, other.z) {
490            self.z = Some(z * other_z);
491        }
492    }
493}
494impl<M: Clone> Mul<f64> for &VectorPoint<M> {
495    type Output = VectorPoint<M>;
496    fn mul(self, other: f64) -> Self::Output {
497        VectorPoint {
498            x: self.x * other,
499            y: self.y * other,
500            z: self.z.map(|z| z * other),
501            m: None,
502            t: None,
503        }
504    }
505}
506impl<M: Clone> MulAssign<f64> for VectorPoint<M> {
507    fn mul_assign(&mut self, other: f64) {
508        self.x *= other;
509        self.y *= other;
510        if let Some(z) = self.z {
511            self.z = Some(z * other);
512        }
513    }
514}
515impl<M: Clone> Rem<f64> for VectorPoint<M> {
516    type Output = VectorPoint<M>;
517
518    fn rem(self, modulus: f64) -> Self::Output {
519        self.modulo(modulus)
520    }
521}
522impl<M: Clone> RemAssign<f64> for VectorPoint<M> {
523    fn rem_assign(&mut self, modulus: f64) {
524        let modulus = fabs(modulus);
525        self.x = fmod(self.x, modulus);
526        self.y = fmod(self.y, modulus);
527        if let Some(z) = self.z {
528            self.z = Some(fmod(z, modulus));
529        }
530    }
531}
532impl<M: Clone> Eq for VectorPoint<M> {}
533impl<M: Clone> Ord for VectorPoint<M> {
534    fn cmp(&self, other: &VectorPoint<M>) -> Ordering {
535        match self.x.partial_cmp(&other.x) {
536            Some(Ordering::Equal) => {}
537            other => return other.unwrap_or(Ordering::Greater), /* Handle cases where `x` comparison is not equal */
538        }
539        match self.y.partial_cmp(&other.y) {
540            Some(Ordering::Equal) => {}
541            other => return other.unwrap_or(Ordering::Greater), /* Handle cases where `y` comparison is not equal */
542        }
543        match self.z.partial_cmp(&other.z) {
544            Some(order) => order,
545            None => Ordering::Equal, // This handles the NaN case safely
546        }
547    }
548}
549impl<M: Clone> PartialEq for VectorPoint<M> {
550    fn eq(&self, other: &VectorPoint<M>) -> bool {
551        self.x == other.x && self.y == other.y && self.z == other.z
552    }
553}
554impl<M: Clone> PartialOrd for VectorPoint<M> {
555    fn partial_cmp(&self, other: &VectorPoint<M>) -> Option<Ordering> {
556        Some(self.cmp(other))
557    }
558}
s2json_core/geometry/vector_point.rs

s2json_core/geometry/
vector_point.rs
