rust_3d/
interpolate_2d.rs

1/*
2Copyright 2016 Martin Buck
3
4Permission is hereby granted, free of charge, to any person obtaining a copy
5of this software and associated documentation files (the "Software"),
6to deal in the Software without restriction, including without limitation the
7rights to use, copy, modify, merge, publish, distribute, sublicense,
8and/or sell copies of the Software, and to permit persons to whom the Software
9is furnished to do so, subject to the following conditions:
10
11The above copyright notice and this permission notice shall
12be included all copies or substantial portions of the Software.
13
14THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
15EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
16MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
17IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM,
18DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
19TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE
20OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
21*/
22
23//! interpolations within 2D space. E.g. bezier, linear, cosine
24
25use std::f64::consts::PI;
26
27use crate::*;
28
29//------------------------------------------------------------------------------
30
31fn factorial(number: usize) -> usize {
32    let mut result = 1;
33
34    for i in 1..number + 1 {
35        result *= i;
36    }
37    result
38}
39
40fn binominal_coefficient(n: usize, k: usize) -> usize {
41    factorial(n) / (factorial(k) * factorial(n - k))
42}
43
44fn bernstein_polynomial(n: usize, i: usize, t: f64) -> f64 {
45    (binominal_coefficient(n, i) as f64) * t.powi(i as i32) * (1.0 - t).powi((n - i) as i32)
46}
47
48fn control_polygon<P>(path: &PointCloud2D<P>, n_points: usize, t: f64) -> P
49where
50    P: IsBuildable2D,
51{
52    let mut x: f64 = 0.0;
53    let mut y: f64 = 0.0;
54
55    for i in 0..n_points + 1 {
56        let bp = bernstein_polynomial(n_points, i, t);
57        x += bp * path.data[i].x();
58        y += bp * path.data[i].y();
59    }
60    P::new(x, y)
61}
62
63/// Returns the Bezier interpolation of the given base points
64pub fn interpolate_bezier<P>(
65    base_points: &PointCloud2D<P>,
66    n_points: usize,
67) -> Result<PointCloud2D<P>>
68where
69    P: IsBuildable2D,
70{
71    if base_points.len() < 2 {
72        return Err(ErrorKind::TooFewPoints);
73    }
74
75    let mut pc = PointCloud2D::with_capacity(n_points);
76    let p_dist = 1.0 / (n_points as f64);
77
78    for i in 0..n_points {
79        pc.push(control_polygon(
80            base_points,
81            base_points.len() - 1,
82            (i as f64) * p_dist,
83        ));
84    }
85    Ok(pc)
86}
87
88/// Returns the Cosine interpolation of the given base points
89pub fn interpolate_cosine<P>(
90    base_points: &PointCloud2D<P>,
91    n_points: usize,
92) -> Result<PointCloud2D<P>>
93where
94    P: IsBuildable2D,
95{
96    if base_points.len() < 2 {
97        return Err(ErrorKind::TooFewPoints);
98    }
99
100    let mut pc = PointCloud2D::with_capacity(n_points);
101    let p_dist = base_points.length() / (n_points - 1) as f64;
102
103    for i in 0..n_points {
104        let mut traveled: f64 = 0.0;
105        let mut traveled_before: f64 = 0.0;
106
107        for j in 1..base_points.len() {
108            let ref p_prev = base_points.data[j - 1];
109            let ref p_now = base_points.data[j];
110
111            traveled +=
112                ((p_now.x() - p_prev.x()).powi(2) + (p_now.y() - p_prev.y()).powi(2)).sqrt();
113
114            if traveled >= p_dist * (i as f64) {
115                let proportion =
116                    ((i as f64) * p_dist - traveled_before) / (traveled - traveled_before);
117                let proportion2 = (1.0 - (proportion * PI).cos()) / 2.0;
118                pc.push(P::new(
119                    p_prev.x() + proportion * (p_now.x() - p_prev.x()),
120                    p_prev.y() * (1.0 - proportion2) + p_now.y() * proportion2,
121                ));
122                break;
123            }
124            traveled_before = traveled;
125        }
126    }
127    Ok(pc)
128}
129
130/// Returns the linear interpolation of the given base points
131pub fn interpolation_linear<P>(
132    base_points: &PointCloud2D<P>,
133    n_points: usize,
134) -> Result<PointCloud2D<P>>
135where
136    P: IsBuildable2D,
137{
138    if base_points.len() < 2 {
139        return Err(ErrorKind::TooFewPoints);
140    }
141
142    let mut pc = PointCloud2D::with_capacity(n_points);
143    let p_dist = base_points.length() / (n_points - 1) as f64;
144
145    for i in 0..n_points {
146        let mut traveled: f64 = 0.0;
147        let mut traveled_before: f64 = 0.0;
148
149        for j in 1..base_points.len() {
150            //@todo fails if path too small, handle this
151            let ref p_prev = base_points.data[j - 1];
152            let ref p_now = base_points.data[j];
153
154            traveled +=
155                ((p_now.x() - p_prev.x()).powi(2) + (p_now.y() - p_prev.y()).powi(2)).sqrt();
156
157            if traveled >= p_dist * (i as f64) {
158                let proportion =
159                    ((i as f64) * p_dist - traveled_before) / (traveled - traveled_before);
160                pc.push(P::new(
161                    p_prev.x() + proportion * (p_now.x() - p_prev.x()),
162                    p_prev.y() + proportion * (p_now.y() - p_prev.y()),
163                ));
164                break;
165            }
166            traveled_before = traveled;
167        }
168    }
169    Ok(pc)
170}