From 5bf7ad1208069d38da5f1b18d5ea69091abafc36 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:27:09 +0800 Subject: [PATCH 01/11] feat(sedona-gdal): add geometry and spatial ref primitives --- c/sedona-gdal/src/errors.rs | 16 +++ c/sedona-gdal/src/geo_transform.rs | 152 ++++++++++++++++++++++ c/sedona-gdal/src/lib.rs | 4 + c/sedona-gdal/src/spatial_ref.rs | 155 ++++++++++++++++++++++ c/sedona-gdal/src/vector/geometry.rs | 162 +++++++++++++++++++++++ c/sedona-gdal/src/vector/mod.rs | 20 +++ c/sedona-gdal/src/vsi.rs | 187 +++++++++++++++++++++++++++ 7 files changed, 696 insertions(+) create mode 100644 c/sedona-gdal/src/geo_transform.rs create mode 100644 c/sedona-gdal/src/spatial_ref.rs create mode 100644 c/sedona-gdal/src/vector/geometry.rs create mode 100644 c/sedona-gdal/src/vector/mod.rs create mode 100644 c/sedona-gdal/src/vsi.rs diff --git a/c/sedona-gdal/src/errors.rs b/c/sedona-gdal/src/errors.rs index 6f2ad1a0f..f344667ec 100644 --- a/c/sedona-gdal/src/errors.rs +++ b/c/sedona-gdal/src/errors.rs @@ -20,6 +20,7 @@ //! Original code is licensed under MIT. use std::ffi::NulError; +use std::num::TryFromIntError; use thiserror::Error; @@ -45,8 +46,23 @@ pub enum GdalError { #[error("Bad argument: {0}")] BadArgument(String), + #[error("GDAL method '{method_name}' returned a NULL pointer. Error msg: '{msg}'")] + NullPointer { + method_name: &'static str, + msg: String, + }, + + #[error("OGR method '{method_name}' returned error: '{err:?}'")] + OgrError { err: i32, method_name: &'static str }, + + #[error("Unable to unlink mem file: {file_name}")] + UnlinkMemFile { file_name: String }, + #[error("FFI NUL error: {0}")] FfiNulError(#[from] NulError), + + #[error(transparent)] + IntConversionError(#[from] TryFromIntError), } pub type Result = std::result::Result; diff --git a/c/sedona-gdal/src/geo_transform.rs b/c/sedona-gdal/src/geo_transform.rs new file mode 100644 index 000000000..5504e8501 --- /dev/null +++ b/c/sedona-gdal/src/geo_transform.rs @@ -0,0 +1,152 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. +//! +//! GeoTransform type and extension trait. +//! +//! The [`apply`](GeoTransformEx::apply) and [`invert`](GeoTransformEx::invert) +//! methods are pure-Rust reimplementations of GDAL's `GDALApplyGeoTransform` +//! and `GDALInvGeoTransform` (from `alg/gdaltransformer.cpp`). No FFI call or +//! thread-local state is needed. + +use crate::errors; +use crate::errors::GdalError; + +/// An affine geo-transform: six coefficients mapping pixel/line to projection coordinates. +/// +/// - `[0]`: x-coordinate of the upper-left corner of the upper-left pixel. +/// - `[1]`: W-E pixel resolution (pixel width). +/// - `[2]`: row rotation (typically zero). +/// - `[3]`: y-coordinate of the upper-left corner of the upper-left pixel. +/// - `[4]`: column rotation (typically zero). +/// - `[5]`: N-S pixel resolution (pixel height, negative for North-up). +pub type GeoTransform = [f64; 6]; + +/// Extension methods on [`GeoTransform`]. +pub trait GeoTransformEx { + /// Apply the geo-transform to a pixel/line coordinate, returning (geo_x, geo_y). + fn apply(&self, x: f64, y: f64) -> (f64, f64); + + /// Invert this geo-transform, returning the inverse coefficients for + /// computing (geo_x, geo_y) -> (x, y) transformations. + fn invert(&self) -> errors::Result; +} + +impl GeoTransformEx for GeoTransform { + /// Pure-Rust equivalent of GDAL's `GDALApplyGeoTransform`. + fn apply(&self, x: f64, y: f64) -> (f64, f64) { + let geo_x = self[0] + x * self[1] + y * self[2]; + let geo_y = self[3] + x * self[4] + y * self[5]; + (geo_x, geo_y) + } + + /// Pure-Rust equivalent of GDAL's `GDALInvGeoTransform`. + fn invert(&self) -> errors::Result { + let gt = self; + + // Fast path: no rotation/skew — avoid determinant and precision issues. + if gt[2] == 0.0 && gt[4] == 0.0 && gt[1] != 0.0 && gt[5] != 0.0 { + return Ok([ + -gt[0] / gt[1], + 1.0 / gt[1], + 0.0, + -gt[3] / gt[5], + 0.0, + 1.0 / gt[5], + ]); + } + + // General case: 2x2 matrix inverse via adjugate / determinant. + let det = gt[1] * gt[5] - gt[2] * gt[4]; + let magnitude = gt[1] + .abs() + .max(gt[2].abs()) + .max(gt[4].abs().max(gt[5].abs())); + + if det.abs() <= 1e-10 * magnitude * magnitude { + return Err(GdalError::BadArgument( + "Geo transform is uninvertible".to_string(), + )); + } + + let inv_det = 1.0 / det; + + Ok([ + (gt[2] * gt[3] - gt[0] * gt[5]) * inv_det, + gt[5] * inv_det, + -gt[2] * inv_det, + (-gt[1] * gt[3] + gt[0] * gt[4]) * inv_det, + -gt[4] * inv_det, + gt[1] * inv_det, + ]) + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn test_apply_no_rotation() { + // Origin at (100, 200), 10m pixels, north-up + let gt: GeoTransform = [100.0, 10.0, 0.0, 200.0, 0.0, -10.0]; + let (x, y) = gt.apply(5.0, 3.0); + assert!((x - 150.0).abs() < 1e-12); + assert!((y - 170.0).abs() < 1e-12); + } + + #[test] + fn test_apply_with_rotation() { + let gt: GeoTransform = [100.0, 10.0, 2.0, 200.0, 3.0, -10.0]; + let (x, y) = gt.apply(5.0, 3.0); + // 100 + 5*10 + 3*2 = 156 + assert!((x - 156.0).abs() < 1e-12); + // 200 + 5*3 + 3*(-10) = 185 + assert!((y - 185.0).abs() < 1e-12); + } + + #[test] + fn test_invert_no_rotation() { + let gt: GeoTransform = [100.0, 10.0, 0.0, 200.0, 0.0, -10.0]; + let inv = gt.invert().unwrap(); + // Round-trip: apply then apply inverse should recover pixel/line. + let (geo_x, geo_y) = gt.apply(7.0, 4.0); + let (px, ln) = inv.apply(geo_x, geo_y); + assert!((px - 7.0).abs() < 1e-10); + assert!((ln - 4.0).abs() < 1e-10); + } + + #[test] + fn test_invert_with_rotation() { + let gt: GeoTransform = [100.0, 10.0, 2.0, 200.0, 3.0, -10.0]; + let inv = gt.invert().unwrap(); + let (geo_x, geo_y) = gt.apply(7.0, 4.0); + let (px, ln) = inv.apply(geo_x, geo_y); + assert!((px - 7.0).abs() < 1e-10); + assert!((ln - 4.0).abs() < 1e-10); + } + + #[test] + fn test_invert_singular() { + // Determinant is zero: both rows are proportional. + let gt: GeoTransform = [0.0, 1.0, 2.0, 0.0, 2.0, 4.0]; + assert!(gt.invert().is_err()); + } +} diff --git a/c/sedona-gdal/src/lib.rs b/c/sedona-gdal/src/lib.rs index b64a2275c..0646a241a 100644 --- a/c/sedona-gdal/src/lib.rs +++ b/c/sedona-gdal/src/lib.rs @@ -29,4 +29,8 @@ pub mod global; // --- High-level wrappers --- pub mod config; pub mod cpl; +pub mod geo_transform; pub mod raster; +pub mod spatial_ref; +pub mod vector; +pub mod vsi; diff --git a/c/sedona-gdal/src/spatial_ref.rs b/c/sedona-gdal/src/spatial_ref.rs new file mode 100644 index 000000000..61aa37d10 --- /dev/null +++ b/c/sedona-gdal/src/spatial_ref.rs @@ -0,0 +1,155 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ffi::{CStr, CString}; +use std::ptr; + +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; + +/// An OGR spatial reference system. +pub struct SpatialRef { + api: &'static GdalApi, + c_srs: OGRSpatialReferenceH, +} + +unsafe impl Send for SpatialRef {} + +impl Drop for SpatialRef { + fn drop(&mut self) { + if !self.c_srs.is_null() { + unsafe { call_gdal_api!(self.api, OSRRelease, self.c_srs) }; + } + } +} + +impl SpatialRef { + /// Create a new SpatialRef from a WKT string. + pub fn from_wkt(api: &'static GdalApi, wkt: &str) -> Result { + let c_wkt = CString::new(wkt)?; + let c_srs = unsafe { call_gdal_api!(api, OSRNewSpatialReference, c_wkt.as_ptr()) }; + if c_srs.is_null() { + return Err(GdalError::NullPointer { + method_name: "OSRNewSpatialReference", + msg: "failed to create spatial reference from WKT".to_string(), + }); + } + Ok(Self { api, c_srs }) + } + + /// Create a SpatialRef by cloning a borrowed C handle via `OSRClone`. + /// + /// # Safety + /// + /// The caller must ensure `c_srs` is a valid `OGRSpatialReferenceH`. + pub unsafe fn from_c_srs_clone( + api: &'static GdalApi, + c_srs: OGRSpatialReferenceH, + ) -> Result { + let cloned = call_gdal_api!(api, OSRClone, c_srs); + if cloned.is_null() { + return Err(GdalError::NullPointer { + method_name: "OSRClone", + msg: "failed to clone spatial reference".to_string(), + }); + } + Ok(Self { api, c_srs: cloned }) + } + + /// Return the raw C handle. + pub fn c_srs(&self) -> OGRSpatialReferenceH { + self.c_srs + } + + /// Export to PROJJSON string. + pub fn to_projjson(&self) -> Result { + unsafe { + let mut ptr: *mut std::os::raw::c_char = ptr::null_mut(); + let rv = call_gdal_api!( + self.api, + OSRExportToPROJJSON, + self.c_srs, + &mut ptr, + ptr::null() + ); + if rv != crate::gdal_dyn_bindgen::OGRERR_NONE || ptr.is_null() { + return Err(GdalError::NullPointer { + method_name: "OSRExportToPROJJSON", + msg: "returned null".to_string(), + }); + } + let result = CStr::from_ptr(ptr).to_string_lossy().into_owned(); + call_gdal_api!(self.api, VSIFree, ptr as *mut std::ffi::c_void); + Ok(result) + } + } +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use crate::errors::GdalError; + use crate::global::with_global_gdal_api; + use crate::spatial_ref::SpatialRef; + + const WGS84_WKT: &str = r#"GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]]"#; + + #[test] + fn test_from_wkt() { + with_global_gdal_api(|api| { + let srs = SpatialRef::from_wkt(api, WGS84_WKT).unwrap(); + assert!(!srs.c_srs().is_null()); + }) + .unwrap(); + } + + #[test] + fn test_from_wkt_invalid() { + with_global_gdal_api(|api| { + let err = SpatialRef::from_wkt(api, "WGS\u{0}84"); + assert!(matches!(err, Err(GdalError::FfiNulError(_)))); + }) + .unwrap(); + } + + #[test] + fn test_to_projjson() { + with_global_gdal_api(|api| { + let srs = SpatialRef::from_wkt(api, WGS84_WKT).unwrap(); + let projjson = srs.to_projjson().unwrap(); + assert!( + projjson.contains("WGS 84"), + "unexpected projjson: {projjson}" + ); + }) + .unwrap(); + } + + #[test] + fn test_from_c_srs_clone() { + with_global_gdal_api(|api| { + let srs = SpatialRef::from_wkt(api, WGS84_WKT).unwrap(); + let cloned = unsafe { SpatialRef::from_c_srs_clone(api, srs.c_srs()) }.unwrap(); + assert_eq!(srs.to_projjson().unwrap(), cloned.to_projjson().unwrap()); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/vector/geometry.rs b/c/sedona-gdal/src/vector/geometry.rs new file mode 100644 index 000000000..aefc153a6 --- /dev/null +++ b/c/sedona-gdal/src/vector/geometry.rs @@ -0,0 +1,162 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ffi::CString; +use std::ptr; + +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; + +/// An OGR geometry. +pub struct Geometry { + api: &'static GdalApi, + c_geom: OGRGeometryH, +} + +unsafe impl Send for Geometry {} + +impl Drop for Geometry { + fn drop(&mut self) { + if !self.c_geom.is_null() { + unsafe { call_gdal_api!(self.api, OGR_G_DestroyGeometry, self.c_geom) }; + } + } +} + +impl Geometry { + /// Create a geometry from WKB bytes. + pub fn from_wkb(api: &'static GdalApi, wkb: &[u8]) -> Result { + let wkb_len: i32 = wkb.len().try_into()?; + let mut c_geom: OGRGeometryH = ptr::null_mut(); + let rv = unsafe { + call_gdal_api!( + api, + OGR_G_CreateFromWkb, + wkb.as_ptr() as *const std::ffi::c_void, + ptr::null_mut(), // hSRS + &mut c_geom, + wkb_len + ) + }; + if rv != OGRERR_NONE { + return Err(GdalError::OgrError { + err: rv, + method_name: "OGR_G_CreateFromWkb", + }); + } + if c_geom.is_null() { + return Err(GdalError::NullPointer { + method_name: "OGR_G_CreateFromWkb", + msg: "returned null geometry".to_string(), + }); + } + Ok(Self { api, c_geom }) + } + + /// Create a geometry from WKT string. + pub fn from_wkt(api: &'static GdalApi, wkt: &str) -> Result { + let c_wkt = CString::new(wkt)?; + let mut wkt_ptr = c_wkt.as_ptr() as *mut std::os::raw::c_char; + let mut c_geom: OGRGeometryH = ptr::null_mut(); + let rv = unsafe { + call_gdal_api!( + api, + OGR_G_CreateFromWkt, + &mut wkt_ptr, + ptr::null_mut(), // hSRS + &mut c_geom + ) + }; + if rv != OGRERR_NONE { + return Err(GdalError::OgrError { + err: rv, + method_name: "OGR_G_CreateFromWkt", + }); + } + if c_geom.is_null() { + return Err(GdalError::NullPointer { + method_name: "OGR_G_CreateFromWkt", + msg: "returned null geometry".to_string(), + }); + } + Ok(Self { api, c_geom }) + } + + /// Return the raw C geometry handle. + pub fn c_geometry(&self) -> OGRGeometryH { + self.c_geom + } + + /// Get the bounding envelope. + pub fn envelope(&self) -> Envelope { + let mut env = OGREnvelope { + MinX: 0.0, + MaxX: 0.0, + MinY: 0.0, + MaxY: 0.0, + }; + unsafe { call_gdal_api!(self.api, OGR_G_GetEnvelope, self.c_geom, &mut env) }; + Envelope { + MinX: env.MinX, + MaxX: env.MaxX, + MinY: env.MinY, + MaxY: env.MaxY, + } + } + + /// Export to ISO WKB. + pub fn wkb(&self) -> Result> { + let size = unsafe { call_gdal_api!(self.api, OGR_G_WkbSize, self.c_geom) }; + if size < 0 { + return Err(GdalError::BadArgument(format!( + "OGR_G_WkbSize returned negative size: {size}" + ))); + } + let mut buf = vec![0u8; size as usize]; + let rv = unsafe { + call_gdal_api!( + self.api, + OGR_G_ExportToIsoWkb, + self.c_geom, + wkbNDR, // little-endian + buf.as_mut_ptr() + ) + }; + if rv != OGRERR_NONE { + return Err(GdalError::OgrError { + err: rv, + method_name: "OGR_G_ExportToIsoWkb", + }); + } + Ok(buf) + } +} + +/// Bounding envelope. +#[derive(Debug, Clone, Copy)] +#[allow(non_snake_case)] +pub struct Envelope { + pub MinX: f64, + pub MaxX: f64, + pub MinY: f64, + pub MaxY: f64, +} diff --git a/c/sedona-gdal/src/vector/mod.rs b/c/sedona-gdal/src/vector/mod.rs new file mode 100644 index 000000000..14c13d8a5 --- /dev/null +++ b/c/sedona-gdal/src/vector/mod.rs @@ -0,0 +1,20 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +pub mod geometry; + +pub use geometry::{Envelope, Geometry}; diff --git a/c/sedona-gdal/src/vsi.rs b/c/sedona-gdal/src/vsi.rs new file mode 100644 index 000000000..7c4e93f2f --- /dev/null +++ b/c/sedona-gdal/src/vsi.rs @@ -0,0 +1,187 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. +//! +//! GDAL Virtual File System (VSI) wrappers. + +use std::ffi::CString; + +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; + +/// Creates a new VSI in-memory file from a given buffer. +/// +/// The data is copied into GDAL-allocated memory (via `VSIMalloc`) so that +/// GDAL can safely free it with `VSIFree` when ownership is taken. +pub fn create_mem_file(api: &'static GdalApi, file_name: &str, data: Vec) -> Result<()> { + let c_file_name = CString::new(file_name)?; + let len = data.len(); + + // Allocate via GDAL's allocator so GDAL can safely free it. + let gdal_buf = unsafe { call_gdal_api!(api, VSIMalloc, len) } as *mut u8; + if gdal_buf.is_null() { + return Err(GdalError::NullPointer { + method_name: "VSIMalloc", + msg: format!("failed to allocate {len} bytes"), + }); + } + + // Copy data into GDAL-allocated buffer + unsafe { + std::ptr::copy_nonoverlapping(data.as_ptr(), gdal_buf, len); + } + // Rust Vec is dropped here, freeing the Rust-allocated memory. + + let handle = unsafe { + call_gdal_api!( + api, + VSIFileFromMemBuffer, + c_file_name.as_ptr(), + gdal_buf, + len as i64, + 1 // bTakeOwnership = true — GDAL will VSIFree gdal_buf + ) + }; + + if handle.is_null() { + // GDAL did not take ownership, so we must free. + unsafe { call_gdal_api!(api, VSIFree, gdal_buf as *mut std::ffi::c_void) }; + return Err(GdalError::NullPointer { + method_name: "VSIFileFromMemBuffer", + msg: String::new(), + }); + } + + unsafe { + call_gdal_api!(api, VSIFCloseL, handle); + } + + Ok(()) +} + +/// Unlink (delete) a VSI in-memory file. +pub fn unlink_mem_file(api: &'static GdalApi, file_name: &str) -> Result<()> { + let c_file_name = CString::new(file_name)?; + + let rv = unsafe { call_gdal_api!(api, VSIUnlink, c_file_name.as_ptr()) }; + + if rv != 0 { + return Err(GdalError::UnlinkMemFile { + file_name: file_name.to_string(), + }); + } + + Ok(()) +} + +/// Copies the bytes of the VSI in-memory file, taking ownership and freeing the GDAL memory. +pub fn get_vsi_mem_file_bytes_owned(api: &'static GdalApi, file_name: &str) -> Result> { + let c_file_name = CString::new(file_name)?; + + let owned_bytes = unsafe { + let mut length: i64 = 0; + let bytes = call_gdal_api!( + api, + VSIGetMemFileBuffer, + c_file_name.as_ptr(), + &mut length, + 1 // bUnlinkAndSeize = true + ); + + if bytes.is_null() { + return Err(GdalError::NullPointer { + method_name: "VSIGetMemFileBuffer", + msg: String::new(), + }); + } + + if length < 0 { + call_gdal_api!(api, VSIFree, bytes.cast::()); + return Err(GdalError::BadArgument(format!( + "VSIGetMemFileBuffer returned negative length: {length}" + ))); + } + + let slice = std::slice::from_raw_parts(bytes, length as usize); + let vec = slice.to_vec(); + + call_gdal_api!(api, VSIFree, bytes.cast::()); + + vec + }; + + Ok(owned_bytes) +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use super::*; + use crate::global::with_global_gdal_api; + + #[test] + fn create_and_retrieve_mem_file() { + let file_name = "/vsimem/525ebf24-a030-4677-bb4e-a921741cabe0"; + + with_global_gdal_api(|api| { + create_mem_file(api, file_name, vec![1_u8, 2, 3, 4]).unwrap(); + + let bytes = get_vsi_mem_file_bytes_owned(api, file_name).unwrap(); + + assert_eq!(bytes, vec![1_u8, 2, 3, 4]); + + // mem file must not be there anymore + assert!(matches!( + unlink_mem_file(api, file_name).unwrap_err(), + GdalError::UnlinkMemFile { + file_name: err_file_name + } + if err_file_name == file_name + )); + }) + .unwrap(); + } + + #[test] + fn create_and_unlink_mem_file() { + let file_name = "/vsimem/bbf5f1d6-c1e9-4469-a33b-02cd9173132d"; + + with_global_gdal_api(|api| { + create_mem_file(api, file_name, vec![1_u8, 2, 3, 4]).unwrap(); + + unlink_mem_file(api, file_name).unwrap(); + }) + .unwrap(); + } + + #[test] + fn no_mem_file() { + with_global_gdal_api(|api| { + assert!(matches!( + get_vsi_mem_file_bytes_owned(api, "foobar").unwrap_err(), + GdalError::NullPointer { + method_name: "VSIGetMemFileBuffer", + msg, + } + if msg.is_empty() + )); + }) + .unwrap(); + } +} From 08b3a49e8422cf5575883e80682143a855557333 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:41:46 +0800 Subject: [PATCH 02/11] refactor(sedona-gdal): use flat vector module file --- c/sedona-gdal/src/spatial_ref.rs | 3 ++- c/sedona-gdal/src/{vector/mod.rs => vector.rs} | 0 2 files changed, 2 insertions(+), 1 deletion(-) rename c/sedona-gdal/src/{vector/mod.rs => vector.rs} (100%) diff --git a/c/sedona-gdal/src/spatial_ref.rs b/c/sedona-gdal/src/spatial_ref.rs index 61aa37d10..77559a0de 100644 --- a/c/sedona-gdal/src/spatial_ref.rs +++ b/c/sedona-gdal/src/spatial_ref.rs @@ -24,6 +24,7 @@ use std::ptr; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::OGRERR_NONE; use crate::gdal_dyn_bindgen::*; /// An OGR spatial reference system. @@ -91,7 +92,7 @@ impl SpatialRef { &mut ptr, ptr::null() ); - if rv != crate::gdal_dyn_bindgen::OGRERR_NONE || ptr.is_null() { + if rv != OGRERR_NONE || ptr.is_null() { return Err(GdalError::NullPointer { method_name: "OSRExportToPROJJSON", msg: "returned null".to_string(), diff --git a/c/sedona-gdal/src/vector/mod.rs b/c/sedona-gdal/src/vector.rs similarity index 100% rename from c/sedona-gdal/src/vector/mod.rs rename to c/sedona-gdal/src/vector.rs From e84c3a0f2440f94f73b28c5a3565eef7a57ca684 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Fri, 13 Mar 2026 14:23:04 +0800 Subject: [PATCH 03/11] refactor(sedona-gdal): remove vector re-export aliases --- c/sedona-gdal/src/vector.rs | 2 -- 1 file changed, 2 deletions(-) diff --git a/c/sedona-gdal/src/vector.rs b/c/sedona-gdal/src/vector.rs index 14c13d8a5..10c038edc 100644 --- a/c/sedona-gdal/src/vector.rs +++ b/c/sedona-gdal/src/vector.rs @@ -16,5 +16,3 @@ // under the License. pub mod geometry; - -pub use geometry::{Envelope, Geometry}; From 5d718aebf2a1694ef6cb270100d8bdcdd22e83a8 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:28:06 +0800 Subject: [PATCH 04/11] feat(sedona-gdal): add dataset and vector/raster wrappers --- c/sedona-gdal/src/dataset.rs | 470 +++++++++++++++++++++++++ c/sedona-gdal/src/driver.rs | 246 +++++++++++++ c/sedona-gdal/src/errors.rs | 3 + c/sedona-gdal/src/lib.rs | 2 + c/sedona-gdal/src/raster.rs | 3 + c/sedona-gdal/src/raster/rasterband.rs | 368 +++++++++++++++++++ c/sedona-gdal/src/vector.rs | 2 + c/sedona-gdal/src/vector/feature.rs | 218 ++++++++++++ c/sedona-gdal/src/vector/layer.rs | 118 +++++++ 9 files changed, 1430 insertions(+) create mode 100644 c/sedona-gdal/src/dataset.rs create mode 100644 c/sedona-gdal/src/driver.rs create mode 100644 c/sedona-gdal/src/raster/rasterband.rs create mode 100644 c/sedona-gdal/src/vector/feature.rs create mode 100644 c/sedona-gdal/src/vector/layer.rs diff --git a/c/sedona-gdal/src/dataset.rs b/c/sedona-gdal/src/dataset.rs new file mode 100644 index 000000000..b286535b0 --- /dev/null +++ b/c/sedona-gdal/src/dataset.rs @@ -0,0 +1,470 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ffi::{CStr, CString}; +use std::ptr; + +use crate::cpl::CslStringList; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; +use crate::raster::types::{DatasetOptions, GdalDataType as RustGdalDataType}; +use crate::raster::RasterBand; +use crate::spatial_ref::SpatialRef; +use crate::vector::Layer; + +/// A GDAL dataset. +pub struct Dataset { + api: &'static GdalApi, + c_dataset: GDALDatasetH, + owned: bool, +} + +unsafe impl Send for Dataset {} + +impl Drop for Dataset { + fn drop(&mut self) { + if self.owned && !self.c_dataset.is_null() { + unsafe { call_gdal_api!(self.api, GDALClose, self.c_dataset) }; + } + } +} + +impl Dataset { + /// Open a dataset with extended options. + pub fn open_ex( + api: &'static GdalApi, + path: &str, + open_flags: GDALOpenFlags, + allowed_drivers: Option<&[&str]>, + open_options: Option<&[&str]>, + sibling_files: Option<&[&str]>, + ) -> Result { + let c_path = CString::new(path)?; + + // Build CslStringLists from Option<&[&str]>. + // None → null pointer (use GDAL default). + // Some(&[]) → pointer to [null] (explicitly empty list). + let drivers_csl = allowed_drivers + .map(|v| CslStringList::try_from_iter(v.iter().copied())) + .transpose()?; + let options_csl = open_options + .map(|v| CslStringList::try_from_iter(v.iter().copied())) + .transpose()?; + let siblings_csl = sibling_files + .map(|v| CslStringList::try_from_iter(v.iter().copied())) + .transpose()?; + + let c_dataset = unsafe { + call_gdal_api!( + api, + GDALOpenEx, + c_path.as_ptr(), + open_flags, + drivers_csl + .as_ref() + .map_or(ptr::null(), |csl| csl.as_ptr() as *const *const _), + options_csl + .as_ref() + .map_or(ptr::null(), |csl| csl.as_ptr() as *const *const _), + siblings_csl + .as_ref() + .map_or(ptr::null(), |csl| csl.as_ptr() as *const *const _) + ) + }; + + if c_dataset.is_null() { + return Err(api.last_cpl_err(CE_Failure as u32)); + } + + Ok(Self { + api, + c_dataset, + owned: true, + }) + } + + /// Create a new owned Dataset from a C handle. + pub(crate) fn new_owned(api: &'static GdalApi, c_dataset: GDALDatasetH) -> Self { + Self { + api, + c_dataset, + owned: true, + } + } + + /// Wrap an existing C dataset handle (non-owning). + /// + /// # Safety + /// + /// The caller must ensure the handle is valid and outlives this `Dataset`. + pub unsafe fn from_c_dataset(api: &'static GdalApi, c_dataset: GDALDatasetH) -> Self { + Self { + api, + c_dataset, + owned: false, + } + } + + /// Return the raw C dataset handle. + pub fn c_dataset(&self) -> GDALDatasetH { + self.c_dataset + } + + /// Return raster size as (x_size, y_size). + pub fn raster_size(&self) -> (usize, usize) { + let x = unsafe { call_gdal_api!(self.api, GDALGetRasterXSize, self.c_dataset) }; + let y = unsafe { call_gdal_api!(self.api, GDALGetRasterYSize, self.c_dataset) }; + (x as usize, y as usize) + } + + /// Return the number of raster bands. + pub fn raster_count(&self) -> usize { + unsafe { call_gdal_api!(self.api, GDALGetRasterCount, self.c_dataset) as usize } + } + + /// Get a raster band (1-indexed). + pub fn rasterband(&self, band_index: usize) -> Result> { + let band_index_i32 = i32::try_from(band_index)?; + let c_band = + unsafe { call_gdal_api!(self.api, GDALGetRasterBand, self.c_dataset, band_index_i32) }; + if c_band.is_null() { + return Err(GdalError::NullPointer { + method_name: "GDALGetRasterBand", + msg: format!("band index {band_index}"), + }); + } + Ok(RasterBand::new(self.api, c_band, self)) + } + + /// Get the geo-transform. + pub fn geo_transform(&self) -> Result<[f64; 6]> { + let mut gt = [0.0f64; 6]; + let rv = unsafe { + call_gdal_api!( + self.api, + GDALGetGeoTransform, + self.c_dataset, + gt.as_mut_ptr() + ) + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(gt) + } + + /// Set the geo-transform. + pub fn set_geo_transform(&self, gt: &[f64; 6]) -> Result<()> { + let rv = unsafe { + call_gdal_api!( + self.api, + GDALSetGeoTransform, + self.c_dataset, + gt.as_ptr() as *mut f64 + ) + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Get the projection string. + pub fn projection(&self) -> String { + unsafe { + let ptr = call_gdal_api!(self.api, GDALGetProjectionRef, self.c_dataset); + if ptr.is_null() { + String::new() + } else { + CStr::from_ptr(ptr).to_string_lossy().into_owned() + } + } + } + + /// Set the projection string. + pub fn set_projection(&self, projection: &str) -> Result<()> { + let c_projection = CString::new(projection)?; + let rv = unsafe { + call_gdal_api!( + self.api, + GDALSetProjection, + self.c_dataset, + c_projection.as_ptr() + ) + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Get the spatial reference. + pub fn spatial_ref(&self) -> Result { + let c_srs = unsafe { call_gdal_api!(self.api, GDALGetSpatialRef, self.c_dataset) }; + if c_srs.is_null() { + return Err(GdalError::NullPointer { + method_name: "GDALGetSpatialRef", + msg: "returned null".to_string(), + }); + } + // GDALGetSpatialRef returns a borrowed reference — clone it via OSRClone. + unsafe { SpatialRef::from_c_srs_clone(self.api, c_srs) } + } + + /// Set the spatial reference. + pub fn set_spatial_ref(&self, srs: &SpatialRef) -> Result<()> { + let rv = + unsafe { call_gdal_api!(self.api, GDALSetSpatialRef, self.c_dataset, srs.c_srs()) }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Create a copy of this dataset to a new file using the given driver. + pub fn create_copy( + &self, + driver: &crate::driver::Driver, + filename: &str, + options: &[&str], + ) -> Result { + let c_filename = CString::new(filename)?; + let csl = CslStringList::try_from_iter(options.iter().copied())?; + + let c_ds = unsafe { + call_gdal_api!( + self.api, + GDALCreateCopy, + driver.c_driver(), + c_filename.as_ptr(), + self.c_dataset, + 0, // bStrict + csl.as_ptr(), + ptr::null_mut(), + ptr::null_mut() + ) + }; + if c_ds.is_null() { + return Err(self.api.last_cpl_err(CE_Failure as u32)); + } + Ok(Dataset { + api: self.api, + c_dataset: c_ds, + owned: true, + }) + } + + /// Create a new vector layer. + pub fn create_layer(&self, options: LayerOptions<'_>) -> Result> { + let c_name = CString::new(options.name)?; + let c_srs = options.srs.map_or(ptr::null_mut(), |s| s.c_srs()); + + let csl = CslStringList::try_from_iter(options.options.unwrap_or(&[]).iter().copied())?; + + let c_layer = unsafe { + call_gdal_api!( + self.api, + GDALDatasetCreateLayer, + self.c_dataset, + c_name.as_ptr(), + c_srs, + options.ty, + csl.as_ptr() + ) + }; + if c_layer.is_null() { + return Err(GdalError::NullPointer { + method_name: "GDALDatasetCreateLayer", + msg: format!("failed to create layer '{}'", options.name), + }); + } + Ok(Layer::new(self.api, c_layer, self)) + } + + /// Get the GDAL API reference. + pub fn api(&self) -> &'static GdalApi { + self.api + } + + /// Open a dataset using a `DatasetOptions` struct (georust-compatible convenience). + pub fn open_ex_with_options( + api: &'static GdalApi, + path: &str, + options: DatasetOptions<'_>, + ) -> Result { + Self::open_ex( + api, + path, + options.open_flags, + options.allowed_drivers, + options.open_options, + options.sibling_files, + ) + } + + /// Add a raster band backed by an existing memory buffer (zero-copy). + /// + /// This wraps `GDALAddBand` with the `DATAPOINTER`, `PIXELOFFSET`, and `LINEOFFSET` + /// options, allowing you to attach existing memory to a MEM dataset without copying. + /// + /// # Arguments + /// * `data_type` - The GDAL data type of the band. + /// * `data_ptr` - Pointer to the band pixel data. + /// * `pixel_offset` - Byte offset between consecutive pixels. `None` defaults to the + /// byte size of `data_type`. + /// * `line_offset` - Byte offset between consecutive lines. `None` defaults to + /// `pixel_offset * width`. + /// + /// # Safety + /// + /// The caller must ensure that `data_ptr` points to a valid buffer of at least + /// `height * line_offset` bytes (or `height * width * data_type.byte_size()` when + /// using defaults), and that the buffer outlives this dataset. + pub unsafe fn add_band_with_data( + &self, + data_type: RustGdalDataType, + data_ptr: *const u8, + pixel_offset: Option, + line_offset: Option, + ) -> Result<()> { + let data_pointer = format!("DATAPOINTER={data_ptr:p}"); + + let mut options = CslStringList::with_capacity(3); + options.add_string(&data_pointer)?; + + if let Some(pixel) = pixel_offset { + options.set_name_value("PIXELOFFSET", &pixel.to_string())?; + } + + if let Some(line) = line_offset { + options.set_name_value("LINEOFFSET", &line.to_string())?; + } + + let err = call_gdal_api!( + self.api, + GDALAddBand, + self.c_dataset, + data_type.to_c(), + options.as_ptr() + ); + if err != CE_None { + return Err(self.api.last_cpl_err(err as u32)); + } + Ok(()) + } + + /// Mark this dataset as owning its handle (for `Drop`). + pub fn set_owned(&mut self, owned: bool) { + self.owned = owned; + } +} + +/// Options for creating a vector layer. +pub struct LayerOptions<'a> { + pub name: &'a str, + pub srs: Option<&'a SpatialRef>, + pub ty: OGRwkbGeometryType, + /// Additional driver-specific options, in the form `"name=value"`. + pub options: Option<&'a [&'a str]>, +} + +impl Default for LayerOptions<'_> { + fn default() -> Self { + Self { + name: "", + srs: None, + ty: OGRwkbGeometryType::wkbUnknown, + options: None, + } + } +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use crate::driver::DriverManager; + use crate::global::with_global_gdal_api; + + #[test] + fn test_geo_transform_roundtrip() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create("", 256, 256, 1).unwrap(); + + let gt = [0.0, 1.0, 0.0, 0.0, 0.0, -1.0]; + ds.set_geo_transform(>).unwrap(); + let got = ds.geo_transform().unwrap(); + assert_eq!(gt, got); + }) + .unwrap(); + } + + #[test] + fn test_geo_transform_unset() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create("", 256, 256, 1).unwrap(); + + // MEM driver without an explicit set_geo_transform returns an error + assert!(ds.geo_transform().is_err()); + }) + .unwrap(); + } + + #[test] + fn test_set_projection_roundtrip() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create("", 256, 256, 1).unwrap(); + + let wkt = r#"GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]]"#; + ds.set_projection(wkt).unwrap(); + let got = ds.projection(); + // The returned WKT may be reformatted by GDAL, so just check it contains WGS 84 + assert!(got.contains("WGS 84"), "Expected WGS 84 in: {got}"); + }) + .unwrap(); + } + + #[test] + fn test_dataset_raster_count() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + + let ds1 = driver.create("", 64, 64, 1).unwrap(); + assert_eq!(ds1.raster_count(), 1); + + let ds3 = driver.create("", 64, 64, 3).unwrap(); + assert_eq!(ds3.raster_count(), 3); + }) + .unwrap(); + } + + #[test] + fn test_dataset_raster_size() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create("", 123, 456, 1).unwrap(); + assert_eq!(ds.raster_size(), (123, 456)); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/driver.rs b/c/sedona-gdal/src/driver.rs new file mode 100644 index 000000000..b3d82379f --- /dev/null +++ b/c/sedona-gdal/src/driver.rs @@ -0,0 +1,246 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ffi::CString; +use std::ptr; + +use crate::dataset::Dataset; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; +use crate::raster::types::GdalDataType as RustGdalDataType; +use crate::raster::GdalType; + +/// A GDAL driver. +pub struct Driver { + api: &'static GdalApi, + c_driver: GDALDriverH, +} + +impl Driver { + /// Wrap an existing C driver handle. + /// + /// # Safety + /// + /// The caller must ensure the handle is valid. + pub unsafe fn from_c_driver(api: &'static GdalApi, c_driver: GDALDriverH) -> Self { + Self { api, c_driver } + } + + /// Return the raw C driver handle. + pub fn c_driver(&self) -> GDALDriverH { + self.c_driver + } + + /// Create a new raster dataset (with u8 band type). + pub fn create( + &self, + filename: &str, + size_x: usize, + size_y: usize, + bands: usize, + ) -> Result { + self.create_with_band_type::(filename, size_x, size_y, bands) + } + + /// Create a new raster dataset with a specific band type. + pub fn create_with_band_type( + &self, + filename: &str, + size_x: usize, + size_y: usize, + bands: usize, + ) -> Result { + let c_filename = CString::new(filename)?; + let x: i32 = size_x.try_into()?; + let y: i32 = size_y.try_into()?; + let b: i32 = bands.try_into()?; + let c_ds = unsafe { + call_gdal_api!( + self.api, + GDALCreate, + self.c_driver, + c_filename.as_ptr(), + x, + y, + b, + T::gdal_ordinal(), + ptr::null_mut() + ) + }; + if c_ds.is_null() { + return Err(self.api.last_cpl_err(CE_Failure as u32)); + } + Ok(Dataset::new_owned(self.api, c_ds)) + } + + /// Create a new raster dataset with a runtime data type. + /// + /// Unlike [`create_with_band_type`](Self::create_with_band_type), this accepts a + /// [`GdalDataType`](RustGdalDataType) enum value instead of a compile-time generic, + /// which is useful when the data type is only known at runtime. + pub fn create_with_data_type( + &self, + filename: &str, + size_x: usize, + size_y: usize, + bands: usize, + data_type: RustGdalDataType, + ) -> Result { + let c_filename = CString::new(filename)?; + let x: i32 = size_x.try_into()?; + let y: i32 = size_y.try_into()?; + let b: i32 = bands.try_into()?; + let c_ds = unsafe { + call_gdal_api!( + self.api, + GDALCreate, + self.c_driver, + c_filename.as_ptr(), + x, + y, + b, + data_type.to_c(), + ptr::null_mut() + ) + }; + if c_ds.is_null() { + return Err(self.api.last_cpl_err(CE_Failure as u32)); + } + Ok(Dataset::new_owned(self.api, c_ds)) + } + + /// Create a new dataset (vector-only, no raster bands). + pub fn create_vector_only(&self, filename: &str) -> Result { + let c_filename = CString::new(filename)?; + let c_ds = unsafe { + call_gdal_api!( + self.api, + GDALCreate, + self.c_driver, + c_filename.as_ptr(), + 1, + 1, + 0, + GDALDataType::GDT_Unknown, + ptr::null_mut() + ) + }; + if c_ds.is_null() { + return Err(self.api.last_cpl_err(CE_Failure as u32)); + } + Ok(Dataset::new_owned(self.api, c_ds)) + } +} + +/// Driver manager for looking up drivers by name. +pub struct DriverManager; + +impl DriverManager { + pub fn get_driver_by_name(api: &'static GdalApi, name: &str) -> Result { + let c_name = CString::new(name)?; + let c_driver = unsafe { call_gdal_api!(api, GDALGetDriverByName, c_name.as_ptr()) }; + if c_driver.is_null() { + // `GDALGetDriverByName` just returns `null` and sets no error message + return Err(GdalError::NullPointer { + method_name: "GDALGetDriverByName", + msg: format!("driver '{name}' not found"), + }); + } + Ok(Driver { api, c_driver }) + } +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use crate::driver::DriverManager; + use crate::errors::GdalError; + use crate::global::with_global_gdal_api; + use crate::raster::types::GdalDataType; + + #[test] + fn test_get_driver_by_name() { + with_global_gdal_api(|api| { + let gtiff = DriverManager::get_driver_by_name(api, "GTiff").unwrap(); + assert!(!gtiff.c_driver().is_null()); + let mem = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + assert!(!mem.c_driver().is_null()); + }) + .unwrap(); + } + + #[test] + fn test_get_driver_by_name_invalid() { + with_global_gdal_api(|api| { + let err = DriverManager::get_driver_by_name(api, "NO_SUCH_DRIVER"); + assert!(matches!(err, Err(GdalError::NullPointer { .. }))); + }) + .unwrap(); + } + + #[test] + fn test_driver_create() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create("", 32, 16, 2).unwrap(); + assert_eq!(ds.raster_size(), (32, 16)); + assert_eq!(ds.raster_count(), 2); + }) + .unwrap(); + } + + #[test] + fn test_driver_create_with_band_type() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create_with_band_type::("", 10, 20, 1).unwrap(); + assert_eq!(ds.raster_count(), 1); + let ds = driver.create_with_band_type::("", 10, 20, 2).unwrap(); + assert_eq!(ds.raster_count(), 2); + let ds = driver.create_with_band_type::("", 10, 20, 3).unwrap(); + assert_eq!(ds.raster_count(), 3); + }) + .unwrap(); + } + + #[test] + fn test_driver_create_with_data_type() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver + .create_with_data_type("", 8, 8, 1, GdalDataType::UInt16) + .unwrap(); + assert_eq!(ds.raster_count(), 1); + }) + .unwrap(); + } + + #[test] + fn test_driver_create_vector_only() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let ds = driver.create_vector_only("").unwrap(); + assert_eq!(ds.raster_count(), 0); + assert_eq!(ds.raster_size(), (1, 1)); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/errors.rs b/c/sedona-gdal/src/errors.rs index f344667ec..3ff5d03db 100644 --- a/c/sedona-gdal/src/errors.rs +++ b/c/sedona-gdal/src/errors.rs @@ -63,6 +63,9 @@ pub enum GdalError { #[error(transparent)] IntConversionError(#[from] TryFromIntError), + + #[error("Buffer length {0} does not match raster size {1:?}")] + BufferSizeMismatch(usize, (usize, usize)), } pub type Result = std::result::Result; diff --git a/c/sedona-gdal/src/lib.rs b/c/sedona-gdal/src/lib.rs index 0646a241a..c331f5dfa 100644 --- a/c/sedona-gdal/src/lib.rs +++ b/c/sedona-gdal/src/lib.rs @@ -29,6 +29,8 @@ pub mod global; // --- High-level wrappers --- pub mod config; pub mod cpl; +pub mod dataset; +pub mod driver; pub mod geo_transform; pub mod raster; pub mod spatial_ref; diff --git a/c/sedona-gdal/src/raster.rs b/c/sedona-gdal/src/raster.rs index 1ddc9b2ed..a0ce55b5f 100644 --- a/c/sedona-gdal/src/raster.rs +++ b/c/sedona-gdal/src/raster.rs @@ -15,4 +15,7 @@ // specific language governing permissions and limitations // under the License. +pub mod rasterband; pub mod types; + +pub use rasterband::{actual_block_size, RasterBand}; diff --git a/c/sedona-gdal/src/raster/rasterband.rs b/c/sedona-gdal/src/raster/rasterband.rs new file mode 100644 index 000000000..12f17616e --- /dev/null +++ b/c/sedona-gdal/src/raster/rasterband.rs @@ -0,0 +1,368 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::marker::PhantomData; + +use crate::dataset::Dataset; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::raster::types::{Buffer, GdalType, ResampleAlg}; +use crate::{gdal_dyn_bindgen::*, raster::types::GdalDataType}; + +/// A raster band of a dataset. +pub struct RasterBand<'a> { + api: &'static GdalApi, + c_rasterband: GDALRasterBandH, + _dataset: PhantomData<&'a Dataset>, +} + +impl<'a> RasterBand<'a> { + pub(crate) fn new( + api: &'static GdalApi, + c_rasterband: GDALRasterBandH, + _dataset: &'a Dataset, + ) -> Self { + Self { + api, + c_rasterband, + _dataset: PhantomData, + } + } + + /// Return the raw C raster band handle. + pub fn c_rasterband(&self) -> GDALRasterBandH { + self.c_rasterband + } + + /// Read a region of the band as a typed buffer. + /// + /// If `e_resample_alg` is `None`, nearest-neighbour resampling is used. + pub fn read_as( + &self, + window: (isize, isize), + window_size: (usize, usize), + size: (usize, usize), + e_resample_alg: Option, + ) -> Result> { + let len = size.0 * size.1; + // Safety: all GdalType implementations are numeric primitives (u8, i8, u16, ..., f64), + // for which zeroed memory is a valid bit pattern. + let mut data: Vec = vec![unsafe { std::mem::zeroed() }; len]; + + let resample_alg = e_resample_alg.unwrap_or(ResampleAlg::NearestNeighbour); + let mut extra_arg = GDALRasterIOExtraArg { + eResampleAlg: resample_alg.to_gdal(), + ..GDALRasterIOExtraArg::default() + }; + + let rv = unsafe { + call_gdal_api!( + self.api, + GDALRasterIOEx, + self.c_rasterband, + GF_Read, + i32::try_from(window.0)?, + i32::try_from(window.1)?, + i32::try_from(window_size.0)?, + i32::try_from(window_size.1)?, + data.as_mut_ptr() as *mut std::ffi::c_void, + i32::try_from(size.0)?, + i32::try_from(size.1)?, + T::gdal_ordinal(), + 0, // nPixelSpace (auto) + 0, // nLineSpace (auto) + &mut extra_arg + ) + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + + Ok(Buffer::new(size, data)) + } + + /// Write a buffer to this raster band. + pub fn write( + &self, + window: (isize, isize), + window_size: (usize, usize), + buffer: &mut Buffer, + ) -> Result<()> { + let expected_len = buffer.shape.0 * buffer.shape.1; + if buffer.data.len() != expected_len { + return Err(GdalError::BufferSizeMismatch( + buffer.data.len(), + buffer.shape, + )); + } + let rv = unsafe { + call_gdal_api!( + self.api, + GDALRasterIO, + self.c_rasterband, + GF_Write, + i32::try_from(window.0)?, + i32::try_from(window.1)?, + i32::try_from(window_size.0)?, + i32::try_from(window_size.1)?, + buffer.data.as_mut_ptr() as *mut std::ffi::c_void, + i32::try_from(buffer.shape.0)?, + i32::try_from(buffer.shape.1)?, + T::gdal_ordinal(), + 0, // nPixelSpace (auto) + 0 // nLineSpace (auto) + ) + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Get the data type of this band. + pub fn band_type(&self) -> GdalDataType { + GdalDataType::from_c(self.c_band_type()).unwrap_or(GdalDataType::Unknown) + } + + /// Get the GDAL data type of this band. + pub fn c_band_type(&self) -> GDALDataType { + unsafe { call_gdal_api!(self.api, GDALGetRasterDataType, self.c_rasterband) } + } + + /// Get band size as (x_size, y_size). + pub fn size(&self) -> (usize, usize) { + let x = unsafe { call_gdal_api!(self.api, GDALGetRasterBandXSize, self.c_rasterband) }; + let y = unsafe { call_gdal_api!(self.api, GDALGetRasterBandYSize, self.c_rasterband) }; + (x as usize, y as usize) + } + + /// Get the block size as (x_size, y_size). + pub fn block_size(&self) -> (usize, usize) { + let mut x: i32 = 0; + let mut y: i32 = 0; + unsafe { + call_gdal_api!( + self.api, + GDALGetBlockSize, + self.c_rasterband, + &mut x, + &mut y + ) + }; + (x as usize, y as usize) + } + + /// Get the no-data value. Returns `Some(value)` if set, `None` otherwise. + pub fn no_data_value(&self) -> Option { + let mut success: i32 = 0; + let value = unsafe { + call_gdal_api!( + self.api, + GDALGetRasterNoDataValue, + self.c_rasterband, + &mut success + ) + }; + if success != 0 { + Some(value) + } else { + None + } + } + + /// Set or clear the no-data value. + pub fn set_no_data_value(&self, value: Option) -> Result<()> { + let rv = if let Some(val) = value { + unsafe { call_gdal_api!(self.api, GDALSetRasterNoDataValue, self.c_rasterband, val) } + } else { + unsafe { call_gdal_api!(self.api, GDALDeleteRasterNoDataValue, self.c_rasterband) } + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Set or clear the no-data value as u64. + pub fn set_no_data_value_u64(&self, value: Option) -> Result<()> { + let rv = if let Some(val) = value { + unsafe { + call_gdal_api!( + self.api, + GDALSetRasterNoDataValueAsUInt64, + self.c_rasterband, + val + ) + } + } else { + unsafe { call_gdal_api!(self.api, GDALDeleteRasterNoDataValue, self.c_rasterband) } + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Set or clear the no-data value as i64. + pub fn set_no_data_value_i64(&self, value: Option) -> Result<()> { + let rv = if let Some(val) = value { + unsafe { + call_gdal_api!( + self.api, + GDALSetRasterNoDataValueAsInt64, + self.c_rasterband, + val + ) + } + } else { + unsafe { call_gdal_api!(self.api, GDALDeleteRasterNoDataValue, self.c_rasterband) } + }; + if rv != CE_None { + return Err(self.api.last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Get the GDAL API reference. + pub fn api(&self) -> &'static GdalApi { + self.api + } +} + +/// Compute the actual block size (clamped to raster extent) for a given block index. +pub fn actual_block_size( + band: &RasterBand<'_>, + block_index: (usize, usize), +) -> Result<(usize, usize)> { + let (block_x, block_y) = band.block_size(); + let (raster_x, raster_y) = band.size(); + let x_off = block_index.0 * block_x; + let y_off = block_index.1 * block_y; + if x_off >= raster_x || y_off >= raster_y { + return Err(GdalError::BadArgument(format!( + "block index ({}, {}) is out of bounds for raster size ({}, {})", + block_index.0, block_index.1, raster_x, raster_y + ))); + } + let actual_x = if x_off + block_x > raster_x { + raster_x - x_off + } else { + block_x + }; + let actual_y = if y_off + block_y > raster_y { + raster_y - y_off + } else { + block_y + }; + Ok((actual_x, actual_y)) +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use crate::dataset::Dataset; + use crate::driver::DriverManager; + use crate::gdal_dyn_bindgen::*; + use crate::global::with_global_gdal_api; + use crate::raster::types::ResampleAlg; + + fn fixture(name: &str) -> String { + sedona_testing::data::test_raster(name).unwrap() + } + + #[test] + fn test_read_raster() { + with_global_gdal_api(|api| { + let path = fixture("tinymarble.tif"); + let dataset = Dataset::open_ex(api, &path, GDAL_OF_READONLY, None, None, None).unwrap(); + let rb = dataset.rasterband(1).unwrap(); + let rv = rb.read_as::((20, 30), (2, 3), (2, 3), None).unwrap(); + assert_eq!(rv.shape, (2, 3)); + assert_eq!(rv.data(), [7, 7, 7, 10, 8, 12]); + }) + .unwrap(); + } + + #[test] + fn test_read_raster_with_default_resample() { + with_global_gdal_api(|api| { + let path = fixture("tinymarble.tif"); + let dataset = Dataset::open_ex(api, &path, GDAL_OF_READONLY, None, None, None).unwrap(); + let rb = dataset.rasterband(1).unwrap(); + let rv = rb.read_as::((20, 30), (4, 4), (2, 2), None).unwrap(); + assert_eq!(rv.shape, (2, 2)); + // Default is NearestNeighbour; exact values are GDAL-version-dependent + // when downsampling from 4x4 to 2x2. Just verify shape and non-emptiness. + assert_eq!(rv.data().len(), 4); + }) + .unwrap(); + } + + #[test] + fn test_read_raster_with_average_resample() { + with_global_gdal_api(|api| { + let path = fixture("tinymarble.tif"); + let dataset = Dataset::open_ex(api, &path, GDAL_OF_READONLY, None, None, None).unwrap(); + let rb = dataset.rasterband(1).unwrap(); + let rv = rb + .read_as::((20, 30), (4, 4), (2, 2), Some(ResampleAlg::Average)) + .unwrap(); + assert_eq!(rv.shape, (2, 2)); + // Average resampling; exact values are GDAL-version-dependent. + // Verify shape and that results differ from the non-resampled full read. + assert_eq!(rv.data().len(), 4); + }) + .unwrap(); + } + + #[test] + fn test_get_no_data_value() { + with_global_gdal_api(|api| { + // tinymarble.tif has no nodata + let path = fixture("tinymarble.tif"); + let dataset = Dataset::open_ex(api, &path, GDAL_OF_READONLY, None, None, None).unwrap(); + let rb = dataset.rasterband(1).unwrap(); + assert!(rb.no_data_value().is_none()); + + // labels.tif has nodata=255 + let path = fixture("labels.tif"); + let dataset = Dataset::open_ex(api, &path, GDAL_OF_READONLY, None, None, None).unwrap(); + let rb = dataset.rasterband(1).unwrap(); + assert_eq!(rb.no_data_value(), Some(255.0)); + }) + .unwrap(); + } + + #[test] + #[allow(clippy::float_cmp)] + fn test_set_no_data_value() { + with_global_gdal_api(|api| { + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let dataset = driver.create("", 20, 10, 1).unwrap(); + let rasterband = dataset.rasterband(1).unwrap(); + assert_eq!(rasterband.no_data_value(), None); + assert!(rasterband.set_no_data_value(Some(1.23)).is_ok()); + assert_eq!(rasterband.no_data_value(), Some(1.23)); + assert!(rasterband.set_no_data_value(None).is_ok()); + assert_eq!(rasterband.no_data_value(), None); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/vector.rs b/c/sedona-gdal/src/vector.rs index 10c038edc..52ff4e1bf 100644 --- a/c/sedona-gdal/src/vector.rs +++ b/c/sedona-gdal/src/vector.rs @@ -15,4 +15,6 @@ // specific language governing permissions and limitations // under the License. +pub mod feature; pub mod geometry; +pub mod layer; diff --git a/c/sedona-gdal/src/vector/feature.rs b/c/sedona-gdal/src/vector/feature.rs new file mode 100644 index 000000000..de3013fca --- /dev/null +++ b/c/sedona-gdal/src/vector/feature.rs @@ -0,0 +1,218 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ffi::CString; +use std::marker::PhantomData; + +use crate::errors::Result; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; + +/// An OGR feature. +pub struct Feature<'a> { + api: &'static GdalApi, + c_feature: OGRFeatureH, + _lifetime: PhantomData<&'a ()>, +} + +impl Drop for Feature<'_> { + fn drop(&mut self) { + if !self.c_feature.is_null() { + unsafe { call_gdal_api!(self.api, OGR_F_Destroy, self.c_feature) }; + } + } +} + +impl<'a> Feature<'a> { + pub(crate) fn new(api: &'static GdalApi, c_feature: OGRFeatureH) -> Self { + Self { + api, + c_feature, + _lifetime: PhantomData, + } + } + + /// Get the geometry reference (borrowed, not owned — do not destroy). + /// + /// Returns None if the feature has no geometry. + pub fn geometry(&self) -> Option> { + let c_geom = unsafe { call_gdal_api!(self.api, OGR_F_GetGeometryRef, self.c_feature) }; + if c_geom.is_null() { + None + } else { + Some(BorrowedGeometry { + api: self.api, + c_geom, + _lifetime: PhantomData, + }) + } + } + + /// Get a field's index by name. Returns an error if the field is not found. + pub fn field_index(&self, name: &str) -> Result { + let c_name = CString::new(name)?; + let idx = unsafe { + call_gdal_api!( + self.api, + OGR_F_GetFieldIndex, + self.c_feature, + c_name.as_ptr() + ) + }; + if idx < 0 { + return Err(crate::errors::GdalError::BadArgument(format!( + "field '{name}' not found" + ))); + } + Ok(idx) + } + + /// Get a field value as f64. + pub fn field_as_double(&self, field_index: i32) -> f64 { + unsafe { + call_gdal_api!( + self.api, + OGR_F_GetFieldAsDouble, + self.c_feature, + field_index + ) + } + } + + /// Get a field value as i32. + /// + /// Returns `Some(value)` if the field is set and not null, `None` otherwise. + pub fn field_as_integer(&self, field_index: i32) -> Option { + let is_set = unsafe { + call_gdal_api!( + self.api, + OGR_F_IsFieldSetAndNotNull, + self.c_feature, + field_index + ) + }; + if is_set != 0 { + Some(unsafe { + call_gdal_api!( + self.api, + OGR_F_GetFieldAsInteger, + self.c_feature, + field_index + ) + }) + } else { + None + } + } +} + +/// A geometry borrowed from a feature (not owned — will NOT be destroyed). +pub struct BorrowedGeometry<'a> { + api: &'static GdalApi, + c_geom: OGRGeometryH, + _lifetime: PhantomData<&'a ()>, +} + +impl<'a> BorrowedGeometry<'a> { + /// Return the raw C geometry handle. + pub fn c_geometry(&self) -> OGRGeometryH { + self.c_geom + } + + /// Export to ISO WKB. + pub fn wkb(&self) -> Result> { + let size = unsafe { call_gdal_api!(self.api, OGR_G_WkbSize, self.c_geom) }; + if size < 0 { + return Err(crate::errors::GdalError::BadArgument(format!( + "OGR_G_WkbSize returned negative size: {size}" + ))); + } + let mut buf = vec![0u8; size as usize]; + let rv = unsafe { + call_gdal_api!( + self.api, + OGR_G_ExportToIsoWkb, + self.c_geom, + wkbNDR, + buf.as_mut_ptr() + ) + }; + if rv != OGRERR_NONE { + return Err(crate::errors::GdalError::OgrError { + err: rv, + method_name: "OGR_G_ExportToIsoWkb", + }); + } + Ok(buf) + } + + /// Get the bounding envelope. + pub fn envelope(&self) -> crate::vector::Envelope { + let mut env = OGREnvelope { + MinX: 0.0, + MaxX: 0.0, + MinY: 0.0, + MaxY: 0.0, + }; + unsafe { call_gdal_api!(self.api, OGR_G_GetEnvelope, self.c_geom, &mut env) }; + crate::vector::Envelope { + MinX: env.MinX, + MaxX: env.MaxX, + MinY: env.MinY, + MaxY: env.MaxY, + } + } +} + +/// An OGR field definition. +pub struct FieldDefn { + api: &'static GdalApi, + c_field_defn: OGRFieldDefnH, +} + +impl Drop for FieldDefn { + fn drop(&mut self) { + if !self.c_field_defn.is_null() { + unsafe { call_gdal_api!(self.api, OGR_Fld_Destroy, self.c_field_defn) }; + } + } +} + +impl FieldDefn { + /// Create a new field definition. + pub fn new(api: &'static GdalApi, name: &str, field_type: OGRFieldType) -> Result { + let c_name = CString::new(name)?; + let c_field_defn = + unsafe { call_gdal_api!(api, OGR_Fld_Create, c_name.as_ptr(), field_type) }; + if c_field_defn.is_null() { + return Err(crate::errors::GdalError::NullPointer { + method_name: "OGR_Fld_Create", + msg: format!("failed to create field definition '{name}'"), + }); + } + Ok(Self { api, c_field_defn }) + } + + /// Return the raw C handle. + pub fn c_field_defn(&self) -> OGRFieldDefnH { + self.c_field_defn + } +} diff --git a/c/sedona-gdal/src/vector/layer.rs b/c/sedona-gdal/src/vector/layer.rs new file mode 100644 index 000000000..0e23a4246 --- /dev/null +++ b/c/sedona-gdal/src/vector/layer.rs @@ -0,0 +1,118 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::marker::PhantomData; + +use crate::dataset::Dataset; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; +use crate::vector::feature::{Feature, FieldDefn}; + +/// An OGR layer (borrowed from a Dataset). +pub struct Layer<'a> { + api: &'static GdalApi, + c_layer: OGRLayerH, + _dataset: PhantomData<&'a Dataset>, +} + +impl<'a> Layer<'a> { + pub(crate) fn new(api: &'static GdalApi, c_layer: OGRLayerH, _dataset: &'a Dataset) -> Self { + Self { + api, + c_layer, + _dataset: PhantomData, + } + } + + /// Return the raw C layer handle. + pub fn c_layer(&self) -> OGRLayerH { + self.c_layer + } + + /// Reset reading to the first feature. + pub fn reset_reading(&self) { + unsafe { call_gdal_api!(self.api, OGR_L_ResetReading, self.c_layer) }; + } + + /// Get the next feature (returns None when exhausted). + pub fn next_feature(&self) -> Option> { + let c_feature = unsafe { call_gdal_api!(self.api, OGR_L_GetNextFeature, self.c_layer) }; + if c_feature.is_null() { + None + } else { + Some(Feature::new(self.api, c_feature)) + } + } + + /// Create a field on this layer. + pub fn create_field(&self, field_defn: &FieldDefn) -> Result<()> { + let rv = unsafe { + call_gdal_api!( + self.api, + OGR_L_CreateField, + self.c_layer, + field_defn.c_field_defn(), + 1 // bApproxOK + ) + }; + if rv != OGRERR_NONE { + return Err(GdalError::OgrError { + err: rv, + method_name: "OGR_L_CreateField", + }); + } + Ok(()) + } + + /// Get the number of features in this layer. + /// + /// If `force` is true, the count will be computed even if it is expensive. + pub fn feature_count(&self, force: bool) -> i64 { + unsafe { + call_gdal_api!( + self.api, + OGR_L_GetFeatureCount, + self.c_layer, + if force { 1 } else { 0 } + ) + } + } + + /// Iterate over all features. + pub fn features(&self) -> FeatureIterator<'_> { + self.reset_reading(); + FeatureIterator { layer: self } + } +} + +/// Iterator over features in a layer. +pub struct FeatureIterator<'a> { + layer: &'a Layer<'a>, +} + +impl<'a> Iterator for FeatureIterator<'a> { + type Item = Feature<'a>; + + fn next(&mut self) -> Option { + self.layer.next_feature() + } +} From d16130f7a10e7310d1ad01c6b78a8f60f7ddb8d0 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:43:28 +0800 Subject: [PATCH 05/11] refactor(sedona-gdal): prefer imported core types --- c/sedona-gdal/src/dataset.rs | 3 ++- c/sedona-gdal/src/vector/feature.rs | 17 ++++++++--------- 2 files changed, 10 insertions(+), 10 deletions(-) diff --git a/c/sedona-gdal/src/dataset.rs b/c/sedona-gdal/src/dataset.rs index b286535b0..111d8f8ef 100644 --- a/c/sedona-gdal/src/dataset.rs +++ b/c/sedona-gdal/src/dataset.rs @@ -23,6 +23,7 @@ use std::ffi::{CStr, CString}; use std::ptr; use crate::cpl::CslStringList; +use crate::driver::Driver; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; @@ -243,7 +244,7 @@ impl Dataset { /// Create a copy of this dataset to a new file using the given driver. pub fn create_copy( &self, - driver: &crate::driver::Driver, + driver: &Driver, filename: &str, options: &[&str], ) -> Result { diff --git a/c/sedona-gdal/src/vector/feature.rs b/c/sedona-gdal/src/vector/feature.rs index de3013fca..a5bee5588 100644 --- a/c/sedona-gdal/src/vector/feature.rs +++ b/c/sedona-gdal/src/vector/feature.rs @@ -22,9 +22,10 @@ use std::ffi::CString; use std::marker::PhantomData; -use crate::errors::Result; +use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; +use crate::vector::Envelope; /// An OGR feature. pub struct Feature<'a> { @@ -78,9 +79,7 @@ impl<'a> Feature<'a> { ) }; if idx < 0 { - return Err(crate::errors::GdalError::BadArgument(format!( - "field '{name}' not found" - ))); + return Err(GdalError::BadArgument(format!("field '{name}' not found"))); } Ok(idx) } @@ -141,7 +140,7 @@ impl<'a> BorrowedGeometry<'a> { pub fn wkb(&self) -> Result> { let size = unsafe { call_gdal_api!(self.api, OGR_G_WkbSize, self.c_geom) }; if size < 0 { - return Err(crate::errors::GdalError::BadArgument(format!( + return Err(GdalError::BadArgument(format!( "OGR_G_WkbSize returned negative size: {size}" ))); } @@ -156,7 +155,7 @@ impl<'a> BorrowedGeometry<'a> { ) }; if rv != OGRERR_NONE { - return Err(crate::errors::GdalError::OgrError { + return Err(GdalError::OgrError { err: rv, method_name: "OGR_G_ExportToIsoWkb", }); @@ -165,7 +164,7 @@ impl<'a> BorrowedGeometry<'a> { } /// Get the bounding envelope. - pub fn envelope(&self) -> crate::vector::Envelope { + pub fn envelope(&self) -> Envelope { let mut env = OGREnvelope { MinX: 0.0, MaxX: 0.0, @@ -173,7 +172,7 @@ impl<'a> BorrowedGeometry<'a> { MaxY: 0.0, }; unsafe { call_gdal_api!(self.api, OGR_G_GetEnvelope, self.c_geom, &mut env) }; - crate::vector::Envelope { + Envelope { MinX: env.MinX, MaxX: env.MaxX, MinY: env.MinY, @@ -203,7 +202,7 @@ impl FieldDefn { let c_field_defn = unsafe { call_gdal_api!(api, OGR_Fld_Create, c_name.as_ptr(), field_type) }; if c_field_defn.is_null() { - return Err(crate::errors::GdalError::NullPointer { + return Err(GdalError::NullPointer { method_name: "OGR_Fld_Create", msg: format!("failed to create field definition '{name}'"), }); From 0e554999b525642f456d0d2753c0bdfd0029a3ce Mon Sep 17 00:00:00 2001 From: Kontinuation Date: Wed, 11 Mar 2026 15:44:16 +0800 Subject: [PATCH 06/11] fix(sedona-gdal): import raster types via module path --- c/sedona-gdal/src/driver.rs | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/c/sedona-gdal/src/driver.rs b/c/sedona-gdal/src/driver.rs index b3d82379f..f0e6c5963 100644 --- a/c/sedona-gdal/src/driver.rs +++ b/c/sedona-gdal/src/driver.rs @@ -27,7 +27,7 @@ use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; use crate::raster::types::GdalDataType as RustGdalDataType; -use crate::raster::GdalType; +use crate::raster::types::GdalType; /// A GDAL driver. pub struct Driver { From efb5606f531fa1ff00d4e08db1370b4a7a90ad75 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Fri, 13 Mar 2026 14:32:55 +0800 Subject: [PATCH 07/11] refactor(sedona-gdal): remove wrapper re-export aliases --- c/sedona-gdal/src/dataset.rs | 4 ++-- c/sedona-gdal/src/raster.rs | 2 -- c/sedona-gdal/src/vector/feature.rs | 2 +- 3 files changed, 3 insertions(+), 5 deletions(-) diff --git a/c/sedona-gdal/src/dataset.rs b/c/sedona-gdal/src/dataset.rs index 111d8f8ef..e7cc6384e 100644 --- a/c/sedona-gdal/src/dataset.rs +++ b/c/sedona-gdal/src/dataset.rs @@ -27,10 +27,10 @@ use crate::driver::Driver; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; +use crate::raster::rasterband::RasterBand; use crate::raster::types::{DatasetOptions, GdalDataType as RustGdalDataType}; -use crate::raster::RasterBand; use crate::spatial_ref::SpatialRef; -use crate::vector::Layer; +use crate::vector::layer::Layer; /// A GDAL dataset. pub struct Dataset { diff --git a/c/sedona-gdal/src/raster.rs b/c/sedona-gdal/src/raster.rs index a0ce55b5f..389d9d73b 100644 --- a/c/sedona-gdal/src/raster.rs +++ b/c/sedona-gdal/src/raster.rs @@ -17,5 +17,3 @@ pub mod rasterband; pub mod types; - -pub use rasterband::{actual_block_size, RasterBand}; diff --git a/c/sedona-gdal/src/vector/feature.rs b/c/sedona-gdal/src/vector/feature.rs index a5bee5588..b29813799 100644 --- a/c/sedona-gdal/src/vector/feature.rs +++ b/c/sedona-gdal/src/vector/feature.rs @@ -25,7 +25,7 @@ use std::marker::PhantomData; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; -use crate::vector::Envelope; +use crate::vector::geometry::Envelope; /// An OGR feature. pub struct Feature<'a> { From d742c27c111b11645e76758382a47b636489e0f4 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:28:50 +0800 Subject: [PATCH 08/11] feat(sedona-gdal): add raster operations and vrt support --- c/sedona-gdal/src/lib.rs | 1 + c/sedona-gdal/src/raster.rs | 3 + c/sedona-gdal/src/raster/polygonize.rs | 274 ++++++++++++++ c/sedona-gdal/src/raster/rasterize.rs | 261 +++++++++++++ c/sedona-gdal/src/raster/rasterize_affine.rs | 368 +++++++++++++++++++ c/sedona-gdal/src/vrt.rs | 336 +++++++++++++++++ 6 files changed, 1243 insertions(+) create mode 100644 c/sedona-gdal/src/raster/polygonize.rs create mode 100644 c/sedona-gdal/src/raster/rasterize.rs create mode 100644 c/sedona-gdal/src/raster/rasterize_affine.rs create mode 100644 c/sedona-gdal/src/vrt.rs diff --git a/c/sedona-gdal/src/lib.rs b/c/sedona-gdal/src/lib.rs index c331f5dfa..f61384515 100644 --- a/c/sedona-gdal/src/lib.rs +++ b/c/sedona-gdal/src/lib.rs @@ -35,4 +35,5 @@ pub mod geo_transform; pub mod raster; pub mod spatial_ref; pub mod vector; +pub mod vrt; pub mod vsi; diff --git a/c/sedona-gdal/src/raster.rs b/c/sedona-gdal/src/raster.rs index 389d9d73b..839674298 100644 --- a/c/sedona-gdal/src/raster.rs +++ b/c/sedona-gdal/src/raster.rs @@ -15,5 +15,8 @@ // specific language governing permissions and limitations // under the License. +pub mod polygonize; pub mod rasterband; +pub mod rasterize; +pub mod rasterize_affine; pub mod types; diff --git a/c/sedona-gdal/src/raster/polygonize.rs b/c/sedona-gdal/src/raster/polygonize.rs new file mode 100644 index 000000000..58cba2736 --- /dev/null +++ b/c/sedona-gdal/src/raster/polygonize.rs @@ -0,0 +1,274 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +use std::ptr; + +use crate::cpl::CslStringList; +use crate::errors::Result; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; +use crate::raster::RasterBand; +use crate::vector::Layer; + +#[derive(Clone, Debug, Default)] +pub struct PolygonizeOptions { + /// Use 8 connectedness (diagonal pixels are considered connected). + /// + /// If `false` (default), 4 connectedness is used. + pub eight_connected: bool, + + /// Name of a dataset from which to read the geotransform. + /// + /// This is useful if the source band has no related dataset, which is typical for mask bands. + /// + /// Corresponds to GDAL's `DATASET_FOR_GEOREF=dataset_name` option. + pub dataset_for_georef: Option, + + /// Interval in number of features at which transactions must be flushed. + /// + /// - `0` means that no transactions are opened. + /// - a negative value means a single transaction. + /// + /// Corresponds to GDAL's `COMMIT_INTERVAL=num` option. + pub commit_interval: Option, +} + +impl PolygonizeOptions { + /// Build a [`CslStringList`] from this options struct. + pub fn to_options_list(&self) -> Result { + let mut options = CslStringList::new(); + + if self.eight_connected { + options.set_name_value("8CONNECTED", "8")?; + } + + if let Some(ref ds) = self.dataset_for_georef { + options.set_name_value("DATASET_FOR_GEOREF", ds)?; + } + + if let Some(interval) = self.commit_interval { + options.set_name_value("COMMIT_INTERVAL", &interval.to_string())?; + } + + Ok(options) + } +} + +/// Polygonize a raster band into a vector layer. +/// +/// Uses `GDALPolygonize` (integer pixel values). +pub fn polygonize( + api: &'static GdalApi, + src_band: &RasterBand<'_>, + mask_band: Option<&RasterBand<'_>>, + out_layer: &Layer<'_>, + pixel_value_field: i32, + options: &PolygonizeOptions, +) -> Result<()> { + let mask = mask_band.map_or(ptr::null_mut(), |b| b.c_rasterband()); + let csl = options.to_options_list()?; + + let rv = unsafe { + call_gdal_api!( + api, + GDALPolygonize, + src_band.c_rasterband(), + mask, + out_layer.c_layer(), + pixel_value_field, + csl.as_ptr(), + ptr::null_mut(), // pfnProgress + ptr::null_mut() // pProgressData + ) + }; + if rv != CE_None { + return Err(api.last_cpl_err(rv as u32)); + } + Ok(()) +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn test_polygonizeoptions_as_ptr() { + let c_options = PolygonizeOptions::default().to_options_list().unwrap(); + assert_eq!(c_options.fetch_name_value("8CONNECTED"), None); + assert_eq!(c_options.fetch_name_value("DATASET_FOR_GEOREF"), None); + assert_eq!(c_options.fetch_name_value("COMMIT_INTERVAL"), None); + + let c_options = PolygonizeOptions { + eight_connected: true, + dataset_for_georef: Some("/vsimem/georef.tif".to_string()), + commit_interval: Some(12345), + } + .to_options_list() + .unwrap(); + assert_eq!(c_options.fetch_name_value("8CONNECTED"), Some("8".into())); + assert_eq!( + c_options.fetch_name_value("DATASET_FOR_GEOREF"), + Some("/vsimem/georef.tif".into()) + ); + assert_eq!( + c_options.fetch_name_value("COMMIT_INTERVAL"), + Some("12345".into()) + ); + } + + #[cfg(feature = "gdal-sys")] + #[test] + fn test_polygonize_connectivity_affects_regions() { + use crate::dataset::LayerOptions; + use crate::driver::DriverManager; + use crate::global::with_global_gdal_api; + use crate::raster::Buffer; + use crate::vector::FieldDefn; + + with_global_gdal_api(|api| { + let mem_driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let raster_ds = mem_driver.create("", 3, 3, 1).unwrap(); + let band = raster_ds.rasterband(1).unwrap(); + + // 3x3 raster with diagonal 1s: + // 1 0 0 + // 0 1 0 + // 0 0 1 + let mut data = Buffer::new((3, 3), vec![1u8, 0, 0, 0, 1, 0, 0, 0, 1]); + band.write((0, 0), (3, 3), &mut data).unwrap(); + + let gpkg_path = "/vsimem/test_polygonize_connectivity.gpkg"; + let gpkg_driver = DriverManager::get_driver_by_name(api, "GPKG").unwrap(); + let vector_ds = gpkg_driver.create_vector_only(gpkg_path).unwrap(); + + // 4-connected output + let layer_4 = vector_ds + .create_layer(LayerOptions { + name: "four", + srs: None, + ty: OGRwkbGeometryType::wkbPolygon, + options: None, + }) + .unwrap(); + let field_defn = FieldDefn::new(api, "val", OGRFieldType::OFTInteger).unwrap(); + layer_4.create_field(&field_defn).unwrap(); + + polygonize(api, &band, None, &layer_4, 0, &PolygonizeOptions::default()).unwrap(); + + let ones_4 = layer_4 + .features() + .filter_map(|f| f.field_as_integer(0)) + .filter(|v| *v == 1) + .count(); + assert_eq!(ones_4, 3); + + // 8-connected output + let layer_8 = vector_ds + .create_layer(LayerOptions { + name: "eight", + srs: None, + ty: OGRwkbGeometryType::wkbPolygon, + options: None, + }) + .unwrap(); + let field_defn = FieldDefn::new(api, "val", OGRFieldType::OFTInteger).unwrap(); + layer_8.create_field(&field_defn).unwrap(); + + polygonize( + api, + &band, + None, + &layer_8, + 0, + &PolygonizeOptions { + eight_connected: true, + dataset_for_georef: None, + commit_interval: None, + }, + ) + .unwrap(); + + let ones_8 = layer_8 + .features() + .filter_map(|f| f.field_as_integer(0)) + .filter(|v| *v == 1) + .count(); + assert_eq!(ones_8, 1); + + crate::vsi::unlink_mem_file(api, gpkg_path).unwrap(); + }) + .unwrap(); + } + + #[cfg(feature = "gdal-sys")] + #[test] + fn test_polygonize_with_mask_band_restricts_output() { + use crate::dataset::LayerOptions; + use crate::driver::DriverManager; + use crate::global::with_global_gdal_api; + use crate::raster::Buffer; + use crate::vector::FieldDefn; + + with_global_gdal_api(|api| { + let mem_driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let raster_ds = mem_driver.create("", 3, 3, 2).unwrap(); + + let value_band = raster_ds.rasterband(1).unwrap(); + let mask_band = raster_ds.rasterband(2).unwrap(); + + // Value band: all 7s. + let mut values = Buffer::new((3, 3), vec![7u8; 9]); + value_band.write((0, 0), (3, 3), &mut values).unwrap(); + + // Mask: only the center pixel is included. + let mut mask = Buffer::new((3, 3), vec![0u8, 0, 0, 0, 1, 0, 0, 0, 0]); + mask_band.write((0, 0), (3, 3), &mut mask).unwrap(); + + let gpkg_path = "/vsimem/test_polygonize_mask.gpkg"; + let gpkg_driver = DriverManager::get_driver_by_name(api, "GPKG").unwrap(); + let vector_ds = gpkg_driver.create_vector_only(gpkg_path).unwrap(); + + let layer = vector_ds + .create_layer(LayerOptions { + name: "masked", + srs: None, + ty: OGRwkbGeometryType::wkbPolygon, + options: None, + }) + .unwrap(); + let field_defn = FieldDefn::new(api, "val", OGRFieldType::OFTInteger).unwrap(); + layer.create_field(&field_defn).unwrap(); + + polygonize( + api, + &value_band, + Some(&mask_band), + &layer, + 0, + &PolygonizeOptions::default(), + ) + .unwrap(); + + assert_eq!(layer.feature_count(true), 1); + let only_val = layer.features().next().unwrap().field_as_integer(0); + assert_eq!(only_val, Some(7)); + + crate::vsi::unlink_mem_file(api, gpkg_path).unwrap(); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/raster/rasterize.rs b/c/sedona-gdal/src/raster/rasterize.rs new file mode 100644 index 000000000..dbcad37ca --- /dev/null +++ b/c/sedona-gdal/src/raster/rasterize.rs @@ -0,0 +1,261 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Ported (and contains copied code) from georust/gdal: +//! . +//! Original code is licensed under MIT. + +use std::ptr; + +use crate::cpl::CslStringList; +use crate::dataset::Dataset; +use crate::errors::Result; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::*; +use crate::vector::Geometry; + +/// Source of burn values. +#[derive(Copy, Clone, Debug)] +pub enum BurnSource { + /// Use whatever `burn_values` argument is supplied to `rasterize`. + UserSupplied, + + /// Add the geometry's Z value to whatever `burn_values` argument + /// is supplied to `rasterize`. + Z, +} + +/// Algorithm for merging new raster values with existing values. +#[derive(Copy, Clone, Debug)] +pub enum MergeAlgorithm { + /// Overwrite existing value (default). + Replace, + /// Add new value to existing value (useful for heatmaps). + Add, +} + +/// Optimization mode for rasterization. +#[derive(Copy, Clone, Debug)] +pub enum OptimizeMode { + /// Let GDAL decide (default). + Automatic, + /// Force raster-based scan (iterates over pixels). + Raster, + /// Force vector-based scan (iterates over geometry edges). + Vector, +} + +/// Options that specify how to rasterize geometries. +#[derive(Copy, Clone, Debug)] +pub struct RasterizeOptions { + /// Set to `true` to set all pixels touched by the line or polygons, + /// not just those whose center is within the polygon or that are + /// selected by Bresenham's line algorithm. Defaults to `false`. + pub all_touched: bool, + + /// May be set to `BurnSource::Z` to use the Z values of the geometries. + /// `burn_value` is added to this before burning. Defaults to + /// `BurnSource::UserSupplied` in which case just the `burn_value` is burned. + pub source: BurnSource, + + /// May be `MergeAlgorithm::Replace` (default) or `MergeAlgorithm::Add`. + /// `Replace` overwrites existing values; `Add` adds to them. + pub merge_algorithm: MergeAlgorithm, + + /// The height in lines of the chunk to operate on. `0` (default) lets GDAL + /// choose based on cache size. Not used in `OPTIM=RASTER` mode. + pub chunk_y_size: usize, + + /// Optimization mode for rasterization. + pub optimize: OptimizeMode, +} + +impl Default for RasterizeOptions { + fn default() -> Self { + RasterizeOptions { + all_touched: false, + source: BurnSource::UserSupplied, + merge_algorithm: MergeAlgorithm::Replace, + chunk_y_size: 0, + optimize: OptimizeMode::Automatic, + } + } +} + +impl RasterizeOptions { + /// Build a [`CslStringList`] from this options struct. + pub fn to_options_list(self) -> Result { + let mut options = CslStringList::with_capacity(5); + + options.set_name_value( + "ALL_TOUCHED", + if self.all_touched { "TRUE" } else { "FALSE" }, + )?; + + options.set_name_value( + "MERGE_ALG", + match self.merge_algorithm { + MergeAlgorithm::Replace => "REPLACE", + MergeAlgorithm::Add => "ADD", + }, + )?; + + options.set_name_value("CHUNKYSIZE", &self.chunk_y_size.to_string())?; + + options.set_name_value( + "OPTIM", + match self.optimize { + OptimizeMode::Automatic => "AUTO", + OptimizeMode::Raster => "RASTER", + OptimizeMode::Vector => "VECTOR", + }, + )?; + + if let BurnSource::Z = self.source { + options.set_name_value("BURN_VALUE_FROM", "Z")?; + } + + Ok(options) + } +} + +/// Rasterize geometries onto a dataset. +/// +/// There must be one burn value for every geometry. Each burn value is +/// replicated across all bands internally, matching the +/// `GDALRasterizeGeometries` contract of `nGeomCount * nBandCount` entries. +pub fn rasterize( + api: &'static GdalApi, + dataset: &Dataset, + band_list: &[i32], + geometries: &[&Geometry], + burn_values: &[f64], + options: Option, +) -> Result<()> { + if band_list.is_empty() { + return Err(crate::errors::GdalError::BadArgument( + "`band_list` must not be empty".to_string(), + )); + } + if burn_values.len() != geometries.len() { + return Err(crate::errors::GdalError::BadArgument(format!( + "burn_values length ({}) must match geometries length ({})", + burn_values.len(), + geometries.len() + ))); + } + let raster_count = dataset.raster_count(); + for &band in band_list { + let is_good = band > 0 && (band as usize) <= raster_count; + if !is_good { + return Err(crate::errors::GdalError::BadArgument(format!( + "Band index {} is out of bounds", + band + ))); + } + } + + let geom_handles: Vec = geometries.iter().map(|g| g.c_geometry()).collect(); + + // Replicate each burn value across all bands, matching the GDAL C API + // contract that expects nGeomCount * nBandCount burn values. + let expanded_burn_values: Vec = burn_values + .iter() + .flat_map(|burn| std::iter::repeat_n(burn, band_list.len())) + .copied() + .collect(); + + let opts = options.unwrap_or_default(); + let csl = opts.to_options_list()?; + + let n_band_count: i32 = band_list.len().try_into()?; + let n_geom_count: i32 = geom_handles.len().try_into()?; + + let rv = unsafe { + call_gdal_api!( + api, + GDALRasterizeGeometries, + dataset.c_dataset(), + n_band_count, + band_list.as_ptr(), + n_geom_count, + geom_handles.as_ptr(), + ptr::null_mut(), // pfnTransformer + ptr::null_mut(), // pTransformArg + expanded_burn_values.as_ptr(), + csl.as_ptr(), + ptr::null_mut(), // pfnProgress + ptr::null_mut() // pProgressData + ) + }; + if rv != CE_None { + return Err(api.last_cpl_err(rv as u32)); + } + Ok(()) +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn test_rasterizeoptions_as_ptr() { + let c_options = RasterizeOptions::default().to_options_list().unwrap(); + assert_eq!( + c_options.fetch_name_value("ALL_TOUCHED"), + Some("FALSE".to_string()) + ); + assert_eq!(c_options.fetch_name_value("BURN_VALUE_FROM"), None); + assert_eq!( + c_options.fetch_name_value("MERGE_ALG"), + Some("REPLACE".to_string()) + ); + assert_eq!( + c_options.fetch_name_value("CHUNKYSIZE"), + Some("0".to_string()) + ); + assert_eq!( + c_options.fetch_name_value("OPTIM"), + Some("AUTO".to_string()) + ); + } + + #[cfg(feature = "gdal-sys")] + #[test] + fn test_rasterize() { + crate::global::with_global_gdal_api(|api| { + let wkt = "POLYGON ((2 2, 2 4.25, 4.25 4.25, 4.25 2, 2 2))"; + let poly = Geometry::from_wkt(api, wkt).unwrap(); + + let driver = crate::driver::DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let dataset = driver.create("", 5, 5, 1).unwrap(); + + let bands = [1]; + let geometries = [&poly]; + let burn_values = [1.0]; + rasterize(api, &dataset, &bands, &geometries, &burn_values, None).unwrap(); + + let rb = dataset.rasterband(1).unwrap(); + let values = rb.read_as::((0, 0), (5, 5), (5, 5), None).unwrap(); + assert_eq!( + values.data(), + [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0] + ); + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/raster/rasterize_affine.rs b/c/sedona-gdal/src/raster/rasterize_affine.rs new file mode 100644 index 000000000..68dad7bbb --- /dev/null +++ b/c/sedona-gdal/src/raster/rasterize_affine.rs @@ -0,0 +1,368 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! Fast affine-transformer rasterize wrapper. +//! +//! GDALRasterizeGeometries() will internally call GDALCreateGenImgProjTransformer2() +//! if pfnTransformer is NULL, even in the common case where only a GeoTransform-based +//! affine conversion from georeferenced coords to pixel/line is needed. +//! +//! This module supplies a minimal GDALTransformerFunc that applies the dataset +//! GeoTransform (and its inverse), avoiding expensive transformer creation. + +use std::ffi::{c_int, c_void}; +use std::ptr; + +use crate::dataset::Dataset; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::gdal_dyn_bindgen::{CE_Failure, CE_None}; +use crate::geo_transform::{GeoTransform, GeoTransformEx}; +use crate::vector::Geometry; + +#[repr(C)] +struct AffineTransformArg { + gt: GeoTransform, + inv_gt: GeoTransform, +} + +unsafe extern "C" fn affine_transformer( + p_transformer_arg: *mut c_void, + b_dst_to_src: c_int, + n_point_count: c_int, + x: *mut f64, + y: *mut f64, + _z: *mut f64, + pan_success: *mut c_int, +) -> c_int { + if p_transformer_arg.is_null() || x.is_null() || y.is_null() || pan_success.is_null() { + return 0; + } + if n_point_count < 0 { + return 0; + } + + // Treat transformer arg as immutable. + let arg = &*(p_transformer_arg as *const AffineTransformArg); + let affine = if b_dst_to_src == 0 { + &arg.inv_gt + } else { + &arg.gt + }; + + let n = n_point_count as usize; + for i in 0..n { + // SAFETY: x/y/pan_success are assumed to point to arrays of length n_point_count. + let xin = unsafe { *x.add(i) }; + let yin = unsafe { *y.add(i) }; + let (xout, yout) = affine.apply(xin, yin); + unsafe { + *x.add(i) = xout; + *y.add(i) = yout; + *pan_success.add(i) = 1; + } + } + + 1 +} + +/// Rasterize geometries with an affine transformer derived from the destination dataset. +/// +/// This mirrors [`rasterize()`](super::rasterize) but avoids GDAL's slow default +/// transformer creation. +/// +/// Assumptions: +/// - Geometry coordinates are already in the destination dataset georeferenced coordinate space. +/// - Only GeoTransform-based affine conversion is supported (no GCP/RPC/geolocs). +pub fn rasterize_affine( + api: &'static GdalApi, + dataset: &Dataset, + bands: &[usize], + geometries: &[Geometry], + burn_values: &[f64], + all_touched: bool, +) -> Result<()> { + if bands.is_empty() { + return Err(GdalError::BadArgument( + "`bands` must not be empty".to_string(), + )); + } + if burn_values.len() != geometries.len() { + return Err(GdalError::BadArgument(format!( + "Burn values length ({}) must match geometries length ({})", + burn_values.len(), + geometries.len() + ))); + } + + let raster_count = dataset.raster_count(); + for band in bands { + let is_good = *band > 0 && *band <= raster_count; + if !is_good { + return Err(GdalError::BadArgument(format!( + "Band index {} is out of bounds", + *band + ))); + } + } + + let bands_i32: Vec = bands.iter().map(|&band| band as c_int).collect(); + + let c_options = if all_touched { + [c"ALL_TOUCHED=TRUE".as_ptr(), ptr::null_mut()] + } else { + [c"ALL_TOUCHED=FALSE".as_ptr(), ptr::null_mut()] + }; + + let geometries_c: Vec<_> = geometries.iter().map(|geo| geo.c_geometry()).collect(); + let burn_values_expanded: Vec = burn_values + .iter() + .flat_map(|burn| std::iter::repeat_n(burn, bands_i32.len())) + .copied() + .collect(); + + let gt = dataset.geo_transform().map_err(|_e| { + GdalError::BadArgument( + "Missing geotransform: only geotransform-based affine rasterize is supported" + .to_string(), + ) + })?; + let inv_gt = gt.invert().map_err(|_e| { + GdalError::BadArgument( + "Non-invertible geotransform: only geotransform-based affine rasterize is supported" + .to_string(), + ) + })?; + let mut arg = AffineTransformArg { gt, inv_gt }; + + unsafe { + let error = call_gdal_api!( + api, + GDALRasterizeGeometries, + dataset.c_dataset(), + bands_i32.len() as c_int, + bands_i32.as_ptr(), + geometries_c.len() as c_int, + geometries_c.as_ptr(), + (affine_transformer as *const ()).cast::() as *mut c_void, + (&mut arg as *mut AffineTransformArg).cast::(), + burn_values_expanded.as_ptr(), + c_options.as_ptr() as *mut *mut i8, + ptr::null_mut(), + ptr::null_mut() + ); + if error != CE_None { + return Err(api.last_cpl_err(CE_Failure as u32)); + } + } + Ok(()) +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use super::*; + + use crate::driver::{Driver, DriverManager}; + use crate::global::with_global_gdal_api; + use crate::raster::rasterize; + use crate::raster::types::Buffer; + + fn mem_driver(api: &'static GdalApi) -> Driver { + DriverManager::get_driver_by_name(api, "MEM").unwrap() + } + + fn make_dataset_u8( + api: &'static GdalApi, + width: usize, + height: usize, + gt: GeoTransform, + ) -> Result { + let driver = mem_driver(api); + let ds = driver.create_with_band_type::("", width, height, 1)?; + ds.set_geo_transform(>)?; + let band = ds.rasterband(1)?; + let mut buf = Buffer::new((width, height), vec![0u8; width * height]); + band.write((0, 0), (width, height), &mut buf)?; + Ok(ds) + } + + fn read_u8(ds: &Dataset, width: usize, height: usize) -> Vec { + let band = ds.rasterband(1).unwrap(); + let buf = band + .read_as::((0, 0), (width, height), (width, height), None) + .unwrap(); + buf.data().to_vec() + } + + fn poly_from_pixel_rect( + api: &'static GdalApi, + gt: &GeoTransform, + x0: f64, + y0: f64, + x1: f64, + y1: f64, + ) -> Geometry { + let (wx0, wy0) = gt.apply(x0, y0); + let (wx1, wy1) = gt.apply(x1, y0); + let (wx2, wy2) = gt.apply(x1, y1); + let (wx3, wy3) = gt.apply(x0, y1); + let wkt = + format!("POLYGON (({wx0} {wy0}, {wx1} {wy1}, {wx2} {wy2}, {wx3} {wy3}, {wx0} {wy0}))"); + Geometry::from_wkt(api, &wkt).unwrap() + } + + fn line_from_pixel_points( + api: &'static GdalApi, + gt: &GeoTransform, + pts: &[(f64, f64)], + ) -> Geometry { + assert!(pts.len() >= 2); + let mut s = String::from("LINESTRING ("); + for (i, (px, py)) in pts.iter().copied().enumerate() { + let (wx, wy) = gt.apply(px, py); + if i > 0 { + s.push_str(", "); + } + s.push_str(&format!("{wx} {wy}")); + } + s.push(')'); + Geometry::from_wkt(api, &s).unwrap() + } + + #[test] + fn test_rasterize_affine_matches_baseline_north_up() { + with_global_gdal_api(|api| { + let (w, h) = (32usize, 24usize); + let gt: GeoTransform = [100.0, 2.0, 0.0, 200.0, 0.0, -2.0]; + + let geom_baseline = poly_from_pixel_rect(api, >, 3.2, 4.7, 20.4, 18.1); + let geom_affine = poly_from_pixel_rect(api, >, 3.2, 4.7, 20.4, 18.1); + + let ds_baseline = make_dataset_u8(api, w, h, gt).unwrap(); + let ds_affine = make_dataset_u8(api, w, h, gt).unwrap(); + + let geom_refs: Vec<&Geometry> = vec![&geom_baseline]; + rasterize(api, &ds_baseline, &[1], &geom_refs, &[1.0], None).unwrap(); + rasterize_affine(api, &ds_affine, &[1], &[geom_affine], &[1.0], false).unwrap(); + + assert_eq!(read_u8(&ds_affine, w, h), read_u8(&ds_baseline, w, h)); + }) + .unwrap(); + } + + #[test] + fn test_rasterize_affine_matches_baseline_rotated_gt_all_touched() { + with_global_gdal_api(|api| { + let (w, h) = (40usize, 28usize); + // Rotated/skewed GeoTransform. + let gt: GeoTransform = [10.0, 1.2, 0.15, 50.0, -0.1, -1.1]; + + let geom_baseline = poly_from_pixel_rect(api, >, 5.25, 4.5, 25.75, 20.25); + let geom_affine = poly_from_pixel_rect(api, >, 5.25, 4.5, 25.75, 20.25); + + let ds_baseline = make_dataset_u8(api, w, h, gt).unwrap(); + let ds_affine = make_dataset_u8(api, w, h, gt).unwrap(); + + let geom_refs: Vec<&Geometry> = vec![&geom_baseline]; + rasterize( + api, + &ds_baseline, + &[1], + &geom_refs, + &[1.0], + Some(crate::raster::RasterizeOptions { + all_touched: true, + ..Default::default() + }), + ) + .unwrap(); + rasterize_affine(api, &ds_affine, &[1], &[geom_affine], &[1.0], true).unwrap(); + + assert_eq!(read_u8(&ds_affine, w, h), read_u8(&ds_baseline, w, h)); + }) + .unwrap(); + } + + #[test] + fn test_rasterize_affine_matches_baseline_linestring() { + with_global_gdal_api(|api| { + let (w, h) = (64usize, 48usize); + // Rotated/skewed GeoTransform. + let gt: GeoTransform = [5.0, 1.0, 0.2, 100.0, -0.15, -1.05]; + + // A polyline with many vertices, defined in pixel/line space. + let mut pts: Vec<(f64, f64)> = Vec::new(); + for i in 0..200 { + let t = i as f64 / 199.0; + let x = 2.625 + t * ((w as f64) - 5.25); + let y = 5.25 + (t * 6.0).sin() * 8.0 + t * ((h as f64) - 12.25); + pts.push((x, y)); + } + let geom_baseline = line_from_pixel_points(api, >, &pts); + let geom_affine = line_from_pixel_points(api, >, &pts); + + let ds_baseline = make_dataset_u8(api, w, h, gt).unwrap(); + let ds_affine = make_dataset_u8(api, w, h, gt).unwrap(); + + let geom_refs: Vec<&Geometry> = vec![&geom_baseline]; + rasterize(api, &ds_baseline, &[1], &geom_refs, &[1.0], None).unwrap(); + rasterize_affine(api, &ds_affine, &[1], &[geom_affine], &[1.0], false).unwrap(); + + let got = read_u8(&ds_affine, w, h); + let expected = read_u8(&ds_baseline, w, h); + if got != expected { + let mut diffs = Vec::new(); + for (i, (a, b)) in got + .iter() + .copied() + .zip(expected.iter().copied()) + .enumerate() + { + if a != b { + let x = i % w; + let y = i / w; + diffs.push((x, y, a, b)); + } + } + panic!( + "raster mismatch: {} differing pixels; first 10: {:?}", + diffs.len(), + &diffs[..diffs.len().min(10)] + ); + } + }) + .unwrap(); + } + + #[test] + fn test_rasterize_affine_fails_on_noninvertible_gt() { + with_global_gdal_api(|api| { + let (w, h) = (8usize, 8usize); + let gt: GeoTransform = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]; + let ds = make_dataset_u8(api, w, h, gt).unwrap(); + let geom = Geometry::from_wkt(api, "POINT (0 0)").unwrap(); + let err = rasterize_affine(api, &ds, &[1], &[geom], &[1.0], true).unwrap_err(); + match err { + GdalError::BadArgument(msg) => { + assert!(msg.contains("Non-invertible geotransform")); + } + other => panic!("Unexpected error: {other:?}"), + } + }) + .unwrap(); + } +} diff --git a/c/sedona-gdal/src/vrt.rs b/c/sedona-gdal/src/vrt.rs new file mode 100644 index 000000000..009b62452 --- /dev/null +++ b/c/sedona-gdal/src/vrt.rs @@ -0,0 +1,336 @@ +// Licensed to the Apache Software Foundation (ASF) under one +// or more contributor license agreements. See the NOTICE file +// distributed with this work for additional information +// regarding copyright ownership. The ASF licenses this file +// to you under the Apache License, Version 2.0 (the +// "License"); you may not use this file except in compliance +// with the License. You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, +// software distributed under the License is distributed on an +// "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY +// KIND, either express or implied. See the License for the +// specific language governing permissions and limitations +// under the License. + +//! GDAL VRT (Virtual Raster) API wrappers. + +use std::ffi::CString; +use std::ops::{Deref, DerefMut}; +use std::ptr::null_mut; + +use crate::cpl::CslStringList; +use crate::dataset::Dataset; +use crate::errors::{GdalError, Result}; +use crate::gdal_api::{call_gdal_api, GdalApi}; +use crate::raster::RasterBand; +use crate::{gdal_dyn_bindgen::*, raster::types::GdalDataType}; + +/// Special value indicating that nodata is not set for a VRT source. +/// Matches `VRT_NODATA_UNSET` from GDAL's `gdal_vrt.h`. +pub const NODATA_UNSET: f64 = -1234.56; + +/// A VRT (Virtual Raster) dataset. +pub struct VrtDataset { + dataset: Dataset, +} + +unsafe impl Send for VrtDataset {} + +impl VrtDataset { + /// Create a new empty VRT dataset with the given dimensions. + pub fn create(api: &'static GdalApi, x_size: usize, y_size: usize) -> Result { + let x: i32 = x_size.try_into()?; + let y: i32 = y_size.try_into()?; + let c_dataset = unsafe { call_gdal_api!(api, VRTCreate, x, y) }; + + if c_dataset.is_null() { + return Err(GdalError::NullPointer { + method_name: "VRTCreate", + msg: String::new(), + }); + } + + Ok(VrtDataset { + dataset: Dataset::new_owned(api, c_dataset), + }) + } + + /// Consume this VRT and return the underlying `Dataset`, transferring ownership. + pub fn as_dataset(self) -> Dataset { + let VrtDataset { dataset } = self; + dataset + } + + /// Add a new band to the VRT dataset. + /// + /// Returns the 1-based index of the newly created band on success. + pub fn add_band(&mut self, data_type: GdalDataType, options: Option<&[&str]>) -> Result { + let csl = CslStringList::try_from_iter(options.unwrap_or(&[]).iter().copied())?; + + // Preserve null semantics: pass null when no options given. + let opts_ptr = if csl.is_empty() { + null_mut() + } else { + csl.as_ptr() + }; + + let rv = unsafe { + call_gdal_api!( + self.dataset.api(), + GDALAddBand, + self.dataset.c_dataset(), + data_type.to_c(), + opts_ptr + ) + }; + + if rv != CE_None { + return Err(self.dataset.api().last_cpl_err(rv as u32)); + } + + Ok(self.raster_count()) + } + + /// Fetch a band from the VRT dataset as a `VrtRasterBand`. + pub fn rasterband(&self, band_index: usize) -> Result> { + let band = self.dataset.rasterband(band_index)?; + Ok(VrtRasterBand { band }) + } +} + +impl Deref for VrtDataset { + type Target = Dataset; + + fn deref(&self) -> &Self::Target { + &self.dataset + } +} + +impl DerefMut for VrtDataset { + fn deref_mut(&mut self) -> &mut Self::Target { + &mut self.dataset + } +} + +impl AsRef for VrtDataset { + fn as_ref(&self) -> &Dataset { + &self.dataset + } +} + +/// A raster band within a VRT dataset. +pub struct VrtRasterBand<'a> { + band: RasterBand<'a>, +} + +impl<'a> VrtRasterBand<'a> { + /// Returns the raw GDAL raster band handle. + pub fn c_rasterband(&self) -> GDALRasterBandH { + self.band.c_rasterband() + } + + /// Adds a simple source to this VRT band. + /// + /// # Arguments + /// * `source_band` - The source raster band to read from + /// * `src_window` - Source window as `(x_offset, y_offset, x_size, y_size)` in pixels + /// * `dst_window` - Destination window as `(x_offset, y_offset, x_size, y_size)` in pixels + /// * `resampling` - Optional resampling method (e.g. "near", "bilinear", "cubic"). + /// * `nodata` - Optional nodata value for the source. If None, uses `NODATA_UNSET`. + pub fn add_simple_source( + &self, + source_band: &RasterBand<'a>, + src_window: (i32, i32, i32, i32), + dst_window: (i32, i32, i32, i32), + resampling: Option<&str>, + nodata: Option, + ) -> Result<()> { + let c_resampling = resampling.map(CString::new).transpose()?; + + let resampling_ptr = c_resampling + .as_ref() + .map(|s| s.as_ptr()) + .unwrap_or(null_mut()); + + let nodata_value = nodata.unwrap_or(NODATA_UNSET); + + let rv = unsafe { + call_gdal_api!( + self.band.api(), + VRTAddSimpleSource, + self.band.c_rasterband(), + source_band.c_rasterband(), + src_window.0, + src_window.1, + src_window.2, + src_window.3, + dst_window.0, + dst_window.1, + dst_window.2, + dst_window.3, + resampling_ptr, + nodata_value + ) + }; + + if rv != CE_None { + return Err(self.band.api().last_cpl_err(rv as u32)); + } + Ok(()) + } + + /// Sets the nodata value for this VRT band. + pub fn set_no_data_value(&self, nodata: f64) -> Result<()> { + let rv = unsafe { + call_gdal_api!( + self.band.api(), + GDALSetRasterNoDataValue, + self.band.c_rasterband(), + nodata + ) + }; + + if rv != CE_None { + return Err(self.band.api().last_cpl_err(rv as u32)); + } + Ok(()) + } +} + +impl<'a> Deref for VrtRasterBand<'a> { + type Target = RasterBand<'a>; + + fn deref(&self) -> &Self::Target { + &self.band + } +} + +impl<'a> DerefMut for VrtRasterBand<'a> { + fn deref_mut(&mut self) -> &mut Self::Target { + &mut self.band + } +} + +impl<'a> AsRef> for VrtRasterBand<'a> { + fn as_ref(&self) -> &RasterBand<'a> { + &self.band + } +} + +#[cfg(all(test, feature = "gdal-sys"))] +mod tests { + use crate::dataset::Dataset; + use crate::gdal_dyn_bindgen::GDAL_OF_READONLY; + use crate::global::with_global_gdal_api; + use crate::raster::types::GdalDataType; + use crate::vrt::{VrtDataset, NODATA_UNSET}; + + fn fixture(name: &str) -> String { + sedona_testing::data::test_raster(name).unwrap() + } + + #[test] + fn test_vrt_create() { + with_global_gdal_api(|api| { + let vrt = VrtDataset::create(api, 100, 100).unwrap(); + assert_eq!(vrt.raster_count(), 0); + assert!(!vrt.c_dataset().is_null()); + }) + .unwrap(); + } + + #[test] + fn test_vrt_add_band() { + with_global_gdal_api(|api| { + let mut vrt = VrtDataset::create(api, 100, 100).unwrap(); + let band_idx = vrt.add_band(GdalDataType::Float32, None).unwrap(); + assert_eq!(band_idx, 1); + assert_eq!(vrt.raster_count(), 1); + + let band_idx = vrt.add_band(GdalDataType::UInt8, None).unwrap(); + assert_eq!(band_idx, 2); + assert_eq!(vrt.raster_count(), 2); + }) + .unwrap(); + } + + #[test] + fn test_vrt_set_geo_transform() { + with_global_gdal_api(|api| { + let vrt = VrtDataset::create(api, 100, 100).unwrap(); + let transform = [0.0, 1.0, 0.0, 100.0, 0.0, -1.0]; + vrt.set_geo_transform(&transform).unwrap(); + assert_eq!(vrt.geo_transform().unwrap(), transform); + }) + .unwrap(); + } + + #[test] + fn test_vrt_set_projection() { + with_global_gdal_api(|api| { + let vrt = VrtDataset::create(api, 100, 100).unwrap(); + vrt.set_projection("EPSG:4326").unwrap(); + assert!(vrt.projection().contains("4326")); + }) + .unwrap(); + } + + #[test] + fn test_vrt_add_simple_source() { + with_global_gdal_api(|api| { + let source = Dataset::open_ex( + api, + &fixture("tinymarble.tif"), + GDAL_OF_READONLY, + None, + None, + None, + ) + .unwrap(); + let source_band_type = source.rasterband(1).unwrap().band_type(); + + let mut vrt = VrtDataset::create(api, 1, 1).unwrap(); + vrt.add_band(source_band_type, None).unwrap(); + + let source_band = source.rasterband(1).unwrap(); + let vrt_band = vrt.rasterband(1).unwrap(); + + vrt_band + .add_simple_source(&source_band, (0, 0, 1, 1), (0, 0, 1, 1), None, None) + .unwrap(); + + let source_px = source_band + .read_as::((0, 0), (1, 1), (1, 1), None) + .unwrap() + .data()[0]; + let vrt_px = vrt_band + .read_as::((0, 0), (1, 1), (1, 1), None) + .unwrap() + .data()[0]; + + assert_eq!(vrt_px, source_px); + }) + .unwrap(); + } + + #[test] + fn test_vrt_nodata_unset() { + assert_eq!(NODATA_UNSET, -1234.56); + } + + #[test] + #[allow(clippy::float_cmp)] + fn test_vrt_set_no_data_value() { + with_global_gdal_api(|api| { + let mut vrt = VrtDataset::create(api, 1, 1).unwrap(); + vrt.add_band(GdalDataType::UInt8, None).unwrap(); + let band = vrt.rasterband(1).unwrap(); + band.set_no_data_value(-9999.0).unwrap(); + assert_eq!(band.no_data_value(), Some(-9999.0)); + }) + .unwrap(); + } +} From 40df37c96c1747c1fded8a8c28e41dc4cdaf2821 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Tue, 10 Mar 2026 01:50:38 +0800 Subject: [PATCH 09/11] refactor(sedona-gdal): prefer imported raster helpers --- c/sedona-gdal/src/raster/polygonize.rs | 6 ++++-- c/sedona-gdal/src/raster/rasterize.rs | 14 ++++++++------ c/sedona-gdal/src/raster/rasterize_affine.rs | 4 ++-- 3 files changed, 14 insertions(+), 10 deletions(-) diff --git a/c/sedona-gdal/src/raster/polygonize.rs b/c/sedona-gdal/src/raster/polygonize.rs index 58cba2736..0315bd656 100644 --- a/c/sedona-gdal/src/raster/polygonize.rs +++ b/c/sedona-gdal/src/raster/polygonize.rs @@ -138,6 +138,7 @@ mod tests { use crate::global::with_global_gdal_api; use crate::raster::Buffer; use crate::vector::FieldDefn; + use crate::vsi::unlink_mem_file; with_global_gdal_api(|api| { let mem_driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); @@ -209,7 +210,7 @@ mod tests { .count(); assert_eq!(ones_8, 1); - crate::vsi::unlink_mem_file(api, gpkg_path).unwrap(); + unlink_mem_file(api, gpkg_path).unwrap(); }) .unwrap(); } @@ -222,6 +223,7 @@ mod tests { use crate::global::with_global_gdal_api; use crate::raster::Buffer; use crate::vector::FieldDefn; + use crate::vsi::unlink_mem_file; with_global_gdal_api(|api| { let mem_driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); @@ -267,7 +269,7 @@ mod tests { let only_val = layer.features().next().unwrap().field_as_integer(0); assert_eq!(only_val, Some(7)); - crate::vsi::unlink_mem_file(api, gpkg_path).unwrap(); + unlink_mem_file(api, gpkg_path).unwrap(); }) .unwrap(); } diff --git a/c/sedona-gdal/src/raster/rasterize.rs b/c/sedona-gdal/src/raster/rasterize.rs index dbcad37ca..f557d10e6 100644 --- a/c/sedona-gdal/src/raster/rasterize.rs +++ b/c/sedona-gdal/src/raster/rasterize.rs @@ -23,7 +23,7 @@ use std::ptr; use crate::cpl::CslStringList; use crate::dataset::Dataset; -use crate::errors::Result; +use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; use crate::vector::Geometry; @@ -147,12 +147,12 @@ pub fn rasterize( options: Option, ) -> Result<()> { if band_list.is_empty() { - return Err(crate::errors::GdalError::BadArgument( + return Err(GdalError::BadArgument( "`band_list` must not be empty".to_string(), )); } if burn_values.len() != geometries.len() { - return Err(crate::errors::GdalError::BadArgument(format!( + return Err(GdalError::BadArgument(format!( "burn_values length ({}) must match geometries length ({})", burn_values.len(), geometries.len() @@ -162,7 +162,7 @@ pub fn rasterize( for &band in band_list { let is_good = band > 0 && (band as usize) <= raster_count; if !is_good { - return Err(crate::errors::GdalError::BadArgument(format!( + return Err(GdalError::BadArgument(format!( "Band index {} is out of bounds", band ))); @@ -211,6 +211,8 @@ pub fn rasterize( #[cfg(test)] mod tests { use super::*; + use crate::driver::DriverManager; + use crate::global::with_global_gdal_api; #[test] fn test_rasterizeoptions_as_ptr() { @@ -237,11 +239,11 @@ mod tests { #[cfg(feature = "gdal-sys")] #[test] fn test_rasterize() { - crate::global::with_global_gdal_api(|api| { + with_global_gdal_api(|api| { let wkt = "POLYGON ((2 2, 2 4.25, 4.25 4.25, 4.25 2, 2 2))"; let poly = Geometry::from_wkt(api, wkt).unwrap(); - let driver = crate::driver::DriverManager::get_driver_by_name(api, "MEM").unwrap(); + let driver = DriverManager::get_driver_by_name(api, "MEM").unwrap(); let dataset = driver.create("", 5, 5, 1).unwrap(); let bands = [1]; diff --git a/c/sedona-gdal/src/raster/rasterize_affine.rs b/c/sedona-gdal/src/raster/rasterize_affine.rs index 68dad7bbb..7a4ae7a2c 100644 --- a/c/sedona-gdal/src/raster/rasterize_affine.rs +++ b/c/sedona-gdal/src/raster/rasterize_affine.rs @@ -178,8 +178,8 @@ mod tests { use crate::driver::{Driver, DriverManager}; use crate::global::with_global_gdal_api; - use crate::raster::rasterize; use crate::raster::types::Buffer; + use crate::raster::{rasterize, RasterizeOptions}; fn mem_driver(api: &'static GdalApi) -> Driver { DriverManager::get_driver_by_name(api, "MEM").unwrap() @@ -284,7 +284,7 @@ mod tests { &[1], &geom_refs, &[1.0], - Some(crate::raster::RasterizeOptions { + Some(RasterizeOptions { all_touched: true, ..Default::default() }), From 4ae2be59a7abcaa3d1824d5019be85a20ccab5b1 Mon Sep 17 00:00:00 2001 From: Kontinuation Date: Wed, 11 Mar 2026 15:45:18 +0800 Subject: [PATCH 10/11] fix(sedona-gdal): update polygonize raster type imports --- c/sedona-gdal/src/raster/polygonize.rs | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/c/sedona-gdal/src/raster/polygonize.rs b/c/sedona-gdal/src/raster/polygonize.rs index 0315bd656..4e5425c1b 100644 --- a/c/sedona-gdal/src/raster/polygonize.rs +++ b/c/sedona-gdal/src/raster/polygonize.rs @@ -136,7 +136,7 @@ mod tests { use crate::dataset::LayerOptions; use crate::driver::DriverManager; use crate::global::with_global_gdal_api; - use crate::raster::Buffer; + use crate::raster::types::Buffer; use crate::vector::FieldDefn; use crate::vsi::unlink_mem_file; @@ -221,7 +221,7 @@ mod tests { use crate::dataset::LayerOptions; use crate::driver::DriverManager; use crate::global::with_global_gdal_api; - use crate::raster::Buffer; + use crate::raster::types::Buffer; use crate::vector::FieldDefn; use crate::vsi::unlink_mem_file; From e1d4617130f6afc6f52e2be200a3554f895d5287 Mon Sep 17 00:00:00 2001 From: Kristin Cowalcijk Date: Fri, 13 Mar 2026 14:36:09 +0800 Subject: [PATCH 11/11] refactor(sedona-gdal): remove raster re-export aliases --- c/sedona-gdal/src/raster/polygonize.rs | 8 ++++---- c/sedona-gdal/src/raster/rasterize.rs | 2 +- c/sedona-gdal/src/raster/rasterize_affine.rs | 4 ++-- c/sedona-gdal/src/vrt.rs | 2 +- 4 files changed, 8 insertions(+), 8 deletions(-) diff --git a/c/sedona-gdal/src/raster/polygonize.rs b/c/sedona-gdal/src/raster/polygonize.rs index 4e5425c1b..6140bebcc 100644 --- a/c/sedona-gdal/src/raster/polygonize.rs +++ b/c/sedona-gdal/src/raster/polygonize.rs @@ -21,8 +21,8 @@ use crate::cpl::CslStringList; use crate::errors::Result; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; -use crate::raster::RasterBand; -use crate::vector::Layer; +use crate::raster::rasterband::RasterBand; +use crate::vector::layer::Layer; #[derive(Clone, Debug, Default)] pub struct PolygonizeOptions { @@ -137,7 +137,7 @@ mod tests { use crate::driver::DriverManager; use crate::global::with_global_gdal_api; use crate::raster::types::Buffer; - use crate::vector::FieldDefn; + use crate::vector::feature::FieldDefn; use crate::vsi::unlink_mem_file; with_global_gdal_api(|api| { @@ -222,7 +222,7 @@ mod tests { use crate::driver::DriverManager; use crate::global::with_global_gdal_api; use crate::raster::types::Buffer; - use crate::vector::FieldDefn; + use crate::vector::feature::FieldDefn; use crate::vsi::unlink_mem_file; with_global_gdal_api(|api| { diff --git a/c/sedona-gdal/src/raster/rasterize.rs b/c/sedona-gdal/src/raster/rasterize.rs index f557d10e6..55e8f919d 100644 --- a/c/sedona-gdal/src/raster/rasterize.rs +++ b/c/sedona-gdal/src/raster/rasterize.rs @@ -26,7 +26,7 @@ use crate::dataset::Dataset; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::*; -use crate::vector::Geometry; +use crate::vector::geometry::Geometry; /// Source of burn values. #[derive(Copy, Clone, Debug)] diff --git a/c/sedona-gdal/src/raster/rasterize_affine.rs b/c/sedona-gdal/src/raster/rasterize_affine.rs index 7a4ae7a2c..431a94f44 100644 --- a/c/sedona-gdal/src/raster/rasterize_affine.rs +++ b/c/sedona-gdal/src/raster/rasterize_affine.rs @@ -32,7 +32,7 @@ use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; use crate::gdal_dyn_bindgen::{CE_Failure, CE_None}; use crate::geo_transform::{GeoTransform, GeoTransformEx}; -use crate::vector::Geometry; +use crate::vector::geometry::Geometry; #[repr(C)] struct AffineTransformArg { @@ -178,8 +178,8 @@ mod tests { use crate::driver::{Driver, DriverManager}; use crate::global::with_global_gdal_api; + use crate::raster::rasterize::{rasterize, RasterizeOptions}; use crate::raster::types::Buffer; - use crate::raster::{rasterize, RasterizeOptions}; fn mem_driver(api: &'static GdalApi) -> Driver { DriverManager::get_driver_by_name(api, "MEM").unwrap() diff --git a/c/sedona-gdal/src/vrt.rs b/c/sedona-gdal/src/vrt.rs index 009b62452..41e960e60 100644 --- a/c/sedona-gdal/src/vrt.rs +++ b/c/sedona-gdal/src/vrt.rs @@ -25,7 +25,7 @@ use crate::cpl::CslStringList; use crate::dataset::Dataset; use crate::errors::{GdalError, Result}; use crate::gdal_api::{call_gdal_api, GdalApi}; -use crate::raster::RasterBand; +use crate::raster::rasterband::RasterBand; use crate::{gdal_dyn_bindgen::*, raster::types::GdalDataType}; /// Special value indicating that nodata is not set for a VRT source.