Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL
4 : * Purpose: "gdal vector clean-coverage" subcommand
5 : * Author: Daniel Baston
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2025, ISciences LLC
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "gdalalg_vector_clean_coverage.h"
14 :
15 : #include "cpl_error.h"
16 : #include "gdal_priv.h"
17 : #include "gdalalg_vector_geom.h"
18 : #include "ogr_geos.h"
19 : #include "ogrsf_frmts.h"
20 :
21 : #include <cinttypes>
22 :
23 : #ifndef _
24 : #define _(x) (x)
25 : #endif
26 :
27 : //! @cond Doxygen_Suppress
28 :
29 59 : GDALVectorCleanCoverageAlgorithm::GDALVectorCleanCoverageAlgorithm(
30 59 : bool standaloneStep)
31 : : GDALVectorPipelineStepAlgorithm(NAME, DESCRIPTION, HELP_URL,
32 59 : standaloneStep)
33 : {
34 59 : AddActiveLayerArg(&m_activeLayer);
35 : AddArg("snapping-distance", 0, _("Distance tolerance for snapping nodes"),
36 118 : &m_opts.snappingTolerance)
37 59 : .SetMinValueIncluded(0);
38 :
39 : AddArg("merge-strategy", 0,
40 : _("Algorithm to assign overlaps to neighboring polygons"),
41 118 : &m_opts.mergeStrategy)
42 295 : .SetChoices({"longest-border", "max-area", "min-area", "min-index"});
43 :
44 : AddArg("maximum-gap-width", 0, _("Maximum width of a gap to be closed"),
45 118 : &m_opts.maximumGapWidth)
46 59 : .SetMinValueIncluded(0);
47 59 : }
48 :
49 : #if defined HAVE_GEOS && \
50 : (GEOS_VERSION_MAJOR > 3 || \
51 : (GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR >= 14))
52 :
53 : class GDALVectorCleanCoverageOutputLayer final
54 : : public GDALGeosNonStreamingAlgorithmLayer
55 : {
56 : public:
57 13 : GDALVectorCleanCoverageOutputLayer(
58 : OGRLayer &srcLayer, int geomFieldIndex,
59 : const GDALVectorCleanCoverageAlgorithm::Options &opts)
60 13 : : GDALGeosNonStreamingAlgorithmLayer(srcLayer, geomFieldIndex),
61 13 : m_opts(opts), m_cleanParams(GetCoverageCleanParams())
62 : {
63 13 : }
64 :
65 : ~GDALVectorCleanCoverageOutputLayer() override;
66 :
67 333 : const OGRFeatureDefn *GetLayerDefn() const override
68 : {
69 333 : return m_srcLayer.GetLayerDefn();
70 : }
71 :
72 37 : int TestCapability(const char *pszCap) const override
73 : {
74 37 : if (EQUAL(pszCap, OLCFastFeatureCount))
75 : {
76 0 : return m_srcLayer.TestCapability(pszCap);
77 : }
78 :
79 37 : return false;
80 : }
81 :
82 13 : GEOSCoverageCleanParams *GetCoverageCleanParams() const
83 : {
84 : GEOSCoverageCleanParams *params =
85 13 : GEOSCoverageCleanParams_create_r(m_poGeosContext);
86 :
87 13 : if (!params)
88 : {
89 0 : CPLError(CE_Failure, CPLE_AppDefined,
90 : "Failed to create coverage clean parameters");
91 0 : return nullptr;
92 : }
93 :
94 13 : if (!GEOSCoverageCleanParams_setSnappingDistance_r(
95 13 : m_poGeosContext, params, m_opts.snappingTolerance))
96 : {
97 0 : CPLError(CE_Failure, CPLE_AppDefined,
98 : "Failed to set snapping tolerance");
99 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
100 0 : return nullptr;
101 : }
102 :
103 13 : if (!GEOSCoverageCleanParams_setGapMaximumWidth_r(
104 13 : m_poGeosContext, params, m_opts.maximumGapWidth))
105 : {
106 0 : CPLError(CE_Failure, CPLE_AppDefined,
107 : "Failed to set maximum gap width");
108 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
109 0 : return nullptr;
110 : }
111 :
112 : int mergeStrategy;
113 13 : if (m_opts.mergeStrategy == "longest-border")
114 : {
115 10 : mergeStrategy = GEOS_MERGE_LONGEST_BORDER;
116 : }
117 3 : else if (m_opts.mergeStrategy == "max-area")
118 : {
119 1 : mergeStrategy = GEOS_MERGE_MAX_AREA;
120 : }
121 2 : else if (m_opts.mergeStrategy == "min-area")
122 : {
123 1 : mergeStrategy = GEOS_MERGE_MIN_AREA;
124 : }
125 1 : else if (m_opts.mergeStrategy == "min-index")
126 : {
127 1 : mergeStrategy = GEOS_MERGE_MIN_INDEX;
128 : }
129 : else
130 : {
131 0 : CPLError(CE_Failure, CPLE_AppDefined,
132 : "Unknown overlap merge strategy: %s",
133 0 : m_opts.mergeStrategy.c_str());
134 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
135 0 : return nullptr;
136 : }
137 :
138 13 : if (!GEOSCoverageCleanParams_setOverlapMergeStrategy_r(
139 13 : m_poGeosContext, params, mergeStrategy))
140 : {
141 0 : CPLError(CE_Failure, CPLE_AppDefined,
142 : "Failed to set overlap merge strategy");
143 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
144 0 : return nullptr;
145 : }
146 :
147 13 : return params;
148 : }
149 :
150 47 : bool PolygonsOnly() const override
151 : {
152 47 : return true;
153 : }
154 :
155 34 : bool SkipEmpty() const override
156 : {
157 34 : return false;
158 : }
159 :
160 10 : bool ProcessGeos() override
161 : {
162 : // Perform coverage cleaning
163 10 : GEOSGeometry *coll = GEOSGeom_createCollection_r(
164 : m_poGeosContext, GEOS_GEOMETRYCOLLECTION, m_apoGeosInputs.data(),
165 10 : static_cast<unsigned int>(m_apoGeosInputs.size()));
166 :
167 10 : if (coll == nullptr)
168 : {
169 0 : return false;
170 : }
171 :
172 10 : m_apoGeosInputs.clear();
173 :
174 10 : m_poGeosResultAsCollection =
175 10 : GEOSCoverageCleanWithParams_r(m_poGeosContext, coll, m_cleanParams);
176 10 : GEOSGeom_destroy_r(m_poGeosContext, coll);
177 :
178 10 : return m_poGeosResultAsCollection != nullptr;
179 : }
180 :
181 : CPL_DISALLOW_COPY_ASSIGN(GDALVectorCleanCoverageOutputLayer)
182 :
183 : private:
184 : const GDALVectorCleanCoverageAlgorithm::Options &m_opts;
185 : GEOSCoverageCleanParams *m_cleanParams;
186 : };
187 :
188 26 : GDALVectorCleanCoverageOutputLayer::~GDALVectorCleanCoverageOutputLayer()
189 : {
190 13 : if (m_poGeosContext != nullptr)
191 : {
192 13 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, m_cleanParams);
193 : }
194 26 : }
195 :
196 14 : bool GDALVectorCleanCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &ctxt)
197 : {
198 14 : auto poSrcDS = m_inputDataset[0].GetDatasetRef();
199 28 : auto poDstDS = std::make_unique<GDALVectorNonStreamingAlgorithmDataset>();
200 :
201 28 : GDALVectorAlgorithmLayerProgressHelper progressHelper(ctxt);
202 :
203 30 : for (auto &&poSrcLayer : poSrcDS->GetLayers())
204 : {
205 20 : if (m_activeLayer.empty() ||
206 4 : m_activeLayer == poSrcLayer->GetDescription())
207 : {
208 13 : progressHelper.AddProcessedLayer(*poSrcLayer);
209 : }
210 : else
211 : {
212 3 : progressHelper.AddPassThroughLayer(*poSrcLayer);
213 : }
214 : }
215 :
216 14 : if (!progressHelper.HasProcessedLayers())
217 : {
218 1 : ReportError(CE_Failure, CPLE_AppDefined,
219 : "Specified layer '%s' was not found",
220 : m_activeLayer.c_str());
221 1 : return false;
222 : }
223 :
224 25 : for (auto [poSrcLayer, bProcessed, layerProgressFunc, layerProgressData] :
225 35 : progressHelper)
226 : {
227 14 : if (bProcessed)
228 : {
229 13 : constexpr int geomFieldIndex = 0; // TODO parametrize
230 :
231 : auto poLayer = std::make_unique<GDALVectorCleanCoverageOutputLayer>(
232 13 : *poSrcLayer, geomFieldIndex, m_opts);
233 :
234 13 : if (!poDstDS->AddProcessedLayer(std::move(poLayer),
235 : layerProgressFunc,
236 : layerProgressData.get()))
237 : {
238 3 : return false;
239 : }
240 : }
241 : else
242 : {
243 1 : poDstDS->AddPassThroughLayer(*poSrcLayer);
244 : }
245 : }
246 :
247 10 : m_outputDataset.Set(std::move(poDstDS));
248 :
249 10 : return true;
250 : }
251 :
252 : #else
253 :
254 : bool GDALVectorCleanCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &)
255 : {
256 : ReportError(CE_Failure, CPLE_AppDefined,
257 : "%s requires GDAL to be built against version 3.14 or later of "
258 : "the GEOS library.",
259 : NAME);
260 : return false;
261 : }
262 : #endif // HAVE_GEOS
263 :
264 : GDALVectorCleanCoverageAlgorithmStandalone::
265 : ~GDALVectorCleanCoverageAlgorithmStandalone() = default;
266 :
267 : //! @endcond
|
