wip

Cargo.toml

@@ -3,6 +3,7 @@ members = [
   "geopath",
   "itm",
   "nasadem",
+  "nasadem/zorderize",
   "propah",
   "terrain",
 ]

nasadem/src/lib.rs

@@ -89,11 +89,11 @@ impl SampleStore {
     fn min(&self) -> i16 {
         match self {
             Self::Tombstone => 0,
-            Self::InMem(samples) => samples.iter().max().copied().unwrap(),
+            Self::InMem(samples) => samples.iter().min().copied().unwrap(),
             Self::MemMap(raw) => (*raw)
                 .chunks_exact(2)
                 .map(|mut bytes| (&mut bytes).read_i16::<BE>().unwrap())
-                .max()
+                .min()
                 .unwrap(),
         }
     }
@@ -224,10 +224,15 @@ impl Tile {
         }
     }
 
+    /// Returns this tile's (x, y) dimensions.
+    pub fn dimensions(&self) -> (usize, usize) {
+        self.dimensions
+    }
+
     /// Returns the number of samples in this tile.
     #[allow(clippy::len_without_is_empty)]
     pub fn len(&self) -> usize {
-        let (x, y) = self.dimensions;
+        let (x, y) = self.dimensions();
         x * y
     }
 
@@ -261,9 +266,9 @@ impl Tile {
         let (idx_x, idx_y) = self.coord_to_xy(coord);
         #[allow(clippy::cast_possible_wrap)]
         if 0 <= idx_x
-            && idx_x < self.dimensions.0 as isize
+            && idx_x < self.dimensions().0 as isize
             && 0 <= idx_y
-            && idx_y < self.dimensions.1 as isize
+            && idx_y < self.dimensions().1 as isize
         {
             #[allow(clippy::cast_sign_loss)]
             let idx_1d = self.xy_to_linear_index((idx_x as usize, idx_y as usize));
@@ -283,7 +288,7 @@ impl Tile {
 
     /// Returns and iterator over `self`'s grid squares.
     pub fn iter(&self) -> impl Iterator<Item = Sample<'_>> + '_ {
-        (0..(self.dimensions.0 * self.dimensions.1)).map(|index| Sample { tile: self, index })
+        (0..(self.dimensions().0 * self.dimensions().1)).map(|index| Sample { tile: self, index })
     }
 
     /// Returns this tile's outline as a polygon.
@@ -326,13 +331,13 @@ impl Tile {
     }
 
     fn linear_index_to_xy(&self, idx: usize) -> (usize, usize) {
-        let y = idx / self.dimensions.0;
-        let x = idx % self.dimensions.1;
-        (x, self.dimensions.1 - 1 - y)
+        let y = idx / self.dimensions().0;
+        let x = idx % self.dimensions().1;
+        (x, self.dimensions().1 - 1 - y)
     }
 
     fn xy_to_linear_index(&self, (x, y): (usize, usize)) -> usize {
-        self.dimensions.0 * (self.dimensions.1 - y - 1) + x
+        self.dimensions().0 * (self.dimensions().1 - y - 1) + x
     }
 
     fn xy_to_polygon(&self, (x, y): (usize, usize)) -> Polygon<C> {
@@ -376,8 +381,11 @@ impl<'a> Sample<'a> {
     }
 
     pub fn polygon(&self) -> Polygon {
-        self.tile
-            .xy_to_polygon(self.tile.linear_index_to_xy(self.index))
+        self.tile.xy_to_polygon(self.index())
+    }
+
+    pub fn index(&self) -> (usize, usize) {
+        self.tile.linear_index_to_xy(self.index)
     }
 }
 

nasadem/zorderize/Cargo.toml

@@ -0,0 +1,10 @@
+[package]
+name = "zorderize"
+version = "0.1.0"
+edition = "2021"
+
+[dependencies]
+clap    = { workspace = true }
+glob    = "0"
+image   = "0"
+nasadem = { path = ".." }

nasadem/zorderize/src/main.rs

@@ -0,0 +1,60 @@
+use clap::Parser;
+use image::{
+    imageops::{resize, FilterType},
+    ImageBuffer, Luma,
+};
+use nasadem::Tile;
+use std::path::PathBuf;
+
+#[derive(Debug, Parser)]
+struct Cmd {
+    /// Source NASADEM hgt file.
+    tile: PathBuf,
+    /// Write output z-order tile to this directory.
+    dest: Option<PathBuf>,
+}
+
+fn tile_to_img(tile: &Tile) -> ImageBuffer<Luma<u16>, Vec<u16>> {
+    let (x_dim, y_dim) = tile.dimensions();
+    let mut img = ImageBuffer::new(x_dim as u32, y_dim as u32);
+    let min_elev = tile.min_elevation() as f32;
+    let max_elev = tile.max_elevation() as f32;
+    let scale = |m: i16| {
+        let m = m as f32;
+        (m - min_elev) / (max_elev - min_elev) * u16::MAX as f32
+    };
+    for sample in tile.iter() {
+        let (x, y) = sample.index();
+        let elev = sample.elevation();
+        let scaled_elev = scale(elev);
+        img.put_pixel(x as u32, (y_dim - 1 - y) as u32, Luma([scaled_elev as u16]))
+    }
+    img
+}
+
+fn scale_to_pow2(img: &ImageBuffer<Luma<u16>, Vec<u16>>) -> ImageBuffer<Luma<u16>, Vec<u16>> {
+    let (scaled_x, scaled_y) = match img.dimensions() {
+        (1201, 1201) => (1024, 1024),
+        (3601, 3601) => (2048, 2048),
+        other => panic!("dimensions {other:?} are not expected for SRTM data."),
+    };
+    resize(img, scaled_x, scaled_y, FilterType::Lanczos3)
+}
+
+fn main() {
+    let cli = Cmd::parse();
+    let tile = Tile::load(&cli.tile).unwrap();
+    let img = tile_to_img(&tile);
+    let scaled = scale_to_pow2(&img);
+    let out = cli.dest.unwrap_or_else(|| {
+        let mut out = cli.tile.clone();
+        if out.is_dir() {
+            let name = cli.tile.file_name().unwrap();
+            out.push(name);
+        } else {
+            out.set_extension("png");
+        }
+        out
+    });
+    scaled.save(out).unwrap();
+}