Skip to content
Merged
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
26 changes: 19 additions & 7 deletions crates/data/tests/planetary_crs.rs
Original file line number Diff line number Diff line change
Expand Up @@ -140,6 +140,14 @@ fn planetographic_systems_without_east_axes_are_refused() {
let dir = tempfile::tempdir().unwrap();
let r = 3_396_190.0f64;
let metres = LONLAT.map(|(lon, lat)| (r * lon.to_radians(), r * lat.to_radians()));
// The planetographic systems are on the flattened Mars (2015) ellipsoid, where equidistant
// cylindrical northings are meridian distances (EPSG method 1028), not r·φ.
let f = 1.0 / 169.894447223612;
let eqc = rustgis_geom::EquidistantCylindrical::new(r, f * (2.0 - f), 0.0, 0.0, 0.0, 0.0, 0.0).unwrap();
let on_ellipsoid = LONLAT.map(|(lon, lat)| {
let c = eqc.forward(Coord::new(lon, lat)).unwrap();
Comment thread
giswqs marked this conversation as resolved.
Resolved
(c.x, c.y)
});
let with_axes = |wkt: &str, axes: &str| format!("{},{axes}]", wkt.strip_suffix(']').unwrap());
let east_north = r#"AXIS["Lon",EAST],AXIS["Lat",NORTH]"#;
let unnamed = MARS_OGRAPHIC_ESRI.replace("GCS_Mars_2015_Ographic", "GCS_Mars");
Expand Down Expand Up @@ -194,14 +202,14 @@ fn planetographic_systems_without_east_axes_are_refused() {
("esri_east", with_axes(MARS_OGRAPHIC_ESRI, east_north), LONLAT),
("wkt2_east", MARS_WEST.replace(",west,", ",east,"), LONLAT),
("ocentric_code", with_axes(&unnamed, r#"AUTHORITY["IAU_2015","49900"]"#), LONLAT),
("eqc_easting", eqc_on(&with_axes(MARS_OGRAPHIC_ESRI, east_north), "Mars_Equirectangular"), metres),
("eqc_easting", eqc_on(&with_axes(MARS_OGRAPHIC_ESRI, east_north), "Mars_Equirectangular"), on_ellipsoid),
(
"eqc_base_easting",
with_axes(&eqc_on(MARS_OGRAPHIC_ESRI, "Mars_Equirectangular"), r#"AXIS["Easting",EAST],AXIS["Northing",NORTH]"#),
metres,
on_ellipsoid,
),
("eqc_ocentric_code", with_axes(MARS_EQC, r#"AUTHORITY["IAU_2015","49910"]"#), metres),
("wkt2_eqc_easting", wkt2_eqc("Mars (2015) / Ographic", "east"), metres),
("wkt2_eqc_easting", wkt2_eqc("Mars (2015) / Ographic", "east"), on_ellipsoid),
];
for (name, prj, native) in &read {
let path = shapefile(dir.path(), name, *native, prj);
Expand Down Expand Up @@ -238,8 +246,12 @@ fn planetographic_systems_without_east_axes_are_refused() {
assert!(planetary_notes(&report.messages).is_empty(), "{:?}", report.messages);
let src = source::open_path(&path, None).unwrap();
assert_eq!(said(&src.info().notes), 1, "{:?}", src.info().notes);
let a = 6_378_137.0f64;
let earth_metres = LONLAT.map(|(lon, lat)| (a * lon.to_radians(), a * lat.to_radians()));
let f = 1.0 / 298.257223563;
let wgs84_eqc = rustgis_geom::EquidistantCylindrical::new(6_378_137.0, f * (2.0 - f), 0.0, 0.0, 0.0, 0.0, 0.0).unwrap();
let earth_metres = LONLAT.map(|(lon, lat)| {
let c = wgs84_eqc.forward(Coord::new(lon, lat)).unwrap();
(c.x, c.y)
});
let earth_eqc = eqc_on(&west, "World_Equidistant_Cylindrical");
assert!(earth_eqc.contains("WGS_1984") && earth_eqc.contains("WEST"), "{earth_eqc}");
let path = shapefile(dir.path(), "earth_eqc_west", earth_metres, &earth_eqc);
Expand Down Expand Up @@ -748,9 +760,9 @@ fn a_mars_projection_in_a_geotiff_is_said_or_refused() {
// where the numbers look like degrees, as they are metres.
let mut unread = eqc.clone();
let at = unread.chunks(4).position(|k| k[0] == 3075).unwrap() * 4 + 3;
unread[at] = 23;
unread[at] = 26;
Comment thread
giswqs marked this conversation as resolved.
Resolved
let e = raster::read_bytes("unread.tif", &geotiff(-180.0, &unread, &[3_396_190.0], "")).unwrap_err().to_string();
assert!(e.contains(ANOTHER_BODY) && e.contains("transformation 23"), "{e}");
assert!(e.contains(ANOTHER_BODY) && e.contains("transformation 26"), "{e}");
// On the Earth it is read as longitude/latitude where the numbers fit, as before.
let earth = raster::read_bytes("earth.tif", &geotiff(-180.0, &unread, &[6_378_137.0], "")).unwrap();
assert_eq!(earth.crs, Crs::LonLat);
Expand Down
Loading