diff --git a/MIGRATION_GUIDE.TXT b/MIGRATION_GUIDE.TXT index c4ad53a63d5e..fd2f57dbe4ff 100644 --- a/MIGRATION_GUIDE.TXT +++ b/MIGRATION_GUIDE.TXT @@ -10,6 +10,19 @@ MIGRATION GUIDE FROM GDAL 3.10 to GDAL 3.11 - If only a specific GDAL Minor version is to be supported, this must now be specified in the find_package call in CMake via a version range specification. +- The following methods + OGRCoordinateTransformation::Transform(size_t nCount, double *x, double *y, + double *z, double *t, int *pabSuccess) and + OGRCoordinateTransformation::TransformWithErrorCodes(size_t nCount, double *x, + double *y, double *z, double *t, int *panErrorCodes) are modified to return + FALSE as soon as at least one point fails to transform (to be consistent with + the other form of Transform() that doesn't take a "t" argument), whereas + previously they would return FALSE only if no transformation was found. + + The GDALTransformerFunc callback and its implementations (GenImgProjTransformer, + RPCTransformer, etc.) are also modified to return FALSE as soon as at least + one point fails to transform. + MIGRATION GUIDE FROM GDAL 3.9 to GDAL 3.10 ------------------------------------------ diff --git a/alg/gdal_alg.h b/alg/gdal_alg.h index e92cd6625313..8bccbb4429d4 100644 --- a/alg/gdal_alg.h +++ b/alg/gdal_alg.h @@ -76,6 +76,22 @@ CPLErr CPL_DLL CPL_STDCALL GDALSieveFilter( * Warp Related. */ +/** + * Callback to transforms points. + * + * @param pTransformerArg return value from a GDALCreateXXXXTransformer() function + * @param bDstToSrc TRUE if transformation is from the destination + * (georeferenced) coordinates to pixel/line or FALSE when transforming + * from pixel/line to georeferenced coordinates. + * @param nPointCount the number of values in the x, y and z arrays. + * @param[in,out] x array containing the X values to be transformed. Must not be NULL. + * @param[in,out] y array containing the Y values to be transformed. Must not be NULL. + * @param[in,out] z array containing the Z values to be transformed. Must not be NULL. + * @param[out] panSuccess array in which a flag indicating success (TRUE) or + * failure (FALSE) of the transformation are placed. Must not be NULL. + * + * @return TRUE if all points have been successfully transformed. + */ typedef int (*GDALTransformerFunc)(void *pTransformerArg, int bDstToSrc, int nPointCount, double *x, double *y, double *z, int *panSuccess); diff --git a/alg/gdal_crs.cpp b/alg/gdal_crs.cpp index 9c134f0b69f5..d58e508f4b4c 100644 --- a/alg/gdal_crs.cpp +++ b/alg/gdal_crs.cpp @@ -422,7 +422,7 @@ void GDALDestroyGCPTransformer(void *pTransformArg) * @param panSuccess array in which a flag indicating success (TRUE) or * failure (FALSE) of the transformation are placed. * - * @return TRUE. + * @return TRUE if all points have been successfully transformed. */ int GDALGCPTransform(void *pTransformArg, int bDstToSrc, int nPointCount, @@ -436,10 +436,12 @@ int GDALGCPTransform(void *pTransformArg, int bDstToSrc, int nPointCount, if (psInfo->bReversed) bDstToSrc = !bDstToSrc; + int bRet = TRUE; for (i = 0; i < nPointCount; i++) { if (x[i] == HUGE_VAL || y[i] == HUGE_VAL) { + bRet = FALSE; panSuccess[i] = FALSE; continue; } @@ -459,7 +461,7 @@ int GDALGCPTransform(void *pTransformArg, int bDstToSrc, int nPointCount, panSuccess[i] = TRUE; } - return TRUE; + return bRet; } /************************************************************************/ diff --git a/alg/gdal_rpc.cpp b/alg/gdal_rpc.cpp index cbff487fc53c..fa342ae613cf 100644 --- a/alg/gdal_rpc.cpp +++ b/alg/gdal_rpc.cpp @@ -1449,10 +1449,15 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, const int nY = static_cast(dfY); const double dfDeltaY = dfY - nY; + int bRet = TRUE; for (int i = 0; i < nPointCount; i++) { if (padfX[i] == HUGE_VAL) + { + bRet = FALSE; + panSuccess[i] = FALSE; continue; + } double dfDEMH = 0.0; const double dfZ_i = padfZ ? padfZ[i] : 0.0; @@ -1502,6 +1507,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, dfDEMH = psTransform->dfDEMMissingValue; else { + bRet = FALSE; panSuccess[i] = FALSE; continue; } @@ -1546,6 +1552,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, { if (!RPCIsValidLongLat(psTransform, padfX[i], padfY[i])) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1565,6 +1572,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, { if (!RPCIsValidLongLat(psTransform, padfX[i], padfY[i])) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1582,6 +1590,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, } else { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1613,6 +1622,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, dfDEMH = psTransform->dfDEMMissingValue; else { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1623,6 +1633,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, if (!RPCIsValidLongLat(psTransform, padfX[i], padfY[i])) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1638,7 +1649,7 @@ GDALRPCTransformWholeLineWithDEM(const GDALRPCTransformInfo *psTransform, VSIFree(padfDEMBuffer); - return TRUE; + return bRet; } /************************************************************************/ @@ -1900,10 +1911,12 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, } } + int bRet = TRUE; for (int i = 0; i < nPointCount; i++) { if (!RPCIsValidLongLat(psTransform, padfX[i], padfY[i])) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1913,6 +1926,7 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, if (!GDALRPCGetHeightAtLongLat(psTransform, padfX[i], padfY[i], &dfHeight)) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1925,13 +1939,15 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, panSuccess[i] = TRUE; } - return TRUE; + return bRet; } if (padfZ == nullptr) { CPLError(CE_Failure, CPLE_NotSupported, "Z array should be provided for reverse RPC computation"); + for (int i = 0; i < nPointCount; i++) + panSuccess[i] = FALSE; return FALSE; } @@ -1940,6 +1956,7 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, /* function uses an iterative method from an initial linear */ /* approximation. */ /* -------------------------------------------------------------------- */ + int bRet = TRUE; for (int i = 0; i < nPointCount; i++) { double dfResultX = 0.0; @@ -1948,6 +1965,7 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, if (!RPCInverseTransformPoint(psTransform, padfX[i], padfY[i], padfZ[i], &dfResultX, &dfResultY)) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1955,6 +1973,7 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, } if (!RPCIsValidLongLat(psTransform, padfX[i], padfY[i])) { + bRet = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -1967,7 +1986,7 @@ int GDALRPCTransform(void *pTransformArg, int bDstToSrc, int nPointCount, panSuccess[i] = TRUE; } - return TRUE; + return bRet; } /************************************************************************/ diff --git a/alg/gdal_tps.cpp b/alg/gdal_tps.cpp index cb259c95e38f..81d9b240c945 100644 --- a/alg/gdal_tps.cpp +++ b/alg/gdal_tps.cpp @@ -324,7 +324,7 @@ void GDALDestroyTPSTransformer(void *pTransformArg) * @param panSuccess array in which a flag indicating success (TRUE) or * failure (FALSE) of the transformation are placed. * - * @return TRUE. + * @return TRUE if all points have been successfully transformed. */ int GDALTPSTransform(void *pTransformArg, int bDstToSrc, int nPointCount, diff --git a/alg/gdalgeoloc.cpp b/alg/gdalgeoloc.cpp index 32f28648cc1f..1cff39151721 100644 --- a/alg/gdalgeoloc.cpp +++ b/alg/gdalgeoloc.cpp @@ -589,6 +589,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, double *padfY, double * /* padfZ */, int *panSuccess) { + int bSuccess = TRUE; GDALGeoLocTransformInfo *psTransform = static_cast(pTransformArg); @@ -607,6 +608,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, { if (padfX[i] == HUGE_VAL || padfY[i] == HUGE_VAL) { + bSuccess = FALSE; panSuccess[i] = FALSE; continue; } @@ -623,6 +625,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, if (!PixelLineToXY(psTransform, dfGeoLocPixel, dfGeoLocLine, padfX[i], padfY[i])) { + bSuccess = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -665,6 +668,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, { if (padfX[i] == HUGE_VAL || padfY[i] == HUGE_VAL) { + bSuccess = FALSE; panSuccess[i] = FALSE; continue; } @@ -688,6 +692,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, dfBMX + 1 < psTransform->nBackMapWidth && dfBMY + 1 < psTransform->nBackMapHeight)) { + bSuccess = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -701,6 +706,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, const auto fBMY_0_0 = pAccessors->backMapYAccessor.Get(iBMX, iBMY); if (fBMX_0_0 == INVALID_BMXY) { + bSuccess = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -914,6 +920,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, } if (!bDone) { + bSuccess = FALSE; panSuccess[i] = FALSE; padfX[i] = HUGE_VAL; padfY[i] = HUGE_VAL; @@ -924,7 +931,7 @@ int GDALGeoLoc::Transform(void *pTransformArg, int bDstToSrc, } } - return TRUE; + return bSuccess; } /*! @endcond */ diff --git a/alg/gdaltransformer.cpp b/alg/gdaltransformer.cpp index 8e5e5e2b02a0..08c059d9149f 100644 --- a/alg/gdaltransformer.cpp +++ b/alg/gdaltransformer.cpp @@ -185,10 +185,9 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput(GDALDatasetH hSrcDS, adfExtent, 0); } -static int GDALSuggestedWarpOutput2_MustAdjustForRightBorder( +static bool GDALSuggestedWarpOutput2_MustAdjustForRightBorder( GDALTransformerFunc pfnTransformer, void *pTransformArg, double *padfExtent, - CPL_UNUSED int nPixels, int nLines, double dfPixelSizeX, - double dfPixelSizeY) + int /* nPixels*/, int nLines, double dfPixelSizeX, double dfPixelSizeY) { double adfX[21] = {}; double adfY[21] = {}; @@ -213,38 +212,35 @@ static int GDALSuggestedWarpOutput2_MustAdjustForRightBorder( int abSuccess[21] = {}; - bool bErr = false; - if (!pfnTransformer(pTransformArg, TRUE, nSamplePoints, adfX, adfY, adfZ, - abSuccess)) - { - bErr = true; - } + pfnTransformer(pTransformArg, TRUE, nSamplePoints, adfX, adfY, adfZ, + abSuccess); - if (!bErr && !pfnTransformer(pTransformArg, FALSE, nSamplePoints, adfX, - adfY, adfZ, abSuccess)) - { - bErr = true; - } + int abSuccess2[21] = {}; + + pfnTransformer(pTransformArg, FALSE, nSamplePoints, adfX, adfY, adfZ, + abSuccess2); nSamplePoints = 0; int nBadCount = 0; - for (double dfRatio = 0.0; !bErr && dfRatio <= 1.01; dfRatio += 0.05) + for (double dfRatio = 0.0; dfRatio <= 1.01; dfRatio += 0.05) { const double expected_x = dfMaxXOut; const double expected_y = dfMaxYOut - dfPixelSizeY * dfRatio * nLines; - if (fabs(adfX[nSamplePoints] - expected_x) > dfPixelSizeX || + if (!abSuccess[nSamplePoints] || !abSuccess2[nSamplePoints] || + fabs(adfX[nSamplePoints] - expected_x) > dfPixelSizeX || fabs(adfY[nSamplePoints] - expected_y) > dfPixelSizeY) + { nBadCount++; + } nSamplePoints++; } return nBadCount == nSamplePoints; } -static int GDALSuggestedWarpOutput2_MustAdjustForBottomBorder( +static bool GDALSuggestedWarpOutput2_MustAdjustForBottomBorder( GDALTransformerFunc pfnTransformer, void *pTransformArg, double *padfExtent, - int nPixels, CPL_UNUSED int nLines, double dfPixelSizeX, - double dfPixelSizeY) + int nPixels, int /* nLines */, double dfPixelSizeX, double dfPixelSizeY) { double adfX[21] = {}; double adfY[21] = {}; @@ -269,28 +265,26 @@ static int GDALSuggestedWarpOutput2_MustAdjustForBottomBorder( int abSuccess[21] = {}; - bool bErr = false; - if (!pfnTransformer(pTransformArg, TRUE, nSamplePoints, adfX, adfY, adfZ, - abSuccess)) - { - bErr = true; - } + pfnTransformer(pTransformArg, TRUE, nSamplePoints, adfX, adfY, adfZ, + abSuccess); - if (!bErr && !pfnTransformer(pTransformArg, FALSE, nSamplePoints, adfX, - adfY, adfZ, abSuccess)) - { - bErr = true; - } + int abSuccess2[21] = {}; + + pfnTransformer(pTransformArg, FALSE, nSamplePoints, adfX, adfY, adfZ, + abSuccess2); nSamplePoints = 0; int nBadCount = 0; - for (double dfRatio = 0.0; !bErr && dfRatio <= 1.01; dfRatio += 0.05) + for (double dfRatio = 0.0; dfRatio <= 1.01; dfRatio += 0.05) { const double expected_x = dfMinXOut + dfPixelSizeX * dfRatio * nPixels; const double expected_y = dfMinYOut; - if (fabs(adfX[nSamplePoints] - expected_x) > dfPixelSizeX || + if (!abSuccess[nSamplePoints] || !abSuccess2[nSamplePoints] || + fabs(adfX[nSamplePoints] - expected_x) > dfPixelSizeX || fabs(adfY[nSamplePoints] - expected_y) > dfPixelSizeY) + { nBadCount++; + } nSamplePoints++; } @@ -530,17 +524,10 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput2(GDALDatasetH hSrcDS, /* -------------------------------------------------------------------- */ /* Transform them to the output coordinate system. */ /* -------------------------------------------------------------------- */ - if (!pfnTransformer(pTransformArg, FALSE, nSamplePoints, padfX, padfY, - padfZ, pabSuccess)) - { - CPLError(CE_Failure, CPLE_AppDefined, - "GDALSuggestedWarpOutput() failed because the passed " - "transformer failed."); - CPLFree(padfX); - CPLFree(padfXRevert); - CPLFree(pabSuccess); - return CE_Failure; - } + CPLTurnFailureIntoWarning(true); + pfnTransformer(pTransformArg, FALSE, nSamplePoints, padfX, padfY, padfZ, + pabSuccess); + CPLTurnFailureIntoWarning(false); constexpr int SIGN_FINAL_UNINIT = -2; constexpr int SIGN_FINAL_INVALID = 0; @@ -601,6 +588,7 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput2(GDALDatasetH hSrcDS, double ayTemp[2] = {padfY[i], padfY[i]}; double azTemp[2] = {padfZ[i], padfZ[i]}; int abSuccess[2] = {FALSE, FALSE}; + CPLTurnFailureIntoWarning(true); if (pfnTransformer(pTransformArg, TRUE, 2, axTemp, ayTemp, azTemp, abSuccess) && fabs(axTemp[0] - axTemp[1]) < 1e-8 && @@ -608,6 +596,7 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput2(GDALDatasetH hSrcDS, { padfX[i] = iSignDiscontinuity * 180.0; } + CPLTurnFailureIntoWarning(false); } } } @@ -626,56 +615,53 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput2(GDALDatasetH hSrcDS, memcpy(padfXRevert, padfX, nSamplePoints * sizeof(double)); memcpy(padfYRevert, padfY, nSamplePoints * sizeof(double)); memcpy(padfZRevert, padfZ, nSamplePoints * sizeof(double)); - if (pfnTransformer(pTransformArg, TRUE, nSamplePoints, padfXRevert, - padfYRevert, padfZRevert, pabSuccess)) + CPLTurnFailureIntoWarning(true); + pfnTransformer(pTransformArg, TRUE, nSamplePoints, padfXRevert, + padfYRevert, padfZRevert, pabSuccess); + CPLTurnFailureIntoWarning(false); + + for (int i = 0; nFailedCount == 0 && i < nSamplePoints; i++) { - for (int i = 0; nFailedCount == 0 && i < nSamplePoints; i++) + if (!pabSuccess[i]) { - if (!pabSuccess[i]) - { - nFailedCount++; - break; - } + nFailedCount++; + break; + } - double dfRatio = (i % nStepsPlusOne) * dfStep; - if (dfRatio > 0.99) - dfRatio = 1.0; + double dfRatio = (i % nStepsPlusOne) * dfStep; + if (dfRatio > 0.99) + dfRatio = 1.0; - double dfExpectedX = 0.0; - double dfExpectedY = 0.0; - if (i < nStepsPlusOne) - { - dfExpectedX = dfRatio * nInXSize; - } - else if (i < 2 * nStepsPlusOne) - { - dfExpectedX = dfRatio * nInXSize; - dfExpectedY = nInYSize; - } - else if (i < 3 * nStepsPlusOne) - { - dfExpectedY = dfRatio * nInYSize; - } - else - { - dfExpectedX = nInXSize; - dfExpectedY = dfRatio * nInYSize; - } - - if (fabs(padfXRevert[i] - dfExpectedX) > - nInXSize / static_cast(nSteps) || - fabs(padfYRevert[i] - dfExpectedY) > - nInYSize / static_cast(nSteps)) - nFailedCount++; + double dfExpectedX = 0.0; + double dfExpectedY = 0.0; + if (i < nStepsPlusOne) + { + dfExpectedX = dfRatio * nInXSize; } - if (nFailedCount != 0) - CPLDebug("WARP", - "At least one point failed after revert transform"); - } - else - { - nFailedCount = 1; + else if (i < 2 * nStepsPlusOne) + { + dfExpectedX = dfRatio * nInXSize; + dfExpectedY = nInYSize; + } + else if (i < 3 * nStepsPlusOne) + { + dfExpectedY = dfRatio * nInYSize; + } + else + { + dfExpectedX = nInXSize; + dfExpectedY = dfRatio * nInYSize; + } + + if (fabs(padfXRevert[i] - dfExpectedX) > + nInXSize / static_cast(nSteps) || + fabs(padfYRevert[i] - dfExpectedY) > + nInYSize / static_cast(nSteps)) + nFailedCount++; } + if (nFailedCount != 0) + CPLDebug("WARP", + "At least one point failed after revert transform"); } /* -------------------------------------------------------------------- */ @@ -707,19 +693,10 @@ CPLErr CPL_STDCALL GDALSuggestedWarpOutput2(GDALDatasetH hSrcDS, CPLAssert(nSamplePoints == nSampleMax); - if (!pfnTransformer(pTransformArg, FALSE, nSamplePoints, padfX, padfY, - padfZ, pabSuccess)) - { - CPLError(CE_Failure, CPLE_AppDefined, - "GDALSuggestedWarpOutput() failed because the passed" - "transformer failed."); - - CPLFree(padfX); - CPLFree(padfXRevert); - CPLFree(pabSuccess); - - return CE_Failure; - } + CPLTurnFailureIntoWarning(true); + pfnTransformer(pTransformArg, FALSE, nSamplePoints, padfX, padfY, padfZ, + pabSuccess); + CPLTurnFailureIntoWarning(false); } /* -------------------------------------------------------------------- */ @@ -2910,11 +2887,12 @@ int GDALGenImgProjTransform(void *pTransformArgIn, int bDstToSrc, pTransformer = psInfo->pSrcTransformer; } + int ret = TRUE; if (pTransformArg != nullptr) { if (!pTransformer(pTransformArg, FALSE, nPointCount, padfX, padfY, padfZ, panSuccess)) - return FALSE; + ret = FALSE; } else { @@ -2942,7 +2920,7 @@ int GDALGenImgProjTransform(void *pTransformArgIn, int bDstToSrc, { if (!psInfo->pReproject(psInfo->pReprojectArg, bDstToSrc, nPointCount, padfX, padfY, padfZ, panSuccess)) - return FALSE; + ret = FALSE; } /* -------------------------------------------------------------------- */ @@ -2965,7 +2943,7 @@ int GDALGenImgProjTransform(void *pTransformArgIn, int bDstToSrc, { if (!pTransformer(pTransformArg, TRUE, nPointCount, padfX, padfY, padfZ, panSuccess)) - return FALSE; + ret = FALSE; } else { @@ -2986,7 +2964,7 @@ int GDALGenImgProjTransform(void *pTransformArgIn, int bDstToSrc, } } - return TRUE; + return ret; } /************************************************************************/ diff --git a/alg/gdalwarpoperation.cpp b/alg/gdalwarpoperation.cpp index f4997c01b16a..fc0bbe9f2956 100644 --- a/alg/gdalwarpoperation.cpp +++ b/alg/gdalwarpoperation.cpp @@ -2579,14 +2579,10 @@ void GDALWarpOperation::ComputeSourceWindowStartingFromSource( /* Transform them to the output pixel coordinate space */ /* -------------------------------------------------------------------- */ - if (!psOptions->pfnTransformer( - psOptions->pTransformerArg, FALSE, nSampleMax, - privateData->adfDstX.data(), privateData->adfDstY.data(), - adfDstZ.data(), privateData->abSuccess.data())) - { - return; - } - + psOptions->pfnTransformer(psOptions->pTransformerArg, FALSE, nSampleMax, + privateData->adfDstX.data(), + privateData->adfDstY.data(), adfDstZ.data(), + privateData->abSuccess.data()); privateData->nStepCount = nStepCount; } @@ -2834,26 +2830,14 @@ bool GDALWarpOperation::ComputeSourceWindowTransformPoints( CPLSetThreadLocalConfigOption("CHECK_WITH_INVERT_PROJ", "YES"); RefreshTransformer(); } - int ret = psOptions->pfnTransformer(psOptions->pTransformerArg, TRUE, - nSamplePoints, padfX, padfY, padfZ, - pabSuccess); + psOptions->pfnTransformer(psOptions->pTransformerArg, TRUE, nSamplePoints, + padfX, padfY, padfZ, pabSuccess); if (bTryWithCheckWithInvertProj) { CPLSetThreadLocalConfigOption("CHECK_WITH_INVERT_PROJ", nullptr); RefreshTransformer(); } - if (!ret) - { - CPLFree(padfX); - CPLFree(pabSuccess); - - CPLError(CE_Failure, CPLE_AppDefined, - "GDALWarperOperation::ComputeSourceWindow() failed because " - "the pfnTransformer failed."); - return false; - } - /* -------------------------------------------------------------------- */ /* Collect the bounds, ignoring any failed points. */ /* -------------------------------------------------------------------- */ diff --git a/apps/gdalwarp_lib.cpp b/apps/gdalwarp_lib.cpp index e0dd6c9b7e19..95bedf8ef0b7 100644 --- a/apps/gdalwarp_lib.cpp +++ b/apps/gdalwarp_lib.cpp @@ -2762,36 +2762,35 @@ static GDALDatasetH GDALWarpDirect(const char *pszDest, GDALDatasetH hDstDS, } } std::vector abSuccess(nPoints); - if (pfnTransformer(hTransformArg, TRUE, nPoints, &adfX[0], - &adfY[0], &adfZ[0], &abSuccess[0])) + pfnTransformer(hTransformArg, TRUE, nPoints, &adfX[0], &adfY[0], + &adfZ[0], &abSuccess[0]); + + double dfMinSrcX = std::numeric_limits::infinity(); + double dfMaxSrcX = -std::numeric_limits::infinity(); + double dfMinSrcY = std::numeric_limits::infinity(); + double dfMaxSrcY = -std::numeric_limits::infinity(); + for (int i = 0; i < nPoints; i++) { - double dfMinSrcX = std::numeric_limits::infinity(); - double dfMaxSrcX = -std::numeric_limits::infinity(); - double dfMinSrcY = std::numeric_limits::infinity(); - double dfMaxSrcY = -std::numeric_limits::infinity(); - for (int i = 0; i < nPoints; i++) + if (abSuccess[i]) { - if (abSuccess[i]) - { - dfMinSrcX = std::min(dfMinSrcX, adfX[i]); - dfMaxSrcX = std::max(dfMaxSrcX, adfX[i]); - dfMinSrcY = std::min(dfMinSrcY, adfY[i]); - dfMaxSrcY = std::max(dfMaxSrcY, adfY[i]); - } - } - if (dfMaxSrcX > dfMinSrcX) - { - dfTargetRatioX = (dfMaxSrcX - dfMinSrcX) / - GDALGetRasterXSize(hDstDS); - } - if (dfMaxSrcY > dfMinSrcY) - { - dfTargetRatioY = (dfMaxSrcY - dfMinSrcY) / - GDALGetRasterYSize(hDstDS); + dfMinSrcX = std::min(dfMinSrcX, adfX[i]); + dfMaxSrcX = std::max(dfMaxSrcX, adfX[i]); + dfMinSrcY = std::min(dfMinSrcY, adfY[i]); + dfMaxSrcY = std::max(dfMaxSrcY, adfY[i]); } - // take the minimum of these ratios #7019 - dfTargetRatio = std::min(dfTargetRatioX, dfTargetRatioY); } + if (dfMaxSrcX > dfMinSrcX) + { + dfTargetRatioX = + (dfMaxSrcX - dfMinSrcX) / GDALGetRasterXSize(hDstDS); + } + if (dfMaxSrcY > dfMinSrcY) + { + dfTargetRatioY = + (dfMaxSrcY - dfMinSrcY) / GDALGetRasterYSize(hDstDS); + } + // take the minimum of these ratios #7019 + dfTargetRatio = std::min(dfTargetRatioX, dfTargetRatioY); } else { diff --git a/autotest/gcore/geoloc.py b/autotest/gcore/geoloc.py index eb5460b2587b..24f613e2f6b0 100755 --- a/autotest/gcore/geoloc.py +++ b/autotest/gcore/geoloc.py @@ -118,10 +118,7 @@ def test_geoloc_fill_line(use_temp_datasets): with gdaltest.config_option("GDAL_GEOLOC_USE_TEMP_DATASETS", use_temp_datasets): warped_ds = gdal.Warp("", ds, format="MEM") assert warped_ds - assert warped_ds.GetRasterBand(1).Checksum() in ( - 22339, - 22336, - ) # 22336 with Intel(R) oneAPI DPC++/C++ Compiler 2022.1.0 + assert warped_ds.GetRasterBand(1).Checksum() == 20177 ############################################################################### diff --git a/ogr/ogr_spatialref.h b/ogr/ogr_spatialref.h index 23b236006fe2..be0fc4a9da67 100644 --- a/ogr/ogr_spatialref.h +++ b/ogr/ogr_spatialref.h @@ -813,8 +813,8 @@ class CPL_DLL OGRCoordinateTransformation * @param pabSuccess array of per-point flags set to TRUE if that point * transforms, or FALSE if it does not. Might be NULL. * - * @return TRUE if a transformation could be found (but not all points may - * have necessarily succeed to transform), otherwise FALSE. + * @return TRUE on success, or FALSE if some or all points fail to + * transform. */ int Transform(size_t nCount, double *x, double *y, double *z = nullptr, int *pabSuccess = nullptr); @@ -835,8 +835,10 @@ class CPL_DLL OGRCoordinateTransformation * @param pabSuccess array of per-point flags set to TRUE if that point * transforms, or FALSE if it does not. Might be NULL. * - * @return TRUE if a transformation could be found (but not all points may - * have necessarily succeed to transform), otherwise FALSE. + * @return TRUE on success, or FALSE if some or all points fail to + * transform. Caution: prior to GDAL 3.11, TRUE could be returned if a + * transformation could be found but not all points may + * have necessarily succeed to transform. */ virtual int Transform(size_t nCount, double *x, double *y, double *z, double *t, int *pabSuccess) = 0; @@ -857,8 +859,10 @@ class CPL_DLL OGRCoordinateTransformation * @param panErrorCodes Output array of nCount value that will be set to 0 * for success, or a non-zero value for failure. Refer to PROJ 8 public * error codes. Might be NULL - * @return TRUE if a transformation could be found (but not all points may - * have necessarily succeed to transform), otherwise FALSE. + * @return TRUE on success, or FALSE if some or all points fail to + * transform. Caution: prior to GDAL 3.11, TRUE could be returned if a + * transformation could be found but not all points may + * have necessarily succeed to transform. * @since GDAL 3.3, and PROJ 8 to be able to use PROJ public error codes */ virtual int TransformWithErrorCodes(size_t nCount, double *x, double *y, diff --git a/ogr/ogrct.cpp b/ogr/ogrct.cpp index 4da6f38a037f..ae93f7163579 100644 --- a/ogr/ogrct.cpp +++ b/ogr/ogrct.cpp @@ -2188,22 +2188,12 @@ int OGRCoordinateTransformation::Transform(size_t nCount, double *x, double *y, if (!pabSuccess) return FALSE; - bool bOverallSuccess = - CPL_TO_BOOL(Transform(nCount, x, y, z, nullptr, pabSuccess)); - - for (size_t i = 0; i < nCount; i++) - { - if (!pabSuccess[i]) - { - bOverallSuccess = false; - break; - } - } + const int bRet = Transform(nCount, x, y, z, nullptr, pabSuccess); if (pabSuccess != pabSuccessIn) CPLFree(pabSuccess); - return bOverallSuccess; + return bRet; } /************************************************************************/ @@ -2219,13 +2209,12 @@ int OGRCoordinateTransformation::TransformWithErrorCodes(size_t nCount, if (nCount == 1) { int nSuccess = 0; - const bool bOverallSuccess = - CPL_TO_BOOL(Transform(nCount, x, y, z, t, &nSuccess)); + const int bRet = Transform(nCount, x, y, z, t, &nSuccess); if (panErrorCodes) { panErrorCodes[0] = nSuccess ? 0 : -1; } - return bOverallSuccess; + return bRet; } std::vector abSuccess; @@ -2240,8 +2229,7 @@ int OGRCoordinateTransformation::TransformWithErrorCodes(size_t nCount, return FALSE; } - const bool bOverallSuccess = - CPL_TO_BOOL(Transform(nCount, x, y, z, t, abSuccess.data())); + const int bRet = Transform(nCount, x, y, z, t, abSuccess.data()); if (panErrorCodes) { @@ -2251,7 +2239,7 @@ int OGRCoordinateTransformation::TransformWithErrorCodes(size_t nCount, } } - return bOverallSuccess; + return bRet; } /************************************************************************/ @@ -2262,8 +2250,7 @@ int OGRProjCT::Transform(size_t nCount, double *x, double *y, double *z, double *t, int *pabSuccess) { - bool bOverallSuccess = - CPL_TO_BOOL(TransformWithErrorCodes(nCount, x, y, z, t, pabSuccess)); + const int bRet = TransformWithErrorCodes(nCount, x, y, z, t, pabSuccess); if (pabSuccess) { @@ -2273,7 +2260,7 @@ int OGRProjCT::Transform(size_t nCount, double *x, double *y, double *z, } } - return bOverallSuccess; + return bRet; } /************************************************************************/ @@ -2405,6 +2392,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, /* Optimized transform from WebMercator to WGS84 */ /* -------------------------------------------------------------------- */ bool bTransformDone = false; + int bRet = TRUE; if (bWebMercatorToWGS84LongLat) { constexpr double REVERSE_SPHERE_RADIUS = 1.0 / 6378137.0; @@ -2417,7 +2405,11 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, double y0 = y[0]; for (size_t i = 0; i < nCount; i++) { - if (x[i] != HUGE_VAL) + if (x[i] == HUGE_VAL) + { + bRet = FALSE; + } + else { x[i] = x[i] * REVERSE_SPHERE_RADIUS; if (x[i] > M_PI) @@ -2697,6 +2689,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, const double yIn = y[i]; if (!std::isfinite(xIn)) { + bRet = FALSE; x[i] = HUGE_VAL; y[i] = HUGE_VAL; if (panErrorCodes) @@ -2725,6 +2718,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, int err = 0; if (std::isnan(coord.xyzt.x)) { + bRet = FALSE; // This shouldn't normally happen if PROJ projections behave // correctly, but e.g inverse laea before PROJ 8.1.1 could // do that for points out of domain. @@ -2743,6 +2737,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, } else if (coord.xyzt.x == HUGE_VAL) { + bRet = FALSE; err = proj_errno(pj); // PROJ should normally emit an error, but in case it does not // (e.g PROJ 6.3 with the +ortho projection), synthetize one @@ -2759,6 +2754,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, if (fabs(coord.xyzt.x - xIn) > dfThreshold || fabs(coord.xyzt.y - yIn) > dfThreshold) { + bRet = FALSE; err = PROJ_ERR_COORD_TRANSFM_OUTSIDE_PROJECTION_DOMAIN; x[i] = HUGE_VAL; y[i] = HUGE_VAL; @@ -2936,7 +2932,7 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, // static_cast(delay * 1000)); #endif - return TRUE; + return bRet; } /************************************************************************/ @@ -2944,39 +2940,42 @@ int OGRProjCT::TransformWithErrorCodes(size_t nCount, double *x, double *y, /************************************************************************/ // --------------------------------------------------------------------------- -static double simple_min(const double *data, const int arr_len) +static double simple_min(const double *data, const int *panErrorCodes, + const int arr_len) { - double min_value = data[0]; - for (int iii = 1; iii < arr_len; iii++) + double min_value = HUGE_VAL; + for (int iii = 0; iii < arr_len; iii++) { - if (data[iii] < min_value) + if ((data[iii] < min_value || min_value == HUGE_VAL) && + panErrorCodes[iii] == 0) min_value = data[iii]; } return min_value; } // --------------------------------------------------------------------------- -static double simple_max(const double *data, const int arr_len) +static double simple_max(const double *data, const int *panErrorCodes, + const int arr_len) { - double max_value = data[0]; - for (int iii = 1; iii < arr_len; iii++) + double max_value = HUGE_VAL; + for (int iii = 0; iii < arr_len; iii++) { if ((data[iii] > max_value || max_value == HUGE_VAL) && - data[iii] != HUGE_VAL) + panErrorCodes[iii] == 0) max_value = data[iii]; } return max_value; } // --------------------------------------------------------------------------- -static int _find_previous_index(const int iii, const double *data, +static int _find_previous_index(const int iii, const int *panErrorCodes, const int arr_len) { // find index of nearest valid previous value if exists int prev_iii = iii - 1; if (prev_iii == -1) // handle wraparound prev_iii = arr_len - 1; - while (data[prev_iii] == HUGE_VAL && prev_iii != iii) + while (panErrorCodes[prev_iii] != 0 && prev_iii != iii) { prev_iii--; if (prev_iii == -1) // handle wraparound @@ -3033,7 +3032,8 @@ but smalller than 240 to account for possible irregularities in distances when re-projecting. Also, 200 ensures latitudes are ignored for axis order handling. ******************************************************************************/ -static double antimeridian_min(const double *data, const int arr_len) +static double antimeridian_min(const double *data, const int *panErrorCodes, + const int arr_len) { double positive_min = HUGE_VAL; double min_value = HUGE_VAL; @@ -3042,9 +3042,9 @@ static double antimeridian_min(const double *data, const int arr_len) for (int iii = 0; iii < arr_len; iii++) { - if (data[iii] == HUGE_VAL) + if (panErrorCodes[iii]) continue; - int prev_iii = _find_previous_index(iii, data, arr_len); + int prev_iii = _find_previous_index(iii, panErrorCodes, arr_len); // check if crossed meridian double delta = data[prev_iii] - data[iii]; // 180 -> -180 @@ -3086,7 +3086,8 @@ static double antimeridian_min(const double *data, const int arr_len) // Note: This requires a densified ring with at least 2 additional // points per edge to correctly handle global extents. // See antimeridian_min docstring for reasoning. -static double antimeridian_max(const double *data, const int arr_len) +static double antimeridian_max(const double *data, const int *panErrorCodes, + const int arr_len) { double negative_max = -HUGE_VAL; double max_value = -HUGE_VAL; @@ -3095,9 +3096,9 @@ static double antimeridian_max(const double *data, const int arr_len) for (int iii = 0; iii < arr_len; iii++) { - if (data[iii] == HUGE_VAL) + if (panErrorCodes[iii]) continue; - int prev_iii = _find_previous_index(iii, data, arr_len); + int prev_iii = _find_previous_index(iii, panErrorCodes, arr_len); // check if crossed meridian double delta = data[prev_iii] - data[iii]; // 180 -> -180 @@ -3119,11 +3120,11 @@ static double antimeridian_max(const double *data, const int arr_len) // negative meridian side max if (negative_meridian && (data[iii] > negative_max || negative_max == HUGE_VAL) && - data[iii] != HUGE_VAL) + panErrorCodes[iii] == 0) negative_max = data[iii]; // track general max value if ((data[iii] > max_value || max_value == HUGE_VAL) && - data[iii] != HUGE_VAL) + panErrorCodes[iii] == 0) max_value = data[iii]; } if (crossed_meridian_count == 2) @@ -3149,17 +3150,14 @@ bool OGRProjCT::ContainsNorthPole(const double xmin, const double ymin, pole_y = 0; pole_x = 90; } - auto inverseCT = GetInverse(); + auto inverseCT = std::unique_ptr(GetInverse()); if (!inverseCT) return false; - bool success = inverseCT->TransformWithErrorCodes( + CPLErrorStateBackuper oBackuper(CPLQuietErrorHandler); + const bool success = inverseCT->TransformWithErrorCodes( 1, &pole_x, &pole_y, nullptr, nullptr, nullptr); - if (success && CPLGetLastErrorType() != CE_None) - CPLErrorReset(); - delete inverseCT; - if (xmin < pole_x && pole_x < xmax && ymax > pole_y && pole_y > ymin) - return true; - return false; + return success && xmin < pole_x && pole_x < xmax && ymax > pole_y && + pole_y > ymin; } // --------------------------------------------------------------------------- @@ -3177,17 +3175,14 @@ bool OGRProjCT::ContainsSouthPole(const double xmin, const double ymin, pole_y = 0; pole_x = -90; } - auto inverseCT = GetInverse(); + auto inverseCT = std::unique_ptr(GetInverse()); if (!inverseCT) return false; - bool success = inverseCT->TransformWithErrorCodes( + CPLErrorStateBackuper oBackuper(CPLQuietErrorHandler); + const bool success = inverseCT->TransformWithErrorCodes( 1, &pole_x, &pole_y, nullptr, nullptr, nullptr); - if (success && CPLGetLastErrorType() != CE_None) - CPLErrorReset(); - delete inverseCT; - if (xmin < pole_x && pole_x < xmax && ymax > pole_y && pole_y > ymin) - return true; - return false; + return success && xmin < pole_x && pole_x < xmax && ymax > pole_y && + pole_y > ymin; } int OGRProjCT::TransformBounds(const double xmin, const double ymin, @@ -3196,8 +3191,6 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, double *out_xmax, double *out_ymax, const int densify_pts) { - CPLErrorReset(); - if (bNoTransform) { *out_xmin = xmin; @@ -3268,10 +3261,12 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, const int boundary_len = side_pts * 4; std::vector x_boundary_array; std::vector y_boundary_array; + std::vector anErrorCodes; try { x_boundary_array.resize(boundary_len); y_boundary_array.resize(boundary_len); + anErrorCodes.resize(boundary_len); } catch (const std::exception &e) // memory allocation failure { @@ -3284,20 +3279,10 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, bool south_pole_in_bounds = false; if (degree_output) { - CPLErrorHandlerPusher oErrorHandlerPusher(CPLQuietErrorHandler); - north_pole_in_bounds = ContainsNorthPole(xmin, ymin, xmax, ymax, output_lon_lat_order); - if (CPLGetLastErrorType() != CE_None) - { - return false; - } south_pole_in_bounds = ContainsSouthPole(xmin, ymin, xmax, ymax, output_lon_lat_order); - if (CPLGetLastErrorType() != CE_None) - { - return false; - } } if (degree_input && xmax < xmin) @@ -3350,26 +3335,35 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, } { - CPLErrorHandlerPusher oErrorHandlerPusher(CPLQuietErrorHandler); + CPLErrorStateBackuper oBackuper(CPLQuietErrorHandler); bool success = TransformWithErrorCodes( boundary_len, &x_boundary_array[0], &y_boundary_array[0], nullptr, - nullptr, nullptr); - if (success && CPLGetLastErrorType() != CE_None) - { - CPLErrorReset(); - } - else if (!success) + nullptr, anErrorCodes.data()); + if (!success) { - return false; + for (int i = 0; i < boundary_len; ++i) + { + if (anErrorCodes[i] == 0) + { + success = true; + break; + } + } + if (!success) + return false; } } if (!degree_output) { - *out_xmin = simple_min(&x_boundary_array[0], boundary_len); - *out_xmax = simple_max(&x_boundary_array[0], boundary_len); - *out_ymin = simple_min(&y_boundary_array[0], boundary_len); - *out_ymax = simple_max(&y_boundary_array[0], boundary_len); + *out_xmin = + simple_min(&x_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_xmax = + simple_max(&x_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_ymin = + simple_min(&y_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_ymax = + simple_max(&y_boundary_array[0], anErrorCodes.data(), boundary_len); if (poSRSTarget->IsProjected()) { @@ -3475,13 +3469,15 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, else if (north_pole_in_bounds && output_lon_lat_order) { *out_xmin = -180; - *out_ymin = simple_min(&y_boundary_array[0], boundary_len); + *out_ymin = + simple_min(&y_boundary_array[0], anErrorCodes.data(), boundary_len); *out_xmax = 180; *out_ymax = 90; } else if (north_pole_in_bounds) { - *out_xmin = simple_min(&x_boundary_array[0], boundary_len); + *out_xmin = + simple_min(&x_boundary_array[0], anErrorCodes.data(), boundary_len); *out_ymin = -180; *out_xmax = 90; *out_ymax = 180; @@ -3491,28 +3487,38 @@ int OGRProjCT::TransformBounds(const double xmin, const double ymin, *out_xmin = -180; *out_ymin = -90; *out_xmax = 180; - *out_ymax = simple_max(&y_boundary_array[0], boundary_len); + *out_ymax = + simple_max(&y_boundary_array[0], anErrorCodes.data(), boundary_len); } else if (south_pole_in_bounds) { *out_xmin = -90; *out_ymin = -180; - *out_xmax = simple_max(&x_boundary_array[0], boundary_len); + *out_xmax = + simple_max(&x_boundary_array[0], anErrorCodes.data(), boundary_len); *out_ymax = 180; } else if (output_lon_lat_order) { - *out_xmin = antimeridian_min(&x_boundary_array[0], boundary_len); - *out_xmax = antimeridian_max(&x_boundary_array[0], boundary_len); - *out_ymin = simple_min(&y_boundary_array[0], boundary_len); - *out_ymax = simple_max(&y_boundary_array[0], boundary_len); + *out_xmin = antimeridian_min(&x_boundary_array[0], anErrorCodes.data(), + boundary_len); + *out_xmax = antimeridian_max(&x_boundary_array[0], anErrorCodes.data(), + boundary_len); + *out_ymin = + simple_min(&y_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_ymax = + simple_max(&y_boundary_array[0], anErrorCodes.data(), boundary_len); } else { - *out_xmin = simple_min(&x_boundary_array[0], boundary_len); - *out_xmax = simple_max(&x_boundary_array[0], boundary_len); - *out_ymin = antimeridian_min(&y_boundary_array[0], boundary_len); - *out_ymax = antimeridian_max(&y_boundary_array[0], boundary_len); + *out_xmin = + simple_min(&x_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_xmax = + simple_max(&x_boundary_array[0], anErrorCodes.data(), boundary_len); + *out_ymin = antimeridian_min(&y_boundary_array[0], anErrorCodes.data(), + boundary_len); + *out_ymax = antimeridian_max(&y_boundary_array[0], anErrorCodes.data(), + boundary_len); } return *out_xmin != HUGE_VAL && *out_ymin != HUGE_VAL &&