This is an automated email from the ASF dual-hosted git repository. jiayu pushed a commit to branch branch-1.8.1 in repository https://gitbox.apache.org/repos/asf/sedona.git
commit a979b5a5c9e31916dbd6d6723034a14c3d2f2626 Author: Jia Yu <[email protected]> AuthorDate: Mon Mar 2 12:23:47 2026 -0700 [GH-2674] Add RS_SetCRS and RS_CRS for custom CRS string support (#2677) --- .../sedona/common/raster/CrsNormalization.java | 243 +++++++ .../sedona/common/raster/RasterAccessors.java | 101 +++ .../apache/sedona/common/raster/RasterEditors.java | 159 +++++ .../common/raster/CrsRoundTripComplianceTest.java | 697 +++++++++++++++++++++ .../sedona/common/raster/RasterAccessorsTest.java | 68 ++ .../sedona/common/raster/RasterEditorsTest.java | 146 +++++ docs/api/sql/Raster-Functions.md | 2 + docs/api/sql/Raster-Operators/RS_CRS.md | 101 +++ docs/api/sql/Raster-Operators/RS_SRID.md | 2 +- docs/api/sql/Raster-Operators/RS_SetCRS.md | 68 ++ pom.xml | 2 +- .../scala/org/apache/sedona/sql/UDF/Catalog.scala | 2 + .../expressions/raster/RasterAccessors.scala | 9 + .../expressions/raster/RasterEditors.scala | 7 + .../org/apache/sedona/sql/rasteralgebraTest.scala | 71 +++ 15 files changed, 1676 insertions(+), 2 deletions(-) diff --git a/common/src/main/java/org/apache/sedona/common/raster/CrsNormalization.java b/common/src/main/java/org/apache/sedona/common/raster/CrsNormalization.java new file mode 100644 index 0000000000..c342360339 --- /dev/null +++ b/common/src/main/java/org/apache/sedona/common/raster/CrsNormalization.java @@ -0,0 +1,243 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one + * or more contributor license agreements. See the NOTICE file + * distributed with this work for additional information + * regarding copyright ownership. The ASF licenses this file + * to you under the Apache License, Version 2.0 (the + * "License"); you may not use this file except in compliance + * with the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, + * software distributed under the License is distributed on an + * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + * KIND, either express or implied. See the License for the + * specific language governing permissions and limitations + * under the License. + */ +package org.apache.sedona.common.raster; + +import java.util.HashMap; +import java.util.HashSet; +import java.util.Locale; +import java.util.Map; +import java.util.Set; +import java.util.concurrent.ConcurrentHashMap; +import java.util.regex.Matcher; +import java.util.regex.Pattern; +import org.geotools.api.referencing.operation.MathTransformFactory; +import org.geotools.api.referencing.operation.OperationMethod; +import org.geotools.api.referencing.operation.Projection; +import org.geotools.referencing.ReferencingFactoryFinder; +import org.geotools.referencing.operation.DefaultMathTransformFactory; + +/** + * Centralized CRS name normalization for bridging GeoTools ↔ proj4sedona naming differences. + * + * <p>GeoTools and proj4sedona use different canonical names for some map projections (e.g. GeoTools + * uses "Mercator_2SP" while proj4sedona uses "Mercator"). This class normalizes names once at the + * import/export boundaries so downstream code operates on compatible names. + * + * <p>All lookups use pre-computed maps; normalization is O(1) after class initialization. + */ +final class CrsNormalization { + + private CrsNormalization() {} + + // Shared regex for extracting projection names from WKT1 strings + static final Pattern PROJECTION_PATTERN = Pattern.compile("PROJECTION\\[\"([^\"]+)\"\\]"); + + // ===================================================================== + // Import direction: proj4sedona → GeoTools + // ===================================================================== + + // Hardcoded fallback for proj4sedona names with no GeoTools alias. + // Keys are pre-normalized (lowercase, no spaces/underscores) for O(1) lookup. + // Verified via exhaustive testing of all 58 proj4sedona registered projection names. + private static final Map<String, String> PROJ4SEDONA_TO_GEOTOOLS; + + static { + Map<String, String> m = new HashMap<>(); + m.put("lambertcylindricalequalarea", "Cylindrical_Equal_Area"); + m.put("extendedtransversemercator", "Transverse_Mercator"); + m.put("lamberttangentialconformalconicprojection", "Lambert_Conformal_Conic"); + m.put("mercatorvarianta", "Mercator_1SP"); + m.put("polarstereographicvarianta", "Polar_Stereographic"); + m.put("polarstereographicvariantb", "Polar_Stereographic"); + m.put("universaltransversemercatorsystem", "Transverse_Mercator"); + m.put("universaltransversemercator", "Transverse_Mercator"); + PROJ4SEDONA_TO_GEOTOOLS = m; + } + + // ===================================================================== + // Export direction: GeoTools → proj4sedona + // ===================================================================== + + // GeoTools canonical names that proj4sedona does not recognize. + private static final Map<String, String> GEOTOOLS_TO_PROJ4SEDONA; + + static { + Map<String, String> m = new HashMap<>(); + m.put("Mercator_2SP", "Mercator"); + GEOTOOLS_TO_PROJ4SEDONA = m; + } + + // ===================================================================== + // GeoTools alias caches (lazy-initialized, thread-safe) + // ===================================================================== + + // aliasCache: exact alias string → canonical OGC name + // normalizedCache: normalized form → set of canonical names (for disambiguation) + private static volatile Map<String, String> aliasCache; + private static volatile Map<String, Set<String>> normalizedCache; + + // ===================================================================== + // Public API + // ===================================================================== + + /** + * Normalize WKT1 projection names from proj4sedona output for GeoTools consumption. Uses a + * three-tier resolution strategy: + * + * <ol> + * <li><b>Exact alias matching</b> — direct lookup against all GeoTools registered aliases (OGC, + * EPSG, GeoTIFF, ESRI, PROJ authorities). + * <li><b>Normalized matching</b> — case-insensitive, ignoring spaces/underscores. Only used + * when the normalized form maps unambiguously to a single canonical name. + * <li><b>Hardcoded fallback</b> — pre-normalized lookup for proj4sedona-specific names that + * have no equivalent in GeoTools' alias database. + * </ol> + * + * @param wkt1 The WKT1 string from proj4sedona. + * @return The WKT1 string with the projection name normalized for GeoTools. + */ + static String normalizeWkt1ForGeoTools(String wkt1) { + return replaceProjectionName(wkt1, CrsNormalization::resolveForGeoTools); + } + + /** + * Normalize WKT1 projection names from GeoTools output for proj4sedona consumption. + * + * @param wkt1 The WKT1 string from GeoTools. + * @return The WKT1 string with the projection name normalized for proj4sedona. + */ + static String normalizeWkt1ForProj4sedona(String wkt1) { + return replaceProjectionName(wkt1, GEOTOOLS_TO_PROJ4SEDONA::get); + } + + // ===================================================================== + // Internal implementation + // ===================================================================== + + /** Functional interface for projection name lookup. */ + @FunctionalInterface + private interface NameResolver { + /** Return the replacement name, or null if no mapping exists. */ + String resolve(String projName); + } + + /** + * Replace the projection name in a WKT1 string using the given resolver. Shared logic for both + * import and export directions. + */ + private static String replaceProjectionName(String wkt1, NameResolver resolver) { + Matcher m = PROJECTION_PATTERN.matcher(wkt1); + if (m.find()) { + String projName = m.group(1); + String resolved = resolver.resolve(projName); + if (resolved != null && !resolved.equals(projName)) { + return wkt1.substring(0, m.start(1)) + resolved + wkt1.substring(m.end(1)); + } + } + return wkt1; + } + + /** + * Three-tier resolution: proj4sedona projection name → canonical GeoTools name. + * + * @return The resolved name, or null if the name is already compatible. + */ + private static String resolveForGeoTools(String projName) { + ensureCachesBuilt(); + + // Tier 1: Exact alias match from GeoTools + String resolved = aliasCache.get(projName); + if (resolved != null) { + return resolved; + } + + // Tier 2: Normalized match (handles space/underscore/case differences) + String normalized = normalizeForMatch(projName); + Set<String> candidates = normalizedCache.get(normalized); + if (candidates != null && candidates.size() == 1) { + String canonical = candidates.iterator().next(); + aliasCache.put(projName, canonical); // cache for next time + return canonical; + } + + // Tier 3: Pre-normalized hardcoded fallback (O(1) lookup) + String fallback = PROJ4SEDONA_TO_GEOTOOLS.get(normalized); + if (fallback != null) { + aliasCache.put(projName, fallback); // cache for next time + return fallback; + } + + return null; // name is already compatible or unknown + } + + /** + * Normalize a string for loose matching: lowercase, remove spaces and underscores. + * + * @param name The name to normalize. + * @return Normalized form (e.g. "Lambert_Conformal_Conic_2SP" → "lambertconformalconic2sp"). + */ + static String normalizeForMatch(String name) { + return name.toLowerCase(Locale.ROOT).replaceAll("[_ ]", ""); + } + + /** + * Build GeoTools alias caches from all registered {@link OperationMethod} objects. Thread-safe + * via double-checked locking. Called at most once. + */ + private static void ensureCachesBuilt() { + if (aliasCache != null) { + return; + } + synchronized (CrsNormalization.class) { + if (aliasCache != null) { + return; + } + MathTransformFactory mtf = ReferencingFactoryFinder.getMathTransformFactory(null); + DefaultMathTransformFactory factory; + if (mtf instanceof DefaultMathTransformFactory) { + factory = (DefaultMathTransformFactory) mtf; + } else { + factory = new DefaultMathTransformFactory(); + } + Set<OperationMethod> methods = factory.getAvailableMethods(Projection.class); + + Map<String, String> aliases = new ConcurrentHashMap<>(); + Map<String, Set<String>> normalized = new HashMap<>(); + + for (OperationMethod method : methods) { + String canonical = method.getName().getCode(); + aliases.put(canonical, canonical); + normalized + .computeIfAbsent(normalizeForMatch(canonical), k -> new HashSet<>()) + .add(canonical); + if (method.getAlias() != null) { + for (Object alias : method.getAlias()) { + String aliasName = alias.toString().replaceAll("^[^:]+:", ""); + aliases.put(aliasName, canonical); + normalized + .computeIfAbsent(normalizeForMatch(aliasName), k -> new HashSet<>()) + .add(canonical); + } + } + } + aliasCache = aliases; + normalizedCache = normalized; + } + } +} diff --git a/common/src/main/java/org/apache/sedona/common/raster/RasterAccessors.java b/common/src/main/java/org/apache/sedona/common/raster/RasterAccessors.java index 49a9223908..4f4a9f5df2 100644 --- a/common/src/main/java/org/apache/sedona/common/raster/RasterAccessors.java +++ b/common/src/main/java/org/apache/sedona/common/raster/RasterAccessors.java @@ -21,8 +21,10 @@ package org.apache.sedona.common.raster; import java.awt.geom.Point2D; import java.awt.image.RenderedImage; import java.util.Arrays; +import java.util.Locale; import java.util.Set; import org.apache.sedona.common.utils.RasterUtils; +import org.datasyslab.proj4sedona.core.Proj; import org.geotools.api.referencing.FactoryException; import org.geotools.api.referencing.ReferenceIdentifier; import org.geotools.api.referencing.crs.CoordinateReferenceSystem; @@ -31,12 +33,14 @@ import org.geotools.coverage.grid.GridCoverage2D; import org.geotools.coverage.grid.GridEnvelope2D; import org.geotools.referencing.crs.DefaultEngineeringCRS; import org.geotools.referencing.operation.transform.AffineTransform2D; +import org.geotools.referencing.wkt.Formattable; import org.locationtech.jts.geom.Coordinate; import org.locationtech.jts.geom.Geometry; import org.locationtech.jts.geom.GeometryFactory; import org.locationtech.jts.geom.Point; public class RasterAccessors { + public static int srid(GridCoverage2D raster) throws FactoryException { CoordinateReferenceSystem crs = raster.getCoordinateReferenceSystem(); if (crs instanceof DefaultEngineeringCRS) { @@ -359,4 +363,101 @@ public class RasterAccessors { (int) meta[10], (int) meta[11]); } + + /** + * Returns the CRS of a raster as PROJJSON string. + * + * @param raster The input raster. + * @return The CRS definition as PROJJSON string, or null if no CRS is set. + */ + public static String crs(GridCoverage2D raster) { + return crs(raster, "projjson"); + } + + /** + * Returns the CRS of a raster in the specified format. + * + * @param raster The input raster. + * @param format The desired output format: "projjson", "wkt2", "wkt1", or "proj". + * @return The CRS definition string in the requested format, or null if no CRS is set. + */ + public static String crs(GridCoverage2D raster, String format) { + String fmt; + if (format == null || format.trim().isEmpty()) { + fmt = "projjson"; + } else { + fmt = format.toLowerCase(Locale.ROOT).trim(); + } + CoordinateReferenceSystem crsDef = raster.getCoordinateReferenceSystem(); + if (crsDef instanceof DefaultEngineeringCRS) { + if (((DefaultEngineeringCRS) crsDef).isWildcard()) { + return null; + } + } + + // Get WKT1 representation from GeoTools (native, no conversion needed) + String wkt1; + if (crsDef instanceof Formattable) { + wkt1 = ((Formattable) crsDef).toWKT(2, false); + } else { + wkt1 = crsDef.toWKT(); + } + + if ("wkt1".equals(fmt) || "wkt".equals(fmt)) { + return wkt1; + } + + // For all other formats, convert through proj4sedona. + // Prefer EPSG SRID when available: GeoTools WKT1 projection names (e.g. Mercator_2SP) + // may not be recognized by proj4sedona, but EPSG codes always work. + try { + Proj proj; + int srid = srid(raster); + if (srid > 0) { + try { + proj = new Proj("EPSG:" + srid); + } catch (Exception e) { + // EPSG code not recognized by proj4sedona, fall back to WKT1 + proj = createProjFromWkt1(wkt1); + } + } else { + proj = createProjFromWkt1(wkt1); + } + switch (fmt) { + case "projjson": + return proj.toProjJson(); + case "wkt2": + return proj.toWkt2(); + case "proj": + case "proj4": + return proj.toProjString(); + default: + throw new IllegalArgumentException( + "Unsupported CRS format: '" + + format + + "'. Supported formats: projjson, wkt2, wkt1, proj"); + } + } catch (IllegalArgumentException e) { + throw e; + } catch (Exception e) { + throw new RuntimeException( + "Failed to convert CRS to format '" + format + "': " + e.getMessage(), e); + } + } + + /** + * Create a Proj object from GeoTools WKT1, with normalization fallback for projection names that + * proj4sedona does not recognize (e.g. Mercator_2SP). + */ + private static Proj createProjFromWkt1(String wkt1) { + try { + return new Proj(wkt1); + } catch (Exception wktError) { + String normalized = CrsNormalization.normalizeWkt1ForProj4sedona(wkt1); + if (!normalized.equals(wkt1)) { + return new Proj(normalized); + } + throw wktError; + } + } } diff --git a/common/src/main/java/org/apache/sedona/common/raster/RasterEditors.java b/common/src/main/java/org/apache/sedona/common/raster/RasterEditors.java index f28fd77c6a..818ba65e83 100644 --- a/common/src/main/java/org/apache/sedona/common/raster/RasterEditors.java +++ b/common/src/main/java/org/apache/sedona/common/raster/RasterEditors.java @@ -26,15 +26,19 @@ import java.awt.image.*; import java.util.Arrays; import java.util.Map; import java.util.Objects; +import java.util.regex.Matcher; +import java.util.regex.Pattern; import javax.media.jai.Interpolation; import javax.media.jai.RasterFactory; import org.apache.sedona.common.FunctionsGeoTools; import org.apache.sedona.common.utils.RasterInterpolate; import org.apache.sedona.common.utils.RasterUtils; +import org.datasyslab.proj4sedona.core.Proj; import org.geotools.api.coverage.grid.GridCoverage; import org.geotools.api.coverage.grid.GridGeometry; import org.geotools.api.metadata.spatial.PixelOrientation; import org.geotools.api.referencing.FactoryException; +import org.geotools.api.referencing.crs.CRSFactory; import org.geotools.api.referencing.crs.CoordinateReferenceSystem; import org.geotools.api.referencing.datum.PixelInCell; import org.geotools.api.referencing.operation.MathTransform; @@ -48,8 +52,11 @@ import org.geotools.coverage.grid.GridEnvelope2D; import org.geotools.coverage.grid.GridGeometry2D; import org.geotools.coverage.processing.Operations; import org.geotools.geometry.jts.ReferencedEnvelope; +import org.geotools.referencing.CRS; +import org.geotools.referencing.ReferencingFactoryFinder; import org.geotools.referencing.crs.DefaultEngineeringCRS; import org.geotools.referencing.operation.transform.AffineTransform2D; +import org.geotools.util.factory.Hints; import org.locationtech.jts.index.strtree.STRtree; public class RasterEditors { @@ -102,7 +109,159 @@ public class RasterEditors { } else { crs = FunctionsGeoTools.sridToCRS(srid); } + return replaceCrs(raster, crs); + } + + /** + * Sets the CRS of a raster using a CRS string. Accepts EPSG codes (e.g. "EPSG:4326"), WKT1, WKT2, + * PROJ strings, and PROJJSON. + * + * @param raster The input raster. + * @param crsString The CRS definition string. + * @return The raster with the new CRS. + */ + public static GridCoverage2D setCrs(GridCoverage2D raster, String crsString) { + CoordinateReferenceSystem crs = parseCrsString(crsString); + return replaceCrs(raster, crs); + } + + /** + * Parse a CRS string in any supported format into a GeoTools CoordinateReferenceSystem. + * + * <p>Parsing priority: + * + * <ol> + * <li>GeoTools CRS.decode — handles authority codes like EPSG:4326 + * <li>GeoTools CRS.parseWKT — handles WKT1 strings + * <li>proj4sedona — handles WKT2, PROJ strings, PROJJSON. If an EPSG authority can be resolved, + * uses CRS.decode for a lossless result. Otherwise falls back to WKT1 conversion. + * </ol> + * + * @param crsString The CRS definition string. + * @return The parsed CoordinateReferenceSystem. + * @throws IllegalArgumentException if the CRS string cannot be parsed. + */ + static CoordinateReferenceSystem parseCrsString(String crsString) { + if (crsString == null || crsString.trim().isEmpty()) { + throw new IllegalArgumentException( + "CRS string must not be null or empty. " + + "Supported formats: EPSG code (e.g. 'EPSG:4326'), WKT1, WKT2, PROJ string, PROJJSON."); + } + + // Step 1: Try GeoTools CRS.decode (handles EPSG:xxxx, AUTO:xxxx, etc.) + try { + return CRS.decode(crsString, true); + } catch (FactoryException e) { + // Not an authority code, continue + } + + // Step 2: Try GeoTools WKT parsing with longitude-first axis order (handles WKT1) + Hints hints = new Hints(Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER, Boolean.TRUE); + CRSFactory crsFactory = ReferencingFactoryFinder.getCRSFactory(hints); + try { + return crsFactory.createFromWKT(crsString); + } catch (FactoryException e) { + // Not WKT1, continue + } + + // Step 3: Use proj4sedona (handles WKT2, PROJ, PROJJSON) + Exception lastError = null; + try { + Proj proj = new Proj(crsString); + + // Try to resolve to an EPSG authority code for a lossless result + String authority = proj.toEpsgCode(); + if (authority != null && !authority.isEmpty()) { + try { + return CRS.decode(authority, true); + } catch (FactoryException ex) { + // Authority code not recognized by GeoTools, fall through to WKT1 + } + } + + // Fallback: convert to WKT1 via proj4sedona and parse with GeoTools. + // proj4sedona may include parameters GeoTools doesn't expect (e.g. standard_parallel_1 + // for projections that don't use it). We handle this by trying several parse strategies: + // 1. Raw WKT1 (proj4sedona's projection names may already be recognized by GeoTools) + // 2. Normalized WKT1 (resolve projection names to canonical OGC names) + // 3. Strip unexpected parameters iteratively + String wkt1 = proj.toWkt1(); + if (wkt1 != null && !wkt1.isEmpty()) { + + // Strategy 1: Try raw WKT1 directly + try { + return crsFactory.createFromWKT(wkt1); + } catch (FactoryException ex) { + // Raw WKT1 failed, continue with normalization + } + + // Strategy 2: Try with normalized projection name + String normalizedWkt = CrsNormalization.normalizeWkt1ForGeoTools(wkt1); + // Strategy 3: If parsing fails due to unexpected parameters, strip them iteratively. + // proj4sedona sometimes includes parameters like standard_parallel_1 for projections + // that don't use it. We parse the error message to identify and remove the offending + // parameter, then retry. + String currentWkt = normalizedWkt; + for (int attempt = 0; attempt < 5; attempt++) { + try { + return crsFactory.createFromWKT(currentWkt); + } catch (FactoryException ex) { + lastError = ex; + String msg = ex.getMessage(); + if (msg != null) { + Matcher paramMatcher = UNEXPECTED_PARAM_PATTERN.matcher(msg); + if (paramMatcher.find()) { + String stripped = stripWktParameter(currentWkt, paramMatcher.group(1)); + if (stripped.equals(currentWkt)) { + break; // Strip was a no-op, give up + } + currentWkt = stripped; + continue; + } + } + break; // Different kind of error, give up + } + } + } + } catch (RuntimeException e) { + lastError = e; + } + + IllegalArgumentException error = + new IllegalArgumentException( + "Cannot parse CRS string. Supported formats: EPSG code (e.g. 'EPSG:4326'), " + + "WKT1, WKT2, PROJ string, PROJJSON. Input: " + + crsString); + if (lastError != null) { + error.addSuppressed(lastError); + } + throw error; + } + + private static final Pattern UNEXPECTED_PARAM_PATTERN = + Pattern.compile("Parameter \"([^\"]+)\" was not expected"); + + /** + * Strip a named PARAMETER from a WKT1 string. Used to remove parameters that proj4sedona includes + * but GeoTools does not expect (e.g. standard_parallel_1 for Transverse Mercator). + * + * @param wkt The WKT1 string. + * @param paramName The parameter name to strip (e.g. "standard_parallel_1"). + * @return The WKT1 string with the parameter removed. + */ + private static String stripWktParameter(String wkt, String paramName) { + // Remove ,PARAMETER["paramName",value] or PARAMETER["paramName",value], + String escaped = Pattern.quote(paramName); + Pattern pattern = Pattern.compile(",\\s*PARAMETER\\[\"" + escaped + "\",[^\\]]*\\]"); + String result = pattern.matcher(wkt).replaceAll(""); + if (result.equals(wkt)) { + Pattern pattern2 = Pattern.compile("PARAMETER\\[\"" + escaped + "\",[^\\]]*\\]\\s*,?"); + result = pattern2.matcher(wkt).replaceAll(""); + } + return result; + } + private static GridCoverage2D replaceCrs(GridCoverage2D raster, CoordinateReferenceSystem crs) { GridCoverageFactory gridCoverageFactory = CoverageFactoryFinder.getGridCoverageFactory(null); MathTransform2D transform = raster.getGridGeometry().getGridToCRS2D(); Map<?, ?> properties = raster.getProperties(); diff --git a/common/src/test/java/org/apache/sedona/common/raster/CrsRoundTripComplianceTest.java b/common/src/test/java/org/apache/sedona/common/raster/CrsRoundTripComplianceTest.java new file mode 100644 index 0000000000..19cc815491 --- /dev/null +++ b/common/src/test/java/org/apache/sedona/common/raster/CrsRoundTripComplianceTest.java @@ -0,0 +1,697 @@ +/* + * Licensed to the Apache Software Foundation (ASF) under one + * or more contributor license agreements. See the NOTICE file + * distributed with this work for additional information + * regarding copyright ownership. The ASF licenses this file + * to you under the Apache License, Version 2.0 (the + * "License"); you may not use this file except in compliance + * with the License. You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, + * software distributed under the License is distributed on an + * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + * KIND, either express or implied. See the License for the + * specific language governing permissions and limitations + * under the License. + */ +package org.apache.sedona.common.raster; + +import static org.junit.Assert.*; + +import java.util.regex.Matcher; +import java.util.regex.Pattern; +import org.geotools.api.referencing.FactoryException; +import org.geotools.coverage.grid.GridCoverage2D; +import org.junit.Test; + +/** + * Round-trip compliance tests for RS_SetCRS and RS_CRS across representative EPSG codes. + * + * <p>For each EPSG code and each format (PROJ, PROJJSON, WKT1, WKT2), this test: + * + * <ol> + * <li>Creates a raster with that CRS via RS_SetCRS("EPSG:xxxx") + * <li>Exports the CRS via RS_CRS(raster, format) + * <li>Re-imports the exported string via RS_SetCRS(exportedString) + * <li>Re-exports via RS_CRS(raster2, format) and verifies the exported string is identical + * </ol> + */ +public class CrsRoundTripComplianceTest extends RasterTestBase { + + private static final Pattern WKT1_PROJECTION_PATTERN = + Pattern.compile("PROJECTION\\[\"([^\"]+)\""); + private static final Pattern WKT1_AUTHORITY_PATTERN = + Pattern.compile("AUTHORITY\\[\"EPSG\",\\s*\"(\\d+)\"\\]"); + + // --------------------------------------------------------------------------- + // PROJ format round-trip tests + // --------------------------------------------------------------------------- + + @Test + public void testProjRoundTrip_Geographic_4326() throws FactoryException { + assertProjRoundTrip(4326); + } + + @Test + public void testProjRoundTrip_Geographic_NAD83_4269() throws FactoryException { + assertProjRoundTrip(4269); + } + + @Test + public void testProjRoundTrip_TransverseMercator_32617() throws FactoryException { + assertProjRoundTrip(32617); + } + + @Test + public void testProjRoundTrip_PseudoMercator_3857() throws FactoryException { + assertProjRoundTrip(3857); + } + + @Test + public void testProjRoundTrip_Mercator1SP_3395() throws FactoryException { + assertProjRoundTrip(3395); + } + + @Test + public void testProjRoundTrip_LambertConformalConic2SP_2154() throws FactoryException { + assertProjRoundTrip(2154); + } + + @Test + public void testProjRoundTrip_LambertAzimuthalEqualArea_Spherical_2163() throws FactoryException { + assertProjRoundTrip(2163); + } + + @Test + public void testProjRoundTrip_AlbersEqualArea_5070() throws FactoryException { + assertProjRoundTrip(5070); + } + + @Test + public void testProjRoundTrip_ObliqueStereographic_28992() throws FactoryException { + assertProjRoundTrip(28992); + } + + @Test + public void testProjRoundTrip_PolarStereographicB_3031() throws FactoryException { + assertProjRoundTrip(3031); + } + + @Test + public void testProjRoundTrip_LambertAzimuthalEqualArea_3035() throws FactoryException { + assertProjRoundTrip(3035); + } + + @Test + public void testProjRoundTrip_Mercator1SP_Spherical_3785() throws FactoryException { + assertProjRoundTrip(3785); + } + + @Test + public void testProjRoundTrip_EquidistantCylindrical_4087() throws FactoryException { + assertProjRoundTrip(4087); + } + + @Test + public void testProjRoundTrip_PolarStereographicA_32661() throws FactoryException { + assertProjRoundTrip(32661); + } + + @Test + public void testProjRoundTrip_TransverseMercator_OSGB_27700() throws FactoryException { + assertProjRoundTrip(27700); + } + + @Test + public void testProjRoundTrip_AlbersEqualArea_Australian_3577() throws FactoryException { + assertProjRoundTrip(3577); + } + + @Test + public void testProjRoundTrip_LambertConformalConic2SP_Vicgrid_3111() throws FactoryException { + assertProjRoundTrip(3111); + } + + @Test + public void testProjRoundTrip_PolarStereographicB_NSIDC_3413() throws FactoryException { + assertProjRoundTrip(3413); + } + + @Test + public void testProjRoundTrip_LambertAzimuthalEqualArea_EASE_6931() throws FactoryException { + assertProjRoundTrip(6931); + } + + // --------------------------------------------------------------------------- + // PROJJSON format round-trip tests + // --------------------------------------------------------------------------- + + @Test + public void testProjJsonRoundTrip_Geographic_4326() throws FactoryException { + assertProjJsonRoundTrip(4326); + } + + @Test + public void testProjJsonRoundTrip_Geographic_NAD83_4269() throws FactoryException { + assertProjJsonRoundTrip(4269); + } + + @Test + public void testProjJsonRoundTrip_TransverseMercator_32617() throws FactoryException { + assertProjJsonRoundTrip(32617); + } + + @Test + public void testProjJsonRoundTrip_PseudoMercator_3857() throws FactoryException { + assertProjJsonRoundTrip(3857); + } + + @Test + public void testProjJsonRoundTrip_Mercator1SP_3395() throws FactoryException { + assertProjJsonRoundTrip(3395); + } + + @Test + public void testProjJsonRoundTrip_LambertConformalConic2SP_2154() throws FactoryException { + assertProjJsonRoundTrip(2154); + } + + @Test + public void testProjJsonRoundTrip_AlbersEqualArea_5070() throws FactoryException { + assertProjJsonRoundTrip(5070); + } + + @Test + public void testProjJsonRoundTrip_ObliqueStereographic_28992() throws FactoryException { + assertProjJsonRoundTrip(28992); + } + + @Test + public void testProjJsonRoundTrip_PolarStereographicB_3031() throws FactoryException { + assertProjJsonRoundTrip(3031); + } + + @Test + public void testProjJsonRoundTrip_LambertAzimuthalEqualArea_3035() throws FactoryException { + assertProjJsonRoundTrip(3035); + } + + @Test + public void testProjJsonRoundTrip_EquidistantCylindrical_4087() throws FactoryException { + assertProjJsonRoundTrip(4087); + } + + @Test + public void testProjJsonRoundTrip_PolarStereographicA_32661() throws FactoryException { + assertProjJsonRoundTrip(32661); + } + + @Test + public void testProjJsonRoundTrip_TransverseMercator_OSGB_27700() throws FactoryException { + assertProjJsonRoundTrip(27700); + } + + @Test + public void testProjJsonRoundTrip_AlbersEqualArea_Australian_3577() throws FactoryException { + assertProjJsonRoundTrip(3577); + } + + @Test + public void testProjJsonRoundTrip_LambertConformalConic2SP_Vicgrid_3111() + throws FactoryException { + assertProjJsonRoundTrip(3111); + } + + @Test + public void testProjJsonRoundTrip_PolarStereographicB_NSIDC_3413() throws FactoryException { + assertProjJsonRoundTrip(3413); + } + + @Test + public void testProjJsonRoundTrip_LambertAzimuthalEqualArea_EASE_6931() throws FactoryException { + assertProjJsonRoundTrip(6931); + } + + // --------------------------------------------------------------------------- + // WKT1 format round-trip tests + // WKT1 includes AUTHORITY["EPSG","xxxx"] so SRID is always preserved. + // --------------------------------------------------------------------------- + + @Test + public void testWkt1RoundTrip_Geographic_4326() throws FactoryException { + assertWkt1RoundTrip(4326); + } + + @Test + public void testWkt1RoundTrip_Geographic_NAD83_4269() throws FactoryException { + assertWkt1RoundTrip(4269); + } + + @Test + public void testWkt1RoundTrip_TransverseMercator_32617() throws FactoryException { + assertWkt1RoundTrip(32617); + } + + @Test + public void testWkt1RoundTrip_PseudoMercator_3857() throws FactoryException { + assertWkt1RoundTrip(3857); + } + + @Test + public void testWkt1RoundTrip_Mercator1SP_3395() throws FactoryException { + assertWkt1RoundTrip(3395); + } + + @Test + public void testWkt1RoundTrip_LambertConformalConic2SP_2154() throws FactoryException { + assertWkt1RoundTrip(2154); + } + + @Test + public void testWkt1RoundTrip_LambertAzimuthalEqualArea_Spherical_2163() throws FactoryException { + assertWkt1RoundTrip(2163); + } + + @Test + public void testWkt1RoundTrip_AlbersEqualArea_5070() throws FactoryException { + assertWkt1RoundTrip(5070); + } + + @Test + public void testWkt1RoundTrip_ObliqueStereographic_28992() throws FactoryException { + assertWkt1RoundTrip(28992); + } + + @Test + public void testWkt1RoundTrip_PolarStereographicB_3031() throws FactoryException { + assertWkt1RoundTrip(3031); + } + + @Test + public void testWkt1RoundTrip_LambertAzimuthalEqualArea_3035() throws FactoryException { + assertWkt1RoundTrip(3035); + } + + @Test + public void testWkt1RoundTrip_Mercator1SP_Spherical_3785() throws FactoryException { + assertWkt1RoundTrip(3785); + } + + @Test + public void testWkt1RoundTrip_EquidistantCylindrical_4087() throws FactoryException { + assertWkt1RoundTrip(4087); + } + + @Test + public void testWkt1RoundTrip_PolarStereographicA_32661() throws FactoryException { + assertWkt1RoundTrip(32661); + } + + @Test + public void testWkt1RoundTrip_TransverseMercator_OSGB_27700() throws FactoryException { + assertWkt1RoundTrip(27700); + } + + @Test + public void testWkt1RoundTrip_AlbersEqualArea_Australian_3577() throws FactoryException { + assertWkt1RoundTrip(3577); + } + + @Test + public void testWkt1RoundTrip_LambertConformalConic2SP_Vicgrid_3111() throws FactoryException { + assertWkt1RoundTrip(3111); + } + + @Test + public void testWkt1RoundTrip_PolarStereographicB_NSIDC_3413() throws FactoryException { + assertWkt1RoundTrip(3413); + } + + @Test + public void testWkt1RoundTrip_LambertAzimuthalEqualArea_EASE_6931() throws FactoryException { + assertWkt1RoundTrip(6931); + } + + @Test + public void testWkt1RoundTrip_Krovak_2065() throws FactoryException { + // Krovak fails for PROJ/PROJJSON export but WKT1 is GeoTools-native, so it works + assertWkt1RoundTrip(2065); + } + + @Test + public void testWkt1RoundTrip_HotineObliqueMercator_2056() throws FactoryException { + // Hotine Oblique Mercator fails for PROJ/PROJJSON export but works for WKT1 + assertWkt1RoundTrip(2056); + } + + // --------------------------------------------------------------------------- + // WKT2 format round-trip tests + // WKT2 goes through proj4sedona for both export and import. + // --------------------------------------------------------------------------- + + @Test + public void testWkt2RoundTrip_Geographic_4326() throws FactoryException { + assertWkt2RoundTrip(4326); + } + + @Test + public void testWkt2RoundTrip_Geographic_NAD83_4269() throws FactoryException { + assertWkt2RoundTrip(4269); + } + + @Test + public void testWkt2RoundTrip_TransverseMercator_32617() throws FactoryException { + assertWkt2RoundTrip(32617); + } + + @Test + public void testWkt2RoundTrip_PseudoMercator_3857() throws FactoryException { + assertWkt2RoundTrip(3857); + } + + @Test + public void testWkt2RoundTrip_Mercator1SP_3395() throws FactoryException { + assertWkt2RoundTrip(3395); + } + + @Test + public void testWkt2RoundTrip_LambertConformalConic2SP_2154() throws FactoryException { + assertWkt2RoundTrip(2154); + } + + @Test + public void testWkt2RoundTrip_LambertAzimuthalEqualArea_Spherical_2163() throws FactoryException { + assertWkt2RoundTrip(2163); + } + + @Test + public void testWkt2RoundTrip_AlbersEqualArea_5070() throws FactoryException { + assertWkt2RoundTrip(5070); + } + + @Test + public void testWkt2RoundTrip_PolarStereographicB_3031() throws FactoryException { + assertWkt2RoundTrip(3031); + } + + @Test + public void testWkt2RoundTrip_LambertAzimuthalEqualArea_3035() throws FactoryException { + assertWkt2RoundTrip(3035); + } + + @Test + public void testWkt2RoundTrip_Mercator1SP_Spherical_3785() throws FactoryException { + assertWkt2RoundTrip(3785); + } + + @Test + public void testWkt2RoundTrip_EquidistantCylindrical_4087() throws FactoryException { + assertWkt2RoundTrip(4087); + } + + @Test + public void testWkt2RoundTrip_TransverseMercator_OSGB_27700() throws FactoryException { + assertWkt2RoundTrip(27700); + } + + @Test + public void testWkt2RoundTrip_AlbersEqualArea_Australian_3577() throws FactoryException { + assertWkt2RoundTrip(3577); + } + + @Test + public void testWkt2RoundTrip_LambertConformalConic2SP_Vicgrid_3111() throws FactoryException { + assertWkt2RoundTrip(3111); + } + + @Test + public void testWkt2RoundTrip_PolarStereographicB_NSIDC_3413() throws FactoryException { + assertWkt2RoundTrip(3413); + } + + @Test + public void testWkt2RoundTrip_PolarStereographicA_32661() throws FactoryException { + assertWkt2RoundTrip(32661); + } + + @Test + public void testWkt2RoundTrip_LambertAzimuthalEqualArea_EASE_6931() throws FactoryException { + assertWkt2RoundTrip(6931); + } + + @Test + public void testWkt2RoundTrip_ObliqueStereographic_28992() throws FactoryException { + assertWkt2RoundTrip(28992); + } + + // --------------------------------------------------------------------------- + // PROJJSON import failures — spherical datums not parseable after round-trip + // --------------------------------------------------------------------------- + + @Test + public void testProjJsonRoundTrip_LambertAzimuthalEqualArea_Spherical_2163_importFails() + throws FactoryException { + assertProjJsonImportFails(2163); + } + + @Test + public void testProjJsonRoundTrip_Mercator1SP_Spherical_3785() throws FactoryException { + assertProjJsonRoundTrip(3785); + } + + // --------------------------------------------------------------------------- + // Export failures — projection types not supported by proj4sedona + // --------------------------------------------------------------------------- + + @Test + public void testSetCrs_LambertCylindricalEqualArea_6933() throws FactoryException { + // EPSG:6933 setCrs works; WKT1 may not contain AUTHORITY so we just verify it's parseable + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D result = RasterEditors.setCrs(baseRaster, "EPSG:6933"); + String wkt1 = RasterAccessors.crs(result, "wkt1"); + assertNotNull("EPSG:6933 should produce valid WKT1", wkt1); + assertTrue("WKT1 should contain PROJCS", wkt1.contains("PROJCS")); + } + + @Test + public void testExportFails_Krovak_2065() throws FactoryException { + assertExportFails(2065); + } + + @Test + public void testExportFails_HotineObliqueMercator_2056() throws FactoryException { + assertExportFails(2056); + } + + // --------------------------------------------------------------------------- + // Helper methods + // --------------------------------------------------------------------------- + + /** + * Assert a full PROJ format round trip: EPSG → RS_CRS("proj") → RS_SetCRS → RS_CRS("proj") → + * RS_SetCRS → RS_CRS("proj"). The first export from EPSG may carry extra metadata, so we verify + * idempotency: the second and third exports (from PROJ string input) must be identical. + */ + private void assertProjRoundTrip(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + // First export from EPSG + String export1 = RasterAccessors.crs(raster1, "proj"); + assertNotNull("EPSG:" + epsg + " export to PROJ should not be null", export1); + + // Re-import from PROJ string and re-export + GridCoverage2D raster2 = RasterEditors.setCrs(baseRaster, export1); + String export2 = RasterAccessors.crs(raster2, "proj"); + assertNotNull("EPSG:" + epsg + " second export to PROJ should not be null", export2); + + // Third round-trip to verify idempotency + GridCoverage2D raster3 = RasterEditors.setCrs(baseRaster, export2); + String export3 = RasterAccessors.crs(raster3, "proj"); + assertNotNull("EPSG:" + epsg + " third export to PROJ should not be null", export3); + + assertEquals( + "EPSG:" + epsg + " PROJ string should be stable after round trip (export2 == export3)", + export2, + export3); + } + + /** + * Assert a full PROJJSON format round trip: EPSG → RS_CRS("projjson") → RS_SetCRS → + * RS_CRS("projjson") → RS_SetCRS → RS_CRS("projjson"). The first export from EPSG may carry extra + * metadata (e.g., datum names), so we verify idempotency: the second and third exports (from + * PROJJSON string input) must be identical. + */ + private void assertProjJsonRoundTrip(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + // First export from EPSG + String export1 = RasterAccessors.crs(raster1, "projjson"); + assertNotNull("EPSG:" + epsg + " export to PROJJSON should not be null", export1); + + // Re-import from PROJJSON string and re-export + GridCoverage2D raster2 = RasterEditors.setCrs(baseRaster, export1); + String export2 = RasterAccessors.crs(raster2, "projjson"); + assertNotNull("EPSG:" + epsg + " second export to PROJJSON should not be null", export2); + + // Third round-trip to verify idempotency + GridCoverage2D raster3 = RasterEditors.setCrs(baseRaster, export2); + String export3 = RasterAccessors.crs(raster3, "projjson"); + assertNotNull("EPSG:" + epsg + " third export to PROJJSON should not be null", export3); + + assertEquals( + "EPSG:" + epsg + " PROJJSON string should be stable after round trip (export2 == export3)", + export2, + export3); + } + + /** + * Assert a full WKT2 format round trip: EPSG → RS_CRS("wkt2") → RS_SetCRS → RS_CRS("wkt2") → + * RS_SetCRS → RS_CRS("wkt2"). The first export from EPSG may carry extra metadata, so we verify + * idempotency: the second and third exports (from WKT2 string input) must be identical. + */ + private void assertWkt2RoundTrip(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + // First export from EPSG + String export1 = RasterAccessors.crs(raster1, "wkt2"); + assertNotNull("EPSG:" + epsg + " export to WKT2 should not be null", export1); + + // Re-import from WKT2 string and re-export + GridCoverage2D raster2 = RasterEditors.setCrs(baseRaster, export1); + String export2 = RasterAccessors.crs(raster2, "wkt2"); + assertNotNull("EPSG:" + epsg + " second export to WKT2 should not be null", export2); + + // Third round-trip to verify idempotency + GridCoverage2D raster3 = RasterEditors.setCrs(baseRaster, export2); + String export3 = RasterAccessors.crs(raster3, "wkt2"); + assertNotNull("EPSG:" + epsg + " third export to WKT2 should not be null", export3); + + assertEquals( + "EPSG:" + epsg + " WKT2 string should be stable after round trip (export2 == export3)", + export2, + export3); + } + + /** + * Assert that PROJJSON export succeeds but re-import fails (spherical datum CRS that proj4sedona + * can export but GeoTools cannot re-parse). + */ + private void assertProjJsonImportFails(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + // Export should succeed + String exported = RasterAccessors.crs(raster1, "projjson"); + assertNotNull("EPSG:" + epsg + " export to PROJJSON should succeed", exported); + + // Re-import should fail + Exception thrown = + assertThrows( + "EPSG:" + epsg + " PROJJSON re-import should fail for spherical datum", + IllegalArgumentException.class, + () -> RasterEditors.setCrs(baseRaster, exported)); + assertTrue( + "Error message should mention CRS parsing", + thrown.getMessage().contains("Cannot parse CRS string")); + } + + /** + * Assert that RS_CRS export fails for projection types not supported by proj4sedona. Tests both + * "proj" and "projjson" formats. + */ + private void assertExportFails(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + assertThrows( + "EPSG:" + epsg + " export to PROJ should fail", + Exception.class, + () -> RasterAccessors.crs(raster1, "proj")); + + assertThrows( + "EPSG:" + epsg + " export to PROJJSON should fail", + Exception.class, + () -> RasterAccessors.crs(raster1, "projjson")); + } + + /** + * Assert a full WKT1 format round trip: EPSG → RS_CRS("wkt1") → RS_SetCRS → RS_CRS("wkt1"). WKT1 + * includes AUTHORITY["EPSG","xxxx"] so SRID is always preserved. + */ + private void assertWkt1RoundTrip(int epsg) throws FactoryException { + GridCoverage2D baseRaster = RasterConstructors.makeEmptyRaster(1, 4, 4, 0, 0, 1); + GridCoverage2D raster1 = RasterEditors.setCrs(baseRaster, "EPSG:" + epsg); + + // Export to WKT1 + String exported = RasterAccessors.crs(raster1, "wkt1"); + assertNotNull("EPSG:" + epsg + " export to WKT1 should not be null", exported); + + // Verify AUTHORITY clause is present with correct EPSG code + String topAuthority = extractTopLevelAuthority(exported); + assertEquals( + "EPSG:" + epsg + " WKT1 should contain top-level AUTHORITY", + String.valueOf(epsg), + topAuthority); + + String projNameBefore = extractWkt1ProjectionName(exported); + + // Re-import and re-export + GridCoverage2D raster2 = RasterEditors.setCrs(baseRaster, exported); + String reExported = RasterAccessors.crs(raster2, "wkt1"); + assertNotNull("EPSG:" + epsg + " re-export to WKT1 should not be null", reExported); + + String projNameAfter = extractWkt1ProjectionName(reExported); + assertEquals( + "EPSG:" + epsg + " projection name should be stable after WKT1 round trip", + projNameBefore, + projNameAfter); + + // WKT1 always preserves SRID via AUTHORITY clause + int sridAfter = RasterAccessors.srid(raster2); + assertEquals( + "EPSG:" + epsg + " SRID should be preserved after WKT1 round trip", epsg, sridAfter); + } + + // --------------------------------------------------------------------------- + // Extraction helpers + // --------------------------------------------------------------------------- + + /** + * Extract PROJECTION name from WKT1, or "Geographic" for GEOGCS without PROJECTION. Handles both + * PROJECTION["name"] and PROJECTION["name", AUTHORITY[...]]. + */ + private String extractWkt1ProjectionName(String wkt1) { + Matcher m = WKT1_PROJECTION_PATTERN.matcher(wkt1); + if (m.find()) { + return m.group(1); + } + // Geographic CRS has no PROJECTION element + if (wkt1.startsWith("GEOGCS[")) { + return "Geographic"; + } + fail( + "WKT1 should contain PROJECTION or be GEOGCS: " + + wkt1.substring(0, Math.min(80, wkt1.length()))); + return null; + } + + /** + * Extract the top-level AUTHORITY EPSG code from WKT1. The top-level AUTHORITY is the last one in + * the string (at the outermost nesting level). + */ + private String extractTopLevelAuthority(String wkt1) { + // Find the last AUTHORITY["EPSG","xxxx"] — that's the top-level one + Matcher m = WKT1_AUTHORITY_PATTERN.matcher(wkt1); + String lastCode = null; + while (m.find()) { + lastCode = m.group(1); + } + assertNotNull("WKT1 should contain AUTHORITY[\"EPSG\",\"xxxx\"]", lastCode); + return lastCode; + } +} diff --git a/common/src/test/java/org/apache/sedona/common/raster/RasterAccessorsTest.java b/common/src/test/java/org/apache/sedona/common/raster/RasterAccessorsTest.java index 347309cf00..3565cac554 100644 --- a/common/src/test/java/org/apache/sedona/common/raster/RasterAccessorsTest.java +++ b/common/src/test/java/org/apache/sedona/common/raster/RasterAccessorsTest.java @@ -470,4 +470,72 @@ public class RasterAccessorsTest extends RasterTestBase { assertEquals(heightInPixel, metadata[11], 1e-9); assertEquals(12, metadata.length); } + + @Test + public void testCrsDefaultFormat() throws FactoryException { + // multiBandRaster has WGS84 CRS + String crs = RasterAccessors.crs(multiBandRaster); + assertNotNull(crs); + // Default format is PROJJSON - should be valid JSON containing CRS info + assertTrue(crs.contains("\"type\"")); + assertTrue(crs.contains("WGS 84") || crs.contains("WGS84")); + } + + @Test + public void testCrsWkt1Format() throws FactoryException { + String crs = RasterAccessors.crs(multiBandRaster, "wkt1"); + assertNotNull(crs); + assertTrue(crs.contains("GEOGCS")); + } + + @Test + public void testCrsWkt2Format() throws FactoryException { + String crs = RasterAccessors.crs(multiBandRaster, "wkt2"); + assertNotNull(crs); + // WKT2 uses GEOGCRS or GEODCRS + assertTrue(crs.contains("GEOGCRS") || crs.contains("GEODCRS") || crs.contains("GeographicCRS")); + } + + @Test + public void testCrsProjFormat() throws FactoryException { + String crs = RasterAccessors.crs(multiBandRaster, "proj"); + assertNotNull(crs); + assertTrue(crs.contains("+proj=")); + } + + @Test + public void testCrsNullForNoCrs() throws FactoryException { + // oneBandRaster has no CRS (SRID=0) + String crs = RasterAccessors.crs(oneBandRaster); + assertNull(crs); + } + + @Test + public void testCrsWithSetCrsRoundTrip() throws FactoryException { + // Set a CRS using a PROJ string, then read it back in various formats + String proj = "+proj=utm +zone=10 +datum=WGS84 +units=m +no_defs"; + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + GridCoverage2D withCrs = RasterEditors.setCrs(raster, proj); + + // Should be able to retrieve CRS in all formats + String projjson = RasterAccessors.crs(withCrs, "projjson"); + assertNotNull(projjson); + assertTrue(projjson.contains("\"type\"")); + + String wkt1 = RasterAccessors.crs(withCrs, "wkt1"); + assertNotNull(wkt1); + assertTrue(wkt1.contains("PROJCS") || wkt1.contains("GEOGCS")); + + String wkt2 = RasterAccessors.crs(withCrs, "wkt2"); + assertNotNull(wkt2); + + String projStr = RasterAccessors.crs(withCrs, "proj"); + assertNotNull(projStr); + assertTrue(projStr.contains("+proj=")); + } + + @Test(expected = IllegalArgumentException.class) + public void testCrsInvalidFormat() throws FactoryException { + RasterAccessors.crs(multiBandRaster, "invalid_format"); + } } diff --git a/common/src/test/java/org/apache/sedona/common/raster/RasterEditorsTest.java b/common/src/test/java/org/apache/sedona/common/raster/RasterEditorsTest.java index 8022ed734c..02f4af1cf0 100644 --- a/common/src/test/java/org/apache/sedona/common/raster/RasterEditorsTest.java +++ b/common/src/test/java/org/apache/sedona/common/raster/RasterEditorsTest.java @@ -4210,4 +4210,150 @@ public class RasterEditorsTest extends RasterTestBase { } } } + + @Test + public void testSetCrsWithEpsgCode() throws FactoryException { + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + assertEquals(0, RasterAccessors.srid(raster)); + + GridCoverage2D result = RasterEditors.setCrs(raster, "EPSG:4326"); + assertEquals(4326, RasterAccessors.srid(result)); + } + + @Test + public void testSetCrsWithWkt1() throws FactoryException { + String wkt1 = + "GEOGCS[\"WGS 84\"," + + "DATUM[\"WGS_1984\"," + + "SPHEROID[\"WGS 84\",6378137,298.257223563]]," + + "PRIMEM[\"Greenwich\",0]," + + "UNIT[\"degree\",0.0174532925199433]," + + "AUTHORITY[\"EPSG\",\"4326\"]]"; + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + GridCoverage2D result = RasterEditors.setCrs(raster, wkt1); + assertEquals(4326, RasterAccessors.srid(result)); + } + + @Test + public void testSetCrsWithProjString() throws FactoryException { + String proj = "+proj=longlat +datum=WGS84 +no_defs"; + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + GridCoverage2D result = RasterEditors.setCrs(raster, proj); + assertEquals(4326, RasterAccessors.srid(result)); + } + + @Test + public void testSetCrsWithProjJson() throws FactoryException { + // PROJJSON for EPSG:3857 + String projjson = + "{\"$schema\":\"https://proj.org/schemas/v0.7/projjson.schema.json\","; + + "\"type\":\"ProjectedCRS\"," + + "\"name\":\"WGS 84 / Pseudo-Mercator\"," + + "\"base_crs\":{\"name\":\"WGS 84\"," + + "\"datum\":{\"type\":\"GeodeticReferenceFrame\"," + + "\"name\":\"World Geodetic System 1984\"," + + "\"ellipsoid\":{\"name\":\"WGS 84\"," + + "\"semi_major_axis\":6378137," + + "\"inverse_flattening\":298.257223563}}," + + "\"coordinate_system\":{\"subtype\":\"ellipsoidal\"," + + "\"axis\":[{\"name\":\"Geodetic latitude\",\"abbreviation\":\"Lat\"," + + "\"direction\":\"north\",\"unit\":\"degree\"}," + + "{\"name\":\"Geodetic longitude\",\"abbreviation\":\"Lon\"," + + "\"direction\":\"east\",\"unit\":\"degree\"}]}}," + + "\"conversion\":{\"name\":\"Popular Visualisation Pseudo-Mercator\"," + + "\"method\":{\"name\":\"Popular Visualisation Pseudo Mercator\"," + + "\"id\":{\"authority\":\"EPSG\",\"code\":1024}}," + + "\"parameters\":[{\"name\":\"Latitude of natural origin\",\"value\":0," + + "\"unit\":\"degree\",\"id\":{\"authority\":\"EPSG\",\"code\":8801}}," + + "{\"name\":\"Longitude of natural origin\",\"value\":0," + + "\"unit\":\"degree\",\"id\":{\"authority\":\"EPSG\",\"code\":8802}}," + + "{\"name\":\"False easting\",\"value\":0," + + "\"unit\":\"metre\",\"id\":{\"authority\":\"EPSG\",\"code\":8806}}," + + "{\"name\":\"False northing\",\"value\":0," + + "\"unit\":\"metre\",\"id\":{\"authority\":\"EPSG\",\"code\":8807}}]}," + + "\"coordinate_system\":{\"subtype\":\"Cartesian\"," + + "\"axis\":[{\"name\":\"Easting\",\"abbreviation\":\"X\"," + + "\"direction\":\"east\",\"unit\":\"metre\"}," + + "{\"name\":\"Northing\",\"abbreviation\":\"Y\"," + + "\"direction\":\"north\",\"unit\":\"metre\"}]}," + + "\"id\":{\"authority\":\"EPSG\",\"code\":3857}}"; + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + GridCoverage2D result = RasterEditors.setCrs(raster, projjson); + assertEquals(3857, RasterAccessors.srid(result)); + } + + @Test + public void testSetCrsWithCustomProj() throws FactoryException { + // Custom Lambert Conformal Conic - no EPSG code + String proj = + "+proj=lcc +lat_1=25 +lat_2=60 +lat_0=42.5 +lon_0=-100 " + + "+x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"; + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + GridCoverage2D result = RasterEditors.setCrs(raster, proj); + // Custom CRS has no EPSG code, SRID should be 0 + assertEquals(0, RasterAccessors.srid(result)); + // But the CRS should be valid and contain the projection info + String crsWkt = RasterAccessors.crs(result, "wkt1"); + Assert.assertTrue(crsWkt.contains("Lambert_Conformal_Conic")); + } + + @Test(expected = IllegalArgumentException.class) + public void testSetCrsWithInvalidString() throws FactoryException { + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + RasterEditors.setCrs(raster, "NOT_A_VALID_CRS"); + } + + /** + * Comprehensive test: verify that RS_SetCRS works with every projection type that proj4sedona + * supports. Each projection short code is tested with appropriate parameters. proj4sedona outputs + * WKT1 with projection names that may differ from GeoTools conventions (e.g. "Azimuthal + * Equidistant" vs "Azimuthal_Equidistant"), and may include parameters not expected by GeoTools + * (e.g. standard_parallel_1 for Transverse Mercator). The normalization and parameter-stripping + * logic in parseCrsString handles both cases. + */ + @Test + public void testSetCrsWithRepresentativeProj4SedonaProjections() throws FactoryException { + GridCoverage2D raster = RasterConstructors.makeEmptyRaster(1, 10, 10, 0, 0, 1); + + // A representative set of projection short codes supported by proj4sedona, each with + // appropriate parameters. + // Format: {shortCode, projString} + String[][] projConfigs = { + { + "aea", + "+proj=aea +lat_0=0 +lon_0=0 +lat_1=30 +lat_2=60 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs" + }, + {"aeqd", "+proj=aeqd +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"cea", "+proj=cea +lat_ts=30 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"eqc", "+proj=eqc +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + { + "eqdc", + "+proj=eqdc +lat_0=0 +lon_0=0 +lat_1=30 +lat_2=60 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs" + }, + { + "etmerc", "+proj=etmerc +lat_0=0 +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs" + }, + {"laea", "+proj=laea +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + { + "lcc", + "+proj=lcc +lat_0=0 +lon_0=0 +lat_1=30 +lat_2=60 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs" + }, + {"merc", "+proj=merc +lat_ts=0 +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"moll", "+proj=moll +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"robin", "+proj=robin +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"sinu", "+proj=sinu +lon_0=0 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"stere", "+proj=stere +lat_0=90 +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"tmerc", "+proj=tmerc +lat_0=0 +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs"}, + {"utm", "+proj=utm +zone=17 +datum=WGS84 +units=m +no_defs"}, + }; + + for (String[] config : projConfigs) { + String code = config[0]; + String projStr = config[1]; + GridCoverage2D result = RasterEditors.setCrs(raster, projStr); + String wkt1 = RasterAccessors.crs(result, "wkt1"); + Assert.assertNotNull("setCrs should produce a valid CRS for +proj=" + code, wkt1); + Assert.assertTrue("WKT1 should contain PROJCS for +proj=" + code, wkt1.contains("PROJCS")); + } + } } diff --git a/docs/api/sql/Raster-Functions.md b/docs/api/sql/Raster-Functions.md index b7ea6e46b9..bddf209a21 100644 --- a/docs/api/sql/Raster-Functions.md +++ b/docs/api/sql/Raster-Functions.md @@ -122,10 +122,12 @@ These functions perform operations on raster objects. | [RS_SetBandNoDataValue](Raster-Operators/RS_SetBandNoDataValue.md) | This sets the no data value for a specified band in the raster. If the band index is not provided, band 1 is assumed by default. Passing a `null` value for `noDataValue` will remove the no data val... | v1.5.0 | | [RS_SetGeoReference](Raster-Operators/RS_SetGeoReference.md) | Sets the Georeference information of an object in a single call. Accepts inputs in `GDAL` and `ESRI` format. Default format is `GDAL`. If all 6 parameters are not provided then will return null. | v1.5.0 | | [RS_SetPixelType](Raster-Operators/RS_SetPixelType.md) | Returns a modified raster with the desired pixel data type. | v1.6.0 | +| [RS_SetCRS](Raster-Operators/RS_SetCRS.md) | Sets the coordinate reference system (CRS) of the raster using a CRS definition string. Accepts EPSG codes, WKT1, WKT2, PROJ strings, and PROJJSON. | v1.9.0 | | [RS_SetSRID](Raster-Operators/RS_SetSRID.md) | Sets the spatial reference system identifier (SRID) of the raster geometry. | v1.4.1 | | [RS_SetValue](Raster-Operators/RS_SetValue.md) | Returns a raster by replacing the value of pixel specified by `colX` and `rowY`. | v1.5.0 | | [RS_SetValues](Raster-Operators/RS_SetValues.md) | Returns a raster by replacing the values of pixels in a specified rectangular region. The top left corner of the region is defined by the `colX` and `rowY` coordinates. The `width` and `height` par... | v1.5.0 | | [RS_SRID](Raster-Operators/RS_SRID.md) | Returns the spatial reference system identifier (SRID) of the raster geometry. | v1.4.1 | +| [RS_CRS](Raster-Operators/RS_CRS.md) | Returns the coordinate reference system (CRS) of the raster as a string in the specified format (projjson, wkt2, wkt1, proj). Defaults to PROJJSON. | v1.9.0 | | [RS_Union](Raster-Operators/RS_Union.md) | Returns a combined multi-band raster from 2 or more input Rasters. The order of bands in the resultant raster will be in the order of the input rasters. For example if `RS_Union` is called on two 2... | v1.6.0 | | [RS_Value](Raster-Operators/RS_Value.md) | Returns the value at the given point in the raster. If no band number is specified it defaults to 1. | v1.4.0 | | [RS_Values](Raster-Operators/RS_Values.md) | Returns the values at the given points or grid coordinates in the raster. If no band number is specified it defaults to 1. | v1.4.0 | diff --git a/docs/api/sql/Raster-Operators/RS_CRS.md b/docs/api/sql/Raster-Operators/RS_CRS.md new file mode 100644 index 0000000000..7e719f2c88 --- /dev/null +++ b/docs/api/sql/Raster-Operators/RS_CRS.md @@ -0,0 +1,101 @@ +<!-- + Licensed to the Apache Software Foundation (ASF) under one + or more contributor license agreements. See the NOTICE file + distributed with this work for additional information + regarding copyright ownership. The ASF licenses this file + to you under the Apache License, Version 2.0 (the + "License"); you may not use this file except in compliance + with the License. You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, + software distributed under the License is distributed on an + "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + KIND, either express or implied. See the License for the + specific language governing permissions and limitations + under the License. + --> + +# RS_CRS + +Introduction: Returns the coordinate reference system (CRS) of a raster as a string in the specified format. If no format is specified, the CRS is returned in PROJJSON format. Returns null if the raster has no CRS defined. + +Format: + +``` +RS_CRS (raster: Raster) +``` + +``` +RS_CRS (raster: Raster, format: String) +``` + +Since: `v1.9.0` + +## Supported output formats + +| Format | Description | +| :--- | :--- | +| `'projjson'` | PROJJSON format (default). Modern, machine-readable JSON representation. | +| `'wkt2'` | Well-Known Text 2 (ISO 19162). Modern standard CRS representation. | +| `'wkt1'` | Well-Known Text 1 (OGC 01-009). Legacy format, widely supported. | +| `'proj'` | PROJ string format. Compact, human-readable representation. | + +## SQL Examples + +Getting CRS in default PROJJSON format: + +```sql +SELECT RS_CRS(raster) FROM raster_table +``` + +Output: + +```json +{ + "$schema": "https://proj.org/schemas/v0.7/projjson.schema.json";, + "type": "GeographicCRS", + "name": "WGS 84", + ... +} +``` + +Getting CRS in WKT1 format: + +```sql +SELECT RS_CRS(raster, 'wkt1') FROM raster_table +``` + +Output: + +``` +GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]] +``` + +Getting CRS in PROJ string format: + +```sql +SELECT RS_CRS(raster, 'proj') FROM raster_table +``` + +Output: + +``` ++proj=longlat +datum=WGS84 +no_defs +type=crs +``` + +Getting CRS in WKT2 format: + +```sql +SELECT RS_CRS(raster, 'wkt2') FROM raster_table +``` + +## Limitations + +The `wkt2`, `proj`, and `projjson` output formats are generated by [proj4sedona](https://github.com/jiayuasu/proj4sedona) from the raster's internal WKT1 CRS. This conversion may cause the following limitations: + +- **Unsupported projection types**: Some projection types (e.g., Krovak, Hotine Oblique Mercator) cannot be exported to `wkt2`, `proj`, or `projjson` formats and will throw an error. Use `'wkt1'` format for these. + +!!!note + `RS_CRS` returns null only when the raster has no CRS defined. Note that `RS_SRID` may return `0` either when no CRS is defined or when a custom (non-EPSG) CRS has been set via `RS_SetCRS`, so `RS_SRID = 0` does not always mean "no CRS". To test for a missing CRS, use `RS_CRS(raster) IS NULL`. The `wkt1` format always produces a lossless representation of the internally stored CRS. diff --git a/docs/api/sql/Raster-Operators/RS_SRID.md b/docs/api/sql/Raster-Operators/RS_SRID.md index 7820828e2e..8d79b94b0c 100644 --- a/docs/api/sql/Raster-Operators/RS_SRID.md +++ b/docs/api/sql/Raster-Operators/RS_SRID.md @@ -19,7 +19,7 @@ # RS_SRID -Introduction: Returns the spatial reference system identifier (SRID) of the raster geometry. +Introduction: Returns the spatial reference system identifier (SRID) of the raster geometry. Returns 0 if the raster has no CRS defined or if the CRS is a custom (non-EPSG) coordinate reference system. To retrieve the full CRS definition for custom CRS, use [RS_CRS](RS_CRS.md). Format: `RS_SRID (raster: Raster)` diff --git a/docs/api/sql/Raster-Operators/RS_SetCRS.md b/docs/api/sql/Raster-Operators/RS_SetCRS.md new file mode 100644 index 0000000000..98d219e02e --- /dev/null +++ b/docs/api/sql/Raster-Operators/RS_SetCRS.md @@ -0,0 +1,68 @@ +<!-- + Licensed to the Apache Software Foundation (ASF) under one + or more contributor license agreements. See the NOTICE file + distributed with this work for additional information + regarding copyright ownership. The ASF licenses this file + to you under the Apache License, Version 2.0 (the + "License"); you may not use this file except in compliance + with the License. You may obtain a copy of the License at + + http://www.apache.org/licenses/LICENSE-2.0 + + Unless required by applicable law or agreed to in writing, + software distributed under the License is distributed on an + "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY + KIND, either express or implied. See the License for the + specific language governing permissions and limitations + under the License. + --> + +# RS_SetCRS + +Introduction: Sets the coordinate reference system (CRS) of a raster using a CRS definition string. Unlike `RS_SetSRID` which only accepts integer EPSG codes, `RS_SetCRS` accepts CRS definitions in multiple formats including EPSG codes, WKT1, WKT2, PROJ strings, and PROJJSON. This function does not reproject/transform the raster data — it only sets the CRS metadata. + +Format: `RS_SetCRS (raster: Raster, crsString: String)` + +Since: `v1.9.0` + +## Supported CRS formats + +| Format | Example | +| :--- | :--- | +| EPSG code | `'EPSG:4326'` | +| WKT1 | `'GEOGCS["WGS 84", DATUM["WGS_1984", ...], ...]'` | +| WKT2 | `'GEOGCRS["WGS 84", DATUM["World Geodetic System 1984", ...], ...]'` | +| PROJ string | `'+proj=longlat +datum=WGS84 +no_defs'` | +| PROJJSON | `'{"type": "GeographicCRS", "name": "WGS 84", ...}'` | + +## SQL Examples + +Setting CRS with an EPSG code: + +```sql +SELECT RS_SetCRS(raster, 'EPSG:4326') FROM raster_table +``` + +Setting CRS with a PROJ string (useful for custom projections): + +```sql +SELECT RS_SetCRS(raster, '+proj=lcc +lat_1=25 +lat_2=60 +lat_0=42.5 +lon_0=-100 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs') +FROM raster_table +``` + +Setting CRS with a WKT1 string: + +```sql +SELECT RS_SetCRS(raster, 'GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563]],PRIMEM["Greenwich",0],UNIT["degree",0.0174532925199433]]') +FROM raster_table +``` + +## Limitations + +Internally, Sedona stores raster CRS in WKT1 format (via GeoTools). When you provide a CRS in WKT2, PROJ, or PROJJSON format, it is converted to WKT1 using [proj4sedona](https://github.com/jiayuasu/proj4sedona). This conversion may cause the following limitations: + +- **SRID not preserved for projected CRS**: When importing PROJ or PROJJSON strings, the EPSG SRID is often lost for projected coordinate systems. Only geographic CRS (e.g., EPSG:4326), Web Mercator (EPSG:3857), and UTM zones reliably preserve their SRID. Use `RS_SetSRID` after `RS_SetCRS` if you need to set a specific SRID. +- **Unsupported projection types**: Some projection types (e.g., Krovak, Hotine Oblique Mercator) are not supported by proj4sedona and will fail for WKT2, PROJ, and PROJJSON formats. Use `'EPSG:xxxx'` or WKT1 for these. + +!!!note + For the most reliable results, use `'EPSG:xxxx'` format when your CRS has a known EPSG code. WKT1 input is also lossless since it is stored natively. WKT2, PROJ, and PROJJSON inputs undergo conversion and may experience the limitations above. diff --git a/pom.xml b/pom.xml index 97e49947f3..9b72df9f7e 100644 --- a/pom.xml +++ b/pom.xml @@ -96,7 +96,7 @@ <scala-collection-compat.version>2.5.0</scala-collection-compat.version> <geoglib.version>1.52</geoglib.version> <caffeine.version>2.9.2</caffeine.version> - <proj4sedona.version>0.0.3</proj4sedona.version> + <proj4sedona.version>0.0.8</proj4sedona.version> <geotools.scope>provided</geotools.scope> <!-- Because it's not in Maven central, make it provided by default --> diff --git a/spark/common/src/main/scala/org/apache/sedona/sql/UDF/Catalog.scala b/spark/common/src/main/scala/org/apache/sedona/sql/UDF/Catalog.scala index a2d483398a..fd10623b18 100644 --- a/spark/common/src/main/scala/org/apache/sedona/sql/UDF/Catalog.scala +++ b/spark/common/src/main/scala/org/apache/sedona/sql/UDF/Catalog.scala @@ -294,12 +294,14 @@ object Catalog extends AbstractCatalog with Logging { function[RS_NumBands](), function[RS_Metadata](), function[RS_SetSRID](), + function[RS_SetCRS](), function[RS_SetGeoReference](), function[RS_SetBandNoDataValue](), function[RS_SetPixelType](), function[RS_SetValues](), function[RS_SetValue](), function[RS_SRID](), + function[RS_CRS](), function[RS_Value](1), function[RS_Values](1), function[RS_Intersects](), diff --git a/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterAccessors.scala b/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterAccessors.scala index ed3cc5327f..7d75bc137a 100644 --- a/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterAccessors.scala +++ b/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterAccessors.scala @@ -38,6 +38,15 @@ private[apache] case class RS_SRID(inputExpressions: Seq[Expression]) } } +private[apache] case class RS_CRS(inputExpressions: Seq[Expression]) + extends InferredExpression( + inferrableFunction1(RasterAccessors.crs), + inferrableFunction2(RasterAccessors.crs)) { + protected def withNewChildrenInternal(newChildren: IndexedSeq[Expression]) = { + copy(inputExpressions = newChildren) + } +} + private[apache] case class RS_Width(inputExpressions: Seq[Expression]) extends InferredExpression(RasterAccessors.getWidth _) { protected def withNewChildrenInternal(newChildren: IndexedSeq[Expression]) = { diff --git a/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterEditors.scala b/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterEditors.scala index 96855548ee..2394a49375 100644 --- a/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterEditors.scala +++ b/spark/common/src/main/scala/org/apache/spark/sql/sedona_sql/expressions/raster/RasterEditors.scala @@ -31,6 +31,13 @@ private[apache] case class RS_SetSRID(inputExpressions: Seq[Expression]) } } +private[apache] case class RS_SetCRS(inputExpressions: Seq[Expression]) + extends InferredExpression(RasterEditors.setCrs _) { + protected def withNewChildrenInternal(newChildren: IndexedSeq[Expression]) = { + copy(inputExpressions = newChildren) + } +} + private[apache] case class RS_SetGeoReference(inputExpressions: Seq[Expression]) extends InferredExpression( inferrableFunction2(RasterEditors.setGeoReference), diff --git a/spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala b/spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala index 1b6490c234..47840f500b 100644 --- a/spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala +++ b/spark/common/src/test/scala/org/apache/sedona/sql/rasteralgebraTest.scala @@ -472,6 +472,77 @@ class rasteralgebraTest extends TestBaseScala with BeforeAndAfter with GivenWhen assert(metadata(9) == metadata_4326(9)) } + it("Passed RS_SetCRS should handle null values") { + val result = sparkSession.sql("select RS_SetCRS(null, null)").first().get(0) + assert(result == null) + } + + it("Passed RS_SetCRS with EPSG code string") { + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + val result = df + .selectExpr("RS_SRID(RS_SetCRS(RS_FromGeoTiff(content), 'EPSG:4326'))") + .first() + .getInt(0) + assert(result == 4326) + } + + it("Passed RS_SetCRS with PROJ string") { + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + val result = df + .selectExpr( + "RS_SRID(RS_SetCRS(RS_FromGeoTiff(content), '+proj=longlat +datum=WGS84 +no_defs'))") + .first() + .getInt(0) + assert(result == 4326) + } + + it("Passed RS_SetCRS with WKT1 string") { + val wkt1 = + "GEOGCS[\"WGS 84\",DATUM[\"WGS_1984\",SPHEROID[\"WGS 84\",6378137,298.257223563]],PRIMEM[\"Greenwich\",0],UNIT[\"degree\",0.0174532925199433]]" + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + // WKT1 without AUTHORITY clause: RS_SRID returns 0, but RS_CRS still works + val srid = + df.selectExpr(s"RS_SRID(RS_SetCRS(RS_FromGeoTiff(content), '${wkt1}'))").first().getInt(0) + assert(srid == 0) + val crs = + df.selectExpr(s"RS_CRS(RS_SetCRS(RS_FromGeoTiff(content), '${wkt1}'), 'wkt1')") + .first() + .getString(0) + assert(crs != null && crs.contains("WGS")) + } + + it("Passed RS_CRS should handle null values") { + val result = sparkSession.sql("select RS_CRS(null)").first().get(0) + assert(result == null) + } + + it("Passed RS_CRS returns PROJJSON by default") { + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + val result = df.selectExpr("RS_CRS(RS_FromGeoTiff(content))").first().getString(0) + assert(result != null) + assert(result.contains("\"type\"")) + } + + it("Passed RS_CRS with wkt1 format") { + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + val result = df.selectExpr("RS_CRS(RS_FromGeoTiff(content), 'wkt1')").first().getString(0) + assert(result != null) + assert(result.contains("PROJCS")) + } + + it("Passed RS_CRS with proj format") { + val df = sparkSession.read.format("binaryFile").load(resourceFolder + "raster/test1.tiff") + val result = df.selectExpr("RS_CRS(RS_FromGeoTiff(content), 'proj')").first().getString(0) + assert(result != null) + assert(result.contains("+proj=")) + } + + it("Passed RS_CRS returns null for raster without CRS") { + val result = + sparkSession.sql("select RS_CRS(RS_MakeEmptyRaster(1, 10, 10, 0, 0, 1))").first().get(0) + assert(result == null) + } + it("Passed RS_SetGeoReference should handle null values") { val result = sparkSession.sql("select RS_SetGeoReference(null, null)").first().get(0) assertNull(result)
