Skip to content

Instantly share code, notes, and snippets.

@JayKickliter
Created January 14, 2024 23:40
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save JayKickliter/8684e27e32521e5e9d5edfebeba54b80 to your computer and use it in GitHub Desktop.
Save JayKickliter/8684e27e32521e5e9d5edfebeba54b80 to your computer and use it in GitHub Desktop.
ESRI Land Cover + GeoTIFF + Rust + GDAL + Projection
use anyhow::Result as AnyRes;
use gdal::{
spatial_ref::{CoordTransform, SpatialRef},
Dataset, GeoTransformEx,
};
fn main() -> AnyRes<()> {
let ds = Dataset::open("lulc2022/17R_20220101-20230101.tif")?;
// WGS84
let src_ref = SpatialRef::from_epsg(4326)?;
// EPSG:32617
let dst_ref = ds.spatial_ref()?;
let mut lat = [31.994247];
let mut lon = [-84.23041];
let transform = CoordTransform::new(&src_ref, &dst_ref)?;
transform.transform_coords(&mut lat, &mut lon, &mut [])?;
println!("{lat:?} {lon:?}");
let transform = ds.geo_transform()?;
let inverse = transform.invert()?;
let (pixel, line) = inverse.apply(lat[0], lon[0]);
println!("(pixel,line): ({pixel},{line})");
Ok(())
}
// 53 │ Upper Left ( 194780.000, 3544360.000) ( 84d13'49.48"W, 31d59'39.29"N)
// 54 │ Lower Left ( 194780.000, 2654230.000) ( 83d59'57.30"W, 23d58'14.51"N)
// 55 │ Upper Right ( 805220.000, 3544360.000) ( 77d46'10.52"W, 31d59'39.29"N)
// 56 │ Lower Right ( 805220.000, 2654230.000) ( 78d 0' 2.70"W, 23d58'14.51"N)
// 57 │ Center ( 500000.000, 3099295.000) ( 81d 0' 0.00"W, 28d 1' 8.01"N)
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment