Skip to content
Merged
Show file tree
Hide file tree
Changes from 6 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions crates/analysis/src/network.rs
Original file line number Diff line number Diff line change
Expand Up @@ -537,15 +537,15 @@ impl Steps {
}

/// The longest step, in degrees of longitude or latitude, that [`polar_reach`] follows an edge by.
const REACH_STEP_DEG: f64 = 0.5;
pub const REACH_STEP_DEG: f64 = 0.5;
/// The most steps it takes along one edge; one needing more (over 512° long, as only a hostile
/// edge is) counts as near as its box.
const REACH_STEPS_MAX: f64 = 1024.0;
pub const REACH_STEPS_MAX: f64 = 1024.0;
/// The most steps one search ([`within`]) takes along all its shapes' edges, a point measured
/// counting as one; past them, a point or edge counts as near as its box
/// ([`least_distance_to_box`]), so a shape may be counted in that's farther than the radius, never
/// left out for want of steps.
const REACH_STEPS_TOTAL: usize = 250_000;
pub const REACH_STEPS_TOTAL: usize = 250_000;
/// More than a degree of longitude or latitude is long anywhere on the ground (meters): a step's
/// length is at most this times the degrees it spans in each.
const DEGREE_M_MAX: f64 = 111_700.0;
Expand Down
63 changes: 52 additions & 11 deletions crates/analysis/src/tools.rs
Original file line number Diff line number Diff line change
Expand Up @@ -98,7 +98,7 @@ pub static TOOLS: &[ToolSpec] = &[
id: "repairGeometry",
name: "Repair Geometry",
toolbox: "Data Management › Features",
summary: "Writes every feature with its geometry made valid: rings closed, repeated vertices and degenerate parts dropped, self-intersections, overlapping parts and holes outside their shell resolved. Says what was repaired in each.",
summary: "Writes every feature with its geometry made valid where it can be: rings closed, repeated vertices and degenerate parts dropped, self-intersections, overlapping parts and holes outside their shell resolved. Features with a missing or non-finite coordinate anywhere, even in a hole, or left with no usable points or, for polygons, no area once repaired, are written without geometry; features too complex to check, or whose repair fails, are written unchanged. Says what was repaired in each.",
Comment thread
coderabbitai[bot] marked this conversation as resolved.
Resolved
params: &[
P { name: "input", label: "Input Features", kind: K::VectorLayer, required: true, default: "", help: "" },
P { name: "output", label: "Output Name", kind: K::Text, required: false, default: "", help: "" },
Expand Down Expand Up @@ -1512,7 +1512,7 @@ fn buffer_with_limits(inp: &dyn Inputs, p: &Json, limit: usize, bytes: usize) ->
}
if coarse > 0 {
r.messages.push(format!(
"{coarse} of {} features' buffers have edges left longer than their outline allows (a ring of thousands of long edges uses up the points added to it): their sides there can stray from the distance by more than a buffer otherwise keeps to (0.012%, 0.05% along a line or outline buffered in stretches, or 1 cm)",
"{coarse} of {} features' buffers have outline edges left unsplit, with no points added along them (a ring of thousands of long edges uses up the points it may add): their sides there can stray from the distance by more than a buffer otherwise keeps to (0.012%, 0.05% along a line or outline buffered in stretches, or 1 cm)",
Comment thread
giswqs marked this conversation as resolved.
Resolved
picked.len()
));
}
Expand Down Expand Up @@ -2342,8 +2342,9 @@ fn size_output_points<'a>(
for f in feats {
let Some(g) = &f.geometry else { continue };
let k = count(g);
// Saturating, as sizes from input are. A row repeated in the selection (a project can
// restore one) is counted again, and its points are written twice.
// Saturating, as sizes from input are. A row repeated in the selection (only a custom
// `Inputs` can give one: a project's selection is normalized on load) is counted again,
// and its points are written twice.
n = n.saturating_add(k);
if estimate {
let own = f.attrs.get(..kept).unwrap_or(&f.attrs);
Expand Down Expand Up @@ -3321,7 +3322,7 @@ fn select_by_location_with_work_limit(inp: &dyn Inputs, p: &Json, limit: usize,
let coarse = coarse.into_inner();
if coarse > 0 {
messages.push(format!(
"{coarse} selecting features' buffers have edges left longer than their outline allows (a ring of thousands of long edges uses up the points added to it): features near their edges there may be selected or not by more than a buffer otherwise keeps to (0.012% of the distance, 0.05% along a line or outline buffered in stretches, or 1 cm)"
"{coarse} selecting features' buffers have outline edges left unsplit, with no points added along them (a ring of thousands of long edges uses up the points it may add): features near their edges there may be selected or not by more than a buffer otherwise keeps to (0.012% of the distance, 0.05% along a line or outline buffered in stretches, or 1 cm)"
));
}
Ok(ToolResult { output: Output::Selection { layer: ds.name.clone(), rows: hits }, messages })
Expand Down Expand Up @@ -3630,6 +3631,9 @@ fn kernel_density_tool(inp: &dyn Inputs, p: &Json) -> Result<ToolResult, String>
// which would land on the map's edge rather than where it is (#1033): its latitude doesn't
// come back.
let mut beyond: usize = 0;
// A multipoint's members with non-finite coordinates, in the features whose other members
// were used (a feature with none usable is counted in `skipped`).
let (mut non_finite, mut non_finite_features) = (0usize, 0usize);
let project = |c: Coord, beyond: &mut usize| {
// Web Mercator's limit is known exactly: the round trip alone would let a latitude a hair
// past it through, onto the edge.
Expand Down Expand Up @@ -3666,17 +3670,20 @@ fn kernel_density_tool(inp: &dyn Inputs, p: &Json) -> Result<ToolResult, String>
}
}
Geometry::MultiPoint(cs) => {
let mut any = false;
let (mut any, mut bad) = (false, 0usize);
for c in cs {
if c.is_finite()
&& let Some(p_proj) = project(*c, &mut beyond)
{
if !c.is_finite() {
bad += 1;
} else if let Some(p_proj) = project(*c, &mut beyond) {
pts.push((p_proj, w));
any = true;
}
}
if !any {
skipped += 1;
} else if bad > 0 {
non_finite = non_finite.saturating_add(bad);
non_finite_features += 1;
}
}
_ => {
Expand All @@ -3688,6 +3695,9 @@ fn kernel_density_tool(inp: &dyn Inputs, p: &Json) -> Result<ToolResult, String>
// Only Web Mercator's limit is known to be a latitude.
let why = if crs == Crs::WebMercator { " (past 85.05° of latitude)" } else { "" };
let beyond_note = (beyond > 0).then(|| format!("{beyond} point(s) beyond what {} can hold{why} were left out", crs.label()));
let non_finite_note = (non_finite > 0).then(|| {
format!("{non_finite} point(s) with non-finite coordinates were left out of {non_finite_features} multipoint feature(s) whose other points were used")
});
if pts.is_empty() {
let none = "input layer has no valid points with positive weights to analyze";
return Err(beyond_note.map_or_else(|| none.to_string(), |note| format!("{none}: {note}")));
Expand Down Expand Up @@ -3753,7 +3763,7 @@ fn kernel_density_tool(inp: &dyn Inputs, p: &Json) -> Result<ToolResult, String>
};
Ok(ToolResult {
output: Output::Layers(vec![OutputLayer::Raster { file_name: format!("{name}.tif"), name, bytes }]),
messages: std::iter::once(msg).chain(beyond_note).collect(),
messages: std::iter::once(msg).chain(beyond_note).chain(non_finite_note).collect(),
})
}

Expand Down Expand Up @@ -8491,7 +8501,7 @@ mod tests {
for dissolve in ["None", "All"] {
let params = serde_json::json!({"input": "l", "distance": 1, "units": "Meters", "dissolve": dissolve});
let res = crate::run("buffer", &m, &params).unwrap_or_else(|e| panic!("{e}"));
let said = res.messages.iter().any(|m| m.starts_with("1 of 1 features' buffers have edges left longer"));
let said = res.messages.iter().any(|m| m.starts_with("1 of 1 features' buffers have outline edges left unsplit"));
assert_eq!(said, coarse, "dissolve {dissolve}: {:?}", res.messages);
}
}
Expand Down Expand Up @@ -8744,6 +8754,37 @@ mod tests {
assert!(note.starts_with("1 point(s) beyond what") && !note.contains("85.05°"), "{:?}", res.messages);
}

/// A multipoint's non-finite members are left out and counted when its other members are
/// used; a clean multipoint says nothing of the kind.
#[test]
fn kernel_density_counts_a_multipoints_non_finite_members() {
let run = |members: Vec<Coord>| {
let mut pts = Dataset::new("pts", GeometryType::Point, vec![]);
pts.push(Some(Geometry::Point(Coord::new(0.0, 60.0))), vec![]);
pts.push(Some(Geometry::MultiPoint(members)), vec![]);
let m = Map(HashMap::from([("pts".to_string(), pts)]));
let params = serde_json::json!({"input": "pts", "crs": "EPSG:3857", "radius": 5000.0, "pixelSize": 1000.0});
crate::run("kernelDensity", &m, &params).map(|r| r.messages).unwrap_or_else(|e| vec![e])
};
let lossy = run(vec![Coord::new(0.01, 60.0), Coord::new(f64::NAN, 60.0), Coord::new(0.0, f64::INFINITY)]);
assert!(lossy.iter().any(|m| m.contains("from 2 points (bandwidth")), "{lossy:?}");
Comment thread
coderabbitai[bot] marked this conversation as resolved.
Resolved
assert!(lossy.iter().all(|m| !m.contains("beyond what")), "{lossy:?}");
let notes: Vec<_> = lossy.iter().filter(|m| m.contains("non-finite")).collect();
assert!(
notes.len() == 1
&& notes.iter().all(|m| m.starts_with("2 point(s) with non-finite coordinates were left out of 1 multipoint feature(s)")),
"{lossy:?}"
Comment thread
giswqs marked this conversation as resolved.
Resolved
);
let unusable = run(vec![Coord::new(f64::NAN, 60.0), Coord::new(0.0, f64::INFINITY)]);
assert!(unusable.iter().any(|m| m.contains("from 1 points (1 skipped, bandwidth")), "{unusable:?}");
Comment thread
coderabbitai[bot] marked this conversation as resolved.
Resolved
assert!(unusable.iter().all(|m| !m.contains("were left out of")), "{unusable:?}");
let clean = run(vec![Coord::new(0.01, 60.0), Coord::new(0.02, 60.0)]);
assert!(
clean.iter().any(|m| m.contains("from 3 points (bandwidth")) && clean.iter().all(|m| !m.contains("non-finite")),
"{clean:?}"
);
}

#[test]
fn kernel_density_analysis_tool_tests() {
// Create point dataset with UTM projection
Expand Down
6 changes: 5 additions & 1 deletion crates/assistant/skills/geoprocessing/SKILL.md
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,11 @@ Based on the RustGIS user guide, Analysis › Geoprocessing tools
`featureEnvelopeToPolygon`, `unsplitLine` (`input`, optional `output`; see the user guide
for envelope axes and endpoint-joining rules); `checkGeometry` (`input`: a layer of the
features with invalid geometry, saying what is wrong with each); `repairGeometry` (`input`:
every feature with its geometry made valid, and a field saying what was repaired)
every feature with its geometry made valid where it can be, and a field saying what was
repaired; features with a missing or non-finite coordinate, even in a hole, or left with no
usable points or, for polygons, no area once repaired are written without geometry; valid
points and lines are kept; features too complex to check, or whose repair fails, are
written unchanged)
- selection: `selectByLocation` (`input`, `select`, `distance`, `invert`), `selectByAttributes`
(`input`, `where`)
- statistics: `countPoints` (Summarize Within: `polygons`, `points`, `statistics`: summary
Expand Down
22 changes: 21 additions & 1 deletion crates/engine/src/cmd/select.rs
Original file line number Diff line number Diff line change
Expand Up @@ -269,8 +269,10 @@ pub(crate) fn warn_bounded(s: &mut Session, bounded: &[(String, usize)]) {
}
let total: usize = bounded.iter().map(|(_, n)| n).sum();
let layers: Vec<String> = bounded.iter().map(|(layer, n)| format!("{layer} {n}")).collect();
let steps = rustgis_map::graphics::thousands(network::REACH_STEPS_TOTAL as f64);
let longest = network::REACH_STEPS_MAX * network::REACH_STEP_DEG;
s.warn(format!(
"{total} feature(s) near a pole ({}) couldn't be measured in full (past a search's 250,000 steps a layer, or along an edge over 512° long) and were counted at the least distance they could be: they may lie farther, even beyond the radius",
"{total} feature(s) near a pole ({}) couldn't be measured in full (past a search's {steps} steps a layer, or along an edge over {longest}° long) and were counted at the least distance they could be: they may lie farther, even beyond the radius",
layers.join(", ")
));
}
Expand Down Expand Up @@ -512,6 +514,24 @@ mod tests {
use rustgis_data::{Dataset, Field, FieldType, Geometry, GeometryType, Value as DataValue};
use rustgis_geom::{lonlat_to_merc, merc_to_lonlat};

/// The pole warning states the search's limits from the constants that set them.
#[test]
fn bounded_warning_states_the_search_limits() {
let mut s = Session::new(Project::default());
let mark = s.notices_given();
warn_bounded(&mut s, &[("Cities".into(), 2)]);
let text = s.notices_since(mark).first().map(|n| n.text.as_str()).unwrap_or_default();
assert!(
text.contains("2 feature(s) near a pole (Cities 2)")
&& text.contains(&format!(
"past a search's {} steps a layer",
rustgis_map::graphics::thousands(network::REACH_STEPS_TOTAL as f64)
))
Comment thread
coderabbitai[bot] marked this conversation as resolved.
Resolved
&& text.contains(&format!("along an edge over {}° long", network::REACH_STEPS_MAX * network::REACH_STEP_DEG)),
"{text}"
);
}

fn fixture() -> (Session, u32, u32) {
let mut s = Session::new(Project::default());
let mut ds =
Expand Down
14 changes: 11 additions & 3 deletions docs/user-guide/analysis.md
Original file line number Diff line number Diff line change
Expand Up @@ -455,7 +455,9 @@ gap between two parts never counts; a line's parts that cross the antimeridian,
where the next begins on its other side, as a GeoJSON file splits a line that crosses it, are one
line. A ring ends where it starts, so a point there is made once, and a line that ends where it
starts (a contour) is walked as a ring too; a part shorter than its offsets gets none. A line that
passes through one place twice gets a point there each time the distance lands on it. Each point
passes through one place twice gets a point there each time the distance lands on it, except at a
closed line's shared start and end, which gets one point (a line through A, B and back to A gets
one at A). Each point
has its feature's attributes and three fields: `orig_fid`, the feature's OID, as Feature Vertices
To Points writes it, and, named as in QGIS's Points along geometry, `distance` (how far along its
part, in `units`, the start offset included) and `angle` (the line's direction there, in degrees
Expand Down Expand Up @@ -936,7 +938,13 @@ CRS metadata. Coordinates are used directly, without projection or snapping.
- **Unsplit Line** (`unsplitLine`) joins line parts at exactly equal endpoints where precisely
two part ends meet and all attribute values agree. Reversed parts and closed cycles are
supported. Branches, attribute changes, and intersections away from endpoints stay separate.
Every output chain retains its common attributes; no field values are discarded.
Every output chain retains its common attributes; no field values are discarded. Each chain of
joined parts is its own output feature, so a multipart input can become several features.
Joins are decided among the selected features only, when there is a selection: an unselected
part meeting an endpoint neither joins it nor stops a join there. Two open parts with
equal attributes that share both endpoints join into one closed line when no other considered
part ends there; a part that is already closed has two ends at its endpoint, so it joins
nothing there.

Missing geometries and empty or too-short parts are skipped with a count. If skipping leaves
no usable output, the tool returns an error. Non-finite coordinates and malformed rows with
Expand Down Expand Up @@ -966,7 +974,7 @@ features were checked.
| Code | What is wrong |
| --- | --- |
| `null_geometry` | No geometry, or one without coordinates |
| `non_finite_coordinate` | A coordinate is not a number or infinite. The row is written without geometry, so it can't spoil the layer's extent |
| `non_finite_coordinate` | A coordinate is not a number or infinite. The row is written without geometry, so it can't spoil the layer's extent. A ring with one is reported for those coordinates and left out of every other check: its points, its edges' crossings and its containment in other rings |
| `too_few_points` | A line part with fewer than 2 distinct points, a ring with fewer than 3, or a polygon part without rings |
| `ring_not_closed` | A ring's last point is not its first |
| `repeated_vertex` | The same vertex twice in a row (a closed ring's last point doesn't count) |
Expand Down
Loading