Skip to content

Commit

Permalink
gdal raster mosaic/stack: turn gdalbuildvrt -strict mode
Browse files Browse the repository at this point in the history
  • Loading branch information
rouault committed Jan 27, 2025
1 parent fcafee9 commit 4f4c0a7
Show file tree
Hide file tree
Showing 4 changed files with 108 additions and 52 deletions.
5 changes: 5 additions & 0 deletions apps/gdalalg_raster_mosaic.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -185,6 +185,11 @@ bool GDALRasterMosaicAlgorithm::RunImpl(GDALProgressFunc pfnProgress,
"VRT");

CPLStringList aosOptions;
aosOptions.push_back("-strict");

aosOptions.push_back("-program_name");
aosOptions.push_back("gdal raster mosaic");

if (!m_resolution.empty())
{
const auto aosTokens =
Expand Down
5 changes: 5 additions & 0 deletions apps/gdalalg_raster_stack.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -181,6 +181,11 @@ bool GDALRasterStackAlgorithm::RunImpl(GDALProgressFunc pfnProgress,

CPLStringList aosOptions;

aosOptions.push_back("-strict");

aosOptions.push_back("-program_name");
aosOptions.push_back("gdal raster stack");

aosOptions.push_back("-separate");

if (!m_resolution.empty())
Expand Down
128 changes: 77 additions & 51 deletions apps/gdalbuildvrt_lib.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -290,6 +290,8 @@ class VRTBuilder
~VRTBuilder();

GDALDataset *Build(GDALProgressFunc pfnProgress, void *pProgressData);

std::string m_osProgramName{};
};

/************************************************************************/
Expand Down Expand Up @@ -514,42 +516,48 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
int bGotGeoTransform = poDS->GetGeoTransform(padfGeoTransform) == CE_None;
if (bSeparate)
{
std::string osProgramName(m_osProgramName);
if (osProgramName == "gdalbuildvrt")
osProgramName += " -separate";

if (bFirst)
{
bHasGeoTransform = bGotGeoTransform;
if (!bHasGeoTransform)
{
if (bUserExtent)
{
CPLError(CE_Warning, CPLE_NotSupported,
"User extent ignored by gdalbuildvrt -separate "
"with ungeoreferenced images.");
CPLError(CE_Warning, CPLE_NotSupported, "%s",
("User extent ignored by " + osProgramName +
"with ungeoreferenced images.")
.c_str());
}
if (resolutionStrategy == USER_RESOLUTION)
{
CPLError(CE_Warning, CPLE_NotSupported,
"User resolution ignored by gdalbuildvrt "
"-separate with ungeoreferenced images.");
CPLError(CE_Warning, CPLE_NotSupported, "%s",
("User resolution ignored by " + osProgramName +
" with ungeoreferenced images.")
.c_str());
}
}
}
else if (bHasGeoTransform != bGotGeoTransform)
{
return "gdalbuildvrt -separate cannot stack ungeoreferenced and "
"georeferenced images.";
return osProgramName + " cannot stack ungeoreferenced and "
"georeferenced images.";
}
else if (!bHasGeoTransform && (nRasterXSize != poDS->GetRasterXSize() ||
nRasterYSize != poDS->GetRasterYSize()))
{
return "gdalbuildvrt -separate cannot stack ungeoreferenced images "
"that have not the same dimensions.";
return osProgramName + " cannot stack ungeoreferenced images "
"that have not the same dimensions.";
}
}
else
{
if (!bGotGeoTransform)
{
return "gdalbuildvrt does not support ungeoreferenced image.";
return m_osProgramName + " does not support ungeoreferenced image.";
}
bHasGeoTransform = TRUE;
}
Expand All @@ -559,11 +567,13 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
if (padfGeoTransform[GEOTRSFRM_ROTATION_PARAM1] != 0 ||
padfGeoTransform[GEOTRSFRM_ROTATION_PARAM2] != 0)
{
return "gdalbuildvrt does not support rotated geo transforms.";
return m_osProgramName +
" does not support rotated geo transforms.";
}
if (padfGeoTransform[GEOTRSFRM_NS_RES] >= 0)
{
return "gdalbuildvrt does not support positive NS resolution.";
return m_osProgramName +
" does not support positive NS resolution.";
}
}

Expand Down Expand Up @@ -800,7 +810,8 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
{
CPLString osExpected = GetProjectionName(pszProjectionRef);
CPLString osGot = GetProjectionName(proj);
return CPLSPrintf("gdalbuildvrt does not support heterogeneous "
return m_osProgramName +
CPLSPrintf(" does not support heterogeneous "
"projection: expected %s, got %s.",
osExpected.c_str(), osGot.c_str());
}
Expand All @@ -816,18 +827,18 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
}
else
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous band "
"numbers: expected %d, got %d.",
nTotalBands, _nBands);
return m_osProgramName +
CPLSPrintf(" does not support heterogeneous band "
"numbers: expected %d, got %d.",
nTotalBands, _nBands);
}
}
else if (bExplicitBandList && _nBands < nMaxSelectedBandNo)
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous band "
"numbers: expected at least %d, got %d.",
nMaxSelectedBandNo, _nBands);
return m_osProgramName +
CPLSPrintf(" does not support heterogeneous band "
"numbers: expected at least %d, got %d.",
nMaxSelectedBandNo, _nBands);
}

for (int j = 0; j < nSelectedBands; j++)
Expand All @@ -838,21 +849,25 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
if (asBandProperties[j].colorInterpretation !=
poBand->GetColorInterpretation())
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous "
"band color interpretation: expected %s, got %s.",
GDALGetColorInterpretationName(
asBandProperties[j].colorInterpretation),
GDALGetColorInterpretationName(
poBand->GetColorInterpretation()));
return m_osProgramName +
CPLSPrintf(
" does not support heterogeneous "
"band color interpretation: expected %s, got "
"%s.",
GDALGetColorInterpretationName(
asBandProperties[j].colorInterpretation),
GDALGetColorInterpretationName(
poBand->GetColorInterpretation()));
}
if (asBandProperties[j].dataType != poBand->GetRasterDataType())
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous "
"band data type: expected %s, got %s.",
GDALGetDataTypeName(asBandProperties[j].dataType),
GDALGetDataTypeName(poBand->GetRasterDataType()));
return m_osProgramName +
CPLSPrintf(" does not support heterogeneous "
"band data type: expected %s, got %s.",
GDALGetDataTypeName(
asBandProperties[j].dataType),
GDALGetDataTypeName(
poBand->GetRasterDataType()));
}
if (asBandProperties[j].colorTable)
{
Expand All @@ -862,7 +877,8 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
if (colorTable == nullptr ||
colorTable->GetColorEntryCount() != nRefColorEntryCount)
{
return "gdalbuildvrt does not support rasters with "
return m_osProgramName +
" does not support rasters with "
"different color tables (different number of "
"color table entries)";
}
Expand Down Expand Up @@ -894,9 +910,9 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
"You're advised to pre-process your "
"rasters with other tools, such as "
"pct2rgb.py or gdal_translate -expand RGB\n"
"to operate gdalbuildvrt on RGB rasters "
"to operate %s on RGB rasters "
"instead",
dsFileName);
dsFileName, m_osProgramName.c_str());
else
CPLError(CE_Warning, CPLE_NotSupported,
"%s has different values than the "
Expand All @@ -915,13 +931,15 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
psDatasetProperties->adfOffset[j] !=
asBandProperties[j].dfOffset))
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous "
"band offset: expected (%d,%f), got (%d,%f).",
static_cast<int>(asBandProperties[j].bHasOffset),
asBandProperties[j].dfOffset,
static_cast<int>(psDatasetProperties->abHasOffset[j]),
psDatasetProperties->adfOffset[j]);
return m_osProgramName +
CPLSPrintf(
" does not support heterogeneous "
"band offset: expected (%d,%f), got (%d,%f).",
static_cast<int>(asBandProperties[j].bHasOffset),
asBandProperties[j].dfOffset,
static_cast<int>(
psDatasetProperties->abHasOffset[j]),
psDatasetProperties->adfOffset[j]);
}

if (psDatasetProperties->abHasScale[j] !=
Expand All @@ -930,13 +948,15 @@ std::string VRTBuilder::AnalyseRaster(GDALDatasetH hDS,
psDatasetProperties->adfScale[j] !=
asBandProperties[j].dfScale))
{
return CPLSPrintf(
"gdalbuildvrt does not support heterogeneous "
"band scale: expected (%d,%f), got (%d,%f).",
static_cast<int>(asBandProperties[j].bHasScale),
asBandProperties[j].dfScale,
static_cast<int>(psDatasetProperties->abHasScale[j]),
psDatasetProperties->adfScale[j]);
return m_osProgramName +
CPLSPrintf(
" does not support heterogeneous "
"band scale: expected (%d,%f), got (%d,%f).",
static_cast<int>(asBandProperties[j].bHasScale),
asBandProperties[j].dfScale,
static_cast<int>(
psDatasetProperties->abHasScale[j]),
psDatasetProperties->adfScale[j]);
}
}
}
Expand Down Expand Up @@ -1752,6 +1772,7 @@ static bool add_file_to_list(const char *filename, const char *tile_index,
*/
struct GDALBuildVRTOptions
{
std::string osProgramName = "gdalbuildvrt";
std::string osTileIndex = "location";
bool bStrict = false;
std::string osResolution{};
Expand Down Expand Up @@ -1925,6 +1946,7 @@ GDALDatasetH GDALBuildVRT(const char *pszDest, int nSrcCount,
sOptions.osOutputSRS.empty() ? nullptr : sOptions.osOutputSRS.c_str(),
sOptions.osResampling.empty() ? nullptr : sOptions.osResampling.c_str(),
sOptions.aosOpenOptions.List(), sOptions.aosCreateOptions);
oBuilder.m_osProgramName = sOptions.osProgramName;

return GDALDataset::ToHandle(
oBuilder.Build(sOptions.pfnProgress, sOptions.pProgressData));
Expand Down Expand Up @@ -2195,6 +2217,10 @@ GDALBuildVRTOptionsGetParser(GDALBuildVRTOptions *psOptions,
"when the value of the mask band of the source is less or "
"equal to the threshold."));

argParser->add_argument("-program_name")
.store_into(psOptions->osProgramName)
.hidden();

if (psOptionsForBinary)
{
if (psOptionsForBinary->osDstFilename.empty())
Expand Down
22 changes: 21 additions & 1 deletion autotest/utilities/test_gdalalg_raster_mosaic.py
Original file line number Diff line number Diff line change
Expand Up @@ -13,7 +13,7 @@

import pytest

from osgeo import gdal
from osgeo import gdal, osr


def get_mosaic_alg():
Expand Down Expand Up @@ -373,3 +373,23 @@ def test_gdalalg_raster_mosaic_tif_creation_options(tmp_vsimem):
with gdal.Open(out_filename) as ds:
assert ds.GetRasterBand(1).Checksum() == 50054
assert ds.GetRasterBand(1).GetBlockSize() == [256, 256]


def test_gdalalg_raster_mosaic_inconsistent_characteristics():

src1_ds = gdal.GetDriverByName("MEM").Create("", 1, 1)
src1_ds.SetGeoTransform([2, 1, 0, 49, 0, -1])
src2_ds = gdal.GetDriverByName("MEM").Create("", 2, 2)
src2_ds.SetGeoTransform([3, 0.5, 0, 49, 0, -0.5])
srs = osr.SpatialReference()
srs.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
srs.SetFromUserInput("WGS84")
src2_ds.SetSpatialRef(srs)

alg = get_mosaic_alg()
alg.GetArg("input").Set([src1_ds, src2_ds])
alg.GetArg("output").Set("")
with pytest.raises(
Exception, match="gdal raster mosaic does not support heterogeneous projection"
):
assert alg.Run()

0 comments on commit 4f4c0a7

Please sign in to comment.