Skip to content

Commit

Permalink
Proper error handling, alter CoordinateComponent to use AoS/SoA
Browse files Browse the repository at this point in the history
terms.
  • Loading branch information
JmsPae committed Jan 13, 2025
1 parent 92901d2 commit bea5eeb
Show file tree
Hide file tree
Showing 4 changed files with 58 additions and 31 deletions.
82 changes: 54 additions & 28 deletions src/vector/geometry.rs
Original file line number Diff line number Diff line change
Expand Up @@ -9,10 +9,10 @@ use std::{

use gdal_sys::{OGRErr, OGRGeometryH, OGRwkbGeometryType};

use crate::errors::*;
use crate::spatial_ref::SpatialRef;
use crate::utils::{_last_null_pointer_err, _string};
use crate::vector::{Envelope, Envelope3D};
use crate::{errors::*, utils::_last_cpl_err};

#[derive(Debug, Clone, Copy)]
#[repr(u8)]
Expand All @@ -21,8 +21,9 @@ pub enum CoordinateComponent {
Y = 0b00010,
Z = 0b00100,
M = 0b01000,
// Whether or not the layout is sequential (eg. XyXy as opposed to XxYy.)
Sequential = 0b10000,
/// Whether the point data layout is of an AoS structure (Array of Structures, e.g. XyzXyz) as
/// opposed to a SoA layout (Structure of Arrays, e.g. XxYyZz.)
AoS = 0b10000,
}

/// Desired coordinates and the output layout.
Expand Down Expand Up @@ -303,7 +304,7 @@ impl Geometry {
let component_size = size_of::<f64>();

let (stride, ptr_offset): (c_int, isize) =
match layout.has_component(CoordinateComponent::Sequential) {
match layout.has_component(CoordinateComponent::AoS) {
true => (
(component_size * coord_size) as c_int,
component_size as isize,
Expand All @@ -315,6 +316,7 @@ impl Geometry {
};

unsafe {
gdal_sys::CPLErrorReset();
let data = in_points.as_ptr() as *mut c_void;

// Per-component ptr offset.
Expand Down Expand Up @@ -353,18 +355,22 @@ impl Geometry {
stride,
);
}

let err = unsafe { gdal_sys::CPLGetLastErrorType() };

if err != 0 {
return Err(_last_cpl_err(err));
}

Ok(())
}

/// Writes all points in the geometry to `out_points` according to the specified `layout`.
///
/// For some geometry types, like polygons, that don't consist of points, 0 will be returned
/// and `out_points` will remain unmodified.
pub fn get_points(&self, out_points: &mut Vec<f64>, layout: CoordinateLayout) -> usize {
/// For some geometry types, like polygons, that don't consist of points, Err will be returned
/// and `out_points` will only be resized.
pub fn get_points(&self, out_points: &mut Vec<f64>, layout: CoordinateLayout) -> Result<usize> {
let num_points = unsafe { gdal_sys::OGR_G_GetPointCount(self.c_geometry()) } as usize;
if num_points == 0 {
return 0;
}

// Number of dims
let coord_size = layout.coordinate_size();
Expand All @@ -377,7 +383,7 @@ impl Geometry {
let component_size = size_of::<f64>();

let (stride, ptr_offset): (c_int, isize) =
match layout.has_component(CoordinateComponent::Sequential) {
match layout.has_component(CoordinateComponent::AoS) {
true => (
(component_size * coord_size) as c_int,
component_size as isize,
Expand All @@ -389,6 +395,8 @@ impl Geometry {
};

unsafe {
gdal_sys::CPLErrorReset();

let data = out_points.as_mut_ptr() as *mut c_void;

// Per-component ptr offset.
Expand All @@ -414,6 +422,8 @@ impl Geometry {
let m_ptr = component_loc(&layout, CoordinateComponent::M, &mut curr_ptr_offset);

// Should be OK to just use the offset even for unused components...
//
// Returns number of points.
gdal_sys::OGR_G_GetPointsZM(
self.c_geometry(),
x_ptr,
Expand All @@ -425,9 +435,15 @@ impl Geometry {
m_ptr,
stride,
);
};

let err = unsafe { gdal_sys::CPLGetLastErrorType() };

if err != 0 {
return Err(_last_cpl_err(err));
}

length
Ok(num_points)
}

/// Get the geometry type ordinal
Expand Down Expand Up @@ -710,9 +726,9 @@ mod tests {
assert!(CoordinateLayout::XyzXyz.has_component(CoordinateComponent::Y));
assert!(CoordinateLayout::XyzXyz.has_component(CoordinateComponent::Z));
assert!(!CoordinateLayout::XyzXyz.has_component(CoordinateComponent::M));
assert!(CoordinateLayout::XyzXyz.has_component(CoordinateComponent::Sequential));
assert!(CoordinateLayout::XyzXyz.has_component(CoordinateComponent::AoS));

assert!(!CoordinateLayout::XxYy.has_component(CoordinateComponent::Sequential));
assert!(!CoordinateLayout::XxYy.has_component(CoordinateComponent::AoS));

assert!(!CoordinateLayout::XymXym.has_component(CoordinateComponent::Z));
assert!(CoordinateLayout::XymXym.has_component(CoordinateComponent::M));
Expand Down Expand Up @@ -861,19 +877,20 @@ mod tests {
ring.add_point_2d((1179091.1646903288, 712782.8838459781));
assert!(!ring.is_empty());
let mut ring_vec: Vec<f64> = Vec::new();
ring.get_points(&mut ring_vec, CoordinateLayout::XyXy);
ring.get_points(&mut ring_vec, CoordinateLayout::XyXy)
.unwrap();
assert_eq!(ring_vec.len(), 6 * 2);
let mut poly = Geometry::empty(wkbPolygon).unwrap();
poly.add_geometry(ring.to_owned()).unwrap();
let mut poly_vec: Vec<f64> = Vec::new();
poly.get_points(&mut poly_vec, CoordinateLayout::XyXy);
// Points are in ring, not containing geometry.
// NB: In Python SWIG bindings, `GetPoints` is fallible.
assert!(poly_vec.is_empty());
let res = poly.get_points(&mut poly_vec, CoordinateLayout::XyXy);
assert!(res.is_err());
assert_eq!(poly.geometry_count(), 1);
let ring_out = poly.get_geometry(0);
let mut ring_out_vec: Vec<f64> = Vec::new();
ring_out.get_points(&mut ring_out_vec, CoordinateLayout::XyXy);
ring_out
.get_points(&mut ring_out_vec, CoordinateLayout::XyXy)
.unwrap();
// NB: `wkb()` shows it to be a `LINEARRING`, but returned type is LineString
assert_eq!(ring_out.geometry_type(), wkbLineString);
assert!(!&ring_out.is_empty());
Expand All @@ -890,7 +907,9 @@ mod tests {
assert!(geom.json().unwrap().contains("Polygon"));
let inner = geom.get_geometry(0);
let mut points: Vec<f64> = Vec::new();
inner.get_points(&mut points, CoordinateLayout::XyXy);
inner
.get_points(&mut points, CoordinateLayout::XyXy)
.unwrap();
assert!(!points.is_empty());
}

Expand Down Expand Up @@ -958,43 +977,50 @@ mod tests {
line.add_point_zm((1.0, 2.0, 0.5, 1.0));
let mut line_points: Vec<f64> = Vec::new();

line.get_points(&mut line_points, CoordinateLayout::XyXy);
line.get_points(&mut line_points, CoordinateLayout::XyXy)
.unwrap();
assert_eq!(line_points.len(), 3 * 2);
assert_eq!(line_points, vec![0.0, 0.0, 1.0, 0.0, 1.0, 2.0]);

line.get_points(&mut line_points, CoordinateLayout::XxYy);
line.get_points(&mut line_points, CoordinateLayout::XxYy)
.unwrap();
assert_eq!(line_points.len(), 3 * 2);
assert_eq!(line_points, vec![0.0, 1.0, 1.0, 0.0, 0.0, 2.0]);

line.get_points(&mut line_points, CoordinateLayout::XyzXyz);
line.get_points(&mut line_points, CoordinateLayout::XyzXyz)
.unwrap();
assert_eq!(line_points.len(), 3 * 3);
assert_eq!(
line_points,
vec![0.0, 0.0, 0.0, 1.0, 0.0, 0.25, 1.0, 2.0, 0.5]
);

line.get_points(&mut line_points, CoordinateLayout::XxYyZz);
line.get_points(&mut line_points, CoordinateLayout::XxYyZz)
.unwrap();
assert_eq!(line_points.len(), 3 * 3);
assert_eq!(
line_points,
vec![0.0, 1.0, 1.0, 0.0, 0.0, 2.0, 0.0, 0.25, 0.5]
);

line.get_points(&mut line_points, CoordinateLayout::XymXym);
line.get_points(&mut line_points, CoordinateLayout::XymXym)
.unwrap();
assert_eq!(line_points.len(), 3 * 3);
assert_eq!(
line_points,
vec![0.0, 0.0, 0.0, 1.0, 0.0, 0.5, 1.0, 2.0, 1.0]
);

line.get_points(&mut line_points, CoordinateLayout::XyzmXyzm);
line.get_points(&mut line_points, CoordinateLayout::XyzmXyzm)
.unwrap();
assert_eq!(line_points.len(), 3 * 4);
assert_eq!(
line_points,
vec![0.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.25, 0.5, 1.0, 2.0, 0.5, 1.0]
);

line.get_points(&mut line_points, CoordinateLayout::XxYyZzMm);
line.get_points(&mut line_points, CoordinateLayout::XxYyZzMm)
.unwrap();
assert_eq!(line_points.len(), 3 * 4);
assert_eq!(
line_points,
Expand Down
3 changes: 2 additions & 1 deletion src/vector/layer.rs
Original file line number Diff line number Diff line change
Expand Up @@ -1272,7 +1272,8 @@ mod tests {
let geom = feature.geometry().unwrap();
assert_eq!(geom.geometry_type(), OGRwkbGeometryType::wkbLineString);
let mut coords: Vec<f64> = Vec::new();
geom.get_points(&mut coords, CoordinateLayout::XyXy);
geom.get_points(&mut coords, CoordinateLayout::XyXy)
.unwrap();
assert_eq!(
coords,
[26.1019276, 44.4302748, 26.1019382, 44.4303191, 26.1020002, 44.4304202]
Expand Down
2 changes: 1 addition & 1 deletion src/vector/mod.rs
Original file line number Diff line number Diff line change
Expand Up @@ -30,7 +30,7 @@
//! // Summarize the geometry
//! let geometry = feature.geometry().unwrap();
//! let geom_type = geometry_type_to_name(geometry.geometry_type());
//! let geom_len = geometry.get_points(&mut Vec::new(), CoordinateLayout::XyXy);
//! let geom_len = geometry.get_points(&mut Vec::new(), CoordinateLayout::XyXy).unwrap();
//! println!(" Feature fid={fid:?}, geometry_type='{geom_type}', geometry_len={geom_len}");
//! // Get all the available fields and print their values
//! for field in feature.fields() {
Expand Down
2 changes: 1 addition & 1 deletion src/vector/ops/conversions/gdal_to_geo.rs
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,7 @@ impl TryFrom<&Geometry> for geo_types::Geometry<f64> {
}
OGRwkbGeometryType::wkbLineString => {
let mut gdal_coords: Vec<f64> = Vec::new();
geo.get_points(&mut gdal_coords, CoordinateLayout::XyXy);
geo.get_points(&mut gdal_coords, CoordinateLayout::XyXy)?;
let coords = gdal_coords
.chunks(2)
.map(|points| geo_types::Coord {
Expand Down

0 comments on commit bea5eeb

Please sign in to comment.