feat: Add geodetic calculation option for Pseudo-Mercator

This commit introduces the useGeodeticCalculations parameter for maps.
When set to true for Pseudo-Mercator projections, it ensures more accurate scale and bounding box computations by accounting for Mercator projection distortions, especially away from the equator.

This change affects map bounds, scale calculation, and WMTS layer scaling. It also includes new unit tests and an example to demonstrate the feature.

Author: fdiazcarsi

core/src/main/java/org/mapfish/print/PseudoMercatorUtils.java

@@ -0,0 +1,33 @@
+package org.mapfish.print;
+
+import org.geotools.api.referencing.crs.CoordinateReferenceSystem;
+
+/**
+ * Utility class for handling WGS 84 with Pseudo-Mercator projection specific calculations.
+ *
+ * @author fdiaz
+ */
+public final class PseudoMercatorUtils {
+
+    private PseudoMercatorUtils() {
+
+    }
+
+    /**
+     * Checks if a given CoordinateReferenceSystem is WGS 84 with Pseudo-Mercator projection.
+     *
+     * @param crs The CoordinateReferenceSystem to check.
+     * @return {@code true} if the CRS is Pseudo-Mercator, {@code false} otherwise.
+     */
+    public static boolean isPseudoMercator(final CoordinateReferenceSystem crs) {
+        String crsNameCode = crs.getName().getCode();
+        String crsId = crs.getIdentifiers().iterator().next().toString().toLowerCase();
+        return (crsNameCode.contains("wgs 84") && (
+                crsNameCode.contains("pseudo-mercator") ||
+                crsNameCode.contains("pseudo mercator") ||
+                crsNameCode.contains("web-mercator") ||
+                crsNameCode.contains("web mercator")
+            )) || "EPSG:3857".equalsIgnoreCase(crsId);
+    }
+
+}

core/src/main/java/org/mapfish/print/attribute/map/BBoxMapBounds.java

@@ -2,13 +2,21 @@
 
 import java.awt.Rectangle;
 import org.geotools.api.referencing.crs.CoordinateReferenceSystem;
+import org.geotools.api.referencing.operation.MathTransform;
+import org.geotools.geometry.jts.JTS;
 import org.geotools.geometry.jts.ReferencedEnvelope;
+import org.geotools.referencing.CRS;
 import org.geotools.referencing.GeodeticCalculator;
 import org.locationtech.jts.geom.Coordinate;
 import org.locationtech.jts.geom.Envelope;
 import org.mapfish.print.FloatingPointUtil;
+import org.mapfish.print.PrintException;
 import org.mapfish.print.map.DistanceUnit;
 import org.mapfish.print.map.Scale;
+import org.mapfish.print.PseudoMercatorUtils;
+
+
+
 
 /**
  * Represent the map bounds with a bounding box.
@@ -23,12 +31,43 @@ public final class BBoxMapBounds extends MapBounds {
    *
    * @param projection the projection these bounds are defined in.
    * @param envelope the bounds
+   * @param useGeodeticCalculations force to use geodetic calculations in PseudoMercator projection
    */
-  public BBoxMapBounds(final CoordinateReferenceSystem projection, final Envelope envelope) {
-    super(projection);
+  public BBoxMapBounds(final CoordinateReferenceSystem projection, final Envelope envelope, final boolean useGeodeticCalculations) {
+    super(projection, useGeodeticCalculations);
     this.bbox = envelope;
   }
 
+  /**
+   * Constructor.
+   *
+   * @param projection the projection these bounds are defined in.
+   * @param envelope the bounds
+   */
+  public BBoxMapBounds(final CoordinateReferenceSystem projection, final Envelope envelope) {
+    this(projection, envelope, false);
+  }
+
+  /**
+   * Constructor.
+   *
+   * @param projection the projection these bounds are defined in.
+   * @param minX min X coordinate for the MapBounds
+   * @param minY min Y coordinate for the MapBounds
+   * @param maxX max X coordinate for the MapBounds
+   * @param maxY max Y coordinate for the MapBounds
+   * @param useGeodeticCalculations force to use geodetic calculations in PseudoMercator projection
+   */
+  public BBoxMapBounds(
+      final CoordinateReferenceSystem projection,
+      final double minX,
+      final double minY,
+      final double maxX,
+      final double maxY,
+      final boolean useGeodeticCalculations) {
+    this(projection, new Envelope(minX, maxX, minY, maxY), useGeodeticCalculations);
+  }
+
   /**
    * Constructor.
    *
@@ -47,6 +86,22 @@ public BBoxMapBounds(
     this(projection, new Envelope(minX, maxX, minY, maxY));
   }
 
+  /**
+   * Create from a bbox.
+   *
+   * @param bbox the bounds.
+   * @param useGeodeticCalculations force to use geodetic calculations in PseudoMercator projection
+   */
+  public BBoxMapBounds(final ReferencedEnvelope bbox, final boolean useGeodeticCalculations) {
+    this(
+        bbox.getCoordinateReferenceSystem(),
+        bbox.getMinX(),
+        bbox.getMinY(),
+        bbox.getMaxX(),
+        bbox.getMaxY(),
+        useGeodeticCalculations);
+  }
+
   /**
    * Create from a bbox.
    *
@@ -80,7 +135,8 @@ public MapBounds adjustedEnvelope(final Rectangle paintArea) {
           centerX - finalDiff,
           this.bbox.getMinY(),
           centerX + finalDiff,
-          this.bbox.getMaxY());
+          this.bbox.getMaxY(),
+          this.useGeodeticCalculations());
     } else {
       double centerY = (this.bbox.getMinY() + this.bbox.getMaxY()) / 2;
       double factor = bboxAspectRatio / paintAreaAspectRatio;
@@ -90,7 +146,8 @@ public MapBounds adjustedEnvelope(final Rectangle paintArea) {
           this.bbox.getMinX(),
           centerY - finalDiff,
           this.bbox.getMaxX(),
-          centerY + finalDiff);
+          centerY + finalDiff,
+          this.useGeodeticCalculations());
     }
   }
 
@@ -107,26 +164,21 @@ public MapBounds adjustBoundsToNearestScale(
         getNearestScale(zoomLevels, tolerance, zoomLevelSnapStrategy, geodetic, paintArea, dpi);
 
     Coordinate center = this.bbox.centre();
-    return new CenterScaleMapBounds(getProjection(), center.x, center.y, newScale);
+    return new CenterScaleMapBounds(getProjection(), center.x, center.y, newScale, this.useGeodeticCalculations());
   }
 
   @Override
   public Scale getScale(final Rectangle paintArea, final double dpi) {
     final ReferencedEnvelope bboxAdjustedToScreen = toReferencedEnvelope(paintArea);
 
-    DistanceUnit projUnit = DistanceUnit.fromProjection(getProjection());
+    CoordinateReferenceSystem crs = getProjection();
+    DistanceUnit projUnit = DistanceUnit.fromProjection(crs);
 
     double geoWidthInInches;
-    if (projUnit == DistanceUnit.DEGREES) {
-      GeodeticCalculator calculator = new GeodeticCalculator(getProjection());
-      final double centerY = bboxAdjustedToScreen.centre().y;
-      calculator.setStartingGeographicPoint(bboxAdjustedToScreen.getMinX(), centerY);
-      calculator.setDestinationGeographicPoint(bboxAdjustedToScreen.getMaxX(), centerY);
-      double geoWidthInEllipsoidUnits = calculator.getOrthodromicDistance();
-      DistanceUnit ellipsoidUnit =
-          DistanceUnit.fromString(calculator.getEllipsoid().getAxisUnit().toString());
-
-      geoWidthInInches = ellipsoidUnit.convertTo(geoWidthInEllipsoidUnits, DistanceUnit.IN);
+    
+    // If it is geodetic/degress OR it is a special case requiring geodetic calculation (PseudoMercator)
+    if (projUnit == DistanceUnit.DEGREES || (this.useGeodeticCalculations() && PseudoMercatorUtils.isPseudoMercator(crs))) {
+        geoWidthInInches = this.computeGeodeticWidthInInches(bboxAdjustedToScreen);    
     } else {
       // (scale * width ) / dpi = geowidith
       geoWidthInInches = projUnit.convertTo(bboxAdjustedToScreen.getWidth(), DistanceUnit.IN);
@@ -135,6 +187,36 @@ public Scale getScale(final Rectangle paintArea, final double dpi) {
     return new Scale(geoWidthInInches * (dpi / paintArea.getWidth()), projUnit, dpi);
   }
 
+  @SuppressWarnings("UseSpecificCatch")
+  private double computeGeodeticWidthInInches(final ReferencedEnvelope bbox) {
+    try {
+        CoordinateReferenceSystem crs = bbox.getCoordinateReferenceSystem();
+        GeodeticCalculator calculator = new GeodeticCalculator(crs); // Use the original CRS for the ellipsoid
+
+        double centerY = bbox.centre().y;
+        Coordinate start = new Coordinate(bbox.getMinX(), centerY);
+        Coordinate end = new Coordinate(bbox.getMaxX(), centerY);
+        
+        if (this.useGeodeticCalculations() && PseudoMercatorUtils.isPseudoMercator(crs)) {
+            // Needs reprojection
+            final MathTransform transform = CRS.findMathTransform(crs, GenericMapAttribute.parseProjection("EPSG:4326", true));
+            start = JTS.transform(start, null, transform);
+            end = JTS.transform(end, null, transform);
+        }
+
+        // --- Common Logic ---
+        calculator.setStartingGeographicPoint(start.x, start.y);
+        calculator.setDestinationGeographicPoint(end.x, end.y);
+        final double orthodromicWidth = calculator.getOrthodromicDistance();
+        final DistanceUnit ellipsoidUnit = DistanceUnit.fromString(calculator.getEllipsoid().getAxisUnit().toString());
+
+        return ellipsoidUnit.convertTo(orthodromicWidth, DistanceUnit.IN);
+
+    } catch (Exception e) {
+        throw new PrintException("Failed to compute geodetic width", e);
+    }
+  }
+
   @Override
   public MapBounds adjustBoundsToRotation(final double rotation) {
     if (FloatingPointUtil.equals(rotation, 0.0)) {
@@ -157,7 +239,7 @@ public MapBounds adjustBoundsToRotation(final double rotation) {
     final double rotatedMinY = this.bbox.getMinY() - heightDifference;
     final double rotatedMaxY = this.bbox.getMaxY() + heightDifference;
 
-    return new BBoxMapBounds(getProjection(), rotatedMinX, rotatedMinY, rotatedMaxX, rotatedMaxY);
+    return new BBoxMapBounds(getProjection(), rotatedMinX, rotatedMinY, rotatedMaxX, rotatedMaxY, this.useGeodeticCalculations());
   }
 
   private double getRotatedWidth(final double rotation) {
@@ -197,7 +279,7 @@ public MapBounds zoomOut(final double factor) {
     double minGeoY = centerY - destHeight / 2.0f;
     double maxGeoY = centerY + destHeight / 2.0f;
 
-    return new BBoxMapBounds(getProjection(), minGeoX, minGeoY, maxGeoX, maxGeoY);
+    return new BBoxMapBounds(getProjection(), minGeoX, minGeoY, maxGeoX, maxGeoY, this.useGeodeticCalculations());
   }
 
   @Override
@@ -241,7 +323,7 @@ public MapBounds expand(final int margin, final Rectangle paintArea) {
     final double minGeoY = centerY - destHeight / 2.0;
     final double maxGeoY = centerY + destHeight / 2.0;
 
-    return new BBoxMapBounds(getProjection(), minGeoX, minGeoY, maxGeoX, maxGeoY);
+    return new BBoxMapBounds(getProjection(), minGeoX, minGeoY, maxGeoX, maxGeoY, this.useGeodeticCalculations());
   }
 
   @Override