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, P: GetXYZ> From<&P> for VectorPoint<M> {
305    fn from(p: &P) -> Self {
306        Self { x: p.x(), y: p.y(), z: p.z(), 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<M1: Clone, M2: Clone> Add<&VectorPoint<M2>> for &VectorPoint<M1> {
320    type Output = VectorPoint<M1>;
321    fn add(self, other: &VectorPoint<M2>) -> Self::Output {
322        VectorPoint {
323            x: self.x + other.x,
324            y: self.y + other.y,
325            // Only add `z` if both `self.z` and `other.z` are `Some`
326            z: match (self.z, other.z) {
327                (Some(z1), Some(z2)) => Some(z1 + z2),
328                _ => None, // If either `z` is None, the result is `None`
329            },
330            m: None, // Combine m values
331            t: None, // Handle `t` as necessary
332        }
333    }
334}
335impl<M1: Clone, M2: Clone> AddAssign<&VectorPoint<M2>> for VectorPoint<M1> {
336    fn add_assign(&mut self, other: &VectorPoint<M2>) {
337        self.x += other.x;
338        self.y += other.y;
339        if let (Some(z), Some(other_z)) = (self.z, other.z) {
340            self.z = Some(z + other_z);
341        }
342    }
343}
344impl<M: Clone> Add<f64> for &VectorPoint<M> {
345    type Output = VectorPoint<M>;
346    fn add(self, other: f64) -> Self::Output {
347        VectorPoint {
348            x: self.x + other,
349            y: self.y + other,
350            z: self.z.map(|z| z + other),
351            m: None,
352            t: None,
353        }
354    }
355}
356impl<M: Clone> AddAssign<f64> for VectorPoint<M> {
357    fn add_assign(&mut self, other: f64) {
358        self.x += other;
359        self.y += other;
360        if let Some(z) = self.z {
361            self.z = Some(z + other);
362        }
363    }
364}
365// Implementing the Sub trait for VectorPoint
366impl<M1: Clone, M2: Clone> Sub<&VectorPoint<M2>> for &VectorPoint<M1> {
367    type Output = VectorPoint<M1>;
368    fn sub(self, other: &VectorPoint<M2>) -> Self::Output {
369        VectorPoint {
370            x: self.x - other.x,
371            y: self.y - other.y,
372            z: match (self.z, other.z) {
373                (Some(z1), Some(z2)) => Some(z1 - z2),
374                _ => None, // If either `z` is None, the result is `None`
375            },
376            m: None, // Combine m values
377            t: None, // Handle `t` as necessary
378        }
379    }
380}
381impl<M1: Clone, M2: Clone> SubAssign<&VectorPoint<M2>> for VectorPoint<M1> {
382    fn sub_assign(&mut self, other: &VectorPoint<M2>) {
383        self.x -= other.x;
384        self.y -= other.y;
385        if let (Some(z), Some(other_z)) = (self.z, other.z) {
386            self.z = Some(z - other_z);
387        }
388    }
389}
390impl<M: Clone> Sub<f64> for &VectorPoint<M> {
391    type Output = VectorPoint<M>;
392    fn sub(self, other: f64) -> Self::Output {
393        VectorPoint {
394            x: self.x - other,
395            y: self.y - other,
396            z: self.z.map(|z| z - other),
397            m: None,
398            t: None,
399        }
400    }
401}
402impl<M: Clone> SubAssign<f64> for VectorPoint<M> {
403    fn sub_assign(&mut self, other: f64) {
404        self.x -= other;
405        self.y -= other;
406        if let Some(z) = self.z {
407            self.z = Some(z - other);
408        }
409    }
410}
411// Implementing the Neg trait for VectorPoint
412impl<M: Clone> Neg for &VectorPoint<M> {
413    type Output = VectorPoint<M>;
414    fn neg(self) -> Self::Output {
415        VectorPoint { x: -self.x, y: -self.y, z: self.z.map(|z| -z), m: None, t: None }
416    }
417}
418// Implementing the Div trait for VectorPoint
419impl<M1: Clone, M2: Clone> Div<&VectorPoint<M2>> for &VectorPoint<M1> {
420    type Output = VectorPoint<M1>;
421    fn div(self, other: &VectorPoint<M2>) -> Self::Output {
422        VectorPoint {
423            x: self.x / other.x,
424            y: self.y / other.y,
425            z: match (self.z, other.z) {
426                (Some(z1), Some(z2)) => Some(z1 / z2),
427                _ => None, // If either `z` is None, the result is `None`
428            },
429            m: None, // Combine m values
430            t: None, // Handle `t` as necessary
431        }
432    }
433}
434impl<M1: Clone, M2: Clone> DivAssign<&VectorPoint<M2>> for VectorPoint<M1> {
435    fn div_assign(&mut self, other: &VectorPoint<M2>) {
436        self.x /= other.x;
437        self.y /= other.y;
438        if let (Some(z), Some(other_z)) = (self.z, other.z) {
439            self.z = Some(z / other_z);
440        }
441    }
442}
443impl<M: Clone> Div<f64> for &VectorPoint<M> {
444    type Output = VectorPoint<M>;
445    fn div(self, other: f64) -> Self::Output {
446        VectorPoint {
447            x: self.x / other,
448            y: self.y / other,
449            z: self.z.map(|z| z / other),
450            m: None,
451            t: None,
452        }
453    }
454}
455impl<M: Clone> DivAssign<f64> for VectorPoint<M> {
456    fn div_assign(&mut self, other: f64) {
457        self.x /= other;
458        self.y /= other;
459        if let Some(z) = self.z {
460            self.z = Some(z / other);
461        }
462    }
463}
464// Implementing the Mul trait for VectorPoint
465impl<M1: Clone, M2: Clone> Mul<&VectorPoint<M2>> for &VectorPoint<M1> {
466    type Output = VectorPoint<M1>;
467    fn mul(self, other: &VectorPoint<M2>) -> Self::Output {
468        VectorPoint {
469            x: self.x * other.x,
470            y: self.y * other.y,
471            z: match (self.z, other.z) {
472                (Some(z1), Some(z2)) => Some(z1 * z2),
473                _ => None, // If either `z` is None, the result is `None`
474            },
475            m: None, // Combine m values
476            t: None, // Handle `t` as necessary
477        }
478    }
479}
480impl<M1: Clone, M2: Clone> MulAssign<&VectorPoint<M2>> for VectorPoint<M1> {
481    fn mul_assign(&mut self, other: &VectorPoint<M2>) {
482        self.x *= other.x;
483        self.y *= other.y;
484        if let (Some(z), Some(other_z)) = (self.z, other.z) {
485            self.z = Some(z * other_z);
486        }
487    }
488}
489impl<M: Clone> Mul<f64> for &VectorPoint<M> {
490    type Output = VectorPoint<M>;
491    fn mul(self, other: f64) -> Self::Output {
492        VectorPoint {
493            x: self.x * other,
494            y: self.y * other,
495            z: self.z.map(|z| z * other),
496            m: None,
497            t: None,
498        }
499    }
500}
501impl<M: Clone> MulAssign<f64> for VectorPoint<M> {
502    fn mul_assign(&mut self, other: f64) {
503        self.x *= other;
504        self.y *= other;
505        if let Some(z) = self.z {
506            self.z = Some(z * other);
507        }
508    }
509}
510impl<M: Clone> Rem<f64> for VectorPoint<M> {
511    type Output = VectorPoint<M>;
512
513    fn rem(self, modulus: f64) -> Self::Output {
514        self.modulo(modulus)
515    }
516}
517impl<M: Clone> RemAssign<f64> for VectorPoint<M> {
518    fn rem_assign(&mut self, modulus: f64) {
519        let modulus = fabs(modulus);
520        self.x = fmod(self.x, modulus);
521        self.y = fmod(self.y, modulus);
522        if let Some(z) = self.z {
523            self.z = Some(fmod(z, modulus));
524        }
525    }
526}
527impl<M: Clone> Eq for VectorPoint<M> {}
528impl<M: Clone> Ord for VectorPoint<M> {
529    fn cmp(&self, other: &VectorPoint<M>) -> Ordering {
530        match self.x.partial_cmp(&other.x) {
531            Some(Ordering::Equal) => {}
532            other => return other.unwrap_or(Ordering::Greater), /* Handle cases where `x` comparison is not equal */
533        }
534        match self.y.partial_cmp(&other.y) {
535            Some(Ordering::Equal) => {}
536            other => return other.unwrap_or(Ordering::Greater), /* Handle cases where `y` comparison is not equal */
537        }
538        match self.z.partial_cmp(&other.z) {
539            Some(order) => order,
540            None => Ordering::Equal, // This handles the NaN case safely
541        }
542    }
543}
544impl<M: Clone> PartialEq for VectorPoint<M> {
545    fn eq(&self, other: &VectorPoint<M>) -> bool {
546        self.x == other.x && self.y == other.y && self.z == other.z
547    }
548}
549impl<M: Clone> PartialOrd for VectorPoint<M> {
550    fn partial_cmp(&self, other: &VectorPoint<M>) -> Option<Ordering> {
551        Some(self.cmp(other))
552    }
553}
s2json_core/geometry/vector_point.rs

s2json_core/geometry/
vector_point.rs
