Always fallback to default coordinate operations when bounding box tr…

…ansfrom fails

Bounding box transforms are inherently approximate, so we can safely
just fallback to the proj default operation if the user-specified
operation fails when transforming a bounding box.

This fixes cases where the user-specified operation involves a
grid shift and the bounds to transform falls outside of the grid.
In this case proj_trans_bounds fails. By falling back to a proj
default operation and re-trying then we get a valid approximate
transformation of the bounding box.

Fixes map rendering issues when users are using a grid shift operation
to transform layers to map CRS

Fixes #60737
Fixes #60753

src/core/proj/qgscoordinatetransform.cpp

@@ -676,9 +676,27 @@ QgsRectangle QgsCoordinateTransform::transformBoundingBox( const QgsRectangle &r
   proj_errno_reset( projData );
   // proj documentation recommends 21 points for densification
   constexpr int DENSIFY_POINTS = 21;
-  const int projResult = proj_trans_bounds( projContext, projData, ( direction == Qgis::TransformDirection::Forward && !d->mIsReversed ) || ( direction == Qgis::TransformDirection::Reverse && d->mIsReversed ) ? PJ_FWD : PJ_INV,
-                         xMin, yMin, xMax, yMax,
-                         &transXMin, &transYMin, &transXMax, &transYMax, DENSIFY_POINTS );
+  int projResult = proj_trans_bounds( projContext, projData, ( direction == Qgis::TransformDirection::Forward && !d->mIsReversed ) || ( direction == Qgis::TransformDirection::Reverse && d->mIsReversed ) ? PJ_FWD : PJ_INV,
+                                      xMin, yMin, xMax, yMax,
+                                      &transXMin, &transYMin, &transXMax, &transYMax, DENSIFY_POINTS );
+
+  if ( ( projResult != 1
+         || !std::isfinite( transXMin )
+         || !std::isfinite( transXMax )
+         || !std::isfinite( transYMin )
+         || !std::isfinite( transYMax ) )
+       && ( d->mAvailableOpCount > 1 || d->mAvailableOpCount == -1 ) // only use fallbacks if more than one operation is possible -- otherwise we've already tried it and it failed
+     )
+  {
+    // fail #1 -- try with getting proj to auto-pick an appropriate coordinate operation for the points
+    if ( PJ *transform = d->threadLocalFallbackProjData() )
+    {
+      projResult = proj_trans_bounds( projContext, transform, ( direction == Qgis::TransformDirection::Forward && !d->mIsReversed ) || ( direction == Qgis::TransformDirection::Reverse && d->mIsReversed ) ? PJ_FWD : PJ_INV,
+                                      xMin, yMin, xMax, yMax,
+                                      &transXMin, &transYMin, &transXMax, &transYMax, DENSIFY_POINTS );
+    }
+  }
+
   if ( projResult != 1
        || !std::isfinite( transXMin )
        || !std::isfinite( transXMax )