1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
use delaunator::{Point, Triangulation, next_halfedge};
use crate::{Voronoi, utils::{dist2, site_of_incoming, self}};

use super::{EMPTY};

/// Iterator that walks through all the edges connected to a provided starting point.
/// The iteration happens in a clock-wise manner.
///
/// Note: this really only returns edges for triangles around the site. On the convex hull, the rightmost edge will not be returned because there is no incoming rightmost edge.
#[derive(Clone, Debug)]
pub struct EdgesAroundSiteIterator<'t> {
    triangulation: &'t Triangulation,
    start: usize,
    next: usize
}

impl<'t> EdgesAroundSiteIterator<'t> {
    /// Creates iterator based on a incoming edge to a site.
    /// This must be the left-most incoming edge to the site to avoid early iteration stop around the convex hull.
    pub fn new(triangulation: &'t Triangulation, incoming_edge: usize) -> Self {
        Self {
            triangulation,
            start: incoming_edge,
            next: incoming_edge
        }
    }
}

impl<'t> Iterator for EdgesAroundSiteIterator<'t> {
    type Item = usize;
    /// Walks all half-edges around the starting point and returning the associated surrounding incoming edges
    fn next(&mut self) -> Option<Self::Item> {
        let incoming = self.next;

        if incoming != EMPTY {
            let outgoing = next_halfedge(incoming);

            // then take the oposite half-edge, it will be the incoming edge for the opposite triangle
            self.next = self.triangulation.halfedges[outgoing];

            // if we are back to the begining, there is nothing else to do
            if self.next == self.start {
                self.next = EMPTY
            }

            Some(incoming)
        } else {
            None
        }
    }
}

/// Iterates over sites that are topologically adjacent.
///
/// Topological neighbors are sites that share a delaunay edge between them.
/// To take into account the effect of voronoi edge clipping use [NeighborSiteIterator]
/// Sites are returned clockwise.
#[derive(Clone, Debug)]
pub struct TopologicalNeighborSiteIterator<'t> {
    iter: EdgesAroundSiteIterator<'t>,
    last_incoming: usize,
}

impl<'t> TopologicalNeighborSiteIterator<'t> {
    /// Creates iterator based on the site.
    pub fn new(voronoi: &'t Voronoi, site: usize) -> Self {
        Self::with_triangulation(voronoi.triangulation(), &voronoi.site_to_incoming_leftmost_halfedge, site)
    }

    /// Creates iterator based on the site.
    pub fn with_triangulation(triangulation: &'t Triangulation, site_to_incoming_leftmost_halfedge: &'t Vec<usize>, site: usize) -> Self {
        let &incoming_leftmost_edge = site_to_incoming_leftmost_halfedge.get(site).expect("Site does not exist");
        Self {
            iter: EdgesAroundSiteIterator::new(&triangulation, incoming_leftmost_edge),
            last_incoming: EMPTY,
        }
    }
}

impl<'t> Iterator for TopologicalNeighborSiteIterator<'t> {
    type Item = usize;
    /// Get the next neighboring site.
    fn next(&mut self) -> Option<Self::Item> {
        if let Some(incoming_edge) = self.iter.next() {
            self.last_incoming = incoming_edge;
            Some(self.iter.triangulation.triangles[incoming_edge])
        } else if self.last_incoming != EMPTY {
            // if we are on hull, the last neighbor is not returned by EdgesAroundSiteIterator because there is no halfedge
            let outgoing = next_halfedge(self.last_incoming);
            self.last_incoming = EMPTY;
            if self.iter.triangulation.halfedges[outgoing] == EMPTY {
                // on hull edge
                Some(site_of_incoming(self.iter.triangulation, outgoing))
            } else {
                // not on hull
                None
            }
        } else {
            None
        }
    }
}

/// Iterates over sites that are adjacent in the voronoi diagram.
///
/// This iterator expands on [TopologicalNeighborSiteIterator] by taking into account the clipping effect of voronoi edges to decide whether two sites are neighbors.
/// In this iterator, two sites are considered neighbors if their voronoi cells share a common edge in the voronoi graph.
/// Sites on the hull may get disconnected by this classification if their common voronoi edge gets clipped away because it lies beyond the bounding box.
#[derive(Clone, Debug)]
pub struct NeighborSiteIterator<'t> {
    voronoi: &'t Voronoi,
    topo_neighbor_iter: TopologicalNeighborSiteIterator<'t>,
    site: usize
}

impl<'t> NeighborSiteIterator<'t> {
    pub fn new(voronoi: &'t Voronoi, site: usize) -> Self {
        Self {
            voronoi,
            topo_neighbor_iter: TopologicalNeighborSiteIterator::new(voronoi, site),
            site
        }
    }
}

impl<'t> Iterator for NeighborSiteIterator<'t> {
    type Item = usize;

    /// Get the next neighboring site.
    fn next(&mut self) -> Option<Self::Item> {
        let prev_last_incoming = self.topo_neighbor_iter.last_incoming;

        if let Some(neighbor) = self.topo_neighbor_iter.next() {
            // if first neighbor and on hull, need special check for clipping
            if prev_last_incoming == EMPTY
                && self.voronoi.triangulation.halfedges[self.topo_neighbor_iter.last_incoming] == EMPTY {
                if utils::has_common_voronoi_edge(self.voronoi, self.site, neighbor) {
                    Some(neighbor)
                } else {
                    // not connected on the voronoi diagram
                    // move to next
                    self.next()
                }
            } else if self.topo_neighbor_iter.last_incoming == EMPTY {
                // last neighbor, need speciail check for clipping
                if utils::has_common_voronoi_edge(self.voronoi, self.site, neighbor) {
                    Some(neighbor)
                } else {
                    None
                }
            } else {
                Some(neighbor)
            }
        } else {
            None
        }
    }
}

/// Iterator that produces a path between two points in the Voronoi diagram that uses a greed approach to minimizes a cost function.
///
/// A cost function is provided that calculates the cost of an edge; edges for all neighbors are evaluated and the least costly is taken.
/// The process is evaluated for the next cell in the path until no edge can be taken that costs less than f64::MAX.
/// If the destionation point is not contained in the Voronoi diagram, the final cell in the path will be the
/// closest to the destination point.
#[derive(Clone, Debug)]
pub struct CellPathIterator<'t, F> {
    site: usize,
    cost_fn: F,
    triangulation: &'t Triangulation,
    site_to_incoming_leftmost_halfedge: &'t Vec<usize>
}

impl<'t, F> CellPathIterator<'t, F> {
    /// Creates iterator based on the starting site and cost function.
    pub fn new(voronoi: &'t Voronoi, site: usize, cost_fn: F) -> Self {
        assert!(site < voronoi.sites.len(), "site {} does not exist", site);

        Self::with_triangulation(voronoi.triangulation(), &voronoi.site_to_incoming_leftmost_halfedge, site, cost_fn)
    }

    pub fn with_triangulation(triangulation: &'t Triangulation, site_to_incoming_leftmost_halfedge: &'t Vec<usize>, site: usize, cost_fn: F) -> Self {
        Self {
            site,
            cost_fn,
            triangulation,
            site_to_incoming_leftmost_halfedge
        }
    }
}

impl<'t, F> Iterator for CellPathIterator<'t, F>
    where F : Fn(usize, usize) -> f64 {
    type Item = usize;

    /// Walks current site neighbor and find the next site in the path
    fn next(&mut self) -> Option<Self::Item> {
        let current_site = self.site;

        if current_site != EMPTY {
            // take the neighbor with least cost
            let next = TopologicalNeighborSiteIterator::with_triangulation(self.triangulation, self.site_to_incoming_leftmost_halfedge, current_site)
                .map(|n| (n, (self.cost_fn)(current_site, n)))
                .min_by(|(_, cost0), (_, cost1)| cost0.partial_cmp(cost1).unwrap());

            // if next neighbor cost is less than f64::MAX, then we can move to it - it is next in the path
            if let Some((n, cost)) = next {
                if cost < f64::MAX {
                    self.site = n;
                } else {
                    // reached end
                    self.site = EMPTY;
                }
            } else {
                // reached end
                self.site = EMPTY;
            }

            Some(current_site)
        } else {
            None
        }
    }
}

/// Produces an iterator that calculates the shortest path from ```start_site``` to a ```dest``` point.
///
/// If destination point is outside voronoi diagram, then the closest point to destination in the voronoi diagram will be returned.
pub fn shortest_path_iter<'v>(voronoi: &'v Voronoi, start_site: usize, dest: Point) -> impl Iterator<Item = usize> + 'v {
    shortest_path_iter_from_triangulation(voronoi.triangulation(), &voronoi.sites(), &voronoi.site_to_incoming_leftmost_halfedge, start_site, dest)
}

pub (crate) fn shortest_path_iter_from_triangulation<'t>(triangulation: &'t Triangulation, sites: &'t Vec<Point>, site_to_incoming_leftmost_halfedge: &'t Vec<usize>, start_site: usize, dest: Point) -> impl Iterator<Item = usize> + 't {
    CellPathIterator::with_triangulation(triangulation, site_to_incoming_leftmost_halfedge, start_site, move |curr, next| {
        // calculate distance
        let dist_to_dest = dist2(&sites[curr], &dest);
        let dist_from_next = dist2(&sites[next], &dest);

        if dist_to_dest <= dist_from_next {
            // if current is closer to dest than next is, make cost to travel to next impossibly high
            f64::MAX
        } else {
            // cost is distance
            dist_from_next
        }
    })
}

#[cfg(test)]
mod test {
    use delaunator::Point;
    use crate::{VoronoiBuilder, utils::test::assert_list_eq};
    use super::*;

    #[test]
    fn iter_neighbors_hull_test() {
        let sites = vec![Point { x: -0.5, y: 0.0 }, Point { x: 0.5, y: 0.0 }, Point { x: 0.0, y: 0.0 }, Point { x: 0.0, y: 0.5 }, Point { x: 0.0, y: -0.5 }];
        let v = VoronoiBuilder::default()
            .set_sites(sites)
            .build()
            .unwrap();

        let neighbors: Vec<usize> = TopologicalNeighborSiteIterator::new(&v, 0).collect();
        assert_eq!(neighbors.len(), 3, "There are 3 neighboring sites");
        assert_eq!(neighbors[0], 4);
        assert_eq!(neighbors[1], 2);
        assert_eq!(neighbors[2], 3);
    }

    #[test]
    fn iter_neighbors_inner_test() {
        let sites = vec![Point { x: -0.5, y: 0.0 }, Point { x: 0.5, y: 0.0 }, Point { x: 0.0, y: 0.0 }, Point { x: 0.0, y: 0.5 }, Point { x: 0.0, y: -0.5 }];
        let v = VoronoiBuilder::default()
            .set_sites(sites)
            .build()
            .unwrap();
        let neighbors: Vec<usize> = TopologicalNeighborSiteIterator::new(&v, 2).collect();
        assert_eq!(neighbors.len(), 4, "There are 4 neighboring sites");
        assert_eq!(neighbors[0], 3);
        assert_eq!(neighbors[1], 0);
        assert_eq!(neighbors[2], 4);
        assert_eq!(neighbors[3], 1);
    }

    #[test]
    fn iter_neighbors_edge_clipped_by_box_test() -> std::io::Result<()> {
        let voronoi = utils::test::new_voronoi_builder_from_asset("degenerated8.json")?
            .build()
            .expect("Some voronoi expected");

        let neighbors = NeighborSiteIterator::new(&voronoi, 0).collect::<Vec<_>>();
        assert_list_eq(&[3, 2], &neighbors, "Visible neighbors of 0");

        let neighbors = NeighborSiteIterator::new(&voronoi, 1).collect::<Vec<_>>();
        assert_list_eq(&[2], &neighbors, "Visible neighbors of 1");

        Ok(())
    }

    #[test]
    fn iter_cell_path_test() {
        let sites = vec![
            Point { x: -0.5, y: 0.0 },
            Point { x: 0.0, y: 0.0 }, Point { x: 0.0, y: 0.5 }, Point { x: 0.0, y: -0.5 },
            Point { x: 0.2, y: 0.0 }, Point { x: 0.2, y: 0.5 }, Point { x: 0.2, y: -0.5 },
            Point { x: 0.5, y: 0.0 },
        ];
        let v = VoronoiBuilder::default()
            .set_sites(sites.clone())
            .build()
            .unwrap();
        let mut path = shortest_path_iter(&v, 0, sites.last().unwrap().clone());
        assert_eq!(Some(0), path.next());
        assert_eq!(Some(1), path.next());
        assert_eq!(Some(4), path.next());
        assert_eq!(Some(7), path.next());
        assert_eq!(None, path.next());
    }
}