diff --git a/doc/syntax/coord/map.qmd b/doc/syntax/coord/map.qmd new file mode 100644 index 000000000..e800aebb3 --- /dev/null +++ b/doc/syntax/coord/map.qmd @@ -0,0 +1,181 @@ +--- +title: Map +--- + +The map coordinate system facilitates the display of geographical data. +It can project data to various coordinate reference systems or projections. +It is used to display choropleths and other types of maps. + +## Default aesthetics + +* **Primary**: `lon` for longitude (horizontal) +* **Secondary**: `lat` for latitude (vertical) + +Users can provide their own aesthetic names if needed. +For example, if using `x` and `y` aesthetics: + +```ggsql +PROJECT x, y TO map +``` + +This maps `x` to longitude and `y` to latitude. +This is a convenience when converting from a Cartesian coordinate system, without having to modify all the mappings. + +## Settings + +* `clip`: Should data be removed if it appears outside the bounds of the coordinate system, or has become invalid due to the projection. + Defaults to `true`. +* `crs` A string giving the target Coordinate Reference System (CRS) holding parameters of a coordinate transformation. + This is only used with unnamed projections and therefore mutually exclusive with the `origin` and `parallel` settings. + By default, matches the Spatial Reference Identifiers (SRIDs) of geometry columns. + Typically takes one of the following forms: + * A [proj string](https://proj.org/en/stable/usage/quickstart.html) like `'+proj=longlat +lon_0=90'` for 2D longitude/latitude with 90° longitude offset. + * An [EPSG code](https://en.wikipedia.org/wiki/EPSG_Geodetic_Parameter_Dataset) like `'EPSG:3395'` for World Mercator projection. +* `source` A string giving the source CRS for the data to use as fallback when not explicit in geometry. + Useful when using non-spatial layers that don't carry SRIDs. + Defaults to `'EPSG:4326'` or longitude/latitude coordinates. +* `bounds` A numeric array of length 4 giving the bounding box of the area to view. + The numbers mean ('xmin', 'ymin', 'xmax', 'ymax') in that order, in the units of the `crs` argument. + Be mindful that your data might be in longitude/latitude degrees, but the `crs` might use meters, or vice versa. + Bounds can have `Inf`/`-Inf` values to indicate the visible area's extent. + Using `null` values indicates to fall back to the data's extent, so only partial bounds need to be indicated. +* `origin` Projection origin location. + This is only used with named projections and therefore mutually exclusive with `crs`. + One of the following: + * A scalar number, setting the origin longitude or 'central meridian' in degrees between −180 and +180. + Corresponds to the `+lon_0` proj string parameter. + * An array of two numbers, where the first number is the scalar option above. + The second number sets the latitude of origin in degrees between −90 and +90. + The latter corresponds to the `+lat_0` proj string parameter. + If unsupported by the projection, the latitude is silently ignored. +* `parallel` Setting for conic or cylindrical projections. + This is only used with named projections and therefore mutually exclusive with `crs`. + If unsupported by the projection, this setting is silently ignored. + One of the following: + * A scalar number, setting the standard parallel as degrees between −90 and +90. + Corresponds to the `+lat_1` proj string parameter. + * An array of two numbers, where the first number is the scalar option above. + The second number sets the second standard parallel as degrees between −90 and +90. + +## Supported projections + +There are two main ways to set the map projection. +One is via a named projection, which allows for the `origin` and `parallel` settings. + +```ggsql +PROJECT TO orthographic SETTING origin => (150, 52) +``` + +The second option is to use the generic `PROJECT TO map` together with the `crs` setting. +The code below is equivalent to the code above. + +```ggsql +PROJECT TO map SETTING crs => '+proj=ortho +lon_0=150 +lat_0=52' +``` + +In this context, 'supported' means that we are reasonably confident that we can draw a serviceable map. +Below follows a table giving an overview of the supported projections. +The projections are recognised by name in the 'string' column or from the 'proj string' CRS setting. +More information about them can be found on the [PRØJ website](https://proj.org/en/stable/operations/projections/index.html). + +| name | string | origin | parallel | proj string | +|---------------------------------------|-------------------------|--------|--------------|------------------| +| Geographic (unprojected) | `geographic` | 0 | NA | `+proj=longlat` | +| Mercator | `mercator` | 0 | NA | `+proj=merc` | +| Orthographic | `orthographic` | (0, 0) | NA | `+proj=ortho` | +| Miller Cylindrical | `miller` | 0 | NA | `+proj=mill` | +| Equidistant Cylindrical (Plate Carée) | `equirectangular` | 0 | NA | `+proj=eqc` | +| Stereographic | `stereographic` | (0, 0) | NA | `+proj=stere` | +| Gnomonic | `gnomonic` | (0, 0) | NA | `+proj=gnom` | +| Cylindrical Equal Area | `equal_area` | 0 | NA | `+proj=cea` | +| Mollweide | `mollweide` | 0 | NA | `+proj=moll` | +| Sanson-Flamsteed Sinusoidal | `sinusoidal` | 0 | NA | `+proj=sinu` | +| Eckert IV | `eckert4` | 0 | NA | `+proj=eck4` | +| Natural Earth | `natural` | 0 | NA | `+proj=natearth` | +| Winkel Tripel | `winkel_tripel` | 0 | NA | `+proj=wintri` | +| Albers Equal Area | `albers` | (0, 0) | (29.5, 45.5) | `+proj=aea` | +| Lambert Conformal Conic | `lambert_conformal` | (0, 0) | (29.5, 45.5) | `+proj=lcc` | +| Lambert Azimuthal | `lambert` | (0, 0) | NA | `+proj=laea` | +| Azimuthal Equidistant | `azimuthal_equidistant` | (0, 0) | NA | `+proj=aeqd` | +| Interrupted Goode Homolosine | `igh` | 0 | NA | `+proj=igh` | +| Robinson | `robinson` | 0 | NA | `+proj=robin` | + +The way to draw an unsupported map is to use the `crs` setting. +Unsupported maps are likely drawn without projection-appropriate clipping. + +## Examples + +Note that depending on your reader, you may need to activate modules for spatial analysis. + +```{ggsql} +-- For example, for DuckDB, one could use: +INSTALL spatial; +LOAD spatial; +``` + +### Via named projection + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial +PROJECT TO robinson +``` + +### Via proj string + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial +PROJECT TO map SETTING crs => '+proj=robin' +``` + +### Setting the origin + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial +PROJECT TO orthographic + SETTING origin => (133.77, -25.27) +``` + +### Zooming in + +The default zoom behaviour is to zoom in where the data is. + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial + FILTER subregion == 'Western Africa' +PROJECT TO orthographic +``` + +An alternative is to zoom using the `bounds` setting. +Note that while `ggsql:world` is defined in degrees longitude/latitude, the `orthographic` projection uses metres as unit. + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial +PROJECT TO orthographic + SETTING bounds => (-1868152, 468481, 1638882, 2917218) +``` + +The `bounds` arguments can take infinites to take the visible area's extreme values, or `null` to take the data's bounds. + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial + FILTER subregion == 'Western Africa' +PROJECT TO orthographic + SETTING bounds => (-1868152, -Inf, 1638882, null) +``` + +### Parallels + +Default conic projections aren't always kind to the southern hemisphere. +This can be changed by tweaking the `parallel` setting. + +```{ggsql} +VISUALISE continent AS fill FROM ggsql:world +DRAW spatial +PROJECT TO albers SETTING parallel => (-29.5, -45.5) +``` \ No newline at end of file diff --git a/doc/syntax/coord/polar.qmd b/doc/syntax/coord/polar.qmd index 24e8af75d..a17e97f30 100644 --- a/doc/syntax/coord/polar.qmd +++ b/doc/syntax/coord/polar.qmd @@ -16,7 +16,7 @@ Users can provide their own aesthetic names if needed. For example, if using `x` PROJECT y, x TO polar ``` -This maps `y` to radius and `x` to angle. This is useful when converting from a cartesian coordinate system without editing all the mappings. +This maps `y` to radius and `x` to angle. This is useful when converting from a Cartesian coordinate system without editing all the mappings. ## Settings * `clip`: Should data be removed if it appears outside the bounds of the coordinate system. Defaults to `true` diff --git a/src/CLAUDE.md b/src/CLAUDE.md index 5b3152792..0e354b6ce 100644 --- a/src/CLAUDE.md +++ b/src/CLAUDE.md @@ -46,11 +46,10 @@ Grammar lives in [`/tree-sitter-ggsql/`](../tree-sitter-ggsql/) — when adding | `duckdb.rs` | DuckDB (in-memory or file) | `duckdb` (default) | | `sqlite.rs` | SQLite | `sqlite` (default) | | `odbc.rs` | ODBC | `odbc` (default) | -| `snowflake.rs` | Snowflake (very small placeholder today) | — | | `connection.rs` | Connection-string parsing for all of the above | — | | `data.rs`, `spec.rs` | `Spec` type returned by `execute()`, plus DataFrame conversion | — | -`SqlDialect` trait in `mod.rs` lets each driver supply its own type names and information-schema queries. PostgreSQL has a feature flag (`postgres`) but no driver file yet. +`SqlDialect` trait in `mod.rs` lets each driver supply its own type names, information-schema queries, and spatial helper methods (`sql_st_transform`, `sql_geometry_to_wkb`, `sql_geometry_bbox`, `sql_spatial_setup`). ### `execute/` @@ -97,15 +96,12 @@ Defined in `Cargo.toml`: | `sqlite` | ✓ | SQLite reader | | `odbc` | ✓ | ODBC reader | | `parquet` | ✓ | Parquet support in readers/data | -| `ipc` | ✓ | Arrow IPC support | +| `spatial` | ✓ | Spatial/geometry support (geozero for WKT↔GeoJSON) | | `vegalite` | ✓ | Vega-Lite writer | | `builtin-data` | ✓ | Bundled penguins/airquality datasets | -| `postgres` | — | PostgreSQL reader (declared, not yet implemented) | -| `ggplot2` | — | ggplot2 writer (declared, not yet implemented) | -| `all-readers` | — | `duckdb` + `postgres` + `sqlite` + `odbc` | -| `all-writers` | — | `vegalite` + `ggplot2` + `plotters` | +| `all-readers` | — | `duckdb` + `sqlite` + `odbc` | -`ggsql-wasm` builds with `default-features = false` plus `vegalite`, `sqlite`, `builtin-data`. `ggsql-jupyter` builds with `duckdb`, `sqlite`, `odbc`, `vegalite`. +`ggsql-wasm` builds with `default-features = false` plus `vegalite`, `sqlite`, `builtin-data`. `ggsql-jupyter` builds with `duckdb`, `vegalite`. ## Testing diff --git a/src/Cargo.toml b/src/Cargo.toml index 264ae7061..401da7134 100644 --- a/src/Cargo.toml +++ b/src/Cargo.toml @@ -36,7 +36,7 @@ bytes = { workspace = true } adbc_core = { version = "0.23", optional = true } # Spatial -geozero = { workspace = true, optional = true, features = ["with-wkb", "with-geojson"] } +geozero = { workspace = true, optional = true, features = ["with-wkb", "with-wkt", "with-geojson"] } # Serialization serde.workspace = true diff --git a/src/execute/mod.rs b/src/execute/mod.rs index c26f7cd0d..af734d3a4 100644 --- a/src/execute/mod.rs +++ b/src/execute/mod.rs @@ -1368,6 +1368,19 @@ pub fn prepare_data_with_reader(query: &str, reader: &dyn Reader) -> Result = HashMap::new(); for (idx, q) in layer_queries.iter().enumerate() { diff --git a/src/execute/scale.rs b/src/execute/scale.rs index 7d1b4614a..18bd0b1c8 100644 --- a/src/execute/scale.rs +++ b/src/execute/scale.rs @@ -1723,6 +1723,8 @@ mod tests { coord, aesthetics, properties: std::collections::HashMap::new(), + map_projection: None, + computed: std::collections::HashMap::new(), }); // Create scale for pos2 (theta in polar) without explicit expand diff --git a/src/lib.rs b/src/lib.rs index a2ec401f6..43d2bf0b0 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -1029,7 +1029,12 @@ mod integration_tests { let json_str = writer.write(&prepared.specs[0], &prepared.data).unwrap(); let vl_spec: serde_json::Value = serde_json::from_str(&json_str).unwrap(); - assert_eq!(vl_spec["layer"][0]["mark"]["type"], "geoshape"); + let layers = vl_spec["layer"].as_array().unwrap(); + let geoshape_layer = layers + .iter() + .find(|l| l["mark"]["type"] == "geoshape") + .expect("should have a geoshape layer"); + assert_eq!(geoshape_layer["mark"]["type"], "geoshape"); let data = vl_spec["data"]["values"].as_array().unwrap(); let layer_key = prepared.specs[0].layers[0].data_key.as_ref().unwrap(); @@ -1044,6 +1049,83 @@ mod integration_tests { assert_eq!(feature["geometry"]["type"], "Polygon"); } + #[cfg(feature = "spatial")] + #[test] + fn test_end_to_end_spatial_world_orthographic() { + let reader = DuckDBReader::from_connection_string("duckdb://memory").unwrap(); + + let query = r#" + VISUALISE FROM ggsql:world + DRAW spatial PROJECT TO orthographic + "#; + + let prepared = execute::prepare_data_with_reader(query, &reader).unwrap(); + + let writer = VegaLiteWriter::new(); + let json_str = writer.write(&prepared.specs[0], &prepared.data).unwrap(); + let vl_spec: serde_json::Value = serde_json::from_str(&json_str).unwrap(); + + assert_eq!(vl_spec["layer"][0]["mark"]["type"], "geoshape"); + + let data = vl_spec["data"]["values"].as_array().unwrap(); + let layer_key = prepared.specs[0].layers[0].data_key.as_ref().unwrap(); + let spatial_rows: Vec<_> = data + .iter() + .filter(|r| r[naming::SOURCE_COLUMN] == layer_key.as_str()) + .collect(); + assert!(!spatial_rows.is_empty()); + // Horizon clipping filters out back-of-globe features; all remaining have geometry + assert!( + spatial_rows.iter().all(|r| !r["geometry"].is_null()), + "all visible features should have valid geometry" + ); + + assert!( + !vl_spec["projection"]["scale"].is_null(), + "projection.scale should be present" + ); + assert!( + !vl_spec["projection"]["translate"].is_null(), + "projection.translate should be present" + ); + } + + #[cfg(feature = "spatial")] + #[test] + fn test_end_to_end_spatial_world_orthographic_antimeridian() { + let reader = DuckDBReader::from_connection_string("duckdb://memory").unwrap(); + + let queries = &[ + "+proj=ortho +lat_0=0 +lon_0=150", + "+proj=ortho +lat_0=52.36 +lon_0=150.90", + ]; + + for crs in queries { + let query = format!( + "VISUALISE FROM ggsql:world \ + DRAW spatial PROJECT TO map SETTING crs => '{crs}'" + ); + + let prepared = execute::prepare_data_with_reader(&query, &reader).unwrap(); + + let writer = VegaLiteWriter::new(); + let json_str = writer.write(&prepared.specs[0], &prepared.data).unwrap(); + let vl_spec: serde_json::Value = serde_json::from_str(&json_str).unwrap(); + + let data = vl_spec["data"]["values"].as_array().unwrap(); + let layer_key = prepared.specs[0].layers[0].data_key.as_ref().unwrap(); + let spatial_rows: Vec<_> = data + .iter() + .filter(|r| r[naming::SOURCE_COLUMN] == layer_key.as_str()) + .collect(); + assert!(!spatial_rows.is_empty(), "no rows for {crs}"); + assert!( + spatial_rows.iter().all(|r| !r["geometry"].is_null()), + "NULL geometry found for {crs}" + ); + } + } + /// Belt-and-braces regression test: a representative basket of error- /// triggering queries must never produce a user-visible message that /// contains an internal aesthetic name (`__ggsql_aes_*`, `pos1`, `pos2`, @@ -1142,4 +1224,73 @@ mod integration_tests { } } } + + #[cfg(feature = "spatial")] + #[test] + fn test_orthographic_amsterdam_center() { + let reader = DuckDBReader::from_connection_string("duckdb://memory").unwrap(); + + let query = r#" + VISUALISE FROM ggsql:world + DRAW spatial PROJECT TO map SETTING crs => '+proj=ortho +lon_0=4.90 +lat_0=52.36' + "#; + + let prepared = execute::prepare_data_with_reader(query, &reader).unwrap(); + let writer = VegaLiteWriter::new(); + let json_str = writer.write(&prepared.specs[0], &prepared.data).unwrap(); + let vl_spec: serde_json::Value = serde_json::from_str(&json_str).unwrap(); + + let data = vl_spec["data"]["values"].as_array().unwrap(); + let layer_key = prepared.specs[0].layers[0].data_key.as_ref().unwrap(); + let spatial_rows: Vec<_> = data + .iter() + .filter(|r| r[naming::SOURCE_COLUMN] == layer_key.as_str()) + .collect(); + assert!(!spatial_rows.is_empty()); + assert!(spatial_rows.iter().any(|r| !r["geometry"].is_null())); + } + + #[cfg(feature = "spatial")] + #[test] + fn test_non_4326_source_to_orthographic() { + let reader = DuckDBReader::from_connection_string("duckdb://memory").unwrap(); + + // Amsterdam in EPSG:3857: (545977.96, 6867712.83) + // Project to orthographic centered on Amsterdam — the point should survive. + let query = r#" + LOAD spatial; + SELECT ST_Point(545977.96, 6867712.83) AS geometry + VISUALISE + DRAW spatial + PROJECT TO map SETTING + source => 'EPSG:3857', + crs => '+proj=ortho +lon_0=4.90 +lat_0=52.36' + "#; + + let prepared = execute::prepare_data_with_reader(query, &reader).unwrap(); + let writer = VegaLiteWriter::new(); + let json_str = writer.write(&prepared.specs[0], &prepared.data).unwrap(); + let vl_spec: serde_json::Value = serde_json::from_str(&json_str).unwrap(); + + let data = vl_spec["data"]["values"].as_array().unwrap(); + let layer_key = prepared.specs[0].layers[0].data_key.as_ref().unwrap(); + let spatial_rows: Vec<_> = data + .iter() + .filter(|r| r[naming::SOURCE_COLUMN] == layer_key.as_str()) + .collect(); + assert_eq!(spatial_rows.len(), 1); + let geom = &spatial_rows[0]["geometry"]; + assert!( + !geom.is_null(), + "Point in source CRS should project successfully" + ); + + // The projected point should be near (0, 0) in orthographic coords + // since the projection is centered on the same location. + let coords = geom["coordinates"].as_array().unwrap(); + let x = coords[0].as_f64().unwrap(); + let y = coords[1].as_f64().unwrap(); + assert!(x.abs() < 2000.0, "Expected x near 0, got {x}"); + assert!(y.abs() < 2000.0, "Expected y near 0, got {y}"); + } } diff --git a/src/parser/builder.rs b/src/parser/builder.rs index 3d5498fa4..28067b652 100644 --- a/src/parser/builder.rs +++ b/src/parser/builder.rs @@ -4,6 +4,7 @@ //! handling all the node types defined in the grammar. use crate::plot::layer::geom::Geom; +use crate::plot::projection::coord::map_projections::MapSpecification; use crate::plot::projection::resolve_coord; use crate::plot::scale::{color_to_hex, is_color_aesthetic, is_user_facet_aesthetic, Transform}; use crate::plot::*; @@ -86,6 +87,15 @@ fn parse_boolean_node(node: &Node, source: &SourceTree) -> bool { text == "true" } +fn parse_infinity_node(node: &Node, source: &SourceTree) -> f64 { + let text = source.get_text(node); + if text.starts_with('-') { + f64::NEG_INFINITY + } else { + f64::INFINITY + } +} + /// Parse an array node into Vec fn parse_array_node(node: &Node, source: &SourceTree) -> Result> { let mut values = Vec::new(); @@ -105,6 +115,7 @@ fn parse_array_node(node: &Node, source: &SourceTree) -> Result ArrayElement::Number(parse_number_node(&elem_child, source)?), "boolean" => ArrayElement::Boolean(parse_boolean_node(&elem_child, source)), "null_literal" => ArrayElement::Null, + "infinity" => ArrayElement::Number(parse_infinity_node(&elem_child, source)), _ => { return Err(GgsqlError::ParseError(format!( "Invalid array element type: {}", @@ -138,6 +149,7 @@ fn parse_value_node(node: &Node, source: &SourceTree, context: &str) -> Result

Ok(ParameterValue::Null), + "infinity" => Ok(ParameterValue::Number(parse_infinity_node(node, source))), _ => Err(GgsqlError::ParseError(format!( "Unexpected {} value type: {}", context, @@ -272,13 +284,15 @@ fn build_visualise_statement(node: &Node, source: &SourceTree) -> Result { } } - // Resolve coord (infer from mappings if not explicit) + // Resolve coord (infer from mappings and layer types if not explicit) // This must happen after parsing but before initialize_aesthetic_context() let layer_mappings: Vec<&Mappings> = spec.layers.iter().map(|l| &l.mappings).collect(); + let layer_geom_types: Vec = spec.layers.iter().map(|l| l.geom.geom_type()).collect(); if let Some(inferred) = resolve_coord( spec.project.as_ref(), &spec.global_mappings, &layer_mappings, + &layer_geom_types, ) .map_err(GgsqlError::ParseError)? { @@ -958,6 +972,7 @@ fn parse_facet_vars(node: &Node, source: &SourceTree) -> Result> { /// Aesthetics are optional and default to the coord's standard names. fn build_project(node: &Node, source: &SourceTree) -> Result { let mut coord = Coord::cartesian(); + let mut coord_type_name: Option = None; let mut properties = HashMap::new(); let mut user_aesthetics: Option> = None; @@ -970,7 +985,9 @@ fn build_project(node: &Node, source: &SourceTree) -> Result { user_aesthetics = Some(source.find_texts(&child, query)); } "project_type" => { - coord = parse_coord_system(&child, source)?; + let text = source.get_text(&child).to_lowercase(); + coord = parse_coord_system(&text)?; + coord_type_name = Some(text); } "project_properties" => { // Find all project_property nodes @@ -1016,10 +1033,23 @@ fn build_project(node: &Node, source: &SourceTree) -> Result { // Validate properties for this coord type validate_project_properties(&coord, &properties)?; + let map_projection = MapSpecification::new(coord_type_name.as_deref(), &properties); + + if let Some(ref spec) = map_projection { + let proj_str = spec.to_proj_str(); + if !proj_str.is_empty() { + properties + .entry("crs".to_string()) + .or_insert_with(|| ParameterValue::String(proj_str)); + } + } + Ok(Projection { coord, aesthetics, properties, + map_projection, + computed: HashMap::new(), }) } @@ -1090,15 +1120,16 @@ fn validate_project_properties( Ok(()) } -/// Parse coord type from a project_type node -fn parse_coord_system(node: &Node, source: &SourceTree) -> Result { - let text = source.get_text(node); - match text.to_lowercase().as_str() { +fn parse_coord_system(name: &str) -> Result { + use crate::plot::projection::coord::map_projections::NAMED_PROJECTIONS; + match name { "cartesian" => Ok(Coord::cartesian()), "polar" => Ok(Coord::polar()), + "map" => Ok(Coord::map(name)), + _ if NAMED_PROJECTIONS.contains(&name) => Ok(Coord::map(name)), _ => Err(GgsqlError::ParseError(format!( "Unknown coord type: {}", - text + name ))), } } @@ -1341,6 +1372,61 @@ mod tests { .contains("conflicts with material aesthetic")); } + #[test] + fn test_project_map_bare() { + let query = r#" + VISUALISE + DRAW point MAPPING lon AS lon, lat AS lat + PROJECT TO map + "#; + + let result = parse_test_query(query); + assert!(result.is_ok()); + let specs = result.unwrap(); + + let project = specs[0].project.as_ref().unwrap(); + assert_eq!(project.coord.coord_kind(), CoordKind::Map); + assert_eq!( + project.aesthetics, + vec!["lon".to_string(), "lat".to_string()] + ); + assert!(!project.properties.contains_key("crs")); + } + + #[test] + fn test_project_mercator_shorthand() { + let query = r#" + VISUALISE + DRAW point MAPPING lon AS lon, lat AS lat + PROJECT TO mercator + "#; + + let result = parse_test_query(query); + assert!(result.is_ok()); + let specs = result.unwrap(); + + let project = specs[0].project.as_ref().unwrap(); + assert_eq!(project.coord.coord_kind(), CoordKind::Map); + assert_eq!( + project.properties.get("crs"), + Some(&ParameterValue::String("+proj=merc +lon_0=0".to_string())) + ); + } + + #[test] + fn test_project_shorthand_crs_override_rejected() { + let query = r#" + VISUALISE + DRAW point MAPPING lon AS lon, lat AS lat + PROJECT TO mercator SETTING crs => '+proj=custom' + "#; + + let result = parse_test_query(query); + assert!(result.is_err()); + let err = result.unwrap_err().to_string(); + assert!(err.contains("Cannot combine a named projection")); + } + // ======================================== // Case Insensitive Keywords Tests // ======================================== diff --git a/src/plot/CLAUDE.md b/src/plot/CLAUDE.md index 6dee23b26..bc4a58d10 100644 --- a/src/plot/CLAUDE.md +++ b/src/plot/CLAUDE.md @@ -35,7 +35,7 @@ plot/ │ ├── orientation.rs Layer transposition (horizontal vs vertical) │ └── position/ identity, stack, dodge, jitter ├── projection/ PROJECT clause -│ └── coord/ cartesian, polar +│ └── coord/ cartesian, polar, map └── scale/ SCALE clause ├── scale_type/ binned, continuous, discrete, identity, ordinal └── transform/ identity, log, sqrt, asinh, exp, square, pseudo_log, @@ -64,7 +64,15 @@ Scale type / transform docs: [`/doc/syntax/scale/type/`](../../doc/syntax/scale/ ### `facet/` and `projection/` -Smaller subsystems. Each has a `types.rs` (data structure) and `resolve.rs` (logic that runs during execution). `projection/coord/` currently has `cartesian` and `polar`. Docs: [`/doc/syntax/clause/facet.qmd`](../../doc/syntax/clause/facet.qmd), [`/doc/syntax/coord/`](../../doc/syntax/coord/). +Each has a `types.rs` (data structure) and `resolve.rs` (logic that runs during execution). `projection/coord/` has three implementations: + +- **`cartesian`** — standard x/y. No special behaviour. +- **`polar`** — radius/angle (for pie/rose plots). +- **`map`** — geographic projections via PROJ strings. Implements `apply_projection_transforms` to: detect source CRS from geometry SRID, make clip boundaries, delegate per-layer spatial transforms, materialize projected layers as temp tables, and resolve frame bbox from user bounds / data extent / world extent. Properties: `crs` (PROJ string), `source` (source EPSG), `clip` (bool), `bounds` ([xmin, ymin, xmax, ymax] with null/Inf fallback semantics). + +`Projection` (in `types.rs`) wraps `Coord` + resolved aesthetics + properties + a `computed` map populated at execution time for the writer (e.g., `panel_boundary`, `bbox`, `graticule_lon`, `graticule_lat`). + +Docs: [`/doc/syntax/clause/facet.qmd`](../../doc/syntax/clause/facet.qmd), [`/doc/syntax/coord/`](../../doc/syntax/coord/). ## Adding a new geom / scale type / coord diff --git a/src/plot/layer/geom/mod.rs b/src/plot/layer/geom/mod.rs index 8a0ab55d8..4141027e7 100644 --- a/src/plot/layer/geom/mod.rs +++ b/src/plot/layer/geom/mod.rs @@ -21,7 +21,7 @@ //! ``` use crate::plot::types::DefaultAestheticValue; -use crate::{DataFrame, Mappings, Result}; +use crate::{naming, DataFrame, Mappings, Result}; use serde::{Deserialize, Serialize}; use std::collections::HashMap; use std::sync::Arc; @@ -77,6 +77,7 @@ pub use tile::Tile; pub use violin::Violin; use crate::plot::aesthetic::AestheticContext; +use crate::plot::projection::Projection; use crate::plot::types::{ParameterValue, Schema}; use crate::reader::SqlDialect; @@ -286,6 +287,24 @@ pub trait GeomTrait: std::fmt::Debug + std::fmt::Display + Send + Sync { Ok(df) } + /// Apply coord-specific projection transformations to a layer query. + /// + /// Called after stat transforms, before data fetch. Each geom decides what + /// projection means for its parameterization: + /// - Spatial: ST_AsWKB (always), plus ST_Transform when Map coord has a CRS + /// - Future geoms: rectangles transform corners, lines segmentize, etc. + /// + /// The default is a no-op (returns query unchanged). + fn apply_projection( + &self, + query: &str, + _projection: &Projection, + _dialect: &dyn SqlDialect, + _clip: bool, + ) -> Result { + Ok(query.to_string()) + } + /// Adjust layer mappings and parameters based on geom-specific logic. /// /// This method is called during layer execution to allow geoms to customize @@ -313,6 +332,45 @@ pub trait GeomTrait: std::fmt::Debug + std::fmt::Display + Send + Sync { } } +/// Project pos1/pos2 columns through the map CRS transform. +/// +/// When the coordinate system is Map with a CRS, wraps the position columns +/// with ST_X/ST_Y(ST_Transform(ST_Point(pos1, pos2), source, crs)). Returns +/// the query unchanged for non-map coords or when source == crs. +pub(crate) fn project_position_columns( + query: &str, + projection: &Projection, + dialect: &dyn SqlDialect, +) -> Result { + use crate::plot::projection::coord::CoordKind; + + if projection.coord.coord_kind() != CoordKind::Map { + return Ok(query.to_string()); + } + let crs = match projection.properties.get("crs") { + Some(ParameterValue::String(s)) => s.as_str(), + _ => return Ok(query.to_string()), + }; + let source = match projection.properties.get("source") { + Some(ParameterValue::String(s)) => s.as_str(), + _ => "EPSG:4326", + }; + if source == crs { + return Ok(query.to_string()); + } + + let pos1 = naming::quote_ident(&naming::aesthetic_column("pos1")); + let pos2 = naming::quote_ident(&naming::aesthetic_column("pos2")); + let point_expr = format!("ST_Point({pos1}, {pos2})"); + let transformed = dialect.sql_st_transform(&point_expr, source, crs); + let proj_col = naming::quote_ident("__ggsql_proj_pt__"); + + Ok(format!( + "SELECT * REPLACE (ST_X({proj_col}) AS {pos1}, ST_Y({proj_col}) AS {pos2}) \ + FROM (SELECT *, {transformed} AS {proj_col} FROM ({query}))" + )) +} + /// True when `parameters["aggregate"]` is set to a non-null string or array. pub(crate) fn has_aggregate_param(parameters: &HashMap) -> bool { matches!( @@ -545,6 +603,17 @@ impl Geom { self.0.post_process(df, parameters) } + /// Apply coord-specific projection transformations + pub fn apply_projection( + &self, + query: &str, + projection: &Projection, + dialect: &dyn SqlDialect, + clip: bool, + ) -> Result { + self.0.apply_projection(query, projection, dialect, clip) + } + /// Adjust layer mappings and parameters based on geom-specific logic pub fn setup_layer( &self, diff --git a/src/plot/layer/geom/point.rs b/src/plot/layer/geom/point.rs index 3e9c55a6c..d812e4dc8 100644 --- a/src/plot/layer/geom/point.rs +++ b/src/plot/layer/geom/point.rs @@ -2,9 +2,13 @@ use super::types::POSITION_VALUES; use super::{ - DefaultAesthetics, DefaultParamValue, GeomTrait, GeomType, ParamConstraint, ParamDefinition, + project_position_columns, DefaultAesthetics, DefaultParamValue, GeomTrait, GeomType, + ParamConstraint, ParamDefinition, }; +use crate::plot::projection::Projection; use crate::plot::types::DefaultAestheticValue; +use crate::reader::SqlDialect; +use crate::Result; /// Point geom - scatter plots and similar #[derive(Debug, Clone, Copy)] @@ -50,6 +54,16 @@ impl GeomTrait for Point { fn aggregate_domain_aesthetics(&self) -> Option<&'static [&'static str]> { Some(&[]) } + + fn apply_projection( + &self, + query: &str, + projection: &Projection, + dialect: &dyn SqlDialect, + _clip: bool, + ) -> Result { + project_position_columns(query, projection, dialect) + } } impl std::fmt::Display for Point { diff --git a/src/plot/layer/geom/spatial.rs b/src/plot/layer/geom/spatial.rs index c73b8beb0..08252efd9 100644 --- a/src/plot/layer/geom/spatial.rs +++ b/src/plot/layer/geom/spatial.rs @@ -1,6 +1,40 @@ use super::{DefaultAesthetics, GeomTrait, GeomType, StatResult}; +use crate::naming; +use crate::plot::projection::coord::map::CLIP_BOUNDARY_TABLE; +use crate::plot::projection::coord::CoordKind; +use crate::plot::projection::Projection; use crate::plot::types::DefaultAestheticValue; -use crate::{naming, Mappings}; +use crate::plot::ParameterValue; +use crate::reader::SqlDialect; +use crate::Mappings; + +fn apply_clip_boundary( + query: &str, + col: &str, + source: &str, + crs: &str, + clip_table: &str, +) -> String { + let clip_geom = format!("(SELECT geom FROM {clip_table})"); + let source_esc = source.replace('\'', "''"); + let crs_esc = crs.replace('\'', "''"); + + let clipped = format!("ST_Intersection({col}, {clip_geom})"); + let geom_expr = format!( + "ST_MakeValid(ST_Transform(\ + {clipped},\ + '{source_esc}', '{crs_esc}', always_xy := true\ + ))", + ); + format!( + "SELECT * REPLACE ({geom_expr} AS {col}) FROM ({query}) \ + WHERE ST_Intersects({col}, {clip_geom})", + col = col, + geom_expr = geom_expr, + query = query, + clip_geom = clip_geom, + ) +} #[derive(Debug, Clone, Copy)] pub struct Spatial; @@ -38,36 +72,63 @@ impl GeomTrait for Spatial { execute_query(&stmt)?; } - // Geometry columns use database-native types that don't have an Arrow equivalent. - // Convert to standard WKB so the writer can parse them with geozero. - let geom_col = naming::aesthetic_column("geometry"); - let col = naming::quote_ident(&geom_col); - - // Skip conversion if the geometry column is already in binary WKB format. - let already_wkb = _schema.iter().any(|c| { - c.name == geom_col - && matches!( - c.dtype, - arrow::datatypes::DataType::Binary | arrow::datatypes::DataType::LargeBinary - ) - }); - - if already_wkb { - Ok(StatResult::Transformed { - query: query.to_string(), - stat_columns: vec![], - dummy_columns: vec![], - consumed_aesthetics: vec![], - }) + Ok(StatResult::Transformed { + query: query.to_string(), + stat_columns: vec![], + dummy_columns: vec![], + consumed_aesthetics: vec![], + }) + } + + fn apply_projection( + &self, + query: &str, + projection: &Projection, + dialect: &dyn SqlDialect, + clip: bool, + ) -> crate::Result { + let col = naming::quote_ident(&naming::aesthetic_column("geometry")); + let is_map = projection.coord.coord_kind() == CoordKind::Map; + + // WORKAROUND(duckdb-rs#714): normalize column to GEOMETRY since it may + // be WKB BLOB from the Arrow export workaround. + let ensure_geom = dialect.sql_ensure_geometry(&col); + let geom_query = format!("SELECT * REPLACE ({ensure_geom} AS {col}) FROM ({query})"); + + let geom_expr = if let (true, Some(ParameterValue::String(crs))) = + (is_map, projection.properties.get("crs")) + { + let source = match projection.properties.get("source") { + Some(ParameterValue::String(s)) => s.as_str(), + _ => "EPSG:4326", + }; + + if clip { + return Ok(apply_clip_boundary( + &geom_query, + &col, + source, + crs, + CLIP_BOUNDARY_TABLE, + )); + } + + dialect.sql_st_transform(&col, source, crs) + } else if is_map { + // Map coord without CRS — keep native geometry (WKB added later by framing) + return Ok(geom_query); } else { + // Non-map coord — convert to WKB directly let wkb_expr = dialect.sql_geometry_to_wkb(&col); - Ok(StatResult::Transformed { - query: format!("SELECT * REPLACE ({wkb_expr} AS {col}) FROM ({query})"), - stat_columns: vec![], - dummy_columns: vec![], - consumed_aesthetics: vec![], - }) - } + return Ok(format!( + "SELECT * REPLACE ({wkb_expr} AS {col}) FROM ({geom_query})" + )); + }; + + // Map coord with CRS — output native projected geometry (WKB added by framing) + Ok(format!( + "SELECT * REPLACE ({geom_expr} AS {col}) FROM ({geom_query})" + )) } } @@ -76,3 +137,107 @@ impl std::fmt::Display for Spatial { write!(f, "spatial") } } + +#[cfg(test)] +mod tests { + use super::*; + use crate::reader::AnsiDialect; + + #[test] + fn test_apply_projection_without_map_coord() { + let spatial = Spatial; + let projection = Projection::cartesian(); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, false) + .unwrap(); + + assert!(result.contains("ST_AsBinary")); + assert!(!result.contains("ST_Transform")); + } + + #[test] + fn test_apply_projection_map_without_crs() { + let spatial = Spatial; + let projection = Projection::map(); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, false) + .unwrap(); + + // Map without CRS ensures GEOMETRY type for the framing step + assert!(result.contains("ST_GeomFromWKB")); + assert!(!result.contains("ST_Transform")); + } + + #[test] + fn test_apply_projection_map_with_crs_no_clip() { + let spatial = Spatial; + let mut projection = Projection::map(); + projection.properties.insert( + "crs".to_string(), + ParameterValue::String("+proj=merc".to_string()), + ); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, false) + .unwrap(); + + // Without clip=true, just ST_Transform + assert!(!result.contains("ST_AsBinary")); + assert!(result.contains("ST_Transform")); + assert!(result.contains("+proj=merc")); + assert!(!result.contains("ST_Intersection")); + } + + #[test] + fn test_apply_projection_mercator_with_clip() { + let spatial = Spatial; + let mut projection = Projection::map(); + projection.properties.insert( + "crs".to_string(), + ParameterValue::String("+proj=merc".to_string()), + ); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, true) + .unwrap(); + + assert!(result.contains("ST_Intersection")); + assert!(result.contains("ST_Intersects")); + assert!(result.contains("__ggsql_clip_boundary__")); + } + + #[test] + fn test_orthographic_with_clip() { + let spatial = Spatial; + let mut projection = Projection::map(); + projection.properties.insert( + "crs".to_string(), + ParameterValue::String("+proj=ortho +lat_0=45 +lon_0=10".to_string()), + ); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, true) + .unwrap(); + + assert!(result.contains("ST_Transform")); + assert!(result.contains("ST_MakeValid")); + assert!(result.contains("ST_Intersection")); + assert!(result.contains("ST_Intersects")); + assert!(result.contains("__ggsql_clip_boundary__")); + } + + #[test] + fn test_gnomonic_with_clip() { + let spatial = Spatial; + let mut projection = Projection::map(); + projection.properties.insert( + "crs".to_string(), + ParameterValue::String("+proj=gnom +lat_0=90 +lon_0=0".to_string()), + ); + let result = spatial + .apply_projection("SELECT * FROM t", &projection, &AnsiDialect, true) + .unwrap(); + + assert!(result.contains("ST_MakeValid")); + assert!(result.contains("ST_Intersection")); + assert!(result.contains("ST_Intersects")); + assert!(result.contains("__ggsql_clip_boundary__")); + } +} diff --git a/src/plot/layer/geom/text.rs b/src/plot/layer/geom/text.rs index ae71ca55e..cd9974b41 100644 --- a/src/plot/layer/geom/text.rs +++ b/src/plot/layer/geom/text.rs @@ -2,11 +2,13 @@ use super::types::POSITION_VALUES; use super::{ - DefaultAesthetics, DefaultParamValue, GeomTrait, GeomType, ParamConstraint, ParamDefinition, - ParameterValue, + project_position_columns, DefaultAesthetics, DefaultParamValue, GeomTrait, GeomType, + ParamConstraint, ParamDefinition, ParameterValue, }; +use crate::plot::projection::Projection; use crate::plot::types::DefaultAestheticValue; use crate::plot::{ArrayConstraint, NumberConstraint}; +use crate::reader::SqlDialect; use crate::{naming, DataFrame, Result}; use std::collections::HashMap; @@ -68,6 +70,16 @@ impl GeomTrait for Text { Some(&[]) } + fn apply_projection( + &self, + query: &str, + projection: &Projection, + dialect: &dyn SqlDialect, + _clip: bool, + ) -> Result { + project_position_columns(query, projection, dialect) + } + fn post_process( &self, df: DataFrame, diff --git a/src/plot/main.rs b/src/plot/main.rs index 4535c0bbf..2fd6502f6 100644 --- a/src/plot/main.rs +++ b/src/plot/main.rs @@ -774,6 +774,8 @@ mod tests { coord: Coord::cartesian(), aesthetics: vec!["y".to_string(), "x".to_string()], properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), }); spec.labels = Some(Labels { labels: HashMap::from([ @@ -804,6 +806,8 @@ mod tests { coord: Coord::polar(), aesthetics: vec!["angle".to_string(), "radius".to_string()], properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), }); spec.labels = Some(Labels { labels: HashMap::from([ diff --git a/src/plot/projection/coord/map.rs b/src/plot/projection/coord/map.rs new file mode 100644 index 000000000..06e4c2e7a --- /dev/null +++ b/src/plot/projection/coord/map.rs @@ -0,0 +1,1166 @@ +//! Map coordinate system implementation + +use std::collections::HashMap; + +use super::{CoordKind, CoordTrait}; +use crate::naming; +use crate::plot::layer::geom::GeomType; +use crate::plot::types::{ + validate_parameter, DefaultParamValue, ParamConstraint, ParamDefinition, TypeConstraint, +}; +use crate::plot::{DataSource, Layer, ParameterValue}; +use crate::reader::SqlDialect; +use crate::DataFrame; + +pub const CLIP_BOUNDARY_TABLE: &str = "__ggsql_clip_boundary__"; + +// --------------------------------------------------------------------------- +// Map coord +// --------------------------------------------------------------------------- + +/// Map coordinate system - for geographic/cartographic projections +#[derive(Debug, Clone, Copy)] +pub struct Map { + coord_type_name: &'static str, +} + +impl Map { + pub fn new(name: &str) -> Self { + use super::map_projections::NAMED_PROJECTIONS; + let coord_type_name = NAMED_PROJECTIONS + .iter() + .find(|&&n| n == name) + .copied() + .unwrap_or("map"); + Self { coord_type_name } + } +} + +impl CoordTrait for Map { + fn coord_kind(&self) -> CoordKind { + CoordKind::Map + } + + fn name(&self) -> &'static str { + self.coord_type_name + } + + fn position_aesthetic_names(&self) -> &'static [&'static str] { + &["lon", "lat"] + } + + fn default_properties(&self) -> &'static [ParamDefinition] { + use crate::plot::types::{ArrayConstraint, NumberConstraint}; + const LON_RANGE: NumberConstraint = NumberConstraint::range(-180.0, 180.0); + const LAT_RANGE: NumberConstraint = NumberConstraint::range(-90.0, 90.0); + const PARAMS: &[ParamDefinition] = &[ + ParamDefinition { + name: "crs", + default: DefaultParamValue::Null, + constraint: ParamConstraint::string(), + }, + ParamDefinition { + name: "source", + default: DefaultParamValue::Null, + constraint: ParamConstraint::string(), + }, + ParamDefinition { + name: "clip", + default: DefaultParamValue::Boolean(true), + constraint: ParamConstraint::boolean(), + }, + // [xmin, ymin, xmax, ymax] in projected coordinates; null uses data bbox, Inf uses world bbox + ParamDefinition { + name: "bounds", + default: DefaultParamValue::Null, + constraint: ParamConstraint { + number: TypeConstraint::Forbidden, + string: TypeConstraint::Forbidden, + boolean: TypeConstraint::Forbidden, + array: TypeConstraint::Constrained( + ArrayConstraint::of_numbers_len(NumberConstraint::unconstrained(), 4) + .with_null_elements(), + ), + allow_null: true, + }, + }, + // origin => 30 (lon only) or origin => (30, 45) (lon, lat) + ParamDefinition { + name: "origin", + default: DefaultParamValue::Null, + constraint: ParamConstraint::number_or_numeric_array( + LON_RANGE, + ArrayConstraint::of_numbers_len(LON_RANGE, 2), + ), + }, + // parallel => 30 (tangent) or parallel => (30, 50) (secant) + ParamDefinition { + name: "parallel", + default: DefaultParamValue::Null, + constraint: ParamConstraint::number_or_numeric_array( + LAT_RANGE, + ArrayConstraint::of_numbers_len(LAT_RANGE, 2), + ), + }, + ]; + PARAMS + } + + fn resolve_properties( + &self, + properties: &HashMap, + ) -> Result, String> { + if self.coord_type_name != "map" && properties.contains_key("crs") { + return Err(format!( + "Cannot combine a named projection ('{}') with a 'crs' string. \ + Use either PROJECT TO {} or PROJECT TO map SETTING crs => '...'", + self.coord_type_name, self.coord_type_name + )); + } + let has_crs = properties.contains_key("crs"); + let has_origin = properties.contains_key("origin"); + let has_parallel = properties.contains_key("parallel"); + if has_crs && (has_origin || has_parallel) { + return Err( + "Cannot combine 'crs' setting with 'origin' or 'parallel'. \ + Use either the CRS string or a named projection with 'origin'/'parallel' settings." + .to_string(), + ); + } + // Delegate to default validation + let defaults = self.default_properties(); + for (key, value) in properties.iter() { + if let Some(param) = defaults.iter().find(|p| p.name == key) { + validate_parameter(key, value, ¶m.constraint)?; + } else { + let allowed: Vec<&str> = defaults.iter().map(|p| p.name).collect(); + return Err(format!( + "{} projection property should be {}, not '{}'", + self.name(), + crate::or_list_quoted(&allowed, '\''), + key + )); + } + } + // ArrayConstraint applies uniform bounds to all elements, so we validate + // the latitude element (index 1) separately against [-90, 90]. + if let Some(ParameterValue::Array(arr)) = properties.get("origin") { + if let Some(crate::plot::types::ArrayElement::Number(lat)) = arr.get(1) { + if *lat < -90.0 || *lat > 90.0 { + return Err(format!( + "origin latitude must be between -90 and 90, got {}", + lat + )); + } + } + } + let mut resolved = properties.clone(); + for param in defaults { + if !resolved.contains_key(param.name) { + if let Some(default) = param.to_parameter_value() { + resolved.insert(param.name.to_string(), default); + } + } + } + Ok(resolved) + } + + fn apply_projection_transforms( + &self, + layers: &[Layer], + layer_queries: &mut [String], + projection: &mut super::super::Projection, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, + ) -> crate::Result<()> { + for stmt in dialect.sql_spatial_setup() { + execute_query(&stmt)?; + } + + // Step 1: Detect source CRS from geometry columns if not explicitly set + if !projection.properties.contains_key("source") { + if let Some(srid) = detect_source_srid(layers, layer_queries, execute_query)? { + projection + .properties + .insert("source".to_string(), ParameterValue::String(srid)); + } + } + + let source = match projection.properties.get("source") { + Some(ParameterValue::String(s)) => s.clone(), + _ => "EPSG:4326".to_string(), + }; + let crs = match projection.properties.get("crs") { + Some(ParameterValue::String(s)) => s.clone(), + _ => { + projection + .properties + .insert("crs".to_string(), ParameterValue::String(source.clone())); + source.clone() + } + }; + + // Validate CRS by attempting a single point transform + let probe = dialect.sql_st_transform("ST_Point(0, 0)", &source, &crs); + if let Err(e) = execute_query(&format!("SELECT {probe}")) { + let msg = e.to_string(); + return Err(crate::GgsqlError::ValidationError(format!( + "Invalid CRS '{}': {}", + crs, + msg.split(':').next_back().unwrap_or(&msg).trim() + ))); + } + + // Step 2: Materialize clip boundary, panel boundary, and world bbox. + let clip_enabled = match projection.properties.get("clip") { + Some(ParameterValue::Boolean(b)) => *b, + _ => true, + }; + let mut world_bbox: Option = None; + let mut boundary_lonlat: Option = None; + + let clip_wkt = clip_enabled + .then_some(projection.map_projection.as_ref()) + .flatten() + .and_then(|map_proj| map_proj.visible_area_wkt().map(|_| map_proj)); + + if let Some(map_proj) = clip_wkt { + let b = materialize_clip_boundary(map_proj, &source, dialect, execute_query)?; + if let Some(wkt) = boundary_to_target_crs(&b, &crs, dialect, execute_query) { + projection + .computed + .insert("panel_boundary".to_string(), ParameterValue::String(wkt)); + } + world_bbox = compute_world_bbox(&source, &crs, dialect, execute_query); + boundary_lonlat = Some(b); + } + let clip = boundary_lonlat.is_some(); + + // Step 3: Apply per-layer projection (ST_Transform, clip to horizon) + for (idx, layer) in layers.iter().enumerate() { + layer_queries[idx] = + layer + .geom + .apply_projection(&layer_queries[idx], projection, dialect, clip)?; + } + + // Step 4: Materialize projected layers as temp tables, compute data bbox, + // then rewrite layer queries to read from those tables. + let user_bbox = projection.properties.get("bounds"); + let needs_data_bbox = needs_data_bbox(user_bbox); + let mut data_bbox: Option = None; + + for (idx, layer) in layers.iter().enumerate() { + let is_annotation = matches!(layer.source, Some(DataSource::Annotation)); + let is_spatial = layer.geom.geom_type() == GeomType::Spatial; + let has_projected_positions = !is_spatial + && source != crs + && layer.mappings.contains_key("pos1") + && layer.mappings.contains_key("pos2"); + + if is_annotation || (!is_spatial && !has_projected_positions) { + continue; + } + + let table_quoted = materialize_layer(idx, &layer_queries[idx], dialect, execute_query)?; + + if needs_data_bbox { + let layer_bbox = + compute_layer_bbox(&table_quoted, is_spatial, &crs, dialect, execute_query); + data_bbox = BBox::merge(data_bbox, layer_bbox)?; + } + + layer_queries[idx] = if is_spatial { + let geom_col_quoted = naming::quote_ident(&naming::aesthetic_column("geometry")); + let wkb_expr = dialect.sql_geometry_to_wkb(&geom_col_quoted); + format!("SELECT * REPLACE ({wkb_expr} AS {geom_col_quoted}) FROM {table_quoted}") + } else { + format!("SELECT * FROM {table_quoted}") + }; + } + + // Step 5: Resolve final frame bbox from user bounds + data bounds + world bounds + let Some(bbox) = resolve_final_bbox(user_bbox, data_bbox, world_bbox) else { + return Ok(()); + }; + projection + .computed + .insert("bbox".to_string(), bbox.as_parameter_value()); + + // Step 6: Generate graticule lines. The graticule is built and clipped + // in EPSG:4326 (independent of source), then projected to target. + let (lon_wkt, lat_wkt) = build_graticule( + &bbox, + boundary_lonlat.as_deref(), + &crs, + dialect, + execute_query, + )?; + if let Some(wkt) = lon_wkt { + projection + .computed + .insert("graticule_lon".to_string(), ParameterValue::String(wkt)); + } + if let Some(wkt) = lat_wkt { + projection + .computed + .insert("graticule_lat".to_string(), ParameterValue::String(wkt)); + } + + Ok(()) + } +} + +impl std::fmt::Display for Map { + fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { + write!(f, "{}", self.name()) + } +} + +// --------------------------------------------------------------------------- +// BBox +// --------------------------------------------------------------------------- + +#[derive(Debug, Clone, PartialEq)] +struct BBox { + xmin: f64, + ymin: f64, + xmax: f64, + ymax: f64, + crs: String, +} + +impl BBox { + fn from_df(df: &DataFrame, crs: &str) -> Option { + use arrow::array::Array; + let batch = df.inner(); + if batch.num_rows() == 0 || batch.num_columns() < 4 { + return None; + } + let get_f64 = |col: usize| -> Option { + batch + .column(col) + .as_any() + .downcast_ref::() + .filter(|a| !a.is_null(0)) + .map(|a| a.value(0)) + }; + match (get_f64(0), get_f64(1), get_f64(2), get_f64(3)) { + (Some(xmin), Some(ymin), Some(xmax), Some(ymax)) => Some(Self { + xmin, + ymin, + xmax, + ymax, + crs: crs.to_string(), + }), + _ => None, + } + } + + fn merge(existing: Option, new: Option) -> crate::Result> { + match (existing, new) { + (Some(a), Some(b)) => { + if a.crs != b.crs { + return Err(crate::GgsqlError::InternalError(format!( + "Cannot merge bounding boxes with different CRS: '{}' vs '{}'", + a.crs, b.crs + ))); + } + Ok(Some(Self { + xmin: a.xmin.min(b.xmin), + ymin: a.ymin.min(b.ymin), + xmax: a.xmax.max(b.xmax), + ymax: a.ymax.max(b.ymax), + crs: a.crs, + })) + } + (Some(b), None) | (None, Some(b)) => Ok(Some(b)), + (None, None) => Ok(None), + } + } + + fn from_array(arr: [f64; 4], crs: &str) -> Self { + Self { + xmin: arr[0].min(arr[2]), + ymin: arr[1].min(arr[3]), + xmax: arr[0].max(arr[2]), + ymax: arr[1].max(arr[3]), + crs: crs.to_string(), + } + } + + fn to_array(&self) -> [f64; 4] { + [self.xmin, self.ymin, self.xmax, self.ymax] + } + + fn clamp(mut self, xmin: f64, ymin: f64, xmax: f64, ymax: f64) -> Self { + self.xmin = self.xmin.clamp(xmin, xmax); + self.ymin = self.ymin.clamp(ymin, ymax); + self.xmax = self.xmax.clamp(xmin, xmax); + self.ymax = self.ymax.clamp(ymin, ymax); + self + } + + fn xrange(&self) -> (f64, f64) { + (self.xmin, self.xmax) + } + + fn yrange(&self) -> (f64, f64) { + (self.ymin, self.ymax) + } + + fn as_parameter_value(&self) -> ParameterValue { + use crate::plot::types::ArrayElement; + ParameterValue::Array(vec![ + ArrayElement::Number(self.xmin), + ArrayElement::Number(self.ymin), + ArrayElement::Number(self.xmax), + ArrayElement::Number(self.ymax), + ]) + } + + fn reproject( + &self, + target_crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, + ) -> Option { + let envelope = format!( + "ST_MakeEnvelope({}, {}, {}, {})", + self.xmin, self.ymin, self.xmax, self.ymax + ); + let transformed = dialect.sql_st_transform(&envelope, &self.crs, target_crs); + let sql = format!( + "SELECT ST_XMin(g) AS xmin, ST_YMin(g) AS ymin, \ + ST_XMax(g) AS xmax, ST_YMax(g) AS ymax \ + FROM (SELECT {transformed} AS g)" + ); + execute_query(&sql) + .ok() + .and_then(|df| Self::from_df(&df, target_crs)) + } +} + +// --------------------------------------------------------------------------- +// Graticule helpers +// --------------------------------------------------------------------------- + +/// Build graticule lines: determine the visible lon/lat extent, generate densified +/// meridians and parallels, clip and project them, and return projected WKT. +fn build_graticule( + bbox: &BBox, + clip_boundary_wkt: Option<&str>, + crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result<(Option, Option)> { + let Some(geo_bbox) = graticule_bbox(bbox, clip_boundary_wkt, dialect, execute_query)? else { + return Ok((None, None)); + }; + + let lon_breaks = graticule_breaks(geo_bbox.xrange()); + let lat_breaks = graticule_breaks(geo_bbox.yrange()); + + if lon_breaks.is_empty() && lat_breaks.is_empty() { + return Ok((None, None)); + } + + // Densification interval based on angular extent + let max_range = (geo_bbox.xmax - geo_bbox.xmin).max(geo_bbox.ymax - geo_bbox.ymin); + let step_deg = if max_range > 90.0 { + 2.0 + } else if max_range > 30.0 { + 1.0 + } else { + 0.5 + }; + + // Clamp meridians away from ±180 to avoid antimeridian issues, and + // deduplicate (e.g. if both -180 and 180 were present, they become the same) + let lon_breaks: Vec = { + let mut clamped: Vec = lon_breaks + .iter() + .map(|&v| { + if v <= -180.0 { + -179.999999 + } else if v >= 180.0 { + 179.999999 + } else { + v + } + }) + .collect(); + clamped.dedup_by(|a, b| (*a - *b).abs() < 0.001); + clamped + }; + + let lon_wkt = if !lon_breaks.is_empty() { + Some(grid_lines_wkt( + &lon_breaks, + geo_bbox.yrange(), + step_deg, + true, + )) + } else { + None + }; + let lat_wkt = if !lat_breaks.is_empty() { + Some(grid_lines_wkt( + &lat_breaks, + geo_bbox.xrange(), + step_deg, + false, + )) + } else { + None + }; + + Ok(( + project_graticule_wkt(lon_wkt, clip_boundary_wkt, crs, dialect, execute_query)?, + project_graticule_wkt(lat_wkt, clip_boundary_wkt, crs, dialect, execute_query)?, + )) +} + +/// Determine the lon/lat bounding box visible in the current frame by inverse-projecting +/// the bbox corners to EPSG:4326. Falls back to the clip boundary extent for azimuthal +/// projections where corners collapse to degenerate values. +fn graticule_bbox( + bbox: &BBox, + clip_boundary_wkt: Option<&str>, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result> { + let mut geo_bbox = match bbox.reproject("EPSG:4326", dialect, execute_query) { + Some(b) => b.clamp(-180.0, -90.0, 180.0, 90.0), + None => return Ok(None), + }; + + // For azimuthal projections the bbox corners often inverse-project to + // degenerate or incomplete values. Use the clip boundary extent which + // correctly represents the visible hemisphere. + if let Some(wkt) = clip_boundary_wkt { + let sql = format!( + "SELECT ST_XMin(g) AS xmin, ST_YMin(g) AS ymin, \ + ST_XMax(g) AS xmax, ST_YMax(g) AS ymax \ + FROM (SELECT ST_GeomFromText('{wkt}') AS g)" + ); + if let Ok(df) = execute_query(&sql) { + if let Some(clip_bbox) = BBox::from_df(&df, "EPSG:4326") { + geo_bbox = clip_bbox; + } + } + } + + // For projections showing the full globe, expand to full range + if geo_bbox.xmax - geo_bbox.xmin > 300.0 { + geo_bbox.xmin = -180.0; + geo_bbox.xmax = 180.0; + } + if geo_bbox.ymax - geo_bbox.ymin > 150.0 { + geo_bbox.ymin = -90.0; + geo_bbox.ymax = 90.0; + } + + Ok(Some(geo_bbox)) +} + +/// Pick pretty graticule break positions for a lon or lat range. +/// Uses standard angular intervals (multiples of 1, 2, 5, 10, 15, 30, 45, 90). +fn graticule_breaks((min, max): (f64, f64)) -> Vec { + let range = max - min; + if range <= 0.0 { + return vec![]; + } + + const STEPS: &[f64] = &[1.0, 2.0, 5.0, 10.0, 15.0, 20.0, 30.0, 45.0, 60.0, 90.0]; + + // Pick the smallest step that gives at most ~7 lines + let step = STEPS + .iter() + .copied() + .find(|&s| range / s <= 8.0) + .unwrap_or(90.0); + + let start = (min / step).ceil() as i64; + let end = (max / step).floor() as i64; + let mut breaks: Vec = (start..=end) + .map(|i| i as f64 * step) + .filter(|&v| v > min && v < max) + .collect(); + + // Include the boundary value when range covers the full extent, + // so the antimeridian/pole gets a line + if min <= -180.0 && !breaks.contains(&-180.0) { + breaks.insert(0, -180.0); + } else if max >= 180.0 && !breaks.contains(&180.0) { + breaks.push(180.0); + } + if min <= -90.0 && !breaks.contains(&-90.0) { + breaks.insert(0, -90.0); + } else if max >= 90.0 && !breaks.contains(&90.0) { + breaks.push(90.0); + } + + breaks +} + +/// Generate a MULTILINESTRING WKT with one line per break value, densified along +/// the varying axis at `step_deg` intervals. +/// - `lon_first = true`: fixed longitude (meridians), varying latitude. +/// - `lon_first = false`: fixed latitude (parallels), varying longitude. +fn grid_lines_wkt( + breaks: &[f64], + (vary_min, vary_max): (f64, f64), + step_deg: f64, + lon_first: bool, +) -> String { + let mut lines: Vec = Vec::with_capacity(breaks.len()); + for &fixed in breaks { + let mut coords = Vec::new(); + let mut v = vary_min; + while v < vary_max { + let (lon, lat) = if lon_first { (fixed, v) } else { (v, fixed) }; + coords.push(format!("{lon:.6} {lat:.6}")); + v += step_deg; + } + let (lon, lat) = if lon_first { + (fixed, vary_max) + } else { + (vary_max, fixed) + }; + coords.push(format!("{lon:.6} {lat:.6}")); + if coords.len() >= 2 { + lines.push(format!("({})", coords.join(", "))); + } + } + format!("MULTILINESTRING({})", lines.join(", ")) +} + +// --------------------------------------------------------------------------- +// Generic helpers +// --------------------------------------------------------------------------- + +/// Execute a query and extract a single string value from the first row, first column. +fn query_scalar_string( + sql: &str, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> Option { + use arrow::array::Array; + let df = execute_query(sql).ok()?; + let batch = df.inner(); + if batch.num_rows() == 0 { + return None; + } + let arr = batch + .column(0) + .as_any() + .downcast_ref::()?; + if arr.is_null(0) { + return None; + } + Some(arr.value(0).to_string()) +} + +/// Compose the clip boundary (visible area minus seam slits), materialize it as a +/// temp table in source CRS for per-layer clipping, and return the WKT in EPSG:4326. +fn materialize_clip_boundary( + map_proj: &super::map_projections::MapSpecification, + source: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result { + let wkt = map_proj.visible_area_wkt().unwrap(); + let half_width = 0.005; + let slit_wkt = map_proj.slit_wkt(half_width); + + let boundary_lonlat = if let Some(slit) = &slit_wkt { + let sql = format!( + "SELECT ST_AsText(ST_Difference(ST_GeomFromText('{wkt}'), ST_GeomFromText('{slit}'))) AS wkt" + ); + query_scalar_string(&sql, execute_query).unwrap_or(wkt) + } else { + wkt + }; + + let source_geom = dialect.sql_st_transform( + &format!("ST_GeomFromText('{boundary_lonlat}')"), + "EPSG:4326", + source, + ); + let body = format!("SELECT {source_geom} AS geom"); + for stmt in dialect.create_or_replace_temp_table_sql(CLIP_BOUNDARY_TABLE, &[], &body) { + execute_query(&stmt)?; + } + + Ok(boundary_lonlat) +} + +/// Project the clip boundary from EPSG:4326 to target CRS, returning the WKT. +fn boundary_to_target_crs( + boundary_lonlat: &str, + crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> Option { + let panel_geom = dialect.sql_st_transform( + &format!("ST_GeomFromText('{boundary_lonlat}')"), + "EPSG:4326", + crs, + ); + let sql = format!("SELECT ST_AsText({panel_geom}) AS wkt"); + query_scalar_string(&sql, execute_query) +} + +/// Materialize a layer query as a temp table, returning the quoted table name. +fn materialize_layer( + idx: usize, + query: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result { + let table_name = format!("{}_proj", naming::layer_key(idx)); + for stmt in dialect.create_or_replace_temp_table_sql(&table_name, &[], query) { + execute_query(&stmt)?; + } + Ok(naming::quote_ident(&table_name)) +} + +/// Compute the bounding box of a single materialized layer table. +fn compute_layer_bbox( + table: &str, + is_spatial: bool, + crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> Option { + let sql = if is_spatial { + let geom_col = naming::quote_ident(&naming::aesthetic_column("geometry")); + dialect.sql_geometry_bbox(&geom_col, table) + } else { + let pos1_col = naming::quote_ident(&naming::aesthetic_column("pos1")); + let pos2_col = naming::quote_ident(&naming::aesthetic_column("pos2")); + format!( + "SELECT MIN({pos1_col}), MIN({pos2_col}), \ + MAX({pos1_col}), MAX({pos2_col}) FROM {table}" + ) + }; + if let Ok(df) = execute_query(&sql) { + BBox::from_df(&df, crs) + } else { + None + } +} + +/// Compute the world bounding box by reading the extent of the materialized +/// clip boundary table, projected to target CRS. +fn compute_world_bbox( + source: &str, + crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> Option { + let projected = dialect.sql_st_transform("geom", source, crs); + let sql = dialect.sql_geometry_bbox(&projected, CLIP_BOUNDARY_TABLE); + if let Ok(df) = execute_query(&sql) { + BBox::from_df(&df, crs) + } else { + None + } +} + +/// Clip (if needed) and project a graticule WKT from EPSG:4326 to the target CRS. +fn project_graticule_wkt( + wkt: Option, + clip_boundary_wkt: Option<&str>, + crs: &str, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result> { + let Some(wkt) = wkt else { return Ok(None) }; + let geom_expr = format!("ST_GeomFromText('{wkt}')"); + let clipped = if let Some(boundary) = clip_boundary_wkt { + // ST_CollectionExtract(..., 2) keeps only linestring components, + // discarding stray points from vertex-on-boundary intersections. + format!( + "ST_CollectionExtract(ST_Intersection({geom_expr}, \ + ST_GeomFromText('{boundary}')), 2)" + ) + } else { + geom_expr + }; + let projected = dialect.sql_st_transform(&clipped, "EPSG:4326", crs); + let sql = format!("SELECT ST_AsText({projected}) AS wkt"); + Ok(query_scalar_string(&sql, execute_query)) +} + +/// Returns true if we need to compute a bbox (bounding box representing the extent of geometry) +/// from the data — i.e. when bounds is absent or has null elements that need filling in. +fn needs_data_bbox(user_bbox: Option<&ParameterValue>) -> bool { + match user_bbox { + Some(ParameterValue::Array(arr)) => { + use crate::plot::types::ArrayElement; + arr.iter().any(|e| !matches!(e, ArrayElement::Number(_))) + } + _ => true, + } +} + +/// Resolve the frame bbox: merge explicit bounds with computed values. +/// - Null elements fall back to the corresponding data-computed bbox. +/// - Inf/-Inf elements fall back to the clip boundary (world) bbox. +fn resolve_final_bbox( + user_bbox: Option<&ParameterValue>, + computed: Option, + world: Option, +) -> Option { + if let Some(ParameterValue::Array(arr)) = user_bbox { + use crate::plot::types::ArrayElement; + let data_fallback = computed.as_ref().map_or([f64::NAN; 4], |b| b.to_array()); + let world_fallback = world.as_ref().map_or([f64::NAN; 4], |b| b.to_array()); + let crs = computed + .as_ref() + .or(world.as_ref()) + .map(|b| b.crs.clone()) + .unwrap_or_default(); + let resolved: Vec = arr + .iter() + .enumerate() + .map(|(i, e)| match e { + ArrayElement::Number(n) if n.is_finite() => *n, + ArrayElement::Number(_) => world_fallback[i], + _ => data_fallback[i], + }) + .collect(); + if resolved.len() == 4 && resolved.iter().all(|v| v.is_finite()) { + return Some(BBox::from_array( + [resolved[0], resolved[1], resolved[2], resolved[3]], + &crs, + )); + } + } + computed +} + +fn detect_source_srid( + layers: &[Layer], + layer_queries: &[String], + execute_query: &dyn Fn(&str) -> crate::Result, +) -> crate::Result> { + let geom_col = naming::quote_ident(&naming::aesthetic_column("geometry")); + let mut detected: Option = None; + + for (idx, layer) in layers.iter().enumerate() { + if layer.geom.geom_type() != GeomType::Spatial { + continue; + } + let sql = format!( + "SELECT ST_SRID({geom_col}) AS srid FROM ({}) WHERE {geom_col} IS NOT NULL LIMIT 1", + layer_queries[idx] + ); + if let Ok(df) = execute_query(&sql) { + let batch = df.inner(); + if batch.num_rows() == 0 { + continue; + } + if let Some(arr) = batch + .column(0) + .as_any() + .downcast_ref::() + { + let srid = arr.value(0); + if srid != 0 { + let crs = format!("EPSG:{srid}"); + if let Some(ref prev) = detected { + if *prev != crs { + return Err(crate::GgsqlError::ValidationError(format!( + "Spatial layers have conflicting CRS: '{}' vs '{}'. \ + Set PROJECT source to specify which CRS the data is in.", + prev, crs + ))); + } + } else { + detected = Some(crs); + } + } + } + } + } + Ok(detected) +} + +// --------------------------------------------------------------------------- +// Tests +// --------------------------------------------------------------------------- + +#[cfg(test)] +mod tests { + use super::*; + use crate::plot::ParameterValue; + use std::collections::HashMap; + + #[test] + fn test_map_properties() { + let map = Map::new("map"); + assert_eq!(map.coord_kind(), CoordKind::Map); + assert_eq!(map.name(), "map"); + assert_eq!(map.position_aesthetic_names(), &["lon", "lat"]); + } + + #[test] + fn test_map_default_properties() { + let map = Map::new("map"); + let defaults = map.default_properties(); + let names: Vec<&str> = defaults.iter().map(|p| p.name).collect(); + assert!(names.contains(&"crs")); + assert!(names.contains(&"source")); + assert!(names.contains(&"clip")); + assert!(names.contains(&"bounds")); + assert!(names.contains(&"origin")); + assert!(names.contains(&"parallel")); + assert_eq!(defaults.len(), 6); + } + + #[test] + fn test_map_accepts_crs_string() { + let map = Map::new("map"); + let mut props = HashMap::new(); + props.insert( + "crs".to_string(), + ParameterValue::String("+proj=merc".to_string()), + ); + + let resolved = map.resolve_properties(&props); + assert!(resolved.is_ok()); + let resolved = resolved.unwrap(); + assert_eq!( + resolved.get("crs").unwrap(), + &ParameterValue::String("+proj=merc".to_string()) + ); + } + + #[test] + fn test_map_rejects_unknown_property() { + let map = Map::new("map"); + let mut props = HashMap::new(); + props.insert( + "unknown".to_string(), + ParameterValue::String("value".to_string()), + ); + + let resolved = map.resolve_properties(&props); + assert!(resolved.is_err()); + let err = resolved.unwrap_err(); + assert!(err.contains("not 'unknown'")); + } + + #[test] + fn test_crs_rejects_origin_and_parallel() { + let map = Map::new("map"); + let mut props = HashMap::new(); + props.insert( + "crs".to_string(), + ParameterValue::String("+proj=ortho".to_string()), + ); + props.insert("origin".to_string(), ParameterValue::Number(30.0)); + + let resolved = map.resolve_properties(&props); + assert!(resolved.is_err()); + let err = resolved.unwrap_err(); + assert!(err.contains("Cannot combine 'crs'")); + } + + #[test] + fn test_origin_rejects_latitude_out_of_range() { + let map = Map::new("albers"); + let mut props = HashMap::new(); + props.insert( + "origin".to_string(), + ParameterValue::Array(vec![ + crate::plot::types::ArrayElement::Number(0.0), + crate::plot::types::ArrayElement::Number(95.0), + ]), + ); + + let resolved = map.resolve_properties(&props); + assert!(resolved.is_err()); + let err = resolved.unwrap_err(); + assert!(err.contains("origin latitude must be between -90 and 90")); + } + + fn bbox(xmin: f64, ymin: f64, xmax: f64, ymax: f64) -> BBox { + BBox::from_array([xmin, ymin, xmax, ymax], "EPSG:4326") + } + + #[test] + fn test_resolve_final_bbox_no_bounds_uses_computed() { + let computed = Some(bbox(0.0, 0.0, 100.0, 200.0)); + assert_eq!(resolve_final_bbox(None, computed.clone(), None), computed); + } + + #[test] + fn test_resolve_final_bbox_no_bounds_no_computed() { + assert_eq!(resolve_final_bbox(None, None, None), None); + } + + #[test] + fn test_resolve_final_bbox_explicit_bounds_override_computed() { + use crate::plot::types::ArrayElement; + let bounds = ParameterValue::Array(vec![ + ArrayElement::Number(10.0), + ArrayElement::Number(20.0), + ArrayElement::Number(30.0), + ArrayElement::Number(40.0), + ]); + let computed = Some(bbox(0.0, 0.0, 100.0, 200.0)); + assert_eq!( + resolve_final_bbox(Some(&bounds), computed, None), + Some(bbox(10.0, 20.0, 30.0, 40.0)) + ); + } + + #[test] + fn test_resolve_final_bbox_null_elements_use_computed() { + use crate::plot::types::ArrayElement; + let bounds = ParameterValue::Array(vec![ + ArrayElement::Null, + ArrayElement::Number(20.0), + ArrayElement::Null, + ArrayElement::Number(40.0), + ]); + let computed = Some(bbox(5.0, 0.0, 95.0, 0.0)); + assert_eq!( + resolve_final_bbox(Some(&bounds), computed, None), + Some(bbox(5.0, 20.0, 95.0, 40.0)) + ); + } + + #[test] + fn test_resolve_final_bbox_inf_elements_use_world() { + use crate::plot::types::ArrayElement; + let bounds = ParameterValue::Array(vec![ + ArrayElement::Number(f64::NEG_INFINITY), + ArrayElement::Number(20.0), + ArrayElement::Number(f64::INFINITY), + ArrayElement::Number(40.0), + ]); + let computed = Some(bbox(5.0, 0.0, 95.0, 0.0)); + let world = Some(bbox(-500.0, -500.0, 500.0, 500.0)); + assert_eq!( + resolve_final_bbox(Some(&bounds), computed, world), + Some(bbox(-500.0, 20.0, 500.0, 40.0)) + ); + } + + #[test] + fn test_resolve_final_bbox_null_without_computed_falls_through() { + use crate::plot::types::ArrayElement; + let bounds = ParameterValue::Array(vec![ + ArrayElement::Null, + ArrayElement::Number(20.0), + ArrayElement::Number(30.0), + ArrayElement::Number(40.0), + ]); + assert_eq!(resolve_final_bbox(Some(&bounds), None, None), None); + } + + #[test] + fn test_resolve_final_bbox_inf_without_world_falls_through() { + use crate::plot::types::ArrayElement; + let bounds = ParameterValue::Array(vec![ + ArrayElement::Number(f64::INFINITY), + ArrayElement::Number(20.0), + ArrayElement::Number(30.0), + ArrayElement::Number(40.0), + ]); + let computed = Some(bbox(5.0, 0.0, 95.0, 200.0)); + assert_eq!( + resolve_final_bbox(Some(&bounds), computed.clone(), None), + computed + ); + } + + #[test] + fn test_merge_bbox() { + let a = Some(bbox(0.0, 10.0, 50.0, 60.0)); + let b = Some(bbox(-5.0, 15.0, 45.0, 70.0)); + assert_eq!( + BBox::merge(a.clone(), b).unwrap(), + Some(bbox(-5.0, 10.0, 50.0, 70.0)) + ); + assert_eq!(BBox::merge(a.clone(), None).unwrap(), a); + assert_eq!(BBox::merge(None, a.clone()).unwrap(), a); + assert_eq!(BBox::merge(None, None).unwrap(), None); + } + + #[test] + fn test_merge_bbox_crs_mismatch() { + let a = Some(BBox::from_array([0.0, 0.0, 1.0, 1.0], "EPSG:4326")); + let b = Some(BBox::from_array([0.0, 0.0, 1.0, 1.0], "EPSG:3857")); + assert!(BBox::merge(a, b).is_err()); + } + + #[test] + fn test_clamp() { + // restricts values that exceed bounds + let b = BBox::from_array([-200.0, -100.0, 200.0, 100.0], "EPSG:4326"); + assert_eq!( + b.clamp(-180.0, -90.0, 180.0, 90.0), + bbox(-180.0, -90.0, 180.0, 90.0) + ); + + // no-op when already within bounds + let b = bbox(10.0, 20.0, 30.0, 40.0); + assert_eq!( + b.clamp(-180.0, -90.0, 180.0, 90.0), + bbox(10.0, 20.0, 30.0, 40.0) + ); + } + + #[test] + fn test_graticule_breaks_world() { + let breaks = graticule_breaks((-180.0, 180.0)); + assert_eq!( + breaks, + vec![-180.0, -135.0, -90.0, -45.0, 0.0, 45.0, 90.0, 135.0] + ); + } + + #[test] + fn test_graticule_breaks_hemisphere() { + let breaks = graticule_breaks((-88.0, 88.0)); + assert_eq!(breaks, vec![-60.0, -30.0, 0.0, 30.0, 60.0]); + } + + #[test] + fn test_graticule_breaks_small_range() { + let breaks = graticule_breaks((10.0, 20.0)); + assert!(!breaks.is_empty()); + for &b in &breaks { + assert!(b > 10.0 && b < 20.0); + } + } + + #[test] + fn test_graticule_breaks_empty_for_zero_range() { + let breaks = graticule_breaks((50.0, 50.0)); + assert!(breaks.is_empty()); + } + + #[test] + fn test_grid_lines_wkt_meridians() { + let wkt = grid_lines_wkt(&[0.0, 30.0], (-90.0, 90.0), 45.0, true); + assert!(wkt.starts_with("MULTILINESTRING("), "{wkt}"); + assert!(wkt.contains("0.000000 -90.000000"), "{wkt}"); + assert!(wkt.contains("30.000000 -90.000000"), "{wkt}"); + assert!(wkt.contains("0.000000 90.000000"), "{wkt}"); + assert!(wkt.contains("30.000000 90.000000"), "{wkt}"); + } + + #[test] + fn test_grid_lines_wkt_parallels() { + let wkt = grid_lines_wkt(&[0.0, 45.0], (-180.0, 180.0), 90.0, false); + assert!(wkt.starts_with("MULTILINESTRING(")); + assert!(wkt.contains("0.000000")); + assert!(wkt.contains("45.000000")); + } +} diff --git a/src/plot/projection/coord/map_projections.rs b/src/plot/projection/coord/map_projections.rs new file mode 100644 index 000000000..aed1b2468 --- /dev/null +++ b/src/plot/projection/coord/map_projections.rs @@ -0,0 +1,1247 @@ +use serde::{Deserialize, Serialize}; +use std::collections::HashMap; +use std::fmt; +use std::sync::Arc; + +use crate::plot::types::ArrayElement; +use crate::plot::ParameterValue; + +// When adding a new named projection: +// 1. Add a struct implementing MapProjectionTrait (in this file) +// 2. Add the name to NAMED_PROJECTIONS below +// 3. Add match arms in MapSpecification::new(): +// - the `crs` branch (maps +proj= code to the new struct) +// - the `else` branch (maps coord type name to the new struct) +pub const NAMED_PROJECTIONS: &[&str] = &[ + "geographic", + "mercator", + "orthographic", + "miller", + "equirectangular", + "stereographic", + "gnomonic", + "equal_area", + "mollweide", + "sinusoidal", + "eckert4", + "natural", + "winkel_tripel", + "albers", + "lambert_conformal", + "lambert", + "azimuthal_equidistant", + "igh", + "robinson", +]; + +// ============================================================================= +// Trait +// ============================================================================= + +pub trait MapProjectionTrait: fmt::Debug + Send + Sync { + fn proj_code(&self) -> &'static str; + fn display_name(&self) -> &'static str; + fn origin(&self) -> (f64, f64); + fn to_proj_str(&self) -> String; + + fn visible_area_wkt(&self) -> Option { + Some(self.clip_shape_wkt()) + } + + fn clip_shape_wkt(&self) -> String { + rectangle_wkt( + -180.0, + self.lat_bounds().0, + 180.0, + self.lat_bounds().1, + self.edge_segments(), + ) + } + + fn lat_bounds(&self) -> (f64, f64) { + (-90.0, 90.0) + } + + fn edge_segments(&self) -> [usize; 4] { + [1, 36, 1, 36] + } + + fn slit_wkt(&self, epsilon: f64) -> Option { + let seam = wrap_lon(self.origin().0 + 180.0); + if (seam - (-180.0)).abs() > epsilon && (seam - 180.0).abs() > epsilon { + let segs = self.edge_segments()[1]; + Some(rectangle_wkt( + seam - epsilon, + -90.0, + seam + epsilon, + 90.0, + [1, segs, 1, segs], + )) + } else { + None + } + } +} + +// ============================================================================= +// Wrapper +// ============================================================================= + +#[derive(Clone)] +pub struct MapSpecification(Arc); + +impl MapSpecification { + pub fn new(name: Option<&str>, properties: &HashMap) -> Option { + let name = name?; + + let obj: Arc = + if let Some(ParameterValue::String(crs)) = properties.get("crs") { + let code = extract_proj_param_str(crs, "+proj=").unwrap_or(""); + let lon_0 = extract_f64_param(crs, "+lon_0=").unwrap_or(0.0); + let lat_0 = extract_f64_param(crs, "+lat_0=").unwrap_or(0.0); + match code { + "longlat" | "latlong" => Arc::new(Geographic { lon_0 }), + "ortho" => Arc::new(Orthographic { lon_0, lat_0 }), + "stere" => Arc::new(Stereographic { lon_0, lat_0 }), + "gnom" => Arc::new(Gnomonic { lon_0, lat_0 }), + "laea" => Arc::new(LambertAzimuthal { lon_0, lat_0 }), + "aeqd" => Arc::new(AzimuthalEquidistant { lon_0, lat_0 }), + "merc" => Arc::new(Mercator { lon_0 }), + "mill" => Arc::new(Miller { lon_0 }), + "eqc" => Arc::new(Equirectangular { lon_0 }), + "cea" => Arc::new(CylindricalEqualArea { lon_0 }), + "robin" => Arc::new(Robinson { lon_0 }), + "moll" => Arc::new(Mollweide { lon_0 }), + "sinu" => Arc::new(Sinusoidal { lon_0 }), + "eck4" => Arc::new(Eckert4 { lon_0 }), + "natearth" => Arc::new(NaturalEarth { lon_0 }), + "igh" => Arc::new(Igh { lon_0 }), + "wintri" => Arc::new(WinkelTripel { lon_0 }), + "aea" => Arc::new(AlbersEqualArea { + lon_0, + lat_0, + lat_1: extract_f64_param(crs, "+lat_1=").unwrap_or(29.5), + lat_2: extract_f64_param(crs, "+lat_2=").unwrap_or(45.5), + }), + "lcc" => Arc::new(LambertConformalConic { + lon_0, + lat_0, + lat_1: extract_f64_param(crs, "+lat_1=").unwrap_or(29.5), + lat_2: extract_f64_param(crs, "+lat_2=").unwrap_or(45.5), + }), + _ => Arc::new(UnknownProj { + raw: crs.to_string(), + lon_0, + lat_0, + }), + } + } else { + // Extract origin: number (lon only) or array (lon, lat) + let (lon_0, lat_0) = match properties.get("origin") { + Some(ParameterValue::Number(lon)) => (*lon, 0.0), + Some(ParameterValue::Array(arr)) => { + let lon = match arr.first() { + Some(ArrayElement::Number(n)) => *n, + _ => 0.0, + }; + let lat = match arr.get(1) { + Some(ArrayElement::Number(n)) => *n, + _ => 0.0, + }; + (lon, lat) + } + _ => (0.0, 0.0), + }; + + // Extract parallel: number (tangent) or array (secant) + let (lat_1, lat_2) = match properties.get("parallel") { + Some(ParameterValue::Number(lat)) => (*lat, *lat), + Some(ParameterValue::Array(arr)) => { + let l1 = match arr.first() { + Some(ArrayElement::Number(n)) => *n, + _ => 29.5, + }; + let l2 = match arr.get(1) { + Some(ArrayElement::Number(n)) => *n, + _ => 45.5, + }; + (l1, l2) + } + _ => (29.5, 45.5), + }; + + match name { + "geographic" => Arc::new(Geographic { lon_0 }), + "mercator" => Arc::new(Mercator { lon_0 }), + "orthographic" => Arc::new(Orthographic { lon_0, lat_0 }), + "miller" => Arc::new(Miller { lon_0 }), + "equirectangular" => Arc::new(Equirectangular { lon_0 }), + "stereographic" => Arc::new(Stereographic { lon_0, lat_0 }), + "gnomonic" => Arc::new(Gnomonic { lon_0, lat_0 }), + "equal_area" => Arc::new(CylindricalEqualArea { lon_0 }), + "mollweide" => Arc::new(Mollweide { lon_0 }), + "sinusoidal" => Arc::new(Sinusoidal { lon_0 }), + "eckert4" => Arc::new(Eckert4 { lon_0 }), + "natural" => Arc::new(NaturalEarth { lon_0 }), + "winkel_tripel" => Arc::new(WinkelTripel { lon_0 }), + "albers" => Arc::new(AlbersEqualArea { + lon_0, + lat_0, + lat_1, + lat_2, + }), + "lambert_conformal" => Arc::new(LambertConformalConic { + lon_0, + lat_0, + lat_1, + lat_2, + }), + "lambert" => Arc::new(LambertAzimuthal { lon_0, lat_0 }), + "azimuthal_equidistant" => Arc::new(AzimuthalEquidistant { lon_0, lat_0 }), + "igh" => Arc::new(Igh { lon_0 }), + "robinson" => Arc::new(Robinson { lon_0 }), + "map" => Arc::new(UnknownProj { + raw: String::new(), + lon_0, + lat_0, + }), + _ => return None, + } + }; + Some(Self(obj)) + } + + pub fn from_proj_str(crs: &str) -> Self { + let mut properties = HashMap::new(); + properties.insert("crs".to_string(), ParameterValue::String(crs.to_string())); + Self::new(Some("map"), &properties).unwrap() + } + + pub fn proj_code(&self) -> &'static str { + self.0.proj_code() + } + + pub fn display_name(&self) -> &'static str { + self.0.display_name() + } + + pub fn origin(&self) -> (f64, f64) { + self.0.origin() + } + + pub fn to_proj_str(&self) -> String { + self.0.to_proj_str() + } + + pub fn visible_area_wkt(&self) -> Option { + self.0.visible_area_wkt() + } + + pub fn slit_wkt(&self, epsilon: f64) -> Option { + self.0.slit_wkt(epsilon) + } +} + +impl fmt::Debug for MapSpecification { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + write!(f, "MapSpecification({})", self.0.to_proj_str()) + } +} + +impl fmt::Display for MapSpecification { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + write!(f, "{}", self.0.display_name()) + } +} + +impl PartialEq for MapSpecification { + fn eq(&self, other: &Self) -> bool { + self.0.to_proj_str() == other.0.to_proj_str() + } +} + +impl Serialize for MapSpecification { + fn serialize(&self, serializer: S) -> Result + where + S: serde::Serializer, + { + self.0.to_proj_str().serialize(serializer) + } +} + +impl<'de> Deserialize<'de> for MapSpecification { + fn deserialize(deserializer: D) -> Result + where + D: serde::Deserializer<'de>, + { + let s = String::deserialize(deserializer)?; + Ok(MapSpecification::from_proj_str(&s)) + } +} + +// ============================================================================= +// Azimuthal projections +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct Orthographic { + pub lon_0: f64, + pub lat_0: f64, +} + +impl MapProjectionTrait for Orthographic { + fn proj_code(&self) -> &'static str { + "ortho" + } + fn display_name(&self) -> &'static str { + "Orthographic" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!("+proj=ortho +lon_0={} +lat_0={}", self.lon_0, self.lat_0) + } + fn clip_shape_wkt(&self) -> String { + hemisphere_polygon_wkt(self.lon_0, self.lat_0, 88.0) + } +} + +#[derive(Debug, Clone)] +pub struct Stereographic { + pub lon_0: f64, + pub lat_0: f64, +} + +impl MapProjectionTrait for Stereographic { + fn proj_code(&self) -> &'static str { + "stere" + } + fn display_name(&self) -> &'static str { + "Stereographic" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!("+proj=stere +lon_0={} +lat_0={}", self.lon_0, self.lat_0) + } + fn clip_shape_wkt(&self) -> String { + hemisphere_polygon_wkt(self.lon_0, self.lat_0, 88.0) + } +} + +#[derive(Debug, Clone)] +pub struct Gnomonic { + pub lon_0: f64, + pub lat_0: f64, +} + +impl MapProjectionTrait for Gnomonic { + fn proj_code(&self) -> &'static str { + "gnom" + } + fn display_name(&self) -> &'static str { + "Gnomonic" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!("+proj=gnom +lon_0={} +lat_0={}", self.lon_0, self.lat_0) + } + fn clip_shape_wkt(&self) -> String { + hemisphere_polygon_wkt(self.lon_0, self.lat_0, 60.0) + } +} + +#[derive(Debug, Clone)] +pub struct LambertAzimuthal { + pub lon_0: f64, + pub lat_0: f64, +} + +impl MapProjectionTrait for LambertAzimuthal { + fn proj_code(&self) -> &'static str { + "laea" + } + fn display_name(&self) -> &'static str { + "Lambert Azimuthal Equal-Area" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!("+proj=laea +lon_0={} +lat_0={}", self.lon_0, self.lat_0) + } + fn visible_area_wkt(&self) -> Option { + todo!("full-globe azimuthal visible area") + } +} + +#[derive(Debug, Clone)] +pub struct AzimuthalEquidistant { + pub lon_0: f64, + pub lat_0: f64, +} + +impl MapProjectionTrait for AzimuthalEquidistant { + fn proj_code(&self) -> &'static str { + "aeqd" + } + fn display_name(&self) -> &'static str { + "Azimuthal Equidistant" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!("+proj=aeqd +lon_0={} +lat_0={}", self.lon_0, self.lat_0) + } + fn visible_area_wkt(&self) -> Option { + todo!("full-globe azimuthal visible area") + } +} + +// ============================================================================= +// Geographic (unprojected) +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct Geographic { + pub lon_0: f64, +} + +impl MapProjectionTrait for Geographic { + fn proj_code(&self) -> &'static str { + "longlat" + } + fn display_name(&self) -> &'static str { + "Geographic" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=longlat +lon_0={}", self.lon_0) + } + fn edge_segments(&self) -> [usize; 4] { + [1, 1, 1, 1] + } +} + +// ============================================================================= +// Cylindrical projections +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct Mercator { + pub lon_0: f64, +} + +impl MapProjectionTrait for Mercator { + fn proj_code(&self) -> &'static str { + "merc" + } + fn display_name(&self) -> &'static str { + "Mercator" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=merc +lon_0={}", self.lon_0) + } + fn lat_bounds(&self) -> (f64, f64) { + (-85.0, 85.0) + } + fn edge_segments(&self) -> [usize; 4] { + [1, 1, 1, 1] + } +} + +#[derive(Debug, Clone)] +pub struct Miller { + pub lon_0: f64, +} + +impl MapProjectionTrait for Miller { + fn proj_code(&self) -> &'static str { + "mill" + } + fn display_name(&self) -> &'static str { + "Miller" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=mill +lon_0={}", self.lon_0) + } + fn edge_segments(&self) -> [usize; 4] { + [1, 1, 1, 1] + } +} + +#[derive(Debug, Clone)] +pub struct Equirectangular { + pub lon_0: f64, +} + +impl MapProjectionTrait for Equirectangular { + fn proj_code(&self) -> &'static str { + "eqc" + } + fn display_name(&self) -> &'static str { + "Equirectangular" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=eqc +lon_0={}", self.lon_0) + } + fn edge_segments(&self) -> [usize; 4] { + [1, 1, 1, 1] + } +} + +#[derive(Debug, Clone)] +pub struct CylindricalEqualArea { + pub lon_0: f64, +} + +impl MapProjectionTrait for CylindricalEqualArea { + fn proj_code(&self) -> &'static str { + "cea" + } + fn display_name(&self) -> &'static str { + "Cylindrical Equal-Area" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=cea +lon_0={}", self.lon_0) + } + fn edge_segments(&self) -> [usize; 4] { + [1, 1, 1, 1] + } +} + +// ============================================================================= +// Pseudocylindrical projections +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct Robinson { + pub lon_0: f64, +} + +impl MapProjectionTrait for Robinson { + fn proj_code(&self) -> &'static str { + "robin" + } + fn display_name(&self) -> &'static str { + "Robinson" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=robin +lon_0={}", self.lon_0) + } +} + +#[derive(Debug, Clone)] +pub struct Mollweide { + pub lon_0: f64, +} + +impl MapProjectionTrait for Mollweide { + fn proj_code(&self) -> &'static str { + "moll" + } + fn display_name(&self) -> &'static str { + "Mollweide" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=moll +lon_0={}", self.lon_0) + } +} + +#[derive(Debug, Clone)] +pub struct Sinusoidal { + pub lon_0: f64, +} + +impl MapProjectionTrait for Sinusoidal { + fn proj_code(&self) -> &'static str { + "sinu" + } + fn display_name(&self) -> &'static str { + "Sinusoidal" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=sinu +lon_0={}", self.lon_0) + } +} + +#[derive(Debug, Clone)] +pub struct Eckert4 { + pub lon_0: f64, +} + +impl MapProjectionTrait for Eckert4 { + fn proj_code(&self) -> &'static str { + "eck4" + } + fn display_name(&self) -> &'static str { + "Eckert IV" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=eck4 +lon_0={}", self.lon_0) + } +} + +#[derive(Debug, Clone)] +pub struct NaturalEarth { + pub lon_0: f64, +} + +impl MapProjectionTrait for NaturalEarth { + fn proj_code(&self) -> &'static str { + "natearth" + } + fn display_name(&self) -> &'static str { + "Natural Earth" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=natearth +lon_0={}", self.lon_0) + } +} + +// ============================================================================= +// Interrupted +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct Igh { + pub lon_0: f64, +} + +impl MapProjectionTrait for Igh { + fn proj_code(&self) -> &'static str { + "igh" + } + fn display_name(&self) -> &'static str { + "Interrupted Goode Homolosine" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=igh +lon_0={}", self.lon_0) + } + fn slit_wkt(&self, epsilon: f64) -> Option { + Some(igh_slit_wkt(self.lon_0, epsilon)) + } +} + +// ============================================================================= +// Conic projections +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct AlbersEqualArea { + pub lon_0: f64, + pub lat_0: f64, + pub lat_1: f64, + pub lat_2: f64, +} + +impl MapProjectionTrait for AlbersEqualArea { + fn proj_code(&self) -> &'static str { + "aea" + } + fn display_name(&self) -> &'static str { + "Albers Equal-Area" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!( + "+proj=aea +lon_0={} +lat_0={} +lat_1={} +lat_2={}", + self.lon_0, self.lat_0, self.lat_1, self.lat_2 + ) + } + fn edge_segments(&self) -> [usize; 4] { + [36, 36, 36, 36] + } +} + +#[derive(Debug, Clone)] +pub struct LambertConformalConic { + pub lon_0: f64, + pub lat_0: f64, + pub lat_1: f64, + pub lat_2: f64, +} + +impl MapProjectionTrait for LambertConformalConic { + fn proj_code(&self) -> &'static str { + "lcc" + } + fn display_name(&self) -> &'static str { + "Lambert Conformal Conic" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + format!( + "+proj=lcc +lon_0={} +lat_0={} +lat_1={} +lat_2={}", + self.lon_0, self.lat_0, self.lat_1, self.lat_2 + ) + } + fn lat_bounds(&self) -> (f64, f64) { + (-80.0, 84.0) + } + fn edge_segments(&self) -> [usize; 4] { + [36, 36, 36, 36] + } +} + +// ============================================================================= +// Standalone +// ============================================================================= + +#[derive(Debug, Clone)] +pub struct WinkelTripel { + pub lon_0: f64, +} + +impl MapProjectionTrait for WinkelTripel { + fn proj_code(&self) -> &'static str { + "wintri" + } + fn display_name(&self) -> &'static str { + "Winkel Tripel" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, 0.0) + } + fn to_proj_str(&self) -> String { + format!("+proj=wintri +lon_0={}", self.lon_0) + } + fn edge_segments(&self) -> [usize; 4] { + [36, 36, 36, 36] + } +} + +// ============================================================================= +// Unknown / fallback +// ============================================================================= + +#[derive(Debug, Clone)] +struct UnknownProj { + raw: String, + lon_0: f64, + lat_0: f64, +} + +impl MapProjectionTrait for UnknownProj { + fn proj_code(&self) -> &'static str { + "unknown" + } + fn display_name(&self) -> &'static str { + "Unknown" + } + fn origin(&self) -> (f64, f64) { + (self.lon_0, self.lat_0) + } + fn to_proj_str(&self) -> String { + self.raw.clone() + } + fn visible_area_wkt(&self) -> Option { + None + } +} + +// ============================================================================= +// Helpers +// ============================================================================= + +pub fn wrap_lon(lon: f64) -> f64 { + ((lon + 180.0) % 360.0 + 360.0) % 360.0 - 180.0 +} + +pub fn extract_proj_param_str<'a>(crs: &'a str, key: &str) -> Option<&'a str> { + let start = crs.find(key)?; + let after = &crs[start + key.len()..]; + let end = after.find([' ', '+']).unwrap_or(after.len()); + Some(&after[..end]) +} + +fn extract_f64_param(crs: &str, key: &str) -> Option { + extract_proj_param_str(crs, key).and_then(|s| s.parse().ok()) +} + +pub fn rectangle_wkt(xmin: f64, ymin: f64, xmax: f64, ymax: f64, segments: [usize; 4]) -> String { + let mut coords: Vec = Vec::new(); + let [top, right, bottom, left] = segments.map(|s| s.max(1)); + for i in 0..top { + let t = i as f64 / top as f64; + coords.push(format!("{:.6} {:.6}", xmin + t * (xmax - xmin), ymax)); + } + for i in 0..right { + let t = i as f64 / right as f64; + coords.push(format!("{:.6} {:.6}", xmax, ymax - t * (ymax - ymin))); + } + for i in 0..bottom { + let t = i as f64 / bottom as f64; + coords.push(format!("{:.6} {:.6}", xmax - t * (xmax - xmin), ymin)); + } + for i in 0..left { + let t = i as f64 / left as f64; + coords.push(format!("{:.6} {:.6}", xmin, ymin + t * (ymax - ymin))); + } + coords.push(format!("{:.6} {:.6}", xmin, ymax)); + format!("POLYGON(({}))", coords.join(", ")) +} + +fn igh_slit_wkt(lon_0: f64, half_width: f64) -> String { + let segs = [1, 36, 1, 36]; + let polygon_ring = |wkt: String| -> String { wkt.strip_prefix("POLYGON").unwrap().to_string() }; + + let mut parts = Vec::new(); + + for offset in [80.0, -20.0, -100.0] { + let lon = wrap_lon(lon_0 + offset); + parts.push(polygon_ring(rectangle_wkt( + lon - half_width, + -90.0, + lon + half_width, + 0.0, + segs, + ))); + } + + let north = wrap_lon(lon_0 - 40.0); + parts.push(polygon_ring(rectangle_wkt( + north - half_width, + 0.0, + north + half_width, + 90.0, + segs, + ))); + + let seam = wrap_lon(lon_0 + 180.0); + if (seam - (-180.0)).abs() > half_width && (seam - 180.0).abs() > half_width { + parts.push(polygon_ring(rectangle_wkt( + seam - half_width, + -90.0, + seam + half_width, + 90.0, + segs, + ))); + } + + format!("MULTIPOLYGON({})", parts.join(", ")) +} + +pub fn hemisphere_polygon_wkt(lon0: f64, lat0: f64, radius_deg: f64) -> String { + let d = radius_deg.to_radians(); + let lat0_r = lat0.to_radians(); + let sin_lat0 = lat0_r.sin(); + let cos_lat0 = lat0_r.cos(); + let sin_d = d.sin(); + let cos_d = d.cos(); + + let n_points = 72; + let mut raw_points: Vec<(f64, f64)> = Vec::with_capacity(n_points); + for i in 0..n_points { + let az = (i as f64 * (360.0 / n_points as f64)).to_radians(); + let lat2 = (sin_lat0 * cos_d + cos_lat0 * sin_d * az.cos()).asin(); + let lon2 = + lon0.to_radians() + (az.sin() * sin_d * cos_lat0).atan2(cos_d - sin_lat0 * lat2.sin()); + let mut lon_deg = lon2.to_degrees(); + lon_deg = ((lon_deg + 180.0) % 360.0 + 360.0) % 360.0 - 180.0; + raw_points.push((lon_deg, lat2.to_degrees())); + } + + let mut points: Vec<(f64, f64)> = Vec::with_capacity(n_points + 4); + for i in 0..raw_points.len() { + points.push(raw_points[i]); + let next = (i + 1) % raw_points.len(); + if (raw_points[next].0 - raw_points[i].0).abs() > 180.0 { + let lat = antimeridian_crossing_lat(raw_points[i], raw_points[next]); + if raw_points[i].0 > 0.0 { + points.push((179.999999, lat)); + points.push((-179.999999, lat)); + } else { + points.push((-179.999999, lat)); + points.push((179.999999, lat)); + } + } + } + + let includes_north_pole = lat0 + radius_deg > 90.0; + let includes_south_pole = lat0 - radius_deg < -90.0; + + if includes_north_pole || includes_south_pole { + build_pole_polygon(&points, includes_north_pole) + } else if find_antimeridian_crossings(&points).len() == 2 { + build_antimeridian_multipolygon(&points) + } else { + build_simple_polygon(&points) + } +} + +fn build_simple_polygon(points: &[(f64, f64)]) -> String { + let mut coords: Vec = points + .iter() + .map(|(lon, lat)| format!("{lon:.6} {lat:.6}")) + .collect(); + coords.push(coords[0].clone()); + format!("POLYGON(({}))", coords.join(", ")) +} + +fn build_pole_polygon(points: &[(f64, f64)], north: bool) -> String { + let mut split_idx = 0; + let mut max_jump = 0.0_f64; + for i in 0..points.len() { + let next = (i + 1) % points.len(); + let jump = (points[next].0 - points[i].0).abs(); + if jump > max_jump { + max_jump = jump; + split_idx = next; + } + } + + let mut ordered: Vec<(f64, f64)> = Vec::with_capacity(points.len()); + for i in 0..points.len() { + ordered.push(points[(split_idx + i) % points.len()]); + } + + let pole_lat = if north { 90.0 } else { -90.0 }; + let first = ordered.first().unwrap(); + let last = ordered.last().unwrap(); + + let mut coords: Vec = Vec::with_capacity(points.len() + 6); + for (lon, lat) in &ordered { + coords.push(format!("{lon:.6} {lat:.6}")); + } + coords.push(format!("{:.6} {pole_lat:.6}", last.0)); + if (last.0 - first.0).abs() > 180.0 { + let mid = (last.0 + first.0) / 2.0; + coords.push(format!("{mid:.6} {pole_lat:.6}")); + } + coords.push(format!("{:.6} {pole_lat:.6}", first.0)); + coords.push(format!("{:.6} {:.6}", first.0, first.1)); + + format!("POLYGON(({}))", coords.join(", ")) +} + +fn find_antimeridian_crossings(points: &[(f64, f64)]) -> Vec { + let n = points.len(); + let mut crossings = Vec::new(); + for i in 0..n { + let next = (i + 1) % n; + if (points[next].0 - points[i].0).abs() > 180.0 { + crossings.push(i); + } + } + crossings +} + +fn build_antimeridian_multipolygon(points: &[(f64, f64)]) -> String { + let n = points.len(); + let crossings = find_antimeridian_crossings(points); + assert_eq!(crossings.len(), 2); + + let c1 = crossings[0]; + let c2 = crossings[1]; + + let lat_c1 = antimeridian_crossing_lat(points[c1], points[(c1 + 1) % n]); + let lat_c2 = antimeridian_crossing_lat(points[c2], points[(c2 + 1) % n]); + + let (east_arc, west_arc, [east_start_lat, east_end_lat, west_start_lat, west_end_lat]) = + split_arcs_at_crossings(points, c1, c2, lat_c1, lat_c2); + + let east_coords = build_side_ring(&east_arc, 180.0, east_start_lat, east_end_lat); + let west_coords = build_side_ring(&west_arc, -180.0, west_start_lat, west_end_lat); + + format!( + "MULTIPOLYGON((({})),(({})))", + east_coords.join(", "), + west_coords.join(", ") + ) +} + +type ArcSplit = (Vec<(f64, f64)>, Vec<(f64, f64)>, [f64; 4]); + +fn split_arcs_at_crossings( + points: &[(f64, f64)], + c1: usize, + c2: usize, + lat_c1: f64, + lat_c2: f64, +) -> ArcSplit { + let n = points.len(); + + let mut arc1: Vec<(f64, f64)> = Vec::new(); + let mut i = (c1 + 1) % n; + loop { + arc1.push(points[i]); + if i == c2 { + break; + } + i = (i + 1) % n; + } + + let mut arc2: Vec<(f64, f64)> = Vec::new(); + i = (c2 + 1) % n; + loop { + arc2.push(points[i]); + if i == c1 { + break; + } + i = (i + 1) % n; + } + + if arc1[0].0 > 0.0 { + (arc1, arc2, [lat_c1, lat_c2, lat_c2, lat_c1]) + } else { + (arc2, arc1, [lat_c2, lat_c1, lat_c1, lat_c2]) + } +} + +fn build_side_ring( + arc: &[(f64, f64)], + meridian_lon: f64, + start_lat: f64, + end_lat: f64, +) -> Vec { + let mut coords: Vec = Vec::with_capacity(arc.len() + 3); + coords.push(format!("{meridian_lon:.6} {start_lat:.6}")); + for (lon, lat) in arc.iter() { + coords.push(format!("{lon:.6} {lat:.6}")); + } + coords.push(format!("{meridian_lon:.6} {end_lat:.6}")); + coords.push(coords[0].clone()); + coords +} + +fn antimeridian_crossing_lat(a: (f64, f64), b: (f64, f64)) -> f64 { + let (lon_a, lat_a) = a; + let (lon_b, lat_b) = b; + let (lon_a_u, lon_b_u) = if lon_a > lon_b { + (lon_a, lon_b + 360.0) + } else { + (lon_a + 360.0, lon_b) + }; + let t = (180.0 - lon_a_u) / (lon_b_u - lon_a_u); + lat_a + t * (lat_b - lat_a) +} + +// ============================================================================= +// Tests +// ============================================================================= + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn from_proj_str_known_projections() { + let proj = MapSpecification::from_proj_str("+proj=ortho +lon_0=10 +lat_0=45"); + assert_eq!(proj.proj_code(), "ortho"); + assert_eq!(proj.origin(), (10.0, 45.0)); + + let proj = MapSpecification::from_proj_str("+proj=merc"); + assert_eq!(proj.proj_code(), "merc"); + assert_eq!(proj.origin(), (0.0, 0.0)); + + let proj = MapSpecification::from_proj_str("+proj=aea +lon_0=5 +lat_1=30 +lat_2=50"); + assert_eq!(proj.proj_code(), "aea"); + assert_eq!( + proj.to_proj_str(), + "+proj=aea +lon_0=5 +lat_0=0 +lat_1=30 +lat_2=50" + ); + } + + #[test] + fn new_from_crs_albers() { + let mut props = HashMap::new(); + props.insert( + "crs".to_string(), + ParameterValue::String("+proj=aea +lon_0=5 +lat_0=10 +lat_1=30 +lat_2=50".to_string()), + ); + let proj = MapSpecification::new(Some("map"), &props).unwrap(); + assert_eq!(proj.proj_code(), "aea"); + assert_eq!(proj.origin(), (5.0, 10.0)); + assert_eq!( + proj.to_proj_str(), + "+proj=aea +lon_0=5 +lat_0=10 +lat_1=30 +lat_2=50" + ); + } + + #[test] + fn new_named_albers_with_settings() { + let mut props = HashMap::new(); + props.insert( + "origin".to_string(), + ParameterValue::Array(vec![ + ArrayElement::Number(-96.0), + ArrayElement::Number(37.5), + ]), + ); + props.insert( + "parallel".to_string(), + ParameterValue::Array(vec![ArrayElement::Number(29.5), ArrayElement::Number(45.5)]), + ); + let proj = MapSpecification::new(Some("albers"), &props).unwrap(); + assert_eq!(proj.proj_code(), "aea"); + assert_eq!(proj.origin(), (-96.0, 37.5)); + assert_eq!( + proj.to_proj_str(), + "+proj=aea +lon_0=-96 +lat_0=37.5 +lat_1=29.5 +lat_2=45.5" + ); + } + + #[test] + fn from_proj_str_unknown_projection() { + let proj = MapSpecification::from_proj_str("+proj=fooproj +lon_0=5"); + assert_eq!(proj.proj_code(), "unknown"); + assert_eq!(proj.origin(), (5.0, 0.0)); + assert_eq!(proj.visible_area_wkt(), None); + } + + #[test] + fn round_trip_serialization() { + let proj = MapSpecification::from_proj_str("+proj=robin +lon_0=10"); + let json = serde_json::to_string(&proj).unwrap(); + let deser: MapSpecification = serde_json::from_str(&json).unwrap(); + assert_eq!(proj, deser); + } + + #[test] + fn visible_area_cylindrical() { + let proj = MapSpecification::from_proj_str("+proj=merc"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!(wkt.starts_with("POLYGON((")); + assert!(wkt.contains("85.000000")); + } + + #[test] + fn visible_area_azimuthal() { + let proj = MapSpecification::from_proj_str("+proj=ortho +lon_0=0 +lat_0=0"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!(wkt.starts_with("POLYGON((") || wkt.starts_with("MULTIPOLYGON(")); + } + + #[test] + fn slit_igh() { + let proj = MapSpecification::from_proj_str("+proj=igh"); + let slit = proj.slit_wkt(0.005).unwrap(); + assert!(slit.starts_with("MULTIPOLYGON(")); + } + + #[test] + fn slit_default_antimeridian() { + let proj = MapSpecification::from_proj_str("+proj=robin +lon_0=10"); + let slit = proj.slit_wkt(0.005).unwrap(); + assert!(slit.starts_with("POLYGON((")); + assert!(slit.contains("-170.")); + } + + #[test] + fn slit_at_dateline_returns_none() { + let proj = MapSpecification::from_proj_str("+proj=robin +lon_0=0"); + assert!(proj.slit_wkt(0.005).is_none()); + } + + #[test] + fn visible_area_gnomonic() { + let proj = MapSpecification::from_proj_str("+proj=gnom +lat_0=90 +lon_0=0"); + assert!(proj.visible_area_wkt().is_some()); + } + + #[test] + fn visible_area_antimeridian_crossing() { + let proj = MapSpecification::from_proj_str("+proj=ortho +lat_0=0 +lon_0=150"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!( + wkt.starts_with("MULTIPOLYGON"), + "lon_0=150 should cross antimeridian: {wkt}" + ); + } + + #[test] + fn visible_area_no_antimeridian_for_centered() { + let proj = MapSpecification::from_proj_str("+proj=ortho +lat_0=0 +lon_0=0"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!( + wkt.starts_with("POLYGON(("), + "lon_0=0 should not cross antimeridian: {wkt}" + ); + } + + #[test] + fn visible_area_pole_routing() { + let proj = MapSpecification::from_proj_str("+proj=ortho +lat_0=52.36 +lon_0=150.90"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!( + wkt.starts_with("POLYGON(("), + "pole case should produce POLYGON: {wkt}" + ); + } + + #[test] + fn visible_area_rectangle_always_polygon() { + let proj = MapSpecification::from_proj_str("+proj=robin +lon_0=-90"); + let wkt = proj.visible_area_wkt().unwrap(); + assert!( + wkt.starts_with("POLYGON(("), + "rectangle projections always produce POLYGON: {wkt}" + ); + } + + #[test] + fn seam_position() { + let proj = MapSpecification::from_proj_str("+proj=robin +lon_0=-90"); + let (lon_0, _) = proj.origin(); + let seam = wrap_lon(lon_0 + 180.0); + assert!((seam - 90.0).abs() < 1e-6, "seam should be at 90°"); + } + + #[test] + fn igh_slit_shift_with_lon_0() { + let igh0 = MapSpecification::from_proj_str("+proj=igh"); + let wkt0 = igh0.slit_wkt(0.005).unwrap(); + assert!(wkt0.starts_with("MULTIPOLYGON("), "{wkt0}"); + assert_eq!(wkt0.matches("((").count(), 4, "{wkt0}"); + + let igh90 = MapSpecification::from_proj_str("+proj=igh +lon_0=90"); + let wkt90 = igh90.slit_wkt(0.005).unwrap(); + assert_eq!(wkt90.matches("((").count(), 5, "{wkt90}"); + assert!(wkt90.contains("170.005"), "south slit near 170°: {wkt90}"); + assert!(wkt90.contains("70.005"), "south slit near 70°: {wkt90}"); + assert!(wkt90.contains("-10.005"), "south slit near -10°: {wkt90}"); + assert!(wkt90.contains("50.005"), "north slit near 50°: {wkt90}"); + } +} diff --git a/src/plot/projection/coord/mod.rs b/src/plot/projection/coord/mod.rs index 05d21b831..308d4193d 100644 --- a/src/plot/projection/coord/mod.rs +++ b/src/plot/projection/coord/mod.rs @@ -25,14 +25,20 @@ use std::collections::HashMap; use std::sync::Arc; use crate::plot::types::{validate_parameter, ParamDefinition}; -use crate::plot::ParameterValue; +use crate::plot::{Layer, ParameterValue}; +use crate::reader::SqlDialect; +use crate::DataFrame; // Coord type implementations mod cartesian; +pub mod map; +pub mod map_projections; mod polar; // Re-export coord type structs pub use cartesian::Cartesian; +pub use map::Map; +pub use map_projections::MapSpecification; pub use polar::Polar; // ============================================================================= @@ -47,6 +53,8 @@ pub enum CoordKind { Cartesian, /// Polar coordinates (for pie charts, rose plots) Polar, + /// Map coordinates (for geographic/cartographic projections) + Map, } // ============================================================================= @@ -122,6 +130,27 @@ pub trait CoordTrait: std::fmt::Debug + std::fmt::Display + Send + Sync { Ok(resolved) } + + /// Orchestrate projection transforms for all layers. + /// + /// Iterates layers and calls each geom's `apply_projection()`. + /// Override to add coord-specific setup (e.g., Map loads the spatial extension). + fn apply_projection_transforms( + &self, + layers: &[Layer], + layer_queries: &mut [String], + projection: &mut super::Projection, + dialect: &dyn SqlDialect, + _execute_query: &dyn Fn(&str) -> crate::Result, + ) -> crate::Result<()> { + for (idx, layer) in layers.iter().enumerate() { + layer_queries[idx] = + layer + .geom + .apply_projection(&layer_queries[idx], projection, dialect, false)?; + } + Ok(()) + } } // ============================================================================= @@ -146,11 +175,17 @@ impl Coord { Self(Arc::new(Polar)) } + /// Create a Map coord type + pub fn map(name: &str) -> Self { + Self(Arc::new(Map::new(name))) + } + /// Create a Coord from a CoordKind pub fn from_kind(kind: CoordKind) -> Self { match kind { CoordKind::Cartesian => Self::cartesian(), CoordKind::Polar => Self::polar(), + CoordKind::Map => Self::map("map"), } } @@ -182,6 +217,24 @@ impl Coord { ) -> Result, String> { self.0.resolve_properties(properties) } + + /// Orchestrate projection transforms for all layers. + pub fn apply_projection_transforms( + &self, + layers: &[Layer], + layer_queries: &mut [String], + projection: &mut super::Projection, + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, + ) -> crate::Result<()> { + self.0.apply_projection_transforms( + layers, + layer_queries, + projection, + dialect, + execute_query, + ) + } } impl std::fmt::Debug for Coord { diff --git a/src/plot/projection/resolve.rs b/src/plot/projection/resolve.rs index e08646e73..0eabc98b5 100644 --- a/src/plot/projection/resolve.rs +++ b/src/plot/projection/resolve.rs @@ -8,6 +8,7 @@ use std::collections::HashMap; use super::coord::{Coord, CoordKind}; use super::Projection; use crate::plot::aesthetic::{MATERIAL_AESTHETICS, POSITION_SUFFIXES}; +use crate::plot::layer::GeomType; use crate::plot::scale::ScaleTypeKind; use crate::plot::{Mappings, ParameterValue, Scale}; use crate::GgsqlError; @@ -18,13 +19,17 @@ const CARTESIAN_PRIMARIES: &[&str] = &["x", "y"]; /// Polar primary aesthetic names const POLAR_PRIMARIES: &[&str] = &["angle", "radius"]; +/// Map primary aesthetic names +const MAP_PRIMARIES: &[&str] = &["lon", "lat"]; + /// Resolve coordinate system for a Plot /// /// If `project` is `Some`, returns `Ok(None)` (keep existing, no changes needed). -/// If `project` is `None`, infers coord from aesthetic mappings: +/// If `project` is `None`, infers coord from aesthetic mappings and layer types: /// - x/y/xmin/xmax/ymin/ymax → Cartesian /// - angle/radius/anglemin/... → Polar -/// - Both → Error +/// - Any Spatial layer → Map +/// - Both cartesian+polar → Error /// - Neither → Ok(None) (caller should use default Cartesian) /// /// Called early in the pipeline, before AestheticContext construction. @@ -32,36 +37,79 @@ pub fn resolve_coord( project: Option<&Projection>, global_mappings: &Mappings, layer_mappings: &[&Mappings], + layer_geom_types: &[GeomType], ) -> Result, String> { // If project is explicitly specified, keep it as-is if project.is_some() { return Ok(None); } + // Check if any layer is spatial + let mut found_map = layer_geom_types.contains(&GeomType::Spatial); + // Collect all explicit aesthetic keys from global and layer mappings let mut found_cartesian = false; let mut found_polar = false; // Check global mappings for aesthetic in global_mappings.aesthetics.keys() { - check_aesthetic(aesthetic, &mut found_cartesian, &mut found_polar); + check_aesthetic( + aesthetic, + &mut found_cartesian, + &mut found_polar, + &mut found_map, + ); } // Check layer mappings for layer_map in layer_mappings { for aesthetic in layer_map.aesthetics.keys() { - check_aesthetic(aesthetic, &mut found_cartesian, &mut found_polar); + check_aesthetic( + aesthetic, + &mut found_cartesian, + &mut found_polar, + &mut found_map, + ); } } - // Determine result - if found_cartesian && found_polar { - return Err( - "Conflicting aesthetics: cannot use both cartesian (x/y) and polar (angle/radius) \ - aesthetics in the same plot. Use PROJECT TO cartesian or PROJECT TO polar to \ - specify the coordinate system explicitly." - .to_string(), - ); + // Determine result — check for conflicts + let conflict_count = [found_cartesian, found_polar, found_map] + .iter() + .filter(|&&v| v) + .count(); + if conflict_count > 1 { + let mut systems = Vec::new(); + if found_cartesian { + systems.push("cartesian (x/y)"); + } + if found_polar { + systems.push("polar (angle/radius)"); + } + if found_map { + systems.push("map (lon/lat)"); + } + return Err(format!( + "Conflicting aesthetics: cannot mix {} aesthetics in the same plot. \ + Use PROJECT TO specify the coordinate system explicitly.", + systems.join(" and "), + )); + } + + if found_map { + let coord = Coord::from_kind(CoordKind::Map); + let aesthetics = coord + .position_aesthetic_names() + .iter() + .map(|s| s.to_string()) + .collect(); + return Ok(Some(Projection { + coord, + aesthetics, + properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), + })); } if found_polar { @@ -76,6 +124,8 @@ pub fn resolve_coord( coord, aesthetics, properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), })); } @@ -91,6 +141,8 @@ pub fn resolve_coord( coord, aesthetics, properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), })); } @@ -98,9 +150,14 @@ pub fn resolve_coord( Ok(None) } -/// Check if an aesthetic name indicates cartesian or polar coordinate system. +/// Check if an aesthetic name indicates a coordinate system. /// Updates the found flags accordingly. -fn check_aesthetic(aesthetic: &str, found_cartesian: &mut bool, found_polar: &mut bool) { +fn check_aesthetic( + aesthetic: &str, + found_cartesian: &mut bool, + found_polar: &mut bool, + found_map: &mut bool, +) { // Skip material aesthetics (color, size, etc.) if MATERIAL_AESTHETICS.contains(&aesthetic) { return; @@ -118,6 +175,11 @@ fn check_aesthetic(aesthetic: &str, found_cartesian: &mut bool, found_polar: &mu if POLAR_PRIMARIES.contains(&primary) { *found_polar = true; } + + // Check against map primaries + if MAP_PRIMARIES.contains(&primary) { + *found_map = true; + } } /// Strip position suffix from an aesthetic name. @@ -137,6 +199,8 @@ fn strip_position_suffix(name: &str) -> &str { /// - **`radar`** (polar only): When null (auto), sets to `true` if the theta /// (pos2) scale is discrete/ordinal. When explicitly `true`, validates that /// the theta scale is indeed discrete. +/// - **clip boundary** (map only): Computes the visible-area WKT polygon for +/// azimuthal projections and stores it in `computed`. pub fn resolve_projection_properties( project: &mut Projection, scales: &[Scale], @@ -214,7 +278,7 @@ mod tests { let global = mappings_with(&["angle", "radius"]); // Would infer polar let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(Some(&project), &global, &layers); + let result = resolve_coord(Some(&project), &global, &layers, &[]); assert!(result.is_ok()); assert!(result.unwrap().is_none()); // None means keep existing } @@ -228,7 +292,7 @@ mod tests { let global = mappings_with(&["x", "y"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -242,7 +306,7 @@ mod tests { let global = mappings_with(&["xmin", "ymax"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -256,7 +320,7 @@ mod tests { let layer = mappings_with(&["x", "y"]); let layers: Vec<&Mappings> = vec![&layer]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -273,7 +337,7 @@ mod tests { let global = mappings_with(&["angle", "radius"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -287,7 +351,7 @@ mod tests { let global = mappings_with(&["anglemin", "radiusmax"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -301,7 +365,7 @@ mod tests { let layer = mappings_with(&["angle", "radius"]); let layers: Vec<&Mappings> = vec![&layer]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -318,7 +382,7 @@ mod tests { let global = mappings_with(&["color", "size", "fill", "opacity"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); assert!(result.unwrap().is_none()); // Neither cartesian nor polar } @@ -328,7 +392,7 @@ mod tests { let global = mappings_with(&["x", "y", "color", "size"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -345,7 +409,7 @@ mod tests { let global = mappings_with(&["x", "angle"]); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_err()); let err = result.unwrap_err(); assert!(err.contains("Conflicting")); @@ -359,10 +423,49 @@ mod tests { let layer = mappings_with(&["angle"]); let layers: Vec<&Mappings> = vec![&layer]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); + assert!(result.is_err()); + let err = result.unwrap_err(); + assert!(err.contains("Conflicting")); + } + + #[test] + fn test_conflict_cartesian_and_map() { + let global = mappings_with(&["x", "lon"]); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_err()); let err = result.unwrap_err(); assert!(err.contains("Conflicting")); + assert!(err.contains("cartesian")); + assert!(err.contains("map")); + } + + #[test] + fn test_conflict_polar_and_map() { + let global = mappings_with(&["angle", "lat"]); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(None, &global, &layers, &[]); + assert!(result.is_err()); + let err = result.unwrap_err(); + assert!(err.contains("Conflicting")); + assert!(err.contains("polar")); + assert!(err.contains("map")); + } + + #[test] + fn test_conflict_cartesian_and_spatial_layer() { + let global = mappings_with(&["x", "y"]); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(None, &global, &layers, &[GeomType::Spatial]); + assert!(result.is_err()); + let err = result.unwrap_err(); + assert!(err.contains("Conflicting")); + assert!(err.contains("cartesian")); + assert!(err.contains("map")); } // ======================================== @@ -374,7 +477,7 @@ mod tests { let global = Mappings::new(); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); assert!(result.unwrap().is_none()); } @@ -389,7 +492,7 @@ mod tests { global.insert("angle", AestheticValue::standard_column("cat")); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); let inferred = result.unwrap(); assert!(inferred.is_some()); @@ -402,11 +505,71 @@ mod tests { let global = Mappings::with_wildcard(); let layers: Vec<&Mappings> = vec![]; - let result = resolve_coord(None, &global, &layers); + let result = resolve_coord(None, &global, &layers, &[]); assert!(result.is_ok()); assert!(result.unwrap().is_none()); // Wildcard alone doesn't infer coord } + // ======================================== + // Test: Spatial layer infers Map + // ======================================== + + #[test] + fn test_infer_map_from_spatial_layer() { + let global = Mappings::new(); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(None, &global, &layers, &[GeomType::Spatial]); + assert!(result.is_ok()); + let inferred = result.unwrap(); + assert!(inferred.is_some()); + let proj = inferred.unwrap(); + assert_eq!(proj.coord.coord_kind(), CoordKind::Map); + assert_eq!(proj.aesthetics, vec!["lon", "lat"]); + } + + #[test] + fn test_infer_map_from_spatial_among_other_layers() { + let global = Mappings::new(); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord( + None, + &global, + &layers, + &[GeomType::Spatial, GeomType::Point], + ); + assert!(result.is_ok()); + let inferred = result.unwrap(); + assert!(inferred.is_some()); + let proj = inferred.unwrap(); + assert_eq!(proj.coord.coord_kind(), CoordKind::Map); + } + + #[test] + fn test_infer_map_from_lon_lat_aesthetics() { + let global = mappings_with(&["lon", "lat"]); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(None, &global, &layers, &[]); + assert!(result.is_ok()); + let inferred = result.unwrap(); + assert!(inferred.is_some()); + let proj = inferred.unwrap(); + assert_eq!(proj.coord.coord_kind(), CoordKind::Map); + } + + #[test] + fn test_explicit_project_overrides_spatial() { + let project = Projection::cartesian(); + let global = Mappings::new(); + let layers: Vec<&Mappings> = vec![]; + + let result = resolve_coord(Some(&project), &global, &layers, &[GeomType::Spatial]); + assert!(result.is_ok()); + assert!(result.unwrap().is_none()); + } + // ======================================== // Test: Helper functions // ======================================== diff --git a/src/plot/projection/types.rs b/src/plot/projection/types.rs index 9edfeea8c..c69401cf7 100644 --- a/src/plot/projection/types.rs +++ b/src/plot/projection/types.rs @@ -5,8 +5,10 @@ use serde::{Deserialize, Serialize}; use std::collections::HashMap; -use super::coord::Coord; -use crate::plot::ParameterValue; +use super::coord::{Coord, MapSpecification}; +use crate::plot::{Layer, ParameterValue}; +use crate::reader::SqlDialect; +use crate::DataFrame; /// Projection (from PROJECT clause) #[derive(Debug, Clone, PartialEq, Serialize, Deserialize)] @@ -20,6 +22,13 @@ pub struct Projection { pub aesthetics: Vec, /// Projection-specific options pub properties: HashMap, + /// Typed map projection struct (present when coord is Map and crs is set). + #[serde(default, skip_serializing_if = "Option::is_none")] + pub map_projection: Option, + /// Values computed at execution time (e.g., clip boundary WKT). + /// Not user-facing; populated by apply_projection_transforms. + #[serde(default, skip_serializing_if = "HashMap::is_empty")] + pub computed: HashMap, } impl Projection { @@ -33,6 +42,11 @@ impl Projection { Self::with_defaults(Coord::polar()) } + /// Create a default Map projection (lon, lat). + pub fn map() -> Self { + Self::with_defaults(Coord::map("map")) + } + fn with_defaults(coord: Coord) -> Self { let aesthetics = coord .position_aesthetic_names() @@ -43,6 +57,8 @@ impl Projection { coord, aesthetics, properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), } } @@ -51,4 +67,16 @@ impl Projection { pub fn position_names(&self) -> Vec<&str> { self.aesthetics.iter().map(|s| s.as_str()).collect() } + + /// Orchestrate projection transforms for all layers. + pub fn apply_projection_transforms( + &mut self, + layers: &[Layer], + layer_queries: &mut [String], + dialect: &dyn SqlDialect, + execute_query: &dyn Fn(&str) -> crate::Result, + ) -> crate::Result<()> { + let coord = self.coord.clone(); + coord.apply_projection_transforms(layers, layer_queries, self, dialect, execute_query) + } } diff --git a/src/reader/duckdb.rs b/src/reader/duckdb.rs index 7e8bf273e..1cfcf0583 100644 --- a/src/reader/duckdb.rs +++ b/src/reader/duckdb.rs @@ -54,8 +54,8 @@ fn register_builtin_datasets_duckdb(sql: &str, conn: &Connection) -> Result<()> })?; } - // Arrow export in duckdb-rs v1.10502.0 aborts on GEOMETRY columns. - // Cast to binary WKB when loading the world dataset. + // WORKAROUND(duckdb-rs#714): Arrow export aborts on GEOMETRY columns. + // Store geometry as WKB so Arrow transport doesn't crash. // https://github.com/duckdb/duckdb-rs/issues/714 let select_expr = if name == "world" { "* REPLACE (ST_AsWKB(geom) AS geom)" @@ -104,6 +104,14 @@ impl super::SqlDialect for DuckDbDialect { format!("ST_AsWKB({column})") } + fn sql_geometry_bbox(&self, column: &str, from: &str) -> String { + format!( + "SELECT ST_XMin(ext) AS xmin, ST_YMin(ext) AS ymin, \ + ST_XMax(ext) AS xmax, ST_YMax(ext) AS ymax \ + FROM (SELECT ST_Extent_Agg({column}) AS ext FROM {from})" + ) + } + fn sql_spatial_setup(&self) -> Vec { vec!["LOAD spatial".into()] } diff --git a/src/reader/mod.rs b/src/reader/mod.rs index a2ca57bcf..e8376bc6b 100644 --- a/src/reader/mod.rs +++ b/src/reader/mod.rs @@ -122,6 +122,40 @@ pub trait SqlDialect { format!("ST_AsBinary({column})") } + /// WORKAROUND(duckdb-rs#714): Ensures a column is native GEOMETRY type. + /// + /// Geometry columns may arrive as WKB BLOB (because Arrow export crashes on + /// native GEOMETRY, forcing pre-conversion). This normalizes both GEOMETRY + /// and BLOB to GEOMETRY so spatial functions work uniformly. + /// When the duckdb-rs bug is fixed, this method can return `column` as-is. + fn sql_ensure_geometry(&self, column: &str) -> String { + format!("ST_GeomFromWKB(CAST({column} AS BLOB))") + } + + /// SQL expression to transform a geometry from one CRS to another. + /// + /// Default uses `ST_Transform(column, source_crs, target_crs)` which works for DuckDB. + fn sql_st_transform(&self, column: &str, source_crs: &str, target_crs: &str) -> String { + format!( + "ST_Transform({}, '{}', '{}', always_xy := true)", + column, + source_crs.replace('\'', "''"), + target_crs.replace('\'', "''") + ) + } + + /// SQL query that computes the bounding box of a geometry column. + /// + /// Must return a single row with columns `xmin`, `ymin`, `xmax`, `ymax` (DOUBLE). + /// `from` is the table or subquery to aggregate over. + fn sql_geometry_bbox(&self, column: &str, from: &str) -> String { + format!( + "SELECT ST_XMin(ext) AS xmin, ST_YMin(ext) AS ymin, \ + ST_XMax(ext) AS xmax, ST_YMax(ext) AS ymax \ + FROM (SELECT ST_Extent({column}) AS ext FROM {from})" + ) + } + /// SQL statements to run before spatial operations. /// /// Override for backends that need an extension loaded (e.g. DuckDB spatial). diff --git a/src/writer/vegalite/CLAUDE.md b/src/writer/vegalite/CLAUDE.md index 0baae99a0..dc1bcf05f 100644 --- a/src/writer/vegalite/CLAUDE.md +++ b/src/writer/vegalite/CLAUDE.md @@ -377,10 +377,41 @@ grep '__ggsql_source__' /tmp/test.vl.json # Source values Paste the spec into [Vega-Lite Editor](https://vega.github.io/editor/#/url/vega-lite/) to visualize. +## Projection Rendering + +The `projection/` subdirectory handles coord-specific VL output. Each coord type implements `ProjectionRenderer` (defined in `projection/mod.rs`): + +``` +writer/vegalite/projection/ +├── mod.rs ProjectionRenderer trait + factory (get_projection_renderer) +├── cartesian.rs Standard x/y axes with domain/breaks from scales +├── polar.rs Arc marks, theta/radius channels, radial axes +└── map.rs Identity projection for pre-projected spatial data +``` + +### ProjectionRenderer trait + +Two concerns per implementation: + +1. **Channel mapping** — `position_channels()` returns the VL encoding names for pos1/pos2 (e.g. `("x", "y")` for cartesian, `("radius", "theta")` for polar). +2. **Spec transformation** — `transform_layers()` modifies the VL spec after layers are built (e.g. polar converts marks to arcs, map adds an identity projection with scale/translate expressions). + +Additional hooks: `background_layers()` / `foreground_layers()` inject layers before/after the data layers (e.g. map renders the projected clip boundary as a geoshape panel background), and `apply_projection()` orchestrates all of these plus clip propagation. + +### Map projection specifics + +`MapProjection` reads computed values from `Projection.computed` (populated at execution time by the `Map` coord): + +- `panel_boundary` (WKT) → converted to GeoJSON for a geoshape background layer. +- `bbox` ([xmin, ymin, xmax, ymax]) → emits VL `projection.scale` and `projection.translate` expressions that frame the data to the viewport. + +The VL projection is always `{"type": "identity", "reflectY": true}` because coordinates arrive pre-projected from the SQL layer. + ## References - **Main implementation**: `src/writer/vegalite/mod.rs` - **Layer rendering**: `src/writer/vegalite/layer.rs` - **Data unification**: `src/writer/vegalite/data.rs` (`unify_datasets()`) - **Renderer trait**: `src/writer/vegalite/layer.rs` (`GeomRenderer` trait) +- **Projection rendering**: `src/writer/vegalite/projection/mod.rs` (`ProjectionRenderer` trait) - **Example renderers**: `LineRenderer`, `BoxplotRenderer`, `ViolinRenderer` in `layer.rs` diff --git a/src/writer/vegalite/encoding.rs b/src/writer/vegalite/encoding.rs index 92c9b4b56..34d98e92f 100644 --- a/src/writer/vegalite/encoding.rs +++ b/src/writer/vegalite/encoding.rs @@ -965,7 +965,7 @@ pub(super) fn map_aesthetic_name( ) -> String { // For internal position aesthetics, map directly to Vega-Lite channel names // based on coord type (ignoring user-facing names) - if let Some(vl_channel) = super::projection::map_position_to_vegalite(aesthetic, renderer) { + if let Some(vl_channel) = renderer.map_position(aesthetic) { return vl_channel; } @@ -1014,10 +1014,10 @@ impl<'a> RenderContext<'a> { renderer: &dyn super::projection::ProjectionRenderer, aesthetic_context: crate::plot::aesthetic::AestheticContext, ) -> Self { - let pos1 = super::projection::map_position_to_vegalite("pos1", renderer).unwrap(); - let pos1_end = super::projection::map_position_to_vegalite("pos1end", renderer).unwrap(); - let pos2 = super::projection::map_position_to_vegalite("pos2", renderer).unwrap(); - let pos2_end = super::projection::map_position_to_vegalite("pos2end", renderer).unwrap(); + let pos1 = renderer.map_position("pos1").unwrap(); + let pos1_end = renderer.map_position("pos1end").unwrap(); + let pos2 = renderer.map_position("pos2").unwrap(); + let pos2_end = renderer.map_position("pos2end").unwrap(); let (pos1_offset, pos2_offset) = renderer.offset_channels(); diff --git a/src/writer/vegalite/layer.rs b/src/writer/vegalite/layer.rs index 2ac08834a..ca1d1f805 100644 --- a/src/writer/vegalite/layer.rs +++ b/src/writer/vegalite/layer.rs @@ -2236,8 +2236,10 @@ impl SpatialRenderer { GgsqlError::WriterError(format!("Failed to convert WKB to GeoJSON: {}", e)) })?; - serde_json::from_slice(&geojson_out) - .map_err(|e| GgsqlError::WriterError(format!("Invalid GeoJSON from WKB: {}", e))) + match serde_json::from_slice(&geojson_out) { + Ok(value) => Ok(value), + Err(_) => Ok(Value::Null), + } } fn parse_geometry_from_array(array: &arrow::array::ArrayRef, idx: usize) -> Result { diff --git a/src/writer/vegalite/mod.rs b/src/writer/vegalite/mod.rs index 727526761..d616fbfe2 100644 --- a/src/writer/vegalite/mod.rs +++ b/src/writer/vegalite/mod.rs @@ -308,9 +308,19 @@ fn build_layer_encoding( channel_name = "fillOpacity".to_string(); } - // Secondary position channels (x2, y2, theta2, radius2) only support - // field/datum/value in Vega-Lite — not type, scale, axis, or title - if matches!(channel_name.as_str(), "x2" | "y2" | "theta2" | "radius2") { + // Secondary position channels (x2, y2, theta2, radius2) and geographic + // channels (longitude, latitude) only support field/datum/value in + // Vega-Lite — not type, scale, axis, or title. + if matches!( + channel_name.as_str(), + "x2" | "y2" + | "theta2" + | "radius2" + | "longitude" + | "latitude" + | "longitude2" + | "latitude2" + ) { let secondary_encoding = match value { AestheticValue::Column { name: col, .. } => json!({"field": col}), AestheticValue::Literal(lit) => json!({"value": lit.to_json()}), @@ -383,9 +393,12 @@ fn build_layer_encoding( // This prevents Vega-Lite from applying its own stack/dodge logic on top of ours. // Only set stack: null on primary position channels (y/radius) — Vega-Lite does // not support 'stack' on secondary channels (y2/radius2) and Altair rejects it. - if let Some(y_enc) = encoding.get_mut(pos2.as_str()) { - if let Some(obj) = y_enc.as_object_mut() { - obj.insert("stack".to_string(), Value::Null); + // Geographic channels (latitude) don't support stack either. + if !matches!(pos2.as_str(), "latitude") { + if let Some(y_enc) = encoding.get_mut(pos2.as_str()) { + if let Some(obj) = y_enc.as_object_mut() { + obj.insert("stack".to_string(), Value::Null); + } } } @@ -3018,6 +3031,8 @@ mod tests { aesthetics: vec!["angle".to_string(), "radius".to_string()], coord: Coord::polar(), properties: HashMap::new(), + map_projection: None, + computed: HashMap::new(), }); let layer = Layer::new(Geom::point()) .with_aesthetic( diff --git a/src/writer/vegalite/projection/cartesian.rs b/src/writer/vegalite/projection/cartesian.rs new file mode 100644 index 000000000..5ff234556 --- /dev/null +++ b/src/writer/vegalite/projection/cartesian.rs @@ -0,0 +1,28 @@ +use super::ProjectionRenderer; + +/// Cartesian projection — standard x/y coordinates. +pub(in crate::writer) struct CartesianProjection { + is_faceted: bool, +} + +impl CartesianProjection { + pub(super) fn new(facet: Option<&crate::plot::Facet>) -> Self { + Self { + is_faceted: facet.is_some_and(|f| !f.get_variables().is_empty()), + } + } +} + +impl ProjectionRenderer for CartesianProjection { + fn is_faceted(&self) -> bool { + self.is_faceted + } + + fn position_channels(&self) -> (&'static str, &'static str) { + ("x", "y") + } + + fn offset_channels(&self) -> (&'static str, &'static str) { + ("xOffset", "yOffset") + } +} diff --git a/src/writer/vegalite/projection/map.rs b/src/writer/vegalite/projection/map.rs new file mode 100644 index 000000000..336bae625 --- /dev/null +++ b/src/writer/vegalite/projection/map.rs @@ -0,0 +1,354 @@ +//! Map projection implementation for Vega-Lite writer +//! +//! For data that has been pre-projected server-side (via ST_Transform), Vega-Lite +//! must use an identity projection so it passes coordinates through without +//! re-projecting via d3-geo. + +use crate::plot::types::ArrayElement; +use crate::plot::{ParameterValue, Projection, Scale}; +use crate::{Plot, Result}; +use serde_json::{json, Value}; + +use super::ProjectionRenderer; + +/// Map projection — pre-projected spatial coordinates. +pub(in crate::writer) struct MapProjection { + is_faceted: bool, + panel_boundary_wkt: Option, + graticule_lon_wkt: Option, + graticule_lat_wkt: Option, + bbox: Option<[f64; 4]>, +} + +impl MapProjection { + pub(super) fn new(project: Option<&Projection>, facet: Option<&crate::plot::Facet>) -> Self { + let get_string = |key: &str| -> Option { + project + .and_then(|p| p.computed.get(key)) + .and_then(|v| match v { + ParameterValue::String(s) => Some(s.clone()), + _ => None, + }) + }; + let panel_boundary_wkt = get_string("panel_boundary"); + let graticule_lon_wkt = get_string("graticule_lon"); + let graticule_lat_wkt = get_string("graticule_lat"); + let bbox = if let Some(ParameterValue::Array(arr)) = + project.and_then(|p| p.computed.get("bbox")) + { + let nums: Vec = arr + .iter() + .filter_map(|e| match e { + ArrayElement::Number(n) => Some(*n), + _ => None, + }) + .collect(); + nums.try_into().ok() + } else { + None + }; + Self { + is_faceted: facet.is_some_and(|f| !f.get_variables().is_empty()), + panel_boundary_wkt, + graticule_lon_wkt, + graticule_lat_wkt, + bbox, + } + } + + fn panel_boundary(&self, theme: &mut Value) -> Vec { + let (fill, stroke) = + if let Some(view) = theme.get_mut("view").and_then(|v| v.as_object_mut()) { + let fill = view.remove("fill").unwrap_or(Value::Null); + let stroke = view.remove("stroke").unwrap_or(Value::Null); + view.insert("stroke".to_string(), Value::Null); + (fill, stroke) + } else { + (Value::Null, Value::Null) + }; + + if let Some(ref wkt) = self.panel_boundary_wkt { + let Some(geojson) = wkt_to_geojson(wkt) else { + return Vec::new(); + }; + vec![geoshape_layer( + geojson, + json!({ "fill": fill, "stroke": stroke }), + )] + } else { + vec![json!({ + "mark": { + "type": "rect", + "fill": fill, + "stroke": stroke, + }, + "encoding": { + "x": {"value": 0}, + "y": {"value": 0}, + "x2": {"value": {"expr": "width"}}, + "y2": {"value": {"expr": "height"}}, + } + })] + } + } + + fn graticule(&self, theme: &Value) -> Vec { + let grid_color = theme + .pointer("/axis/gridColor") + .cloned() + .unwrap_or(json!("#cccccc")); + let grid_width = theme + .pointer("/axis/gridWidth") + .cloned() + .unwrap_or(json!(0.5)); + + let mark = json!({ + "filled": false, + "stroke": grid_color, + "strokeWidth": grid_width, + }); + + [&self.graticule_lon_wkt, &self.graticule_lat_wkt] + .into_iter() + .flatten() + .filter_map(|wkt| wkt_to_geojson(wkt)) + .map(|geojson| geoshape_layer(geojson, mark.clone())) + .collect() + } +} + +impl ProjectionRenderer for MapProjection { + fn is_faceted(&self) -> bool { + self.is_faceted + } + + fn position_channels(&self) -> (&'static str, &'static str) { + ("longitude", "latitude") + } + + fn offset_channels(&self) -> (&'static str, &'static str) { + ("longitude", "latitude") + } + + fn transform_layers(&self, _spec: &Plot, vl_spec: &mut Value) -> Result<()> { + let mut proj = json!({ + "type": "identity", + "reflectY": true + }); + if let Some([xmin, ymin, xmax, ymax]) = self.bbox { + let dx = (xmax - xmin) * 1.1; + let dy = (ymax - ymin) * 1.1; + if dx.is_finite() && dy.is_finite() && dx > 0.0 && dy > 0.0 { + let cx = (xmin + xmax) / 2.0; + let cy = (ymin + ymax) / 2.0; + proj["scale"] = json!({"expr": format!( + "min(width / {dx}, height / {dy})" + )}); + proj["translate"] = json!({"expr": format!( + "[width / 2 - min(width / {dx}, height / {dy}) * {cx}, \ + height / 2 + min(width / {dx}, height / {dy}) * {cy}]" + )}); + } + } + vl_spec["projection"] = proj; + Ok(()) + } + + fn background_layers(&self, _scales: &[Scale], theme: &mut Value) -> Vec { + let mut layers = Vec::new(); + layers.extend(self.panel_boundary(theme)); + layers.extend(self.graticule(theme)); + layers + } +} + +fn geoshape_layer(geojson: Value, mut mark: Value) -> Value { + mark["type"] = json!("geoshape"); + json!({ + "data": {"values": [{"type": "Feature", "geometry": geojson}]}, + "mark": mark + }) +} + +#[cfg(feature = "spatial")] +fn wkt_to_geojson(wkt: &str) -> Option { + use geozero::geojson::GeoJsonWriter; + use geozero::wkt::WktReader; + use geozero::GeozeroDatasource; + use std::io::Cursor; + + let mut reader = WktReader(wkt.as_bytes()); + let mut geojson_out = Vec::new(); + reader + .process_geom(&mut GeoJsonWriter::new(Cursor::new(&mut geojson_out))) + .ok()?; + serde_json::from_slice(&geojson_out).ok() +} + +#[cfg(not(feature = "spatial"))] +fn wkt_to_geojson(_wkt: &str) -> Option { + None +} + +#[cfg(test)] +mod tests { + use super::*; + use crate::plot::{Facet, FacetLayout, Projection}; + + #[test] + fn test_map_projection_identity() { + let renderer = MapProjection::new(None, None); + let mut vl_spec = json!({"layer": []}); + let spec = Plot::default(); + + renderer.transform_layers(&spec, &mut vl_spec).unwrap(); + + assert_eq!(vl_spec["projection"]["type"], "identity"); + assert_eq!(vl_spec["projection"]["reflectY"], true); + } + + #[test] + fn test_map_projection_channels() { + let renderer = MapProjection::new(None, None); + assert_eq!(renderer.position_channels(), ("longitude", "latitude")); + assert_eq!(renderer.offset_channels(), ("longitude", "latitude")); + assert_eq!(renderer.map_position("pos1"), Some("longitude".to_string())); + assert_eq!(renderer.map_position("pos2"), Some("latitude".to_string())); + } + + #[test] + fn test_map_projection_faceted() { + let facet = Facet::new(FacetLayout::Wrap { + variables: vec!["region".to_string()], + }); + let renderer = MapProjection::new(None, Some(&facet)); + assert!(renderer.is_faceted()); + assert_eq!(renderer.panel_size(), None); + } + + #[test] + fn test_map_projection_not_faceted() { + let renderer = MapProjection::new(None, None); + assert!(!renderer.is_faceted()); + assert_eq!( + renderer.panel_size(), + Some((json!("container"), json!("container"))) + ); + } + + #[test] + fn test_background_layer_without_boundary() { + let renderer = MapProjection::new(None, None); + let mut theme = json!({"view": {"fill": "white", "stroke": "gray"}}); + let layers = renderer.background_layers(&[], &mut theme); + // Fallback: rect mark spanning the full view + assert_eq!(layers.len(), 1); + assert_eq!(layers[0]["mark"]["type"], "rect"); + assert_eq!(layers[0]["mark"]["fill"], "white"); + assert_eq!(layers[0]["mark"]["stroke"], "gray"); + } + + #[test] + fn test_background_layer_with_boundary() { + let mut proj = Projection::map(); + proj.computed.insert( + "panel_boundary".to_string(), + ParameterValue::String("POLYGON ((0 0, 1 0, 1 1, 0 1, 0 0))".to_string()), + ); + let renderer = MapProjection::new(Some(&proj), None); + let mut theme = json!({"view": {"fill": "white", "stroke": "gray"}}); + let layers = renderer.background_layers(&[], &mut theme); + + assert_eq!(layers.len(), 1); + let layer = &layers[0]; + assert_eq!(layer["mark"]["type"], "geoshape"); + assert_eq!(layer["mark"]["fill"], "white"); + assert_eq!(layer["mark"]["stroke"], "gray"); + + let geom = &layer["data"]["values"][0]["geometry"]; + assert_eq!(geom["type"], "Polygon"); + assert!(!geom["coordinates"].is_null()); + + // view stroke should be nulled out + assert_eq!(theme["view"]["stroke"], Value::Null); + } + + #[test] + fn test_bbox_emits_scale_translate_exprs() { + use crate::plot::types::ArrayElement; + + let mut proj = Projection::map(); + proj.computed.insert( + "bbox".to_string(), + ParameterValue::Array(vec![ + ArrayElement::Number(0.0), + ArrayElement::Number(0.0), + ArrayElement::Number(100.0), + ArrayElement::Number(200.0), + ]), + ); + let renderer = MapProjection::new(Some(&proj), None); + let mut vl_spec = json!({"layer": []}); + let spec = Plot::default(); + + renderer.transform_layers(&spec, &mut vl_spec).unwrap(); + + let scale = &vl_spec["projection"]["scale"]["expr"]; + let translate = &vl_spec["projection"]["translate"]["expr"]; + assert!(scale.is_string(), "scale should be an expr"); + assert!(translate.is_string(), "translate should be an expr"); + assert!(scale.as_str().unwrap().contains("width")); + assert!(translate.as_str().unwrap().contains("height")); + } + + #[test] + fn test_graticule_layers_rendered() { + let mut proj = Projection::map(); + proj.computed.insert( + "graticule_lon".to_string(), + ParameterValue::String("MULTILINESTRING ((0 -90, 0 90), (30 -90, 30 90))".to_string()), + ); + proj.computed.insert( + "graticule_lat".to_string(), + ParameterValue::String( + "MULTILINESTRING ((-180 0, 180 0), (-180 45, 180 45))".to_string(), + ), + ); + let renderer = MapProjection::new(Some(&proj), None); + let mut theme = json!({"axis": {"gridColor": "#dddddd", "gridWidth": 1}}); + let layers = renderer.background_layers(&[], &mut theme); + + // 1 rect fallback + 2 graticule layers + assert_eq!(layers.len(), 3); + assert_eq!(layers[0]["mark"]["type"], "rect"); + for layer in &layers[1..] { + assert_eq!(layer["mark"]["type"], "geoshape"); + assert_eq!(layer["mark"]["filled"], false); + assert_eq!(layer["mark"]["stroke"], "#dddddd"); + assert_eq!(layer["mark"]["strokeWidth"], 1); + let geom = &layer["data"]["values"][0]["geometry"]; + assert_eq!(geom["type"], "MultiLineString"); + } + } + + #[test] + fn test_graticule_with_panel_boundary() { + let mut proj = Projection::map(); + proj.computed.insert( + "panel_boundary".to_string(), + ParameterValue::String("POLYGON ((0 0, 1 0, 1 1, 0 1, 0 0))".to_string()), + ); + proj.computed.insert( + "graticule_lon".to_string(), + ParameterValue::String("MULTILINESTRING ((0.5 0, 0.5 1))".to_string()), + ); + let renderer = MapProjection::new(Some(&proj), None); + let mut theme = json!({"view": {"fill": "white", "stroke": "gray"}}); + let layers = renderer.background_layers(&[], &mut theme); + + // Panel boundary first, then graticule + assert_eq!(layers.len(), 2); + assert_eq!(layers[0]["mark"]["fill"], "white"); + assert_eq!(layers[1]["mark"]["filled"], false); + } +} diff --git a/src/writer/vegalite/projection/mod.rs b/src/writer/vegalite/projection/mod.rs new file mode 100644 index 000000000..cad3afbcb --- /dev/null +++ b/src/writer/vegalite/projection/mod.rs @@ -0,0 +1,259 @@ +//! Projection rendering for Vega-Lite writer +//! +//! This module provides a trait-based design for projection rendering. +//! Each projection type (cartesian, polar, and future map projections) +//! implements `ProjectionRenderer`, which owns both the VL channel mapping +//! and the spec-level transformations for that projection. + +mod cartesian; +mod map; +mod polar; + +use crate::plot::{CoordKind, ParameterValue, Projection, Scale}; +use crate::{Plot, Result}; +use serde_json::{json, Value}; + +use cartesian::CartesianProjection; +use map::MapProjection; +use polar::PolarProjection; + +const ANGLE_TOLERANCE: f64 = 1.49011611938476e-08; // f64::EPSILON.sqrt() + +// ============================================================================= +// ProjectionRenderer trait +// ============================================================================= + +/// Trait defining how a projection type maps to Vega-Lite. +/// +/// Each implementation owns two concerns: +/// 1. **Channel mapping** — translating internal position aesthetics (pos1, pos2, …) +/// to Vega-Lite encoding channel names. +/// 2. **Spec transformation** — modifying the Vega-Lite spec after layers are built +/// (e.g., converting marks to arcs for polar). +pub(super) trait ProjectionRenderer { + /// Whether the spec uses faceting. + fn is_faceted(&self) -> bool; + + /// Primary and secondary VL channel names for this projection. + /// + /// Returns `(pos1_channel, pos2_channel)`, e.g. `("x", "y")` for cartesian, + /// `("radius", "theta")` for polar. + fn position_channels(&self) -> (&'static str, &'static str); + + /// Offset channel names for this projection. + /// + /// Returns `(pos1_offset, pos2_offset)`, e.g. `("xOffset", "yOffset")`. + fn offset_channels(&self) -> (&'static str, &'static str); + + /// Map internal position aesthetic to Vega-Lite channel name. + /// + /// Returns `Some(channel_name)` for internal position aesthetics (pos1, pos2, etc.), + /// or `None` for material aesthetics. + fn map_position(&self, aesthetic: &str) -> Option { + let (primary, secondary) = self.position_channels(); + match aesthetic { + "pos1" | "pos1min" => Some(primary.to_string()), + "pos2" | "pos2min" => Some(secondary.to_string()), + "pos1end" | "pos1max" => Some(format!("{}2", primary)), + "pos2end" | "pos2max" => Some(format!("{}2", secondary)), + _ => None, + } + } + + /// Panel dimensions as VL values (`"container"` or explicit pixels). + /// + /// Returns `None` for faceted cartesian (VL handles sizing). + fn panel_size(&self) -> Option<(Value, Value)> { + if self.is_faceted() { + None + } else { + Some((json!("container"), json!("container"))) + } + } + + /// Apply projection-specific transformations to the VL spec. + /// + /// Called after layers are built but before faceting. + fn transform_layers(&self, _spec: &Plot, _vl_spec: &mut Value) -> Result<()> { + Ok(()) + } + + /// Vega-Lite layers to prepend before the data layers. + fn background_layers(&self, _scales: &[Scale], _theme: &mut Value) -> Vec { + Vec::new() + } + + /// Vega-Lite layers to append after the data layers. + fn foreground_layers(&self, _scales: &[Scale], _theme: &mut Value) -> Vec { + Vec::new() + } + + /// Apply all projection-specific work: transforms, clip, and panel decoration. + fn apply_projection(&self, spec: &Plot, theme: &mut Value, vl_spec: &mut Value) -> Result<()> { + self.transform_layers(spec, vl_spec)?; + + if let Some(ref project) = spec.project { + if let Some(ParameterValue::Boolean(clip)) = project.properties.get("clip") { + apply_clip_to_layers(vl_spec, *clip); + } + } + + let mut bg = self.background_layers(&spec.scales, theme); + let mut fg = self.foreground_layers(&spec.scales, theme); + if !(bg.is_empty() && fg.is_empty()) { + for layer in &mut bg { + layer["description"] = json!("background"); + } + for layer in &mut fg { + layer["description"] = json!("foreground"); + } + if let Some(layers) = get_layers_mut(vl_spec) { + let data_layers = std::mem::take(layers); + layers.reserve(bg.len() + data_layers.len() + fg.len()); + layers.extend(bg); + layers.extend(data_layers); + layers.extend(fg); + } + } + + Ok(()) + } +} + +// ============================================================================= +// Factory +// ============================================================================= + +/// Get the projection renderer for a projection spec. +/// +/// Returns the appropriate renderer based on the projection's coord kind, +/// or a Cartesian renderer if no projection is specified. +pub(super) fn get_projection_renderer( + project: Option<&Projection>, + facet: Option<&crate::plot::Facet>, + scales: &[Scale], +) -> Box { + match project.map(|p| p.coord.coord_kind()) { + Some(CoordKind::Polar) => Box::new(PolarProjection::new(project, facet, scales)), + Some(CoordKind::Map) => Box::new(MapProjection::new(project, facet)), + Some(CoordKind::Cartesian) | None => Box::new(CartesianProjection::new(facet)), + } +} + +// ============================================================================= +// AxisInfo — reusable across projection types +// ============================================================================= + +pub(in crate::writer) struct AxisInfo { + pub domain: Option<(f64, f64)>, + pub breaks: Vec, + pub labels: Vec<(f64, String)>, + pub suppress: bool, +} + +impl AxisInfo { + pub fn new(aesthetic: &str, scales: &[Scale], facet: Option<&crate::plot::Facet>) -> Self { + let scale = scales.iter().find(|s| s.aesthetic == aesthetic); + let (domain, labels) = match scale { + Some(s) => (s.numeric_domain(), s.break_labels()), + None => (None, Vec::new()), + }; + let domain = domain.filter(|(min, max)| (max - min).abs() > f64::EPSILON); + let breaks = labels.iter().map(|(v, _)| *v).collect(); + let suppress = + facet.is_some_and(|f| f.is_free(aesthetic)) || scale.is_some_and(|s| s.is_dummy()); + Self { + domain, + breaks, + labels, + suppress, + } + } +} + +// ============================================================================= +// Shared helpers +// ============================================================================= + +/// Get mutable reference to the layers array, handling both flat and faceted specs. +/// +/// In a flat spec: `vl_spec["layer"]` +/// In a faceted spec: `vl_spec["spec"]["layer"]` +fn get_layers_mut(vl_spec: &mut Value) -> Option<&mut Vec> { + if vl_spec.get("layer").is_some() { + vl_spec.get_mut("layer").and_then(|l| l.as_array_mut()) + } else { + vl_spec + .get_mut("spec") + .and_then(|s| s.get_mut("layer")) + .and_then(|l| l.as_array_mut()) + } +} + +/// Apply clip setting to all layers +fn apply_clip_to_layers(vl_spec: &mut Value, clip: bool) { + if let Some(layers_arr) = get_layers_mut(vl_spec) { + for layer in layers_arr { + if let Some(mark) = layer.get_mut("mark") { + if mark.is_string() { + let mark_type = mark.as_str().unwrap().to_string(); + *mark = json!({"type": mark_type, "clip": clip}); + } else if let Some(obj) = mark.as_object_mut() { + obj.insert("clip".to_string(), json!(clip)); + } + } + } + } +} + +// ============================================================================= +// Tests +// ============================================================================= + +#[cfg(test)] +mod tests { + use super::*; + use crate::plot::Projection; + + #[test] + fn test_map_position_cartesian() { + let renderer = CartesianProjection::new(None); + assert_eq!(renderer.map_position("pos1"), Some("x".to_string())); + assert_eq!(renderer.map_position("pos2"), Some("y".to_string())); + assert_eq!(renderer.map_position("pos1end"), Some("x2".to_string())); + assert_eq!(renderer.map_position("pos2end"), Some("y2".to_string())); + assert_eq!(renderer.map_position("color"), None); + assert_eq!(renderer.offset_channels(), ("xOffset", "yOffset")); + assert_eq!( + renderer.panel_size(), + Some((json!("container"), json!("container"))) + ); + } + + #[test] + fn test_map_position_polar() { + let renderer = PolarProjection::new(None, None, &[]); + assert_eq!(renderer.map_position("pos1"), Some("radius".to_string())); + assert_eq!(renderer.map_position("pos2"), Some("theta".to_string())); + assert_eq!( + renderer.map_position("pos1end"), + Some("radius2".to_string()) + ); + assert_eq!(renderer.map_position("pos2end"), Some("theta2".to_string())); + assert_eq!(renderer.offset_channels(), ("radiusOffset", "thetaOffset")); + assert_eq!( + renderer.panel_size(), + Some((json!("container"), json!("container"))) + ); + } + + #[test] + fn test_get_projection_renderer() { + let cartesian = get_projection_renderer(None, None, &[]); + assert_eq!(cartesian.position_channels(), ("x", "y")); + + let polar_proj = Projection::polar(); + let polar = get_projection_renderer(Some(&polar_proj), None, &[]); + assert_eq!(polar.position_channels(), ("radius", "theta")); + } +} diff --git a/src/writer/vegalite/projection.rs b/src/writer/vegalite/projection/polar.rs similarity index 86% rename from src/writer/vegalite/projection.rs rename to src/writer/vegalite/projection/polar.rs index 16ec4e61f..328202b88 100644 --- a/src/writer/vegalite/projection.rs +++ b/src/writer/vegalite/projection/polar.rs @@ -1,179 +1,15 @@ -//! Projection rendering for Vega-Lite writer +//! Polar projection implementation for Vega-Lite writer //! -//! This module provides a trait-based design for projection rendering. -//! Each projection type (cartesian, polar, and future map projections) -//! implements `ProjectionRenderer`, which owns both the VL channel mapping -//! and the spec-level transformations for that projection. +//! Handles radius/theta coordinate transformations for pie charts, rose plots, +//! and other circular visualizations. -use crate::plot::{CoordKind, ParameterValue, Projection, Scale}; +use crate::plot::{ParameterValue, Projection, Scale}; use crate::{GgsqlError, Plot, Result}; use serde_json::{json, Value}; -use super::DEFAULT_POLAR_SIZE; - -const ANGLE_TOLERANCE: f64 = 1.49011611938476e-08; // f64::EPSILON.sqrt() - -// ============================================================================= -// ProjectionRenderer trait -// ============================================================================= - -/// Trait defining how a projection type maps to Vega-Lite. -/// -/// Each implementation owns two concerns: -/// 1. **Channel mapping** — translating internal position aesthetics (pos1, pos2, …) -/// to Vega-Lite encoding channel names. -/// 2. **Spec transformation** — modifying the Vega-Lite spec after layers are built -/// (e.g., converting marks to arcs for polar). -pub(super) trait ProjectionRenderer { - /// Whether the spec uses faceting. - fn is_faceted(&self) -> bool; - - /// Primary and secondary VL channel names for this projection. - /// - /// Returns `(pos1_channel, pos2_channel)`, e.g. `("x", "y")` for cartesian, - /// `("radius", "theta")` for polar. - fn position_channels(&self) -> (&'static str, &'static str); - - /// Offset channel names for this projection. - /// - /// Returns `(pos1_offset, pos2_offset)`, e.g. `("xOffset", "yOffset")`. - fn offset_channels(&self) -> (&'static str, &'static str); - - /// Panel dimensions as VL values (`"container"` or explicit pixels). - /// - /// Returns `None` for faceted cartesian (VL handles sizing). - fn panel_size(&self) -> Option<(Value, Value)> { - if self.is_faceted() { - None - } else { - Some((json!("container"), json!("container"))) - } - } - - /// Apply projection-specific transformations to the VL spec. - /// - /// Called after layers are built but before faceting. - fn transform_layers(&self, _spec: &Plot, _vl_spec: &mut Value) -> Result<()> { - Ok(()) - } - - /// Vega-Lite layers to prepend before the data layers. - fn background_layers(&self, _scales: &[Scale], _theme: &mut Value) -> Vec { - Vec::new() - } - - /// Vega-Lite layers to append after the data layers. - fn foreground_layers(&self, _scales: &[Scale], _theme: &mut Value) -> Vec { - Vec::new() - } - - /// Apply all projection-specific work: transforms, clip, and panel decoration. - fn apply_projection(&self, spec: &Plot, theme: &mut Value, vl_spec: &mut Value) -> Result<()> { - self.transform_layers(spec, vl_spec)?; - - if let Some(ref project) = spec.project { - if let Some(ParameterValue::Boolean(clip)) = project.properties.get("clip") { - apply_clip_to_layers(vl_spec, *clip); - } - } - - let mut bg = self.background_layers(&spec.scales, theme); - let mut fg = self.foreground_layers(&spec.scales, theme); - if !(bg.is_empty() && fg.is_empty()) { - for layer in &mut bg { - layer["description"] = json!("background"); - } - for layer in &mut fg { - layer["description"] = json!("foreground"); - } - if let Some(layers) = get_layers_mut(vl_spec) { - let data_layers = std::mem::take(layers); - layers.reserve(bg.len() + data_layers.len() + fg.len()); - layers.extend(bg); - layers.extend(data_layers); - layers.extend(fg); - } - } - - Ok(()) - } -} - -// ============================================================================= -// Factory -// ============================================================================= - -/// Get the projection renderer for a projection spec. -/// -/// Returns the appropriate renderer based on the projection's coord kind, -/// or a Cartesian renderer if no projection is specified. -pub(super) fn get_projection_renderer( - project: Option<&Projection>, - facet: Option<&crate::plot::Facet>, - scales: &[Scale], -) -> Box { - let is_faceted = facet.is_some_and(|f| !f.get_variables().is_empty()); - match project.map(|p| p.coord.coord_kind()) { - Some(CoordKind::Polar) => Box::new(PolarProjection { - panel: PolarContext::new(project, facet, scales), - }), - Some(CoordKind::Cartesian) | None => Box::new(CartesianProjection { is_faceted }), - } -} - -// ============================================================================= -// Channel mapping helpers (used by encoding.rs via the trait) -// ============================================================================= - -/// Map internal position aesthetic to Vega-Lite channel name using the renderer. -/// -/// Returns `Some(channel_name)` for internal position aesthetics (pos1, pos2, etc.), -/// or `None` for material aesthetics. -pub(super) fn map_position_to_vegalite( - aesthetic: &str, - renderer: &dyn ProjectionRenderer, -) -> Option { - let (primary, secondary) = renderer.position_channels(); - - // Match internal position aesthetic patterns - // Convention: min → primary channel (x/y), max → secondary channel (x2/y2) - match aesthetic { - // Primary position and min variants - "pos1" | "pos1min" => Some(primary.to_string()), - "pos2" | "pos2min" => Some(secondary.to_string()), - // End and max variants (Vega-Lite uses x2/y2/theta2/radius2) - "pos1end" | "pos1max" => Some(format!("{}2", primary)), - "pos2end" | "pos2max" => Some(format!("{}2", secondary)), - _ => None, - } -} - -// ============================================================================= -// CartesianProjection -// ============================================================================= - -/// Cartesian projection — standard x/y coordinates. -struct CartesianProjection { - is_faceted: bool, -} - -impl ProjectionRenderer for CartesianProjection { - fn is_faceted(&self) -> bool { - self.is_faceted - } - - fn position_channels(&self) -> (&'static str, &'static str) { - ("x", "y") - } - - fn offset_channels(&self) -> (&'static str, &'static str) { - ("xOffset", "yOffset") - } -} - -// ============================================================================= -// PolarProjection -// ============================================================================= +use super::super::escape_vega_string; +use super::super::DEFAULT_POLAR_SIZE; +use super::{get_layers_mut, AxisInfo, ProjectionRenderer, ANGLE_TOLERANCE}; /// Normalized outer radius (proportion of `min(width, height) / 2`). const POLAR_OUTER: f64 = 1.0; @@ -182,37 +18,6 @@ const POLAR_OUTER: f64 = 1.0; /// `1 - paddingInner` for band scales, which is ~0.9). const POLAR_BAND_FRACTION: f64 = 0.9; -struct AxisInfo { - domain: Option<(f64, f64)>, - breaks: Vec, - labels: Vec<(f64, String)>, - suppress: bool, -} - -impl AxisInfo { - fn new(aesthetic: &str, scales: &[Scale], facet: Option<&crate::plot::Facet>) -> Self { - let scale = scales.iter().find(|s| s.aesthetic == aesthetic); - let (domain, labels) = match scale { - Some(s) => (s.numeric_domain(), s.break_labels()), - None => (None, Vec::new()), - }; - // Set domain to None if zero-range - let domain = domain.filter(|(min, max)| (max - min).abs() > f64::EPSILON); - let breaks = labels.iter().map(|(v, _)| *v).collect(); - // Free facet scales have per-panel domains that don't match the global - // positions used for decoration; dummy scales are stat-injected placeholders - // with no meaningful domain to label. - let suppress = - facet.is_some_and(|f| f.is_free(aesthetic)) || scale.is_some_and(|s| s.is_dummy()); - Self { - domain, - breaks, - labels, - suppress, - } - } -} - /// Resolved geometry and scale context for polar specs. /// /// Holds angular range, radius bounds, VL expression strings for the panel @@ -355,10 +160,22 @@ impl PolarContext { } /// Polar projection — radius/theta coordinates for pie charts, rose plots, etc. -struct PolarProjection { +pub(in crate::writer) struct PolarProjection { panel: PolarContext, } +impl PolarProjection { + pub(super) fn new( + project: Option<&Projection>, + facet: Option<&crate::plot::Facet>, + scales: &[Scale], + ) -> Self { + Self { + panel: PolarContext::new(project, facet, scales), + } + } +} + impl ProjectionRenderer for PolarProjection { fn is_faceted(&self) -> bool { self.panel.is_faceted @@ -996,43 +813,6 @@ fn polygon_ring( }) } -// ============================================================================= -// Shared helpers -// ============================================================================= - -/// Get mutable reference to the layers array, handling both flat and faceted specs. -/// -/// In a flat spec: `vl_spec["layer"]` -/// In a faceted spec: `vl_spec["spec"]["layer"]` -fn get_layers_mut(vl_spec: &mut Value) -> Option<&mut Vec> { - // Try flat structure first, then faceted - if vl_spec.get("layer").is_some() { - vl_spec.get_mut("layer").and_then(|l| l.as_array_mut()) - } else { - vl_spec - .get_mut("spec") - .and_then(|s| s.get_mut("layer")) - .and_then(|l| l.as_array_mut()) - } -} - -/// Apply clip setting to all layers -fn apply_clip_to_layers(vl_spec: &mut Value, clip: bool) { - if let Some(layers_arr) = get_layers_mut(vl_spec) { - for layer in layers_arr { - if let Some(mark) = layer.get_mut("mark") { - if mark.is_string() { - // Convert "point" to {"type": "point", "clip": ...} - let mark_type = mark.as_str().unwrap().to_string(); - *mark = json!({"type": mark_type, "clip": clip}); - } else if let Some(obj) = mark.as_object_mut() { - obj.insert("clip".to_string(), json!(clip)); - } - } - } - } -} - // ============================================================================= // Polar projection transformation // ============================================================================= @@ -1396,7 +1176,7 @@ fn extract_polar_channel( if !strings.is_empty() { let literal: String = strings .iter() - .map(|s| format!("'{}'", super::escape_vega_string(s))) + .map(|s| format!("'{}'", escape_vega_string(s))) .collect::>() .join(","); let arr_expr = format!("[{}]", literal); @@ -1540,7 +1320,7 @@ fn apply_polar_radius_range(encoding: &mut Value, panel: &PolarContext) -> Resul #[cfg(test)] mod tests { use super::*; - use crate::plot::{Facet, FacetLayout}; + use crate::plot::{Facet, FacetLayout, ParameterValue, Projection}; fn faceted() -> Facet { Facet::new(FacetLayout::Wrap { @@ -1550,7 +1330,6 @@ mod tests { #[test] fn test_polar_inner_radius_non_faceted() { - // Non-faceted donut should use dynamic min(width,height) expressions let mut encoding = json!({ "radius": { "field": "dummy", @@ -1579,7 +1358,6 @@ mod tests { #[test] fn test_polar_inner_radius_faceted() { - // Faceted donut should use explicit size calculation let mut encoding = json!({ "radius": { "field": "dummy", @@ -1605,7 +1383,6 @@ mod tests { #[test] fn test_polar_inner_radius_zero() { - // inner = 0 should still apply range (full pie, no donut hole) let mut encoding = json!({ "radius": { "field": "dummy", @@ -1621,68 +1398,12 @@ mod tests { let panel = PolarContext::new(Some(&proj), Some(&f), &[]); apply_polar_radius_range(&mut encoding, &panel).unwrap(); - // Range should be [0, 350/2] for full pie let range = encoding["radius"]["scale"]["range"].as_array().unwrap(); assert_eq!(range.len(), 2); assert_eq!(range[0]["expr"].as_str().unwrap(), "175 * (0)"); assert_eq!(range[1]["expr"].as_str().unwrap(), "175 * (1)"); } - #[test] - fn test_map_position_to_vegalite_cartesian() { - let renderer = CartesianProjection { is_faceted: false }; - assert_eq!( - map_position_to_vegalite("pos1", &renderer), - Some("x".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos2", &renderer), - Some("y".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos1end", &renderer), - Some("x2".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos2end", &renderer), - Some("y2".to_string()) - ); - assert_eq!(map_position_to_vegalite("color", &renderer), None); - assert_eq!(renderer.offset_channels(), ("xOffset", "yOffset")); - assert_eq!( - renderer.panel_size(), - Some((json!("container"), json!("container"))) - ); - } - - #[test] - fn test_map_position_to_vegalite_polar() { - let renderer = PolarProjection { - panel: PolarContext::new(None, None, &[]), - }; - assert_eq!( - map_position_to_vegalite("pos1", &renderer), - Some("radius".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos2", &renderer), - Some("theta".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos1end", &renderer), - Some("radius2".to_string()) - ); - assert_eq!( - map_position_to_vegalite("pos2end", &renderer), - Some("theta2".to_string()) - ); - assert_eq!(renderer.offset_channels(), ("radiusOffset", "thetaOffset")); - assert_eq!( - renderer.panel_size(), - Some((json!("container"), json!("container"))) - ); - } - fn continuous_panel() -> PolarContext { let mut panel = PolarContext::new(None, None, &[]); panel.radial.domain = Some((0.0, 10.0)); @@ -1723,7 +1444,6 @@ mod tests { let transforms = layer["transform"].as_array().unwrap(); - // Should contain pixel-coordinate expressions using width/height signals let x_calc = transforms .iter() .find(|t| t["as"] == "__polar_x__") @@ -1744,11 +1464,9 @@ mod tests { "y should use pixel coordinates, got: {y_expr}" ); - // Encoding should use scale:null (raw pixel positions) assert_eq!(layer["encoding"]["x"]["scale"], json!(null)); assert_eq!(layer["encoding"]["y"]["scale"], json!(null)); - // Original polar channels should be removed assert!(layer["encoding"].get("radius").is_none()); assert!(layer["encoding"].get("theta").is_none()); } @@ -1773,24 +1491,12 @@ mod tests { ); } - #[test] - fn test_get_projection_renderer() { - let cartesian = get_projection_renderer(None, None, &[]); - assert_eq!(cartesian.position_channels(), ("x", "y")); - - let polar_proj = Projection::polar(); - let polar = get_projection_renderer(Some(&polar_proj), None, &[]); - assert_eq!(polar.position_channels(), ("radius", "theta")); - } - #[test] fn test_expr_normalize_radius() { let mut p = PolarContext::new(None, None, &[]); - // domain [0, 10], inner 0.2 p.inner = 0.2; p.radial.domain = Some((0.0, 10.0)); - // scale = (1.0 - 0.2) / (10 - 0) = 0.08 let expr = p.expr_normalize_radius("datum.v"); assert!( expr.contains("0.08"), @@ -1801,7 +1507,6 @@ mod tests { "should reference value, got: {expr}" ); - // domain [5, 15], inner 0 → scale = 1.0 / 10 = 0.1 p.inner = 0.0; p.radial.domain = Some((5.0, 15.0)); let expr = p.expr_normalize_radius("datum.x"); @@ -1810,7 +1515,6 @@ mod tests { "scale factor should be 0.1, got: {expr}" ); - // None domain → fallback to midpoint p.radial.domain = None; let expr = p.expr_normalize_radius("datum.v"); assert!( @@ -1823,13 +1527,11 @@ mod tests { fn test_expr_normalize_theta() { use std::f64::consts::PI; - // domain [0, 100], partial circle 90°–270° (π/2 to 3π/2) let mut panel = PolarContext::new(None, None, &[]); panel.start = PI / 2.0; panel.end = 3.0 * PI / 2.0; panel.angle.domain = Some((0.0, 100.0)); let expr = panel.expr_normalize_theta("datum.v"); - // scale = (3π/2 - π/2) / (100 - 0) = π / 100 ≈ 0.031416 let expected_scale = PI / 100.0; assert!( expr.contains(&format!("{expected_scale}")), @@ -1866,20 +1568,17 @@ mod tests { let layer = &layers[0]; - // Data should contain the break values let values = layer["data"]["values"].as_array().unwrap(); assert_eq!(values.len(), 3); assert_eq!(values[0]["v"], json!(25.0)); assert_eq!(values[1]["v"], json!(50.0)); assert_eq!(values[2]["v"], json!(75.0)); - // Mark should be a stroke-only arc assert_eq!(layer["mark"]["type"], "arc"); assert_eq!(layer["mark"]["fill"], json!(null)); assert_eq!(layer["mark"]["stroke"], "#FFF"); assert_eq!(layer["mark"]["strokeWidth"], 2.0); - // Radius encoding should use an expression let radius_expr = layer["encoding"]["radius"]["value"]["expr"] .as_str() .unwrap(); @@ -1900,21 +1599,17 @@ mod tests { let layer = &layers[0]; - // Data should contain the break values let values = layer["data"]["values"].as_array().unwrap(); assert_eq!(values.len(), 2); - // Mark should be a rule assert_eq!(layer["mark"]["type"], "rule"); assert_eq!(layer["mark"]["stroke"], "#CCC"); - // Should have calculate transforms for x, y, x2, y2 let transforms = layer["transform"].as_array().unwrap(); assert_eq!(transforms.len(), 4); let field_names: Vec<&str> = transforms.iter().filter_map(|t| t["as"].as_str()).collect(); assert_eq!(field_names, vec!["x", "y", "x2", "y2"]); - // Encoding should use scale:null for pixel positions assert_eq!(layer["encoding"]["x"]["scale"], json!(null)); assert_eq!(layer["encoding"]["y"]["scale"], json!(null)); } @@ -1943,7 +1638,6 @@ mod tests { "should produce axis line, ticks, and labels" ); - // Layer 0: axis line (single rule from inner to outer) let line = &layers[0]; assert_eq!(line["mark"]["type"], "rule"); assert_eq!(line["data"]["values"].as_array().unwrap().len(), 1); @@ -1951,7 +1645,6 @@ mod tests { let fields: Vec<&str> = transforms.iter().filter_map(|t| t["as"].as_str()).collect(); assert_eq!(fields, vec!["x", "y", "x2", "y2"]); - // Layer 1: ticks (one per break) let ticks = &layers[1]; assert_eq!(ticks["mark"]["type"], "rule"); assert_eq!(ticks["data"]["values"].as_array().unwrap().len(), 3); @@ -1962,7 +1655,6 @@ mod tests { .collect(); assert_eq!(tick_fields, vec!["cx", "cy", "x", "y", "x2", "y2"]); - // Layer 2: labels (one per break) let labels = &layers[2]; assert_eq!(labels["mark"]["type"], "text"); assert_eq!(labels["data"]["values"].as_array().unwrap().len(), 3); @@ -2009,12 +1701,10 @@ mod tests { "should produce axis arc, ticks, and labels" ); - // Layer 0: axis arc along outer edge let arc = &layers[0]; assert_eq!(arc["mark"]["type"], "arc"); assert_eq!(arc["mark"]["fill"], json!(null)); - // Layer 1: ticks (one per break) let ticks = &layers[1]; assert_eq!(ticks["mark"]["type"], "rule"); assert_eq!(ticks["data"]["values"].as_array().unwrap().len(), 3); @@ -2025,7 +1715,6 @@ mod tests { .collect(); assert_eq!(tick_fields, vec!["theta", "cx", "cy", "x", "y", "x2", "y2"]); - // Layer 2: nested label layer with shared data/transforms/encoding let labels = &layers[2]; assert_eq!(labels["encoding"]["text"]["field"], "label"); assert_eq!(labels["data"]["values"].as_array().unwrap().len(), 3); @@ -2038,7 +1727,6 @@ mod tests { assert_eq!(sub["mark"]["type"], "text"); assert!(sub["mark"]["align"].is_string()); assert!(sub["mark"]["baseline"].is_string()); - // Each sub-layer filters by alignment tag assert!(sub["transform"] .as_array() .unwrap() @@ -2064,14 +1752,6 @@ mod tests { // ========================================================================= // Free scales suppress polar decorations - // - // Polar grid lines and axes are drawn as manual VL layers whose positions - // are computed from the global scale domain. With free scales each facet - // panel has its own domain, so the global positions would be wrong. - // Rather than rendering misleading decorations we suppress them entirely. - // Proper per-panel decorations would require computing per-group domains - // and threading them into the decoration data — a significant lift that - // is not yet implemented. // ========================================================================= fn facet_with_free(free: Vec) -> Facet { @@ -2250,7 +1930,6 @@ mod tests { convert_polar_to_cartesian(&mut layer, &panel).unwrap(); - // 3 categories → domain (0.5, 3.5), full circle → scale = 2π / 3.0 let transforms = layer["transform"].as_array().unwrap(); let theta_calc = transforms .iter() @@ -2454,8 +2133,6 @@ mod tests { convert_polar_to_cartesian(&mut layer, &panel).unwrap(); - // 3 categories → domain (0.5, 3.5), scale = 2π/3 - // With band fraction 0.9: effective scale = 2π/3 * 0.9 let expected = 2.0 * std::f64::consts::PI / 3.0 * POLAR_BAND_FRACTION; let transforms = layer["transform"].as_array().unwrap(); let x_calc = transforms @@ -2494,7 +2171,6 @@ mod tests { convert_polar_to_cartesian(&mut layer, &panel).unwrap(); - // Continuous → full scale = 2π/100, no band fraction let full_scale = 2.0 * std::f64::consts::PI / 100.0; let with_band = full_scale * POLAR_BAND_FRACTION; let transforms = layer["transform"].as_array().unwrap(); @@ -2545,7 +2221,6 @@ mod tests { "should produce axis line, ticks, and labels" ); - // Label data should carry category names, not numeric positions let labels = &layers[2]; let values = labels["data"]["values"].as_array().unwrap(); assert_eq!(values.len(), 3); @@ -2553,7 +2228,6 @@ mod tests { assert_eq!(values[1]["label"], "mid"); assert_eq!(values[2]["label"], "high"); - // Numeric positions should be 1, 2, 3 assert_eq!(values[0]["v"], 1.0); assert_eq!(values[1]["v"], 2.0); assert_eq!(values[2]["v"], 3.0); @@ -2576,7 +2250,6 @@ mod tests { "should produce axis arc, ticks, and labels" ); - // Label data should carry category names let labels = &layers[2]; let values = labels["data"]["values"].as_array().unwrap(); assert_eq!(values.len(), 3); @@ -2630,7 +2303,6 @@ mod tests { let f = faceted(); let panel = PolarContext::new(Some(&proj), Some(&f), &[]); - // Faceted panel should use literal pixel values, not width/height signals assert_eq!(panel.cx, "150"); assert_eq!(panel.cy, "150"); assert_eq!(panel.radius, "150"); @@ -2652,7 +2324,6 @@ mod tests { let layers = proj.grid_rings(&theme); assert_eq!(layers.len(), 1); - // Radius expression should use literal pixels (150), not signals let radius_expr = layers[0]["encoding"]["radius"]["value"]["expr"] .as_str() .unwrap(); @@ -2686,8 +2357,6 @@ mod tests { let panel = PolarContext::new(None, None, &scales); let thetas = &panel.angle_breaks_radians; assert_eq!(thetas.len(), 3); - // 3 categories → domain (0.5, 3.5), breaks at 1, 2, 3 - // theta = 0 + 2π/3 * (break - 0.5) let scale = 2.0 * PI / 3.0; assert!((thetas[0] - scale * 0.5).abs() < 1e-10); assert!((thetas[1] - scale * 1.5).abs() < 1e-10); @@ -2707,7 +2376,6 @@ mod tests { panel.angle_breaks_radians = vec![1.0, 2.0, 3.0]; let layer = polygon_ring(&panel, POLAR_OUTER, None, Value::Null, json!("red")); let values = layer["data"]["values"].as_array().unwrap(); - // 3 thetas + 1 closing vertex = 4 assert_eq!(values.len(), 4); assert_eq!(values[0]["theta"], values[3]["theta"]); } @@ -2723,9 +2391,7 @@ mod tests { panel.angle_breaks_radians = vec![0.5, 1.0, 1.5]; let layer = polygon_ring(&panel, POLAR_OUTER, None, Value::Null, json!("red")); let values = layer["data"]["values"].as_array().unwrap(); - // start + 3 breaks + end + centre(end) + centre(start) + close = 8 assert_eq!(values.len(), 8); - // First and last vertex should be at the same position (closed path) assert_eq!(values[0]["theta"], values[7]["theta"]); assert_eq!(values[0]["r"], values[7]["r"]); } @@ -2745,9 +2411,7 @@ mod tests { let r_start = values[0]["r"].as_f64().unwrap(); let r_break = values[1]["r"].as_f64().unwrap(); let r_end = values[2]["r"].as_f64().unwrap(); - // Break vertex at full radius assert!((r_break - POLAR_OUTER).abs() < 1e-10); - // Start/end vertices corrected inward by cos(π/2) let expected = POLAR_OUTER * (PI / 2.0).cos(); assert!((r_start - expected).abs() < 1e-10); assert!((r_end - expected).abs() < 1e-10); @@ -2759,7 +2423,6 @@ mod tests { panel.angle_breaks_radians = vec![1.0, 2.0, 3.0]; let layer = polygon_ring(&panel, POLAR_OUTER, Some(0.3), json!("white"), Value::Null); let values = layer["data"]["values"].as_array().unwrap(); - // Outer: 3 + 1 closing, Inner: 3 + 1 closing = 8 assert_eq!(values.len(), 8); assert_eq!(layer["mark"]["type"], "line"); assert_eq!(layer["mark"]["fill"], "white"); @@ -2832,7 +2495,6 @@ mod tests { let theme = json!({"axis": {"domainColor": "#333"}}); let layers = proj.angular_axis(&theme); assert!(!layers.is_empty()); - // First layer should be the polygon outline, not an arc assert_eq!(layers[0]["mark"]["type"], "line"); } diff --git a/tree-sitter-ggsql/grammar.js b/tree-sitter-ggsql/grammar.js index b6645f6df..add90a275 100644 --- a/tree-sitter-ggsql/grammar.js +++ b/tree-sitter-ggsql/grammar.js @@ -650,6 +650,7 @@ module.exports = grammar({ $.number, $.boolean, $.null_literal, + $.infinity, $.array ), @@ -975,6 +976,8 @@ module.exports = grammar({ boolean: $ => choice('true', 'false'), + infinity: $ => token(seq(optional(choice('-', '+')), caseInsensitive('Inf'))), + array: $ => choice( seq( '[', @@ -998,7 +1001,8 @@ module.exports = grammar({ $.string, $.number, $.boolean, - $.null_literal + $.null_literal, + $.infinity ), null_literal: $ => caseInsensitive('NULL'), diff --git a/tree-sitter-ggsql/test/corpus/basic.txt b/tree-sitter-ggsql/test/corpus/basic.txt index 72b173176..36f9b15ea 100644 --- a/tree-sitter-ggsql/test/corpus/basic.txt +++ b/tree-sitter-ggsql/test/corpus/basic.txt @@ -3715,3 +3715,40 @@ DRAW point MAPPING x AS x FROM {{ ref('fct_orders') }} (bare_identifier)))) name: (aesthetic_name)))) layer_source: (jinja_template)))))) + +================================================================================ +Infinity literals in SETTING +================================================================================ + +VISUALISE +DRAW point SETTING bounds => [null, 1, Inf, -Inf], limit => inf + +-------------------------------------------------------------------------------- + +(query + (visualise_statement + (visualise_keyword) + (viz_clause + (draw_clause + (geom_type) + (setting_clause + (parameter_assignment + name: (parameter_name + (identifier + (bare_identifier))) + value: (parameter_value + (array + (array_element + (null_literal)) + (array_element + (number)) + (array_element + (infinity)) + (array_element + (infinity))))) + (parameter_assignment + name: (parameter_name + (identifier + (bare_identifier))) + value: (parameter_value + (infinity))))))))