dynamic_poisson_sampling/
lib.rs

1
2//! `dynamic_poisson_sampling` is a very simple crate 
3//! that allows using poisson sampling with dynamic distance. 
4
5
6use rand;
7
8#[doc(hidden)]
9fn distance<const N: usize>(lhs: &[f64; N], rhs: &[f64; N]) -> f64{
10    lhs.iter()
11    .zip(rhs.iter())
12    .map(|(&l, &r)| {
13        let diff = l - r;
14        diff * diff
15    })
16    .sum::<f64>().sqrt()
17}
18
19/// Returns the points generated with a dynamic poisson sampling.  
20/// This function is **SLOW** but very permissive. It doesn't need to know the minimum and maximum radius in advance.  
21/// 
22/// The provided density function `density_func` takes the position of the current checked point as argument
23/// and returns an optional radius. If it is None, the point is discarded.
24/// 
25/// `k` represents the number of generated and checked points for each new one.  
26/// High value will slow down a lot the algorithm but will give more accurate results.  
27/// If you need fast results or accuracy doesn't matter, lower the value.  
28/// Typically, `k` varies from 10 to 30.
29/// 
30/// # Example
31/// ```rust
32/// use dynamic_poisson_sampling::get_points;
33/// let mut rng = rand::thread_rng();
34/// let points = get_points(12, (0.5, 0.5).into(), &mut rng, |pos|{
35///     // bound check to avoid infinite loop
36///     if pos[0] < 0.0 || pos[0] >= 1.0 || pos[1] < 0.0 || pos[1] >= 1.0{
37///         None
38///     }else{
39///         Some(0.05)
40///     }
41/// });
42/// ```
43pub fn get_points<const N: usize, T>(k: u32, first_pos: [f64; N], rng: &mut impl rand::Rng, density_func:T)
44-> Vec<[f64; N]>
45where
46    T: Fn(&[f64; N]) -> Option<f64>
47{
48    struct Point<const N: usize>{
49        pos: [f64; N],
50        min_rad: f64,
51    }
52
53    let mut points: Vec<Point<N>> = Vec::new();
54    
55    //active list
56    let mut active_list: Vec<usize> = Vec::new();
57
58
59
60
61
62    //first point
63    let r_rad = density_func(&first_pos);
64    if r_rad.is_none() { return vec![]; }
65    let r_rad = r_rad.unwrap();
66
67    points.push(Point{
68        pos: first_pos,
69        min_rad: r_rad,
70    });
71    active_list.push(0);
72
73
74
75    loop {
76        let r_ind = rng.gen_range(0..active_list.len());
77        let curr_ind = active_list[r_ind];
78        let (curr_pos, curr_rad) = (points[curr_ind].pos, points[curr_ind].min_rad);
79        
80        // let mut to_remove = true;
81        'k_l: for _ in 0..k {
82            let r_distance = rng.gen_range(curr_rad .. 3.*curr_rad);
83
84
85            let mut deltas: [f64; N] = [0.0; N];
86
87            // Generate random direction
88            for delta in &mut deltas {
89                *delta = rng.gen_range(-1.0..=1.0);
90            }
91            
92            // Add deltas to current position
93            let norm = (deltas.iter().map(|&x| x.powi(2)).sum::<f64>()).sqrt();
94            let mut r_pos: [f64; N] = [0.0; N];
95            for i in 0..N{
96                r_pos[i] = curr_pos[i]+(deltas[i]/norm)*r_distance;
97            }
98
99
100            let r_rad = density_func(&r_pos);
101            if r_rad.is_none() { continue 'k_l; }
102            let r_rad = r_rad.unwrap();
103
104            for point in points.iter() {
105                let distance = distance(&point.pos, &r_pos);
106                if distance < point.min_rad+r_rad {
107                    continue 'k_l;
108                }
109            }
110
111            points.push(Point{
112                pos: r_pos,
113                min_rad: r_rad,
114            });
115        
116            active_list.push(points.len()-1);
117        }
118        
119        active_list.remove(r_ind);
120
121        if active_list.is_empty() {
122            break;
123        }
124    }
125
126    points.iter().map(|p| p.pos).collect()
127}
dynamic_poisson_sampling/lib.rs

dynamic_poisson_sampling/
lib.rs
