Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL
4 : * Purpose: "gdal vector check-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_check_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 :
21 : #include <cinttypes>
22 :
23 : #ifndef _
24 : #define _(x) (x)
25 : #endif
26 :
27 : //! @cond Doxygen_Suppress
28 :
29 54 : GDALVectorCheckCoverageAlgorithm::GDALVectorCheckCoverageAlgorithm(
30 54 : bool standaloneStep)
31 : : GDALVectorPipelineStepAlgorithm(
32 : NAME, DESCRIPTION, HELP_URL,
33 0 : ConstructorOptions()
34 54 : .SetStandaloneStep(standaloneStep)
35 108 : .SetNoCreateEmptyLayersArgument(standaloneStep))
36 : {
37 : AddArg("include-valid", 0,
38 : _("Include valid inputs in output, with empty geometry"),
39 54 : &m_includeValid);
40 :
41 : AddArg("geometry-field", 0, _("Name of geometry field to check"),
42 54 : &m_geomField);
43 :
44 : AddArg("maximum-gap-width", 0, _("Maximum width of a gap to be flagged"),
45 108 : &m_maximumGapWidth)
46 54 : .SetMinValueIncluded(0);
47 54 : }
48 :
49 : #if defined HAVE_GEOS && \
50 : (GEOS_VERSION_MAJOR > 3 || \
51 : (GEOS_VERSION_MAJOR == 3 && GEOS_VERSION_MINOR >= 12))
52 :
53 : class GDALVectorCheckCoverageOutputLayer final
54 : : public GDALGeosNonStreamingAlgorithmLayer
55 : {
56 : public:
57 14 : explicit GDALVectorCheckCoverageOutputLayer(OGRLayer &srcLayer,
58 : int geomFieldIndex,
59 : const std::string &name,
60 : double maximumGapWidth,
61 : bool includeValid)
62 14 : : GDALGeosNonStreamingAlgorithmLayer(srcLayer, geomFieldIndex),
63 14 : m_defn(OGRFeatureDefn::CreateFeatureDefn(name.c_str())),
64 14 : m_maximumGapWidth(maximumGapWidth), m_includeValid(includeValid)
65 : {
66 14 : m_defn->Reference();
67 :
68 14 : const OGRFeatureDefn *poSrcLayerDefn = srcLayer.GetLayerDefn();
69 14 : m_defn->SetGeomType(wkbMultiLineString);
70 28 : m_defn->GetGeomFieldDefn(0)->SetSpatialRef(
71 14 : poSrcLayerDefn->GetGeomFieldDefn(geomFieldIndex)->GetSpatialRef());
72 14 : }
73 :
74 : ~GDALVectorCheckCoverageOutputLayer() override;
75 :
76 297 : const OGRFeatureDefn *GetLayerDefn() const override
77 : {
78 297 : return m_defn;
79 : }
80 :
81 12 : int TestCapability(const char *) const override
82 : {
83 12 : return false;
84 : }
85 :
86 84 : bool PolygonsOnly() const override
87 : {
88 84 : return true;
89 : }
90 :
91 67 : bool SkipEmpty() const override
92 : {
93 67 : return !m_includeValid;
94 : }
95 :
96 14 : bool ProcessGeos() override
97 : {
98 : // Perform coverage checking
99 14 : GEOSGeometry *coll = GEOSGeom_createCollection_r(
100 : m_poGeosContext, GEOS_GEOMETRYCOLLECTION, m_apoGeosInputs.data(),
101 14 : static_cast<unsigned int>(m_apoGeosInputs.size()));
102 :
103 14 : if (coll == nullptr)
104 : {
105 0 : return false;
106 : }
107 :
108 14 : m_apoGeosInputs.clear();
109 :
110 : int geos_result =
111 14 : GEOSCoverageIsValid_r(m_poGeosContext, coll, m_maximumGapWidth,
112 : &m_poGeosResultAsCollection);
113 14 : GEOSGeom_destroy_r(m_poGeosContext, coll);
114 :
115 14 : CPLDebug("CoverageIsValid", "%d", geos_result);
116 :
117 14 : return geos_result != 2;
118 : }
119 :
120 : private:
121 : OGRFeatureDefn *m_defn;
122 : const double m_maximumGapWidth;
123 : const bool m_includeValid;
124 :
125 : CPL_DISALLOW_COPY_ASSIGN(GDALVectorCheckCoverageOutputLayer)
126 : };
127 :
128 28 : GDALVectorCheckCoverageOutputLayer::~GDALVectorCheckCoverageOutputLayer()
129 : {
130 14 : if (m_defn != nullptr)
131 : {
132 14 : m_defn->Release();
133 : }
134 28 : }
135 :
136 12 : bool GDALVectorCheckCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &ctxt)
137 : {
138 12 : auto poSrcDS = m_inputDataset[0].GetDatasetRef();
139 24 : auto poDstDS = std::make_unique<GDALVectorNonStreamingAlgorithmDataset>();
140 :
141 12 : const bool bSingleLayerOutput = m_inputLayerNames.empty()
142 15 : ? poSrcDS->GetLayerCount() == 1
143 3 : : m_inputLayerNames.size() == 1;
144 :
145 24 : GDALVectorAlgorithmLayerProgressHelper progressHelper(ctxt);
146 :
147 28 : for (auto &&poSrcLayer : poSrcDS->GetLayers())
148 : {
149 21 : if (m_inputLayerNames.empty() ||
150 0 : std::find(m_inputLayerNames.begin(), m_inputLayerNames.end(),
151 21 : poSrcLayer->GetDescription()) != m_inputLayerNames.end())
152 : {
153 17 : const auto poSrcLayerDefn = poSrcLayer->GetLayerDefn();
154 17 : if (poSrcLayerDefn->GetGeomFieldCount() == 0)
155 : {
156 2 : if (m_inputLayerNames.empty())
157 1 : continue;
158 1 : ReportError(CE_Failure, CPLE_AppDefined,
159 : "Specified layer '%s' has no geometry field",
160 1 : poSrcLayer->GetDescription());
161 1 : return false;
162 : }
163 :
164 15 : progressHelper.AddProcessedLayer(*poSrcLayer);
165 : }
166 : }
167 :
168 29 : for ([[maybe_unused]] auto [poSrcLayer, bProcessed, layerProgressFunc,
169 54 : layerProgressData] : progressHelper)
170 : {
171 15 : const auto poSrcLayerDefn = poSrcLayer->GetLayerDefn();
172 : const int geomFieldIndex =
173 15 : m_geomField.empty()
174 15 : ? 0
175 2 : : poSrcLayerDefn->GetGeomFieldIndex(m_geomField.c_str());
176 :
177 15 : if (geomFieldIndex == -1)
178 : {
179 1 : ReportError(CE_Failure, CPLE_AppDefined,
180 : "Specified geometry field '%s' does not exist in "
181 : "layer '%s'",
182 1 : m_geomField.c_str(), poSrcLayer->GetDescription());
183 1 : return false;
184 : }
185 :
186 : std::string layerName = bSingleLayerOutput
187 : ? "invalid_edge"
188 34 : : std::string("invalid_edge_")
189 38 : .append(poSrcLayer->GetDescription());
190 :
191 : auto poLayer = std::make_unique<GDALVectorCheckCoverageOutputLayer>(
192 14 : *poSrcLayer, geomFieldIndex, layerName, m_maximumGapWidth,
193 14 : m_includeValid);
194 :
195 14 : if (!poDstDS->AddProcessedLayer(std::move(poLayer), layerProgressFunc,
196 : layerProgressData.get()))
197 : {
198 0 : return false;
199 : }
200 : }
201 :
202 10 : m_outputDataset.Set(std::move(poDstDS));
203 :
204 10 : return true;
205 : }
206 :
207 : #else
208 :
209 : bool GDALVectorCheckCoverageAlgorithm::RunStep(GDALPipelineStepRunContext &)
210 : {
211 : ReportError(CE_Failure, CPLE_AppDefined,
212 : "%s requires GDAL to be built against version 3.12 or later of "
213 : "the GEOS library.",
214 : NAME);
215 : return false;
216 : }
217 : #endif // HAVE_GEOS
218 :
219 : GDALVectorCheckCoverageAlgorithmStandalone::
220 : ~GDALVectorCheckCoverageAlgorithmStandalone() = default;
221 :
222 : //! @endcond
|