Skip to content

Commit 44bb31a

Browse files
rouaultrouault
committed
gdal raster color-map: make it work better with pipelines
Fixes OSGeo#13740
1 parent 195f891 commit 44bb31a

File tree

4 files changed

+93
-144
lines changed

4 files changed

+93
-144
lines changed

‎apps/gdaldem_lib.cpp‎

Lines changed: 67 additions & 144 deletions
Original file line numberDiff line numberDiff line change
@@ -104,6 +104,7 @@
104104
#include "cpl_vsi_virtual.h"
105105
#include "gdal.h"
106106
#include "gdal_priv.h"
107+
#include "vrtdataset.h"
107108

108109
#if defined(__x86_64__) || defined(_M_X64)
109110
#define HAVE_16_SSE_REG
@@ -2204,102 +2205,55 @@ GDALColorRelief(GDALRasterBandH hSrcBand, GDALRasterBandH hDstBand1,
22042205
/* GDALGenerateVRTColorRelief() */
22052206
/************************************************************************/
22062207

2207-
static bool GDALGenerateVRTColorRelief(const char *pszDstFilename,
2208-
GDALDatasetH hSrcDataset,
2209-
GDALRasterBandH hSrcBand,
2210-
const char *pszColorFilename,
2211-
ColorSelectionMode eColorSelectionMode,
2212-
bool bAddAlpha)
2208+
static std::unique_ptr<GDALDataset> GDALGenerateVRTColorRelief(
2209+
const char *pszDest, GDALDatasetH hSrcDataset, GDALRasterBandH hSrcBand,
2210+
const char *pszColorFilename, ColorSelectionMode eColorSelectionMode,
2211+
bool bAddAlpha)
22132212
{
22142213
const auto asColorAssociation = GDALColorReliefParseColorFile(
22152214
hSrcBand, pszColorFilename, eColorSelectionMode);
22162215
if (asColorAssociation.empty())
2217-
return false;
2218-
2219-
VSILFILE *fp = VSIFOpenL(pszDstFilename, "wt");
2220-
if (fp == nullptr)
2221-
{
2222-
return false;
2223-
}
2216+
return nullptr;
22242217

2218+
GDALDataset *poSrcDS = GDALDataset::FromHandle(hSrcDataset);
22252219
const int nXSize = GDALGetRasterBandXSize(hSrcBand);
22262220
const int nYSize = GDALGetRasterBandYSize(hSrcBand);
22272221

2228-
bool bOK =
2229-
VSIFPrintfL(fp, "<VRTDataset rasterXSize=\"%d\" rasterYSize=\"%d\">\n",
2230-
nXSize, nYSize) > 0;
2231-
const char *pszProjectionRef = GDALGetProjectionRef(hSrcDataset);
2232-
if (pszProjectionRef && pszProjectionRef[0] != '\0')
2233-
{
2234-
char *pszEscapedString =
2235-
CPLEscapeString(pszProjectionRef, -1, CPLES_XML);
2236-
bOK &= VSIFPrintfL(fp, " <SRS>%s</SRS>\n", pszEscapedString) > 0;
2237-
VSIFree(pszEscapedString);
2238-
}
2239-
GDALGeoTransform gt;
2240-
if (GDALDataset::FromHandle(hSrcDataset)->GetGeoTransform(gt) == CE_None)
2241-
{
2242-
bOK &= VSIFPrintfL(fp,
2243-
" <GeoTransform> %.16g, %.16g, %.16g, "
2244-
"%.16g, %.16g, %.16g</GeoTransform>\n",
2245-
gt[0], gt[1], gt[2], gt[3], gt[4], gt[5]) > 0;
2246-
}
2247-
22482222
int nBlockXSize = 0;
22492223
int nBlockYSize = 0;
22502224
GDALGetBlockSize(hSrcBand, &nBlockXSize, &nBlockYSize);
22512225

2252-
int bRelativeToVRT = FALSE;
2253-
const CPLString osPath = CPLGetPathSafe(pszDstFilename);
2254-
char *pszSourceFilename = CPLStrdup(CPLExtractRelativePath(
2255-
osPath.c_str(), GDALGetDescription(hSrcDataset), &bRelativeToVRT));
2226+
auto poVRTDS =
2227+
std::make_unique<VRTDataset>(nXSize, nYSize, nBlockXSize, nBlockYSize);
2228+
poVRTDS->SetDescription(pszDest);
2229+
poVRTDS->SetSpatialRef(poSrcDS->GetSpatialRef());
2230+
GDALGeoTransform gt;
2231+
if (poSrcDS->GetGeoTransform(gt) == CE_None)
2232+
{
2233+
poVRTDS->SetGeoTransform(gt);
2234+
}
22562235

22572236
const int nBands = 3 + (bAddAlpha ? 1 : 0);
22582237

22592238
for (int iBand = 0; iBand < nBands; iBand++)
22602239
{
2261-
bOK &=
2262-
VSIFPrintfL(fp, " <VRTRasterBand dataType=\"Byte\" band=\"%d\">\n",
2263-
iBand + 1) > 0;
2264-
bOK &= VSIFPrintfL(
2265-
fp, " <ColorInterp>%s</ColorInterp>\n",
2266-
GDALGetColorInterpretationName(
2267-
static_cast<GDALColorInterp>(GCI_RedBand + iBand))) > 0;
2268-
bOK &= VSIFPrintfL(fp, " <ComplexSource>\n") > 0;
2269-
bOK &= VSIFPrintfL(fp,
2270-
" <SourceFilename "
2271-
"relativeToVRT=\"%d\">%s</SourceFilename>\n",
2272-
bRelativeToVRT, pszSourceFilename) > 0;
2273-
bOK &= VSIFPrintfL(fp, " <SourceBand>%d</SourceBand>\n",
2274-
GDALGetBandNumber(hSrcBand)) > 0;
2275-
bOK &= VSIFPrintfL(fp,
2276-
" <SourceProperties RasterXSize=\"%d\" "
2277-
"RasterYSize=\"%d\" DataType=\"%s\" "
2278-
"BlockXSize=\"%d\" BlockYSize=\"%d\"/>\n",
2279-
nXSize, nYSize,
2280-
GDALGetDataTypeName(GDALGetRasterDataType(hSrcBand)),
2281-
nBlockXSize, nBlockYSize) > 0;
2282-
bOK &=
2283-
VSIFPrintfL(
2284-
fp,
2285-
" "
2286-
"<SrcRect xOff=\"0\" yOff=\"0\" xSize=\"%d\" ySize=\"%d\"/>\n",
2287-
nXSize, nYSize) > 0;
2288-
bOK &=
2289-
VSIFPrintfL(
2290-
fp,
2291-
" "
2292-
"<DstRect xOff=\"0\" yOff=\"0\" xSize=\"%d\" ySize=\"%d\"/>\n",
2293-
nXSize, nYSize) > 0;
2294-
2295-
bOK &= VSIFPrintfL(fp, " <LUT>") > 0;
2240+
poVRTDS->AddBand(GDT_Byte, nullptr);
2241+
auto poVRTBand = cpl::down_cast<VRTSourcedRasterBand *>(
2242+
poVRTDS->GetRasterBand(iBand + 1));
2243+
poVRTBand->SetColorInterpretation(
2244+
static_cast<GDALColorInterp>(GCI_RedBand + iBand));
2245+
2246+
auto poComplexSource = std::make_unique<VRTComplexSource>();
2247+
poVRTBand->ConfigureSource(poComplexSource.get(),
2248+
GDALRasterBand::FromHandle(hSrcBand), FALSE,
2249+
0, 0, nXSize, nYSize, 0, 0, nXSize, nYSize);
2250+
2251+
std::vector<double> adfInputLUT;
2252+
std::vector<double> adfOutputLUT;
22962253

22972254
for (size_t iColor = 0; iColor < asColorAssociation.size(); iColor++)
22982255
{
22992256
const double dfVal = asColorAssociation[iColor].dfVal;
2300-
if (iColor > 0)
2301-
bOK &= VSIFPrintfL(fp, ",") > 0;
2302-
23032257
if (iColor > 0 &&
23042258
eColorSelectionMode == COLOR_SELECTION_NEAREST_ENTRY &&
23052259
dfVal !=
@@ -2308,67 +2262,50 @@ static bool GDALGenerateVRTColorRelief(const char *pszDstFilename,
23082262
{
23092263
const double dfMidVal =
23102264
(dfVal + asColorAssociation[iColor - 1].dfVal) / 2.0;
2311-
bOK &=
2312-
VSIFPrintfL(
2313-
fp, "%.18g:%d",
2314-
std::nextafter(
2315-
dfMidVal, -std::numeric_limits<double>::infinity()),
2316-
(iBand == 0) ? asColorAssociation[iColor - 1].nR
2317-
: (iBand == 1) ? asColorAssociation[iColor - 1].nG
2318-
: (iBand == 2) ? asColorAssociation[iColor - 1].nB
2319-
: asColorAssociation[iColor - 1].nA) > 0;
2320-
bOK &= VSIFPrintfL(
2321-
fp, ",%.18g:%d", dfMidVal,
2322-
(iBand == 0) ? asColorAssociation[iColor].nR
2323-
: (iBand == 1) ? asColorAssociation[iColor].nG
2324-
: (iBand == 2) ? asColorAssociation[iColor].nB
2325-
: asColorAssociation[iColor].nA) > 0;
2265+
adfInputLUT.push_back(std::nextafter(
2266+
dfMidVal, -std::numeric_limits<double>::infinity()));
2267+
adfOutputLUT.push_back(
2268+
(iBand == 0) ? asColorAssociation[iColor - 1].nR
2269+
: (iBand == 1) ? asColorAssociation[iColor - 1].nG
2270+
: (iBand == 2) ? asColorAssociation[iColor - 1].nB
2271+
: asColorAssociation[iColor - 1].nA);
2272+
adfInputLUT.push_back(dfMidVal);
2273+
adfOutputLUT.push_back(
2274+
(iBand == 0) ? asColorAssociation[iColor].nR
2275+
: (iBand == 1) ? asColorAssociation[iColor].nG
2276+
: (iBand == 2) ? asColorAssociation[iColor].nB
2277+
: asColorAssociation[iColor].nA);
23262278
}
23272279
else
23282280
{
23292281
if (eColorSelectionMode == COLOR_SELECTION_EXACT_ENTRY)
23302282
{
2331-
bOK &=
2332-
VSIFPrintfL(
2333-
fp, "%.18g:0,",
2334-
std::nextafter(
2335-
dfVal,
2336-
-std::numeric_limits<double>::infinity())) > 0;
2283+
adfInputLUT.push_back(std::nextafter(
2284+
dfVal, -std::numeric_limits<double>::infinity()));
2285+
adfOutputLUT.push_back(0);
23372286
}
2338-
if (dfVal != static_cast<double>(static_cast<int>(dfVal)))
2339-
bOK &= VSIFPrintfL(fp, "%.18g", dfVal) > 0;
2340-
else
2341-
bOK &= VSIFPrintfL(fp, "%d", static_cast<int>(dfVal)) > 0;
2342-
bOK &= VSIFPrintfL(
2343-
fp, ":%d",
2344-
(iBand == 0) ? asColorAssociation[iColor].nR
2345-
: (iBand == 1) ? asColorAssociation[iColor].nG
2346-
: (iBand == 2) ? asColorAssociation[iColor].nB
2347-
: asColorAssociation[iColor].nA) > 0;
2287+
adfInputLUT.push_back(dfVal);
2288+
adfOutputLUT.push_back(
2289+
(iBand == 0) ? asColorAssociation[iColor].nR
2290+
: (iBand == 1) ? asColorAssociation[iColor].nG
2291+
: (iBand == 2) ? asColorAssociation[iColor].nB
2292+
: asColorAssociation[iColor].nA);
23482293
}
23492294

23502295
if (eColorSelectionMode == COLOR_SELECTION_EXACT_ENTRY)
23512296
{
2352-
bOK &= VSIFPrintfL(
2353-
fp, ",%.18g:0",
2354-
std::nextafter(
2355-
dfVal,
2356-
std::numeric_limits<double>::infinity())) > 0;
2297+
adfInputLUT.push_back(std::nextafter(
2298+
dfVal, std::numeric_limits<double>::infinity()));
2299+
adfOutputLUT.push_back(0);
23572300
}
23582301
}
2359-
bOK &= VSIFPrintfL(fp, "</LUT>\n") > 0;
23602302

2361-
bOK &= VSIFPrintfL(fp, " </ComplexSource>\n") > 0;
2362-
bOK &= VSIFPrintfL(fp, " </VRTRasterBand>\n") > 0;
2363-
}
2303+
poComplexSource->SetLUT(adfInputLUT, adfOutputLUT);
23642304

2365-
CPLFree(pszSourceFilename);
2366-
2367-
bOK &= VSIFPrintfL(fp, "</VRTDataset>\n") > 0;
2368-
2369-
VSIFCloseL(fp);
2305+
poVRTBand->AddSource(std::move(poComplexSource));
2306+
}
23702307

2371-
return bOK;
2308+
return poVRTDS;
23722309
}
23732310

23742311
/************************************************************************/
@@ -4038,31 +3975,17 @@ GDALDatasetH GDALDEMProcessing(const char *pszDest, GDALDatasetH hSrcDataset,
40383975
{
40393976
if (eUtilityMode == COLOR_RELIEF)
40403977
{
4041-
const bool bTmpFile = pszDest[0] == 0;
4042-
const std::string osTmpFile =
4043-
VSIMemGenerateHiddenFilename("tmp.vrt");
4044-
if (bTmpFile)
4045-
pszDest = osTmpFile.c_str();
4046-
GDALDatasetH hDS = nullptr;
4047-
if (GDALGenerateVRTColorRelief(
4048-
pszDest, hSrcDataset, hSrcBand, pszColorFilename,
4049-
psOptions->eColorSelectionMode, psOptions->bAddAlpha))
3978+
auto poDS = GDALGenerateVRTColorRelief(
3979+
pszDest, hSrcDataset, hSrcBand, pszColorFilename,
3980+
psOptions->eColorSelectionMode, psOptions->bAddAlpha);
3981+
if (poDS && pszDest[0] != 0)
40503982
{
4051-
if (bTmpFile)
4052-
{
4053-
const GByte *pabyData =
4054-
VSIGetMemFileBuffer(pszDest, nullptr, false);
4055-
hDS = GDALOpen(reinterpret_cast<const char *>(pabyData),
4056-
GA_Update);
4057-
}
4058-
else
4059-
{
4060-
hDS = GDALOpen(pszDest, GA_Update);
4061-
}
3983+
poDS.reset();
3984+
poDS.reset(GDALDataset::Open(
3985+
pszDest,
3986+
GDAL_OF_UPDATE | GDAL_OF_VERBOSE_ERROR | GDAL_OF_RASTER));
40623987
}
4063-
if (bTmpFile)
4064-
VSIUnlink(pszDest);
4065-
return hDS;
3988+
return GDALDataset::ToHandle(poDS.release());
40663989
}
40673990
else
40683991
{

‎autotest/utilities/test_gdalalg_raster_color_map.py‎

Lines changed: 11 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -138,3 +138,14 @@ def test_gdalalg_raster_empty_color_map(tmp_vsimem, output_format):
138138
alg["output-format"] = output_format
139139
with pytest.raises(Exception, match="No color association found"):
140140
alg.Run()
141+
142+
143+
def test_gdalalg_raster_color_map_in_pipeline(tmp_vsimem):
144+
"""Test fix for https://github.com/OSGeo/gdal/issues/13740"""
145+
146+
with gdal.alg.raster.pipeline(
147+
input="../gcore/data/byte.tif",
148+
pipeline="read ! clip --bbox 440720.000,3750120.000,441920.000,3751320.000 ! color-map --color-map data/color_file.txt ! clip --bbox 440720.000,3750120.000,441920.000,3751320.000 ! write",
149+
output_format="MEM",
150+
) as alg:
151+
assert alg.Output().GetRasterBand(1).Checksum() == 4688

‎frmts/vrt/vrtdataset.h‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1742,6 +1742,9 @@ class CPL_DLL VRTComplexSource CPL_NON_FINAL : public VRTSimpleSource
17421742
void SetPowerScaling(double dfExponent, double dfSrcMin, double dfSrcMax,
17431743
double dfDstMin, double dfDstMax, bool bClip = true);
17441744
void SetColorTableComponent(int nComponent);
1745+
1746+
void SetLUT(const std::vector<double> &adfLUTInputs,
1747+
const std::vector<double> &adfLUTOutputs);
17451748
};
17461749

17471750
/************************************************************************/

‎frmts/vrt/vrtsources.cpp‎

Lines changed: 12 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2947,6 +2947,18 @@ CPLErr VRTComplexSource::XMLInit(const CPLXMLNode *psSrc,
29472947
return CE_None;
29482948
}
29492949

2950+
/************************************************************************/
2951+
/* SetLUT() */
2952+
/************************************************************************/
2953+
2954+
void VRTComplexSource::SetLUT(const std::vector<double> &adfLUTInputs,
2955+
const std::vector<double> &adfLUTOutputs)
2956+
{
2957+
m_adfLUTInputs = adfLUTInputs;
2958+
m_adfLUTOutputs = adfLUTOutputs;
2959+
m_nProcessingFlags |= PROCESSING_FLAG_LUT;
2960+
}
2961+
29502962
/************************************************************************/
29512963
/* LookupValue() */
29522964
/************************************************************************/

0 commit comments

Comments
 (0)