diff --git a/ROADMAP.md b/ROADMAP.md index e7291bf16..483933bd6 100644 --- a/ROADMAP.md +++ b/ROADMAP.md @@ -57,7 +57,7 @@ The bugs users hit, the docs and a few quick wins: good first pull requests. - A pass over the user guide, the comparison pages and the README against what the app does today (vertex editing, spatial join and kriging are listed as missing or undocumented; some pages still say 16 languages). -- Styling: standard deviation, geometric interval and manual breaks _(open)_. +- Styling: manual breaks _(open)_. - Selection by polygon and lasso, and a radius tool on the map (Neighborhood Explorer already selects by radius) _(open)_. @@ -187,7 +187,8 @@ Open questions rather than plans; discuss before building. projections), a search box for layers, places and coordinates, bookmarks, web tile layers (XYZ/TMS/WMS, Mapbox Vector Tiles, MapLibre styles, PMTiles and MBTiles archives), a basemap gallery (OpenFreeMap, USGS Imagery), worldviews. -- Cartography: single, unique values, graduated colours and proportional symbols; layer groups; +- Cartography: single, unique values, graduated colours (natural breaks, equal interval, quantile, + geometric interval, standard deviation) and proportional symbols; layer groups; label classes with SQL queries, placement and annotation; custom colour ramps (and QGIS, GDAL and GMT palette files); masking, display filters, blend modes, aggregation, extrusion; graphics. - Editing: create, move, rotate, scale, merge, split, explode, vertex editing; snapping, templates, elevation, diff --git a/apps/rustgis-python/python/rustgis/__init__.py b/apps/rustgis-python/python/rustgis/__init__.py index f0b4796a3..4e68bfe08 100644 --- a/apps/rustgis-python/python/rustgis/__init__.py +++ b/apps/rustgis-python/python/rustgis/__init__.py @@ -157,7 +157,8 @@ def style(self, type: str = "single", field: Optional[str] = None, **params: Any (colours by class) or ``"proportional"`` (symbol sizes by value). field: The field the colours or sizes come from (not for ``"single"``). **params: More ``layer.symbology`` parameters: ``classes`` (1 to 12), ``method`` - (``"naturalBreaks"``, ``"equalInterval"``, ``"quantile"``), ``ramp``. + (``"naturalBreaks"``, ``"equalInterval"``, ``"quantile"``, + ``"geometricInterval"``, ``"standardDeviation"``), ``ramp``. Returns: The layer, so calls can be chained. diff --git a/crates/assistant/skills/thematic-map/SKILL.md b/crates/assistant/skills/thematic-map/SKILL.md index 71ab8f218..8ffd91f3c 100644 --- a/crates/assistant/skills/thematic-map/SKILL.md +++ b/crates/assistant/skills/thematic-map/SKILL.md @@ -9,7 +9,7 @@ description: Style a layer by its data (graduated colours, categories, proportio field suits graduated or proportional symbols; text or a few distinct values suit unique values. 2. Style it with `run_command` `layer.symbology`: - Numbers: `{"layer": L, "type": "graduated", "field": F, "classes": 5, "method": "naturalBreaks"}` (or - `quantile`, `equalInterval`), with an optional `"ramp"`. + `quantile`, `equalInterval`, `geometricInterval`, `standardDeviation`), with an optional `"ramp"`. For counts or sizes, `"type": "proportional"` usually reads better than colour. - Categories: `{"layer": L, "type": "unique", "field": F}`. More than about 12 categories is hard to read, so suggest grouping first. diff --git a/crates/engine/src/cmd/layer.rs b/crates/engine/src/cmd/layer.rs index 7d63561ac..2c6e058d0 100644 --- a/crates/engine/src/cmd/layer.rs +++ b/crates/engine/src/cmd/layer.rs @@ -109,7 +109,7 @@ pub fn commands() -> Vec { .when(need_layer), C::new("layer.symbology", "Symbology", "Style › Drawing", symbology) .when(need_feature_layer) - .params(r#"{"layer"?, "type": "single"|"unique"|"graduated"|"proportional", "field"?, "classes"?: 1..12, "method"?: "naturalBreaks"|"equalInterval"|"quantile", "ramp"?: string}"#), + .params(r#"{"layer"?, "type": "single"|"unique"|"graduated"|"proportional", "field"?, "classes"?: 1..12, "method"?: "naturalBreaks"|"equalInterval"|"quantile"|"geometricInterval"|"standardDeviation", "ramp"?: string}"#), C::new("layer.ramp", "Color Scheme", "Style › Drawing", |s, p| { let r = s.ramp(s_req(p, "ramp")?)?; with(s, p, |l| match &mut l.raster { diff --git a/crates/map/src/symbology.rs b/crates/map/src/symbology.rs index c5cb4dc6a..ae1d78d2a 100644 --- a/crates/map/src/symbology.rs +++ b/crates/map/src/symbology.rs @@ -311,15 +311,27 @@ pub enum ClassMethod { NaturalBreaks, EqualInterval, Quantile, + /// Class widths that grow by one ratio, for values that pile up at the low end. + GeometricInterval, + /// Classes one standard deviation wide about the mean. + StandardDeviation, } impl ClassMethod { - pub const ALL: [ClassMethod; 3] = [ClassMethod::NaturalBreaks, ClassMethod::EqualInterval, ClassMethod::Quantile]; + pub const ALL: [ClassMethod; 5] = [ + ClassMethod::NaturalBreaks, + ClassMethod::EqualInterval, + ClassMethod::Quantile, + ClassMethod::GeometricInterval, + ClassMethod::StandardDeviation, + ]; pub fn label(self) -> &'static str { match self { ClassMethod::NaturalBreaks => "Natural Breaks (Jenks)", ClassMethod::EqualInterval => "Equal Interval", ClassMethod::Quantile => "Quantile", + ClassMethod::GeometricInterval => "Geometric Interval", + ClassMethod::StandardDeviation => "Standard Deviation", } } pub fn from_name(s: &str) -> Option { @@ -328,6 +340,8 @@ impl ClassMethod { "naturalbreaks" | "natural_breaks" | "jenks" | "natural breaks" => ClassMethod::NaturalBreaks, "equalinterval" | "equal_interval" | "equal interval" => ClassMethod::EqualInterval, "quantile" => ClassMethod::Quantile, + "geometricinterval" | "geometric_interval" | "geometric interval" | "geometric" => ClassMethod::GeometricInterval, + "standarddeviation" | "standard_deviation" | "standard deviation" | "stddev" => ClassMethod::StandardDeviation, _ => return None, }) } @@ -1048,8 +1062,50 @@ pub fn class_breaks(sorted: &[f64], k: usize, method: ClassMethod) -> Vec { ClassMethod::EqualInterval => (1..=k).map(|i| lo + (hi - lo) * i as f64 / k as f64).collect::>(), ClassMethod::Quantile => (1..=k).map(|i| sorted.get(((n * i) / k).saturating_sub(1).min(n - 1)).copied().unwrap_or(hi)).collect(), ClassMethod::NaturalBreaks => jenks(sorted, k), + ClassMethod::GeometricInterval => { + // Positive values: each bound is the one before times (hi / lo)^(1/k), worked out on + // logarithms. Values from zero or below are taken as if the minimum sat at a hundredth + // of the range, which makes the ratio 101^(1/k); each bound is then a share of the way + // from lo to hi. Neither adds large values, so the bounds stay finite and in order. + let range = hi - lo; + (1..=k) + .map(|i| { + let f = i as f64 / k as f64; + let v = if lo > 0.0 { + (lo.ln() + (hi.ln() - lo.ln()) * f).exp() + } else { + let t = (101f64.powf(f) - 1.0) / 100.0; + if range.is_finite() { lo + range * t } else { lo * (1.0 - t) + hi * t } + }; + v.clamp(lo, hi) + }) + .collect() + } + ClassMethod::StandardDeviation => { + // On values scaled by the largest magnitude, so the squares can't overflow. + let s = lo.abs().max(hi.abs()); + let mean_s = sorted.iter().map(|v| v / s).sum::() / n as f64; + let mean = mean_s * s; + let sd = s * (sorted.iter().map(|v| (v / s - mean_s).powi(2)).sum::() / n as f64).sqrt(); + // k − 1 inner bounds one deviation apart, centred on the mean. A bound with no value + // above the one before it is dropped, so no class is empty. + let inner = (1..k).map(|i| mean + sd * (i as f64 - k as f64 / 2.0)).filter(|&v| v > lo && v < hi); + let mut b = Vec::new(); + let mut below = 0; + for v in inner.chain(std::iter::once(hi)) { + let upto = sorted.partition_point(|x| *x <= v); + if upto > below { + b.push(v); + below = upto; + } + } + b + } }; - b.dedup_by(|a, b| (*a - *b).abs() < f64::EPSILON); + // Only bounds that are the same value: any tolerance merges distinct values somewhere (tiny + // ones under an absolute one, large ones under a relative one). Repeated bounds come from + // repeated values (quantiles, Jenks), so they are equal exactly. + b.dedup(); if let Some(l) = b.last_mut() { *l = hi; } @@ -1748,6 +1804,86 @@ mod tests { assert_eq!(class_breaks(&[5.0, 5.0], 4, ClassMethod::NaturalBreaks), vec![5.0]); } + fn close(a: &[f64], b: &[f64]) -> bool { + a.len() == b.len() && a.iter().zip(b).all(|(x, y)| (x - y).abs() < 1e-9 * y.abs().max(1.0)) + } + + #[test] + fn geometric_interval_breaks_grow_by_one_ratio() { + // Positive values: each bound is the one before times (hi / lo)^(1/k). + let v = [1.0, 5.0, 40.0, 300.0, 1000.0]; + let b = class_breaks(&v, 3, ClassMethod::GeometricInterval); + assert!(close(&b, &[10.0, 100.0, 1000.0]), "{b:?}"); + // From zero (or below), the values are shifted so the minimum sits at a hundredth of the + // range: the widths still grow by one ratio, and the last bound is the maximum. + let b = class_breaks(&[-50.0, 0.0, 50.0], 2, ClassMethod::GeometricInterval); + assert!(close(&b, &[-51.0 + 101f64.sqrt(), 50.0]), "{b:?}"); + let b = class_breaks(&[0.0, 100.0], 4, ClassMethod::GeometricInterval); + let widths: Vec = std::iter::once(b[0]).chain(b.windows(2).map(|w| w[1] - w[0])).collect(); + assert!(widths.windows(2).all(|w| w[1] > w[0]), "{b:?}"); + assert_eq!(b.len(), 4); + assert_eq!(b.last(), Some(&100.0)); + // Far from zero the shift is not added and taken away again, so it can't round off. + let b = class_breaks(&[-1e16, -9_999_999_999_999_998.0], 2, ClassMethod::GeometricInterval); + assert!(b.iter().all(|v| v.is_finite()), "{b:?}"); + assert_eq!(b.last(), Some(&-9_999_999_999_999_998.0)); + // Near the largest float nothing overflows: the bounds stay finite and in order. + for v in [[1e308, 1.1e308], [-1e308, 1e308], [1e-308, 1e308]] { + let b = class_breaks(&v, 2, ClassMethod::GeometricInterval); + assert!(b.len() == 2 && b.iter().all(|x| x.is_finite()) && b[0] >= v[0] && b[0] < b[1], "{v:?}: {b:?}"); + } + } + + #[test] + fn standard_deviation_breaks_are_one_deviation_wide_about_the_mean() { + // 0..=100: mean 50, population deviation sqrt(850). + let v: Vec = (0..=100).map(f64::from).collect(); + let sd = 850f64.sqrt(); + // An even count puts a bound on the mean; an odd one centres a class on it. + let b = class_breaks(&v, 4, ClassMethod::StandardDeviation); + assert!(close(&b, &[50.0 - sd, 50.0, 50.0 + sd, 100.0]), "{b:?}"); + let b = class_breaks(&v, 3, ClassMethod::StandardDeviation); + assert!(close(&b, &[50.0 - sd / 2.0, 50.0 + sd / 2.0, 100.0]), "{b:?}"); + // Bounds past the data are dropped rather than leaving empty classes: 8 classes of one + // deviation would reach ±4 deviations, but the data spans less than ±2. + let b = class_breaks(&v, 8, ClassMethod::StandardDeviation); + assert!(close(&b, &[50.0 - sd, 50.0, 50.0 + sd, 100.0]), "{b:?}"); + // A bound with no value above the one before it goes too: mean 4, deviation sqrt(24), so + // (4, 4 + sqrt(24)] would hold nothing between the zeros and the tens. + let b = class_breaks(&[0.0, 0.0, 0.0, 10.0, 10.0], 4, ClassMethod::StandardDeviation); + assert!(close(&b, &[4.0, 10.0]), "{b:?}"); + // Large values: the squared deviations would overflow without scaling. + let b = class_breaks(&[1e160, 2e160], 2, ClassMethod::StandardDeviation); + assert!(close(&b, &[1.5e160, 2e160]), "{b:?}"); + let b = class_breaks(&[-1e308, 1e308], 2, ClassMethod::StandardDeviation); + assert!(close(&b, &[0.0, 1e308]), "{b:?}"); + } + + #[test] + fn close_values_keep_their_classes() { + // Bounds closer than f64::EPSILON in absolute terms are still distinct when the values + // themselves are that small. + for m in ClassMethod::ALL { + let b = class_breaks(&[0.0, 1e-18], 2, m); + assert!(b.len() == 2 && b[0] >= 0.0 && b[0] < b[1], "{m:?}: {b:?}"); + } + // And close values far from zero keep theirs too. + let b = class_breaks(&[1e16, 1e16 + 2.0], 2, ClassMethod::Quantile); + assert_eq!(b, vec![1e16, 1e16 + 2.0]); + } + + #[test] + fn every_method_is_listed_named_and_read_back() { + for m in ClassMethod::ALL { + let name = format!("{m:?}"); + let camel = format!("{}{}", name[..1].to_ascii_lowercase(), &name[1..]); + assert_eq!(ClassMethod::from_name(&camel), Some(m), "{camel}"); + assert_eq!(ClassMethod::from_name(&m.label().replace(" (Jenks)", "")), Some(m), "{}", m.label()); + } + assert_eq!(ClassMethod::from_name("stddev"), Some(ClassMethod::StandardDeviation)); + assert_eq!(ClassMethod::from_name("geometric"), Some(ClassMethod::GeometricInterval)); + } + #[test] fn graduated_symbols_by_value() { let mut d = Dataset::new("t", GeometryType::Polygon, vec![Field::new("v", FieldType::Double)]); diff --git a/docs/user-guide/cartography.md b/docs/user-guide/cartography.md index b2388998b..b0ff86a11 100644 --- a/docs/user-guide/cartography.md +++ b/docs/user-guide/cartography.md @@ -18,8 +18,13 @@ Select a layer in Layers and open **Style › Symbology** (or the Symbology pane is not fully transparent. A layout legend can leave it out with **Hide other values** (see [Layouts](layouts.md)). - **Graduated Colors**: 1 to 12 classes of a numeric field by natural breaks (Jenks), equal - interval or quantile. If no field is chosen, the current renderer's field is kept when it is - numeric; otherwise the first numeric field is used. If the layer has no numeric field, or the + interval, quantile, geometric interval or standard deviation. Geometric interval makes class + widths grow by one ratio, for values piled up at the low end; for a field that reaches zero or + below, the minimum is taken as if it sat at a hundredth of the range, so the widths grow by + 101^(1/classes). Standard deviation makes classes one deviation wide about the mean, and leaves + out a class no value falls in, so a skewed field can get fewer classes than asked. If no field is + chosen, the current renderer's field is kept when it is numeric; otherwise the first numeric + field is used. If the layer has no numeric field, or the switch fails for another reason, the pane shows why and keeps the current renderer. When numeric fields exist, a field picker lets you try another one; the error stays visible if a retry fails. A problem no field can fix, such as a definition query that keeps no features, needs fixing