Merge pull request #10898 from OSGeo/backport-10895-to-release/3.9

[Backport release/3.9] VRTProcessedDataset: fix issue when computing RasterIO window on auxiliary datasets on right-most/bottom-most tiles
OSGeo · Sep 30, 2024 · 88cedd8 · 88cedd8
2 parents a2f4c2d + 9c5f69b
commit 88cedd8
Show file tree

Hide file tree

Showing 2 changed files with 56 additions and 2 deletions.
diff --git a/autotest/gdrivers/vrtprocesseddataset.py b/autotest/gdrivers/vrtprocesseddataset.py
@@ -719,6 +719,58 @@ def test_vrtprocesseddataset_dehazing_different_resolution(tmp_vsimem):
     )
 
 
+###############################################################################
+# Test we properly request auxiliary datasets on the right-most/bottom-most
+# truncated tile
+
+
+def test_vrtprocesseddataset_dehazing_edge_effects(tmp_vsimem):
+
+    src_filename = str(tmp_vsimem / "src.tif")
+    src_ds = gdal.GetDriverByName("GTiff").Create(
+        src_filename,
+        257,
+        257,
+        1,
+        gdal.GDT_Byte,
+        ["TILED=YES", "BLOCKXSIZE=256", "BLOCKYSIZE=256"],
+    )
+    src_ds.GetRasterBand(1).Fill(10)
+    src_ds.SetGeoTransform([0, 1, 0, 0, 0, -1])
+    src_ds.Close()
+
+    gain_filename = str(tmp_vsimem / "gain.tif")
+    gain_ds = gdal.GetDriverByName("GTiff").Create(gain_filename, 1, 1)
+    gain_ds.GetRasterBand(1).Fill(2)
+    gain_ds.SetGeoTransform([0, 257, 0, 0, 0, -257])
+    gain_ds.Close()
+
+    offset_filename = str(tmp_vsimem / "offset.tif")
+    offset_ds = gdal.GetDriverByName("GTiff").Create(offset_filename, 1, 1)
+    offset_ds.GetRasterBand(1).Fill(3)
+    offset_ds.SetGeoTransform([0, 257, 0, 0, 0, -257])
+    offset_ds.Close()
+
+    ds = gdal.Open(
+        f"""<VRTDataset subclass='VRTProcessedDataset'>
+    <Input>
+        <SourceFilename>{src_filename}</SourceFilename>
+    </Input>
+    <ProcessingSteps>
+        <Step>
+            <Algorithm>LocalScaleOffset</Algorithm>
+            <Argument name="gain_dataset_filename_1">{gain_filename}</Argument>
+            <Argument name="gain_dataset_band_1">1</Argument>
+            <Argument name="offset_dataset_filename_1">{offset_filename}</Argument>
+            <Argument name="offset_dataset_band_1">1</Argument>
+        </Step>
+    </ProcessingSteps>
+    </VRTDataset>
+        """
+    )
+    assert ds.GetRasterBand(1).ComputeRasterMinMax() == (17, 17)
+
+
 ###############################################################################
 # Test error cases of LocalScaleOffset algorithm
 

diff --git a/frmts/vrt/vrtprocesseddatasetfunctions.cpp b/frmts/vrt/vrtprocesseddatasetfunctions.cpp
@@ -963,8 +963,10 @@ static bool LoadAuxData(double dfULX, double dfULY, double dfLRX, double dfLRY,
         return false;
     }
 
-    const int nAuxXOff = std::max(0, static_cast<int>(std::round(dfULPixel)));
-    const int nAuxYOff = std::max(0, static_cast<int>(std::round(dfULLine)));
+    const int nAuxXOff = std::clamp(static_cast<int>(std::round(dfULPixel)), 0,
+                                    poAuxBand->GetXSize() - 1);
+    const int nAuxYOff = std::clamp(static_cast<int>(std::round(dfULLine)), 0,
+                                    poAuxBand->GetYSize() - 1);
     const int nAuxX2Off = std::min(poAuxBand->GetXSize(),
                                    static_cast<int>(std::round(dfLRPixel)));
     const int nAuxY2Off =