Merge pull request #2018 from rouault/fix_vrt_mask_implicit_overviews

VRT: fix requesting a downsampling version of the mask band of a source that has masks and overviews
OSGeo · Nov 15, 2019 · e9627bf · e9627bf
2 parents 8f549e5 + a01e314
commit e9627bf
Show file tree

Hide file tree

Showing 3 changed files with 73 additions and 21 deletions.
diff --git a/autotest/gcore/vrtmisc.py b/autotest/gcore/vrtmisc.py
@@ -29,6 +29,7 @@
 # DEALINGS IN THE SOFTWARE.
 ###############################################################################
 
+import gdaltest
 from osgeo import gdal
 from osgeo import osr
 
@@ -515,6 +516,27 @@ def test_vrtmisc_write_srs():
 
     gdal.Unlink(tmpfile)
 
+###############################################################################
+# complex scenario involving masks and implicit overviews
+
+def test_vrtmisc_mask_implicit_overviews():
+
+    with gdaltest.config_option('GDAL_TIFF_INTERNAL_MASK', 'YES'):
+        ds = gdal.Translate('/vsimem/cog.tif', 'data/stefan_full_rgba.tif', options = '-outsize 2048 0 -b 1 -b 2 -b 3 -mask 4')
+        ds.BuildOverviews('NEAR', [2, 4])
+        ds = None
+    gdal.Translate('/vsimem/cog.vrt', '/vsimem/cog.tif')
+    gdal.Translate('/vsimem/out.tif', '/vsimem/cog.vrt', options = '-b mask -outsize 10% 0')
+    gdal.GetDriverByName('GTiff').Delete('/vsimem/cog.tif')
+    gdal.Unlink('/vsimem/cog.vrt')
+    ds = gdal.Open('/vsimem/out.tif')
+    histo = ds.GetRasterBand(1).GetDefaultHistogram()[3]
+    # Check that there are only 0 and 255 in the histogram
+    assert histo[0] + histo[255] == ds.RasterXSize * ds.RasterYSize, histo
+    assert ds.GetRasterBand(1).Checksum() == 46885
+    ds = None
+    gdal.Unlink('/vsimem/out.tif')
+
 ###############################################################################
 # Cleanup.
 

diff --git a/gdal/frmts/vrt/vrtdataset.cpp b/gdal/frmts/vrt/vrtdataset.cpp
@@ -1848,44 +1848,58 @@ void VRTDataset::BuildVirtualOverviews()
 
     int nOverviews = 0;
     GDALRasterBand* poFirstBand = nullptr;
-    for( int iBand = 0; iBand < nBands; iBand++ )
+
+    const auto CheckBandForOverview =
+        [&nOverviews, &poFirstBand, this](GDALRasterBand* poBand)
     {
-        if( !reinterpret_cast<VRTRasterBand *>(
-               papoBands[iBand] )->IsSourcedRasterBand())
-            return;
+        if( !cpl::down_cast<VRTRasterBand *>(poBand)->IsSourcedRasterBand())
+            return false;
 
         VRTSourcedRasterBand* poVRTBand
-            = reinterpret_cast<VRTSourcedRasterBand *>( papoBands[iBand] );
+            = cpl::down_cast<VRTSourcedRasterBand *>(poBand);
         if( poVRTBand->nSources != 1 )
-            return;
+            return false;
         if( !poVRTBand->papoSources[0]->IsSimpleSource() )
-            return;
+            return false;
 
         VRTSimpleSource* poSource
-            = reinterpret_cast<VRTSimpleSource *>( poVRTBand->papoSources[0] );
+            = cpl::down_cast<VRTSimpleSource *>( poVRTBand->papoSources[0] );
         if( !EQUAL(poSource->GetType(), "SimpleSource") &&
             !EQUAL(poSource->GetType(), "ComplexSource") )
-            return;
+            return false;
         GDALRasterBand* poSrcBand = poSource->GetBand();
         if( poSrcBand == nullptr )
-            return;
+            return false;
 
         // To prevent recursion
         m_apoOverviewsBak.push_back(nullptr);
         const int nOvrCount = poSrcBand->GetOverviewCount();
         m_apoOverviewsBak.resize(0);
 
         if( nOvrCount == 0 )
-            return;
-        if( iBand == 0 )
+            return false;
+        if( poFirstBand == nullptr )
         {
             if( poSrcBand->GetXSize() == 0 || poSrcBand->GetYSize() == 0 )
-                return;
+                return false;
             poFirstBand = poSrcBand;
             nOverviews = nOvrCount;
         }
         else if( nOvrCount < nOverviews )
             nOverviews = nOvrCount;
+        return true;
+    };
+
+    for( int iBand = 0; iBand < nBands; iBand++ )
+    {
+        if( !CheckBandForOverview(papoBands[iBand]) )
+            return;
+    }
+
+    if( m_poMaskBand )
+    {
+        if( !CheckBandForOverview(m_poMaskBand) )
+            return;
     }
 
     for( int j = 0; j < nOverviews; j++)
@@ -1904,19 +1918,19 @@ void VRTDataset::BuildVirtualOverviews()
         VRTDataset* poOvrVDS = new VRTDataset(nOvrXSize, nOvrYSize);
         m_apoOverviews.push_back(poOvrVDS);
 
-        for( int i = 0; i < nBands; i++ )
+        const auto CreateOverviewBand =
+            [&poOvrVDS, nOvrXSize, nOvrYSize, dfXRatio, dfYRatio]
+            (GDALRasterBand* poBand)
         {
             VRTSourcedRasterBand* poVRTBand
-                = reinterpret_cast<VRTSourcedRasterBand *>(
-                    GetRasterBand(i+1) );
+                = cpl::down_cast<VRTSourcedRasterBand *>(poBand);
             VRTSourcedRasterBand* poOvrVRTBand = new VRTSourcedRasterBand(
                 poOvrVDS,
-                poOvrVDS->GetRasterCount() + 1,
+                poBand->GetBand(),
                 poVRTBand->GetRasterDataType(),
                 nOvrXSize, nOvrYSize);
-            poOvrVDS->SetBand( poOvrVDS->GetRasterCount() + 1, poOvrVRTBand );
 
-            VRTSimpleSource* poSrcSource = reinterpret_cast<VRTSimpleSource *>(
+            VRTSimpleSource* poSrcSource = cpl::down_cast<VRTSimpleSource *>(
                 poVRTBand->papoSources[0] );
             VRTSimpleSource* poNewSource = nullptr;
             if( EQUAL(poSrcSource->GetType(), "SimpleSource") )
@@ -1927,7 +1941,7 @@ void VRTDataset::BuildVirtualOverviews()
             else if( EQUAL(poSrcSource->GetType(), "ComplexSource") )
             {
               poNewSource = new VRTComplexSource(
-                  reinterpret_cast<VRTComplexSource *>( poSrcSource ),
+                  cpl::down_cast<VRTComplexSource *>( poSrcSource ),
                   dfXRatio, dfYRatio );
             }
             else
@@ -1940,6 +1954,19 @@ void VRTDataset::BuildVirtualOverviews()
                     poNewSource->GetBand()->GetDataset()->Reference();
                 poOvrVRTBand->AddSource(poNewSource);
             }
+
+            return poOvrVRTBand;
+        };
+
+        for( int i = 0; i < nBands; i++ )
+        {
+            poOvrVDS->SetBand( poOvrVDS->GetRasterCount() + 1,
+                               CreateOverviewBand(GetRasterBand(i+1)) );
+        }
+
+        if( m_poMaskBand )
+        {
+            poOvrVDS->SetMaskBand( CreateOverviewBand(m_poMaskBand) );
         }
     }
 }

diff --git a/gdal/frmts/vrt/vrtrasterband.cpp b/gdal/frmts/vrt/vrtrasterband.cpp
@@ -1200,7 +1200,10 @@ GDALRasterBand *VRTRasterBand::GetOverview( int iOverview )
             || iOverview >= static_cast<int>( poVRTDS->m_apoOverviews.size() ) )
             return nullptr;
 
-        return poVRTDS->m_apoOverviews[iOverview]->GetRasterBand(nBand ? nBand : 1);
+        auto poOvrBand = poVRTDS->m_apoOverviews[iOverview]->GetRasterBand(nBand ? nBand : 1);
+        if( m_bIsMaskBand )
+            return poOvrBand->GetMaskBand();
+        return poOvrBand;
     }
 
     return nullptr;