surtgis 0.10.0

High-performance geospatial analysis CLI
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166

//! Handler for landscape ecology subcommands.

use anyhow::{Context, Result};
use std::collections::HashSet;
use std::time::Instant;

use surtgis_algorithms::landscape::{
    class_metrics, label_patches, landscape_metrics, patch_metrics, patches_to_csv,
};

use crate::commands::LandscapeCommands;
use crate::helpers::{done, parse_connectivity, read_dem, write_result_i32};

pub fn handle(algorithm: LandscapeCommands, compress: bool) -> Result<()> {
    match algorithm {
        LandscapeCommands::LabelPatches {
            input,
            output,
            connectivity,
        } => {
            let conn = parse_connectivity(connectivity)?;
            let raster = read_dem(&input)?;
            let start = Instant::now();
            let (labels, num_patches) =
                label_patches(&raster, conn).context("Failed to label patches")?;
            let elapsed = start.elapsed();
            write_result_i32(&labels, &output, compress)?;
            println!("{} patches found", num_patches);
            done("Label patches", &output, elapsed);
        }

        LandscapeCommands::PatchMetrics {
            input,
            output,
            connectivity,
        } => {
            let conn = parse_connectivity(connectivity)?;
            let raster = read_dem(&input)?;
            let start = Instant::now();
            let (labels, num_patches) =
                label_patches(&raster, conn).context("Failed to label patches")?;
            let patches = patch_metrics(&raster, &labels, num_patches)
                .context("Failed to compute patch metrics")?;
            let elapsed = start.elapsed();
            let csv = patches_to_csv(&patches);
            std::fs::write(&output, &csv).context("Failed to write CSV")?;
            println!(
                "{} patches, {} classes",
                patches.len(),
                patches
                    .iter()
                    .map(|p| p.class)
                    .collect::<HashSet<_>>()
                    .len()
            );
            done("Patch metrics", &output, elapsed);
        }

        LandscapeCommands::ClassMetrics {
            input,
            connectivity,
        } => {
            let conn = parse_connectivity(connectivity)?;
            let raster = read_dem(&input)?;
            let start = Instant::now();
            let (labels, num_patches) =
                label_patches(&raster, conn).context("Failed to label patches")?;
            let patches = patch_metrics(&raster, &labels, num_patches)
                .context("Failed to compute patch metrics")?;
            let cm = class_metrics(&raster, &patches).context("Failed to compute class metrics")?;
            let elapsed = start.elapsed();

            println!(
                "{:<10} {:>10} {:>10} {:>8} {:>12} {:>8} {:>10}",
                "Class", "Area(m\u{b2})", "Proportion", "Patches", "MeanArea", "AI", "Cohesion"
            );
            println!("{}", "-".repeat(78));
            for c in &cm {
                println!(
                    "{:<10} {:>10.1} {:>10.4} {:>8} {:>12.1} {:>8.1} {:>10.1}",
                    c.class,
                    c.area_m2,
                    c.proportion,
                    c.num_patches,
                    c.mean_patch_area_m2,
                    c.ai,
                    c.cohesion
                );
            }
            println!("\n  Processing time: {:.2?}", elapsed);
        }

        LandscapeCommands::LandscapeMetrics { input } => {
            let raster = read_dem(&input)?;
            let start = Instant::now();
            let lm = landscape_metrics(&raster).context("Failed to compute landscape metrics")?;
            let elapsed = start.elapsed();

            println!("Landscape Metrics:");
            println!("  SHDI (Shannon):  {:.4}", lm.shdi);
            println!("  SIDI (Simpson):  {:.4}", lm.sidi);
            println!("  Classes:         {}", lm.num_classes);
            println!("  Total cells:     {}", lm.total_cells);
            println!("  Total area:      {:.1} m\u{b2}", lm.total_area_m2);
            println!("  Processing time: {:.2?}", elapsed);
        }

        LandscapeCommands::Analyze {
            input,
            output_labels,
            output_csv,
            connectivity,
        } => {
            let conn = parse_connectivity(connectivity)?;
            let raster = read_dem(&input)?;
            let start = Instant::now();

            // 1. Label patches
            let (labels, num_patches) =
                label_patches(&raster, conn).context("Failed to label patches")?;
            println!("Patches: {}", num_patches);

            // 2. Patch metrics
            let patches = patch_metrics(&raster, &labels, num_patches)
                .context("Failed to compute patch metrics")?;

            // 3. Class metrics
            let cm = class_metrics(&raster, &patches).context("Failed to compute class metrics")?;

            println!(
                "\n{:<10} {:>10} {:>10} {:>8} {:>8} {:>10}",
                "Class", "Proportion", "Patches", "AI", "Cohesion", "MeanArea"
            );
            println!("{}", "-".repeat(66));
            for c in &cm {
                println!(
                    "{:<10} {:>10.4} {:>8} {:>8.1} {:>10.1} {:>10.1}",
                    c.class, c.proportion, c.num_patches, c.ai, c.cohesion, c.mean_patch_area_m2
                );
            }

            // 4. Landscape metrics
            let lm = landscape_metrics(&raster).context("Failed to compute landscape metrics")?;
            println!(
                "\nLandscape: SHDI={:.4}  SIDI={:.4}  classes={}",
                lm.shdi, lm.sidi, lm.num_classes
            );

            // 5. Write outputs
            if let Some(ref lp) = output_labels {
                write_result_i32(&labels, lp, compress)?;
                println!("Labels saved to: {}", lp.display());
            }
            if let Some(ref cp) = output_csv {
                let csv = patches_to_csv(&patches);
                std::fs::write(cp, &csv).context("Failed to write CSV")?;
                println!("Patch CSV saved to: {}", cp.display());
            }

            let elapsed = start.elapsed();
            println!("\n  Processing time: {:.2?}", elapsed);
        }
    }

    Ok(())
}