h3ron_ndarray/
resolution.rs

1use geo_types::{Coord, Rect};
2
3use h3ron::{H3Cell, ToPolygon, H3_MAX_RESOLUTION, H3_MIN_RESOLUTION};
4
5use crate::{
6    error::Error,
7    sphere::{area_squaremeters_linearring, area_squaremeters_rect},
8    transform::Transform,
9    AxisOrder,
10};
11
12pub enum ResolutionSearchMode {
13    /// Chose the h3 resolution where the difference in the area of a pixel and the h3index is
14    /// as small as possible.
15    MinDiff,
16
17    /// Chose the h3 resolution where the area of the h3index is smaller than the area of a pixel.
18    SmallerThanPixel,
19}
20
21/// Find the h3 resolution closed to the size of a pixel in an array
22/// of the given shape with the given transform.
23pub fn nearest_h3_resolution(
24    shape: &[usize],
25    transform: &Transform,
26    axis_order: &AxisOrder,
27    search_mode: ResolutionSearchMode,
28) -> Result<u8, Error> {
29    if shape.len() != 2 {
30        return Err(Error::UnsupportedArrayShape);
31    }
32    if shape[0] == 0 || shape[1] == 0 {
33        return Err(Error::EmptyArray);
34    }
35    let bbox_array = Rect::new(
36        transform * Coord::from((0.0_f64, 0.0_f64)),
37        transform
38            * Coord::from((
39                (shape[axis_order.x_axis()] - 1) as f64,
40                (shape[axis_order.y_axis()] - 1) as f64,
41            )),
42    );
43    let area_pixel = area_squaremeters_rect(&bbox_array)
44        / (shape[axis_order.x_axis()] * shape[axis_order.y_axis()]) as f64;
45    let center_of_array = bbox_array.center();
46
47    let mut nearest_h3_res = 0;
48    let mut area_difference = None;
49    for h3_res in H3_MIN_RESOLUTION..=H3_MAX_RESOLUTION {
50        // calculate the area of the center index to avoid using the approximate values
51        // of the h3ron hexArea functions
52        let area_h3_index = area_squaremeters_linearring(
53            H3Cell::from_coordinate(center_of_array, h3_res)?
54                .to_polygon()?
55                .exterior(),
56        );
57
58        match search_mode {
59            ResolutionSearchMode::SmallerThanPixel => {
60                if area_h3_index <= area_pixel {
61                    nearest_h3_res = h3_res;
62                    break;
63                }
64            }
65
66            ResolutionSearchMode::MinDiff => {
67                let new_area_difference = if area_h3_index > area_pixel {
68                    area_h3_index - area_pixel
69                } else {
70                    area_pixel - area_h3_index
71                };
72                if let Some(old_area_difference) = area_difference {
73                    if old_area_difference < new_area_difference {
74                        nearest_h3_res = h3_res - 1;
75                        break;
76                    } else {
77                        area_difference = Some(new_area_difference);
78                    }
79                } else {
80                    area_difference = Some(new_area_difference);
81                }
82            }
83        }
84    }
85
86    Ok(nearest_h3_res)
87}
88
89#[cfg(test)]
90mod tests {
91    use crate::resolution::{nearest_h3_resolution, ResolutionSearchMode};
92    use crate::transform::Transform;
93    use crate::AxisOrder;
94
95    #[test]
96    fn test_nearest_h3_resolution() {
97        // transform of the included r.tiff
98        let gt = Transform::from_rasterio(&[
99            0.0011965049999999992,
100            0.0,
101            8.11377,
102            0.0,
103            -0.001215135,
104            49.40792,
105        ]);
106        let h3_res1 = nearest_h3_resolution(
107            &[2000_usize, 2000_usize],
108            &gt,
109            &AxisOrder::YX,
110            ResolutionSearchMode::MinDiff,
111        )
112        .unwrap();
113        assert_eq!(h3_res1, 10); // TODO: validate
114
115        let h3_res2 = nearest_h3_resolution(
116            &[2000_usize, 2000_usize],
117            &gt,
118            &AxisOrder::YX,
119            ResolutionSearchMode::SmallerThanPixel,
120        )
121        .unwrap();
122        assert_eq!(h3_res2, 11); // TODO: validate
123    }
124}
h3ron_ndarray/resolution.rs

h3ron_ndarray/
resolution.rs
