Skip to content

Commit 25b46ab

Browse files
authored
Merge pull request OSGeo#11634 from rouault/gdal_raster_clip
Add 'gdal raster clip' and 'gdal raster pipeline read ... ! clip ... ! write ...'
2 parents a66da76 + 0f31531 commit 25b46ab

11 files changed

+508
-0
lines changed

apps/CMakeLists.txt

+1
Original file line numberDiff line numberDiff line change
@@ -13,6 +13,7 @@ add_library(
1313
gdalalg_raster.cpp
1414
gdalalg_raster_buildvrt.cpp
1515
gdalalg_raster_info.cpp
16+
gdalalg_raster_clip.cpp
1617
gdalalg_raster_convert.cpp
1718
gdalalg_raster_edit.cpp
1819
gdalalg_raster_pipeline.cpp

apps/gdalalg_raster.cpp

+2
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,7 @@
1414

1515
#include "gdalalg_raster_buildvrt.h"
1616
#include "gdalalg_raster_info.h"
17+
#include "gdalalg_raster_clip.h"
1718
#include "gdalalg_raster_convert.h"
1819
#include "gdalalg_raster_edit.h"
1920
#include "gdalalg_raster_overview.h"
@@ -40,6 +41,7 @@ class GDALRasterAlgorithm final : public GDALAlgorithm
4041
{
4142
RegisterSubAlgorithm<GDALRasterInfoAlgorithm>();
4243
RegisterSubAlgorithm<GDALRasterConvertAlgorithm>();
44+
RegisterSubAlgorithm<GDALRasterClipAlgorithmStandalone>();
4345
RegisterSubAlgorithm<GDALRasterEditAlgorithmStandalone>();
4446
RegisterSubAlgorithm<GDALRasterOverviewAlgorithm>();
4547
RegisterSubAlgorithm<GDALRasterPipelineAlgorithm>();

apps/gdalalg_raster_clip.cpp

+145
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,145 @@
1+
/******************************************************************************
2+
*
3+
* Project: GDAL
4+
* Purpose: "clip" step of "raster pipeline", or "gdal raster clip" standalone
5+
* Author: Even Rouault <even dot rouault at spatialys.com>
6+
*
7+
******************************************************************************
8+
* Copyright (c) 2025, Even Rouault <even dot rouault at spatialys.com>
9+
*
10+
* SPDX-License-Identifier: MIT
11+
****************************************************************************/
12+
13+
#include "gdalalg_raster_clip.h"
14+
15+
#include "gdal_priv.h"
16+
#include "gdal_utils.h"
17+
18+
#include <algorithm>
19+
20+
//! @cond Doxygen_Suppress
21+
22+
#ifndef _
23+
#define _(x) (x)
24+
#endif
25+
26+
/************************************************************************/
27+
/* GDALRasterClipAlgorithm::GDALRasterClipAlgorithm() */
28+
/************************************************************************/
29+
30+
GDALRasterClipAlgorithm::GDALRasterClipAlgorithm(bool standaloneStep)
31+
: GDALRasterPipelineStepAlgorithm(NAME, DESCRIPTION, HELP_URL,
32+
standaloneStep)
33+
{
34+
AddBBOXArg(&m_bbox, _("Clipping bounding box as xmin,ymin,xmax,ymax"))
35+
.SetMutualExclusionGroup("exclusion-group");
36+
AddArg("bbox-crs", 0, _("CRS of clipping bounding box"), &m_bboxCrs)
37+
.SetIsCRSArg()
38+
.AddHiddenAlias("bbox_srs");
39+
AddArg("like", 0, _("Raster dataset to use as a template for bounds"),
40+
&m_likeDataset, GDAL_OF_RASTER)
41+
.SetMetaVar("DATASET")
42+
.SetMutualExclusionGroup("exclusion-group");
43+
}
44+
45+
/************************************************************************/
46+
/* GDALRasterClipAlgorithm::RunStep() */
47+
/************************************************************************/
48+
49+
bool GDALRasterClipAlgorithm::RunStep(GDALProgressFunc, void *)
50+
{
51+
CPLAssert(m_inputDataset.GetDatasetRef());
52+
CPLAssert(m_outputDataset.GetName().empty());
53+
CPLAssert(!m_outputDataset.GetDatasetRef());
54+
55+
CPLStringList aosOptions;
56+
aosOptions.AddString("-of");
57+
aosOptions.AddString("VRT");
58+
if (!m_bbox.empty())
59+
{
60+
aosOptions.AddString("-projwin");
61+
aosOptions.AddString(CPLSPrintf("%.17g", m_bbox[0])); // minx
62+
aosOptions.AddString(CPLSPrintf("%.17g", m_bbox[3])); // maxy
63+
aosOptions.AddString(CPLSPrintf("%.17g", m_bbox[2])); // maxx
64+
aosOptions.AddString(CPLSPrintf("%.17g", m_bbox[1])); // miny
65+
66+
if (!m_bboxCrs.empty())
67+
{
68+
aosOptions.AddString("-projwin_srs");
69+
aosOptions.AddString(m_bboxCrs.c_str());
70+
}
71+
}
72+
else if (auto poLikeDS = m_likeDataset.GetDatasetRef())
73+
{
74+
double adfGT[6];
75+
if (poLikeDS->GetGeoTransform(adfGT) != CE_None)
76+
{
77+
ReportError(CE_Failure, CPLE_AppDefined,
78+
"Dataset '%s' has no geotransform matrix. Its bounds "
79+
"cannot be established.",
80+
poLikeDS->GetDescription());
81+
return false;
82+
}
83+
auto poLikeSRS = poLikeDS->GetSpatialRef();
84+
if (!poLikeSRS)
85+
{
86+
ReportError(CE_Failure, CPLE_AppDefined,
87+
"Dataset '%s' has no SRS. Its bounds cannot be used.",
88+
poLikeDS->GetDescription());
89+
return false;
90+
}
91+
const double dfTLX = adfGT[0];
92+
const double dfTLY = adfGT[3];
93+
const double dfTRX = adfGT[0] + poLikeDS->GetRasterXSize() * adfGT[1];
94+
const double dfTRY = adfGT[3] + poLikeDS->GetRasterXSize() * adfGT[4];
95+
const double dfBLX = adfGT[0] + poLikeDS->GetRasterYSize() * adfGT[2];
96+
const double dfBLY = adfGT[3] + poLikeDS->GetRasterYSize() * adfGT[5];
97+
const double dfBRX = adfGT[0] + poLikeDS->GetRasterXSize() * adfGT[1] +
98+
poLikeDS->GetRasterYSize() * adfGT[2];
99+
const double dfBRY = adfGT[3] + poLikeDS->GetRasterXSize() * adfGT[4] +
100+
poLikeDS->GetRasterYSize() * adfGT[5];
101+
const double dfMinX =
102+
std::min(std::min(dfTLX, dfTRX), std::min(dfBLX, dfBRX));
103+
const double dfMinY =
104+
std::min(std::min(dfTLY, dfTRY), std::min(dfBLY, dfBRY));
105+
const double dfMaxX =
106+
std::max(std::max(dfTLX, dfTRX), std::max(dfBLX, dfBRX));
107+
const double dfMaxY =
108+
std::max(std::max(dfTLY, dfTRY), std::max(dfBLY, dfBRY));
109+
110+
aosOptions.AddString("-projwin");
111+
aosOptions.AddString(CPLSPrintf("%.17g", dfMinX));
112+
aosOptions.AddString(CPLSPrintf("%.17g", dfMaxY));
113+
aosOptions.AddString(CPLSPrintf("%.17g", dfMaxX));
114+
aosOptions.AddString(CPLSPrintf("%.17g", dfMinY));
115+
116+
const char *const apszOptions[] = {"FORMAT=WKT2", nullptr};
117+
const std::string osWKT = poLikeSRS->exportToWkt(apszOptions);
118+
aosOptions.AddString("-projwin_srs");
119+
aosOptions.AddString(osWKT.c_str());
120+
}
121+
else
122+
{
123+
ReportError(CE_Failure, CPLE_AppDefined,
124+
"Either --bbox or --like must be specified");
125+
return false;
126+
}
127+
128+
GDALTranslateOptions *psOptions =
129+
GDALTranslateOptionsNew(aosOptions.List(), nullptr);
130+
131+
GDALDatasetH hSrcDS = GDALDataset::ToHandle(m_inputDataset.GetDatasetRef());
132+
auto poRetDS =
133+
GDALDataset::FromHandle(GDALTranslate("", hSrcDS, psOptions, nullptr));
134+
GDALTranslateOptionsFree(psOptions);
135+
136+
const bool bOK = poRetDS != nullptr;
137+
if (bOK)
138+
{
139+
m_outputDataset.Set(std::unique_ptr<GDALDataset>(poRetDS));
140+
}
141+
142+
return bOK;
143+
}
144+
145+
//! @endcond

apps/gdalalg_raster_clip.h

+62
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,62 @@
1+
/******************************************************************************
2+
*
3+
* Project: GDAL
4+
* Purpose: "clip" step of "raster pipeline", or "gdal raster clip" standalone
5+
* Author: Even Rouault <even dot rouault at spatialys.com>
6+
*
7+
******************************************************************************
8+
* Copyright (c) 2025, Even Rouault <even dot rouault at spatialys.com>
9+
*
10+
* SPDX-License-Identifier: MIT
11+
****************************************************************************/
12+
13+
#ifndef GDALALG_RASTER_CLIP_INCLUDED
14+
#define GDALALG_RASTER_CLIP_INCLUDED
15+
16+
#include "gdalalg_raster_pipeline.h"
17+
18+
//! @cond Doxygen_Suppress
19+
20+
/************************************************************************/
21+
/* GDALRasterClipAlgorithm */
22+
/************************************************************************/
23+
24+
class GDALRasterClipAlgorithm /* non final */
25+
: public GDALRasterPipelineStepAlgorithm
26+
{
27+
public:
28+
static constexpr const char *NAME = "clip";
29+
static constexpr const char *DESCRIPTION = "Clip a raster dataset.";
30+
static constexpr const char *HELP_URL = "/programs/gdal_raster_clip.html";
31+
32+
static std::vector<std::string> GetAliases()
33+
{
34+
return {};
35+
}
36+
37+
explicit GDALRasterClipAlgorithm(bool standaloneStep = false);
38+
39+
private:
40+
bool RunStep(GDALProgressFunc pfnProgress, void *pProgressData) override;
41+
42+
std::vector<double> m_bbox{};
43+
std::string m_bboxCrs{};
44+
GDALArgDatasetValue m_likeDataset{};
45+
};
46+
47+
/************************************************************************/
48+
/* GDALRasterClipAlgorithmStandalone */
49+
/************************************************************************/
50+
51+
class GDALRasterClipAlgorithmStandalone final : public GDALRasterClipAlgorithm
52+
{
53+
public:
54+
GDALRasterClipAlgorithmStandalone()
55+
: GDALRasterClipAlgorithm(/* standaloneStep = */ true)
56+
{
57+
}
58+
};
59+
60+
//! @endcond
61+
62+
#endif /* GDALALG_RASTER_CLIP_INCLUDED */

apps/gdalalg_raster_pipeline.cpp

+2
Original file line numberDiff line numberDiff line change
@@ -12,6 +12,7 @@
1212

1313
#include "gdalalg_raster_pipeline.h"
1414
#include "gdalalg_raster_read.h"
15+
#include "gdalalg_raster_clip.h"
1516
#include "gdalalg_raster_edit.h"
1617
#include "gdalalg_raster_reproject.h"
1718
#include "gdalalg_raster_write.h"
@@ -161,6 +162,7 @@ GDALRasterPipelineAlgorithm::GDALRasterPipelineAlgorithm(
161162

162163
m_stepRegistry.Register<GDALRasterReadAlgorithm>();
163164
m_stepRegistry.Register<GDALRasterWriteAlgorithm>();
165+
m_stepRegistry.Register<GDALRasterClipAlgorithm>();
164166
m_stepRegistry.Register<GDALRasterEditAlgorithm>();
165167
m_stepRegistry.Register<GDALRasterReprojectAlgorithm>();
166168
}

0 commit comments

Comments
 (0)