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_geometry.h"
19 : #include "ogr_geos.h"
20 : #include "ogrsf_frmts.h"
21 :
22 : #include <cinttypes>
23 :
24 : #ifndef _
25 : #define _(x) (x)
26 : #endif
27 :
28 : //! @cond Doxygen_Suppress
29 :
30 63 : GDALVectorCleanCoverageAlgorithm::GDALVectorCleanCoverageAlgorithm(
31 63 : bool standaloneStep)
32 : : GDALVectorPipelineStepAlgorithm(NAME, DESCRIPTION, HELP_URL,
33 63 : standaloneStep)
34 : {
35 63 : AddActiveLayerArg(&m_activeLayer);
36 : AddArg("snapping-distance", 0, _("Distance tolerance for snapping nodes"),
37 126 : &m_opts.snappingTolerance)
38 63 : .SetMinValueIncluded(0);
39 :
40 : AddArg("merge-strategy", 0,
41 : _("Algorithm to assign overlaps to neighboring polygons"),
42 126 : &m_opts.mergeStrategy)
43 315 : .SetChoices({"longest-border", "max-area", "min-area", "min-index"});
44 :
45 : AddArg("maximum-gap-width", 0, _("Maximum width of a gap to be closed"),
46 126 : &m_opts.maximumGapWidth)
47 63 : .SetMinValueIncluded(0);
48 63 : }
49 :
50 : #if defined HAVE_GEOS && \
51 : (GEOS_VERSION_MAJOR > 3 || \
52 : (GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR >= 14))
53 :
54 : class GDALVectorCleanCoverageOutputDataset final
55 : : public GDALGeosNonStreamingAlgorithmDataset
56 : {
57 : public:
58 14 : GDALVectorCleanCoverageOutputDataset(
59 : const GDALVectorCleanCoverageAlgorithm::Options &opts)
60 14 : : m_opts(opts), m_cleanParams(GetCoverageCleanParams())
61 : {
62 14 : }
63 :
64 : ~GDALVectorCleanCoverageOutputDataset() override;
65 :
66 14 : GEOSCoverageCleanParams *GetCoverageCleanParams() const
67 : {
68 : GEOSCoverageCleanParams *params =
69 14 : GEOSCoverageCleanParams_create_r(m_poGeosContext);
70 :
71 14 : if (!params)
72 : {
73 0 : CPLError(CE_Failure, CPLE_AppDefined,
74 : "Failed to create coverage clean parameters");
75 0 : return nullptr;
76 : }
77 :
78 14 : if (!GEOSCoverageCleanParams_setSnappingDistance_r(
79 14 : m_poGeosContext, params, m_opts.snappingTolerance))
80 : {
81 0 : CPLError(CE_Failure, CPLE_AppDefined,
82 : "Failed to set snapping tolerance");
83 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
84 0 : return nullptr;
85 : }
86 :
87 14 : if (!GEOSCoverageCleanParams_setGapMaximumWidth_r(
88 14 : m_poGeosContext, params, m_opts.maximumGapWidth))
89 : {
90 0 : CPLError(CE_Failure, CPLE_AppDefined,
91 : "Failed to set maximum gap width");
92 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
93 0 : return nullptr;
94 : }
95 :
96 : int mergeStrategy;
97 14 : if (m_opts.mergeStrategy == "longest-border")
98 : {
99 11 : mergeStrategy = GEOS_MERGE_LONGEST_BORDER;
100 : }
101 3 : else if (m_opts.mergeStrategy == "max-area")
102 : {
103 1 : mergeStrategy = GEOS_MERGE_MAX_AREA;
104 : }
105 2 : else if (m_opts.mergeStrategy == "min-area")
106 : {
107 1 : mergeStrategy = GEOS_MERGE_MIN_AREA;
108 : }
109 1 : else if (m_opts.mergeStrategy == "min-index")
110 : {
111 1 : mergeStrategy = GEOS_MERGE_MIN_INDEX;
112 : }
113 : else
114 : {
115 0 : CPLError(CE_Failure, CPLE_AppDefined,
116 : "Unknown overlap merge strategy: %s",
117 0 : m_opts.mergeStrategy.c_str());
118 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
119 0 : return nullptr;
120 : }
121 :
122 14 : if (!GEOSCoverageCleanParams_setOverlapMergeStrategy_r(
123 14 : m_poGeosContext, params, mergeStrategy))
124 : {
125 0 : CPLError(CE_Failure, CPLE_AppDefined,
126 : "Failed to set overlap merge strategy");
127 0 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, params);
128 0 : return nullptr;
129 : }
130 :
131 14 : return params;
132 : }
133 :
134 47 : bool PolygonsOnly() const override
135 : {
136 47 : return true;
137 : }
138 :
139 44 : bool SkipEmpty() const override
140 : {
141 44 : return false;
142 : }
143 :
144 10 : bool ProcessGeos() override
145 : {
146 : // Perform coverage cleaning
147 10 : GEOSGeometry *coll = GEOSGeom_createCollection_r(
148 : m_poGeosContext, GEOS_GEOMETRYCOLLECTION, m_apoGeosInputs.data(),
149 10 : static_cast<unsigned int>(m_apoGeosInputs.size()));
150 :
151 10 : if (coll == nullptr)
152 : {
153 0 : return false;
154 : }
155 :
156 10 : m_apoGeosInputs.clear();
157 :
158 10 : m_poGeosResultAsCollection =
159 10 : GEOSCoverageCleanWithParams_r(m_poGeosContext, coll, m_cleanParams);
160 10 : GEOSGeom_destroy_r(m_poGeosContext, coll);
161 :
162 10 : return m_poGeosResultAsCollection != nullptr;
163 : }
164 :
165 : CPL_DISALLOW_COPY_ASSIGN(GDALVectorCleanCoverageOutputDataset)
166 :
167 : private:
168 : const GDALVectorCleanCoverageAlgorithm::Options &m_opts;
169 : GEOSCoverageCleanParams *m_cleanParams;
170 : };
171 :
172 28 : GDALVectorCleanCoverageOutputDataset::~GDALVectorCleanCoverageOutputDataset()
173 : {
174 14 : if (m_poGeosContext != nullptr)
175 : {
176 14 : GEOSCoverageCleanParams_destroy_r(m_poGeosContext, m_cleanParams);
177 : }
178 28 : }
179 :
180 14 : bool GDALVectorCleanCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &)
181 : {
182 14 : auto poSrcDS = m_inputDataset[0].GetDatasetRef();
183 : auto poDstDS =
184 28 : std::make_unique<GDALVectorCleanCoverageOutputDataset>(m_opts);
185 :
186 14 : bool bFoundActiveLayer = false;
187 :
188 27 : for (auto &&poSrcLayer : poSrcDS->GetLayers())
189 : {
190 20 : if (m_activeLayer.empty() ||
191 4 : m_activeLayer == poSrcLayer->GetDescription())
192 : {
193 13 : if (!poDstDS->AddProcessedLayer(*poSrcLayer))
194 : {
195 3 : return false;
196 : }
197 10 : bFoundActiveLayer = true;
198 : }
199 : else
200 : {
201 3 : poDstDS->AddPassThroughLayer(*poSrcLayer);
202 : }
203 : }
204 :
205 11 : if (!bFoundActiveLayer)
206 : {
207 1 : ReportError(CE_Failure, CPLE_AppDefined,
208 : "Specified layer '%s' was not found",
209 : m_activeLayer.c_str());
210 1 : return false;
211 : }
212 :
213 10 : m_outputDataset.Set(std::move(poDstDS));
214 :
215 10 : return true;
216 : }
217 :
218 : #else
219 :
220 : bool GDALVectorCleanCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &)
221 : {
222 : ReportError(CE_Failure, CPLE_AppDefined,
223 : "%s requires GDAL to be built against version 3.14 or later of "
224 : "the GEOS library.",
225 : NAME);
226 : return false;
227 : }
228 : #endif // HAVE_GEOS
229 :
230 : GDALVectorCleanCoverageAlgorithmStandalone::
231 : ~GDALVectorCleanCoverageAlgorithmStandalone() = default;
232 :
233 : //! @endcond
|