Geometry_impl.hpp
Go to the documentation of this file.
1 //=============================================================================
2 // Copyright (C) 2021-2024 Wageningen University - All Rights Reserved
3 // Author: Gonzalo Mier
4 // BSD-3 License
5 //=============================================================================
6 
7 #pragma once
8 #ifndef FIELDS2COVER_TYPES_GEOMETRY_IMPL_HPP_
9 #define FIELDS2COVER_TYPES_GEOMETRY_IMPL_HPP_
10 
11 #include <gdal/cpl_conv.h>
12 #include <geos_c.h>
13 #include <memory>
14 #include <vector>
15 #include <string>
16 
17 
18 namespace f2c::types {
19 
20 template <class T, OGRwkbGeometryType R>
21 Geometry<T, R>::Geometry() : data_(
22  downCast<T*>(OGRGeometryFactory::createGeometry(R)),
23  [](T* f) {OGRGeometryFactory::destroyGeometry(f);}) {}
24 
25 template <class T, OGRwkbGeometryType R>
26 Geometry<T, R>::Geometry(const T& g) : data_(downCast<T*>(g.clone()),
27  [](T* f) {OGRGeometryFactory::destroyGeometry(f);}) {}
28 
29 template <class T, OGRwkbGeometryType R>
30 Geometry<T, R>::Geometry(std::shared_ptr<T> g) : data_(g) {}
31 
32 template <class T, OGRwkbGeometryType R>
33 Geometry<T, R>::Geometry(T* g, EmptyDestructor) : data_(g, [](T* f) {}) {}
34 
35 template <class T, OGRwkbGeometryType R>
36 Geometry<T, R>::Geometry(const T* g) : data_(downCast<T*>(g->clone()),
37  [](T* f) {OGRGeometryFactory::destroyGeometry(f);}) {}
38 
39 template <class T, OGRwkbGeometryType R>
40 Geometry<T, R>::Geometry(OGRGeometry* g, EmptyDestructor) :
41  data_(downCast<T*>(g), [](T* f) {}) {}
42 
43 template <class T, OGRwkbGeometryType R>
44 Geometry<T, R>::Geometry(const OGRGeometry* g) : data_(downCast<T*>(g->clone()),
45  [](T* f) {OGRGeometryFactory::destroyGeometry(f);}) {}
46 
47 template <class T, OGRwkbGeometryType R>
48 Geometry<T, R>::~Geometry() = default;
49 
50 template <class T, OGRwkbGeometryType R>
51 Geometry<T, R>::Geometry(const Geometry& g) = default;
52 
53 template <class T, OGRwkbGeometryType R>
54 Geometry<T, R>::Geometry(Geometry&& g) = default;
55 
56 template <class T, OGRwkbGeometryType R>
57 typename Geometry<T, R>::Geometry& Geometry<T, R>::operator=(
58  Geometry&& g) = default;
59 
60 template <class T, OGRwkbGeometryType R>
61 typename Geometry<T, R>::Geometry& Geometry<T, R>::operator=(
62  const Geometry& g) = default;
63 
64 template <class T, OGRwkbGeometryType R>
65 std::shared_ptr<T> Geometry<T, R>::operator->() {return data_;}
66 
67 template <class T, OGRwkbGeometryType R>
68 std::shared_ptr<const T> Geometry<T, R>::operator->() const {return data_;}
69 
70 template <class T, OGRwkbGeometryType R>
71 T* Geometry<T, R>::get() {return data_.get();}
72 
73 template <class T, OGRwkbGeometryType R>
74 const T* Geometry<T, R>::get() const {return data_.get();}
75 
76 template <class T, OGRwkbGeometryType R>
77 bool Geometry<T, R>::operator!=(const Geometry<T, R>& geom2) const {
78  return !geom2->Equals(data_.get());
79 }
80 
81 template <class T, OGRwkbGeometryType R>
82 bool Geometry<T, R>::operator==(const Geometry<T, R>& geom2) const {
83  return geom2->Equals(data_.get());
84 }
85 
86 
87 template <class T, OGRwkbGeometryType R>
88 double Geometry<T, R>::getDimMinX() const {
89  OGREnvelope envelope;
90  data_->getEnvelope(&envelope);
91  return envelope.MinX;
92 }
93 
94 template <class T, OGRwkbGeometryType R>
95 double Geometry<T, R>::getDimMaxX() const {
96  OGREnvelope envelope;
97  data_->getEnvelope(&envelope);
98  return envelope.MaxX;
99 }
100 
101 template <class T, OGRwkbGeometryType R>
102 double Geometry<T, R>::getDimMinY() const {
103  OGREnvelope envelope;
104  data_->getEnvelope(&envelope);
105  return envelope.MinY;
106 }
107 
108 template <class T, OGRwkbGeometryType R>
109 double Geometry<T, R>::getDimMaxY() const {
110  OGREnvelope envelope;
111  data_->getEnvelope(&envelope);
112  return envelope.MaxY;
113 }
114 
115 template <class T, OGRwkbGeometryType R>
116 double Geometry<T, R>::getHeight() const {
117  OGREnvelope envelope;
118  data_->getEnvelope(&envelope);
119  return envelope.MaxY - envelope.MinY;
120 }
121 
122 template <class T, OGRwkbGeometryType R>
123 double Geometry<T, R>::getWidth() const {
124  OGREnvelope envelope;
125  data_->getEnvelope(&envelope);
126  return envelope.MaxX - envelope.MinX;
127 }
128 
129 template <class T, OGRwkbGeometryType R>
130 double Geometry<T, R>::getMinSafeLength() const {
131  OGREnvelope envelope;
132  data_->getEnvelope(&envelope);
133  return (envelope.MaxX - envelope.MinX) +
134  (envelope.MaxY - envelope.MinY);
135 }
136 
137 template <class T, OGRwkbGeometryType R>
138 template <class T2, OGRwkbGeometryType R2>
139 double Geometry<T, R>::distance(const Geometry<T2, R2>& p) const {
140  return data_->Distance(p.get());
141 }
142 
143 template <class T, OGRwkbGeometryType R>
144 template <class T2, OGRwkbGeometryType R2>
145 bool Geometry<T, R>::disjoint(const Geometry<T2, R2>& geom) const {
146  return data_->Disjoint(geom.get());
147 }
148 
149 template <class T, OGRwkbGeometryType R>
150 template <class T2, OGRwkbGeometryType R2>
151 bool Geometry<T, R>::crosses(const Geometry<T2, R2>& geom) const {
152  return data_->Crosses(geom.get());
153 }
154 
155 template <class T, OGRwkbGeometryType R>
156 template <class T2, OGRwkbGeometryType R2>
157 bool Geometry<T, R>::touches(const Geometry<T2, R2>& geom) const {
158  return data_->Touches(geom.get());
159 }
160 
161 template <class T, OGRwkbGeometryType R>
162 template <class T2, OGRwkbGeometryType R2>
163 bool Geometry<T, R>::within(const Geometry<T2, R2>& geom) const {
164  return data_->Within(geom.get());
165 }
166 
167 template <class T, OGRwkbGeometryType R>
168 template <class T2, OGRwkbGeometryType R2>
169 bool Geometry<T, R>::intersects(const Geometry<T2, R2>& geom) const {
170  return data_->Intersects(geom.get());
171 }
172 
173 template <class T, OGRwkbGeometryType R>
174 double Geometry<T, R>::mod_2pi(double val) {
175  return mod(val, boost::math::constants::two_pi<double>());
176 }
177 
178 template <class T, OGRwkbGeometryType R>
179 double Geometry<T, R>::getAngContinuity(double prev_val, double val) {
180  const auto two_pi = boost::math::constants::two_pi<double>();
181  auto val_id = round(prev_val / two_pi);
182  auto getNearVal =
183  [two_pi] (int i, double value) {return i * two_pi + mod_2pi(value);};
184  for (int i = -2; i <= 0; ++i) {
185  if (fabs(getNearVal(val_id + i + 1, val) - prev_val) >
186  fabs(getNearVal(val_id + i, val) - prev_val)) {
187  return getNearVal(val_id + i, val);
188  }
189  }
190  return getNearVal(val_id + 1, val);
191 }
192 
193 template <class T, OGRwkbGeometryType R>
194 std::vector<double> Geometry<T, R>::getAngContinuity(
195  const std::vector<double>& val) {
196  std::vector<double> res {val[0]};
197  for (size_t i = 1; i < val.size(); ++i) {
198  res.emplace_back(getAngContinuity(res[i - 1], val[i]));
199  }
200  return res;
201 }
202 
203 template <class T, OGRwkbGeometryType R>
204 double Geometry<T, R>::getAngleDiffAbs(double a, double b) {
205  const auto pi = boost::math::constants::pi<double>();
206  return pi - fabs(mod_2pi(b - a) - pi);
207 }
208 
209 template <class T, OGRwkbGeometryType R>
210 double Geometry<T, R>::getAngleAvg(double a, double b) {
211  const auto pi = boost::math::constants::pi<double>();
212  double ma = mod_2pi(a);
213  double mb = mod_2pi(b);
214  if (abs(mb - ma) < pi) {
215  return 0.5 * (ma + mb);
216  } else {
217  return 0.5 * (ma + mb) + pi;
218  }
219 }
220 
221 template <class T, OGRwkbGeometryType R>
222 bool Geometry<T, R>::isEmpty() const {
223  return (!data_ || data_->IsEmpty());
224 }
225 
226 template <class T, OGRwkbGeometryType R>
227 std::string Geometry<T, R>::exportToWkt() const {
228  char *pszWKT = nullptr;
229  data_->exportToWkt(&pszWKT);
230  std::string result{pszWKT};
231  CPLFree(pszWKT);
232  return result;
233 }
234 
235 template <class T, OGRwkbGeometryType R>
236 std::string Geometry<T, R>::exportToGML() const {
237  char *s = data_->exportToGML();
238  std::string str(s);
239  CPLFree(s);
240  return str;
241 }
242 
243 template <class T, OGRwkbGeometryType R>
244 std::string Geometry<T, R>::exportToKML() const {
245  char *s = data_->exportToKML();
246  std::string str(s);
247  CPLFree(s);
248  return str;
249 }
250 
251 template <class T, OGRwkbGeometryType R>
252 std::string Geometry<T, R>::exportToJson() const {
253  char *s = data_->exportToJson();
254  std::string str(s);
255  CPLFree(s);
256  return str;
257 }
258 
259 template <class T, OGRwkbGeometryType R>
260 void Geometry<T, R>::importFromWkt(const std::string& text) {
261  auto char_text = text.c_str();
262  data_->importFromWkt(&char_text);
263 }
264 
265 // ###############################
266 // Code extracted from:
267 // https://github.com/OSGeo/gdal/blob/b0aa6065a39b252cb8306e9c2e2535d6dda0fb55/port/cpl_conv.h#L397
268 template <class T, OGRwkbGeometryType R>
269 template <typename To, typename From>
270 inline To Geometry<T, R>::downCast(From *f) const {
271  static_assert(
272  (std::is_base_of<From, typename std::remove_pointer<To>::type>::value),
273  "target type not derived from source type");
274  CPLAssert(f == nullptr || dynamic_cast<To>(f) != nullptr);
275  return static_cast<To>(f);
276 }
277 // ###############################
278 
279 template <class T, OGRwkbGeometryType R>
280 template<typename RetType>
281 RetType Geometry<T, R>::destroyResGeom(OGRGeometry* g) {
282  RetType res(g);
283  OGRGeometryFactory::destroyGeometry(g);
284  return res;
285 }
286 
287 // ###############################
288 // Code extracted from:
289 // https://github.com/OSGeo/gdal/blob/717dcc0eed252e2f78c142b1f7866e49c5511224/ogr/ogrgeometry.cpp#L4309
290 template <class T, OGRwkbGeometryType R>
291 OGRGeometry* Geometry<T, R>::OGRBuffer(double dfDist, int side) const {
292  OGRGeometry *poOGRProduct = nullptr;
293 
294  GEOSContextHandle_t hGEOSCtxt = OGRGeometry::createGEOSContext();
295  GEOSGeom hGeosGeom = this->data_->exportToGEOS(hGEOSCtxt);
296  if (hGeosGeom != nullptr) {
297  GEOSBufferParams* hBufParams = GEOSBufferParams_create_r(hGEOSCtxt);
298  GEOSBufferParams_setEndCapStyle_r(hGEOSCtxt, hBufParams,
299  GEOSBufCapStyles::GEOSBUF_CAP_FLAT);
300  GEOSBufferParams_setJoinStyle_r(hGEOSCtxt, hBufParams,
301  GEOSBufJoinStyles::GEOSBUF_JOIN_MITRE);
302  GEOSBufferParams_setSingleSided_r(hGEOSCtxt, hBufParams, side != 0);
303 
304  GEOSGeom hGeosProduct = GEOSBufferWithParams_r(hGEOSCtxt, hGeosGeom,
305  hBufParams, side >= 0 ? dfDist : -dfDist);
306  GEOSGeom_destroy_r(hGEOSCtxt, hGeosGeom);
307  GEOSBufferParams_destroy_r(hGEOSCtxt, hBufParams);
308 
309  poOGRProduct = buildGeometryFromGEOS(hGEOSCtxt, hGeosProduct,
310  this->data_.get(), nullptr);
311  }
312  OGRGeometry::freeGEOSContext(hGEOSCtxt);
313  return poOGRProduct;
314 }
316 
317 // ###############################
318 // Code extracted from:
319 // https://github.com/OSGeo/gdal/blob/717dcc0eed252e2f78c142b1f7866e49c5511224/ogr/ogrgeometry.cpp#L4309
320 template <class T, OGRwkbGeometryType R>
321 OGRGeometry* Geometry<T, R>::OGRGeometryRebuildCurves(
322  const OGRGeometry *poGeom, const OGRGeometry *poOtherGeom,
323  OGRGeometry *poOGRProduct) const {
324  if (poOGRProduct != nullptr &&
325  wkbFlatten(poOGRProduct->getGeometryType()) != wkbPoint &&
326  (poGeom->hasCurveGeometry(true) ||
327  (poOtherGeom && poOtherGeom->hasCurveGeometry(true)))) {
328  OGRGeometry *poCurveGeom = poOGRProduct->getCurveGeometry();
329  delete poOGRProduct;
330  return poCurveGeom;
331  }
332  return poOGRProduct;
333 }
334 template <class T, OGRwkbGeometryType R>
335 OGRGeometry* Geometry<T, R>::buildGeometryFromGEOS(
336  GEOSContextHandle_t hGEOSCtxt, GEOSGeom hGeosProduct,
337  const OGRGeometry *poSelf, const OGRGeometry *poOtherGeom) const {
338  OGRGeometry *poOGRProduct = nullptr;
339  if (hGeosProduct != nullptr) {
340  poOGRProduct = OGRGeometryFactory::createFromGEOS(hGEOSCtxt, hGeosProduct);
341  if (poOGRProduct != nullptr &&
342  poSelf->getSpatialReference() != nullptr &&
343  (poOtherGeom == nullptr ||
344  (poOtherGeom->getSpatialReference() != nullptr &&
345  poOtherGeom->getSpatialReference()->IsSame(
346  poSelf->getSpatialReference())))) {
347  poOGRProduct->assignSpatialReference(poSelf->getSpatialReference());
348  }
349  poOGRProduct = OGRGeometryRebuildCurves(poSelf, poOtherGeom, poOGRProduct);
350  GEOSGeom_destroy_r(hGEOSCtxt, hGeosProduct);
351  }
352  return poOGRProduct;
353 }
354 // ###############################
355 
356 
357 template <class T, OGRwkbGeometryType R>
358 double Geometry<T, R>::mod(double a, double b) {
359  return fmod(fmod(a, b) + b, b);
360 }
361 
362 
363 } // namespace f2c::types
364 
365 #endif // FIELDS2COVER_TYPES_GEOMETRY_IMPL_HPP_
f2c::types
Types used by fields2cover library.
Definition: Cell.h:20
f2c::types::Geometry
Definition: Geometry.h:26
f2c::types::Geometry::exportToJson
std::string exportToJson() const
Definition: Geometry_impl.hpp:252
f2c::types::Geometry::exportToKML
std::string exportToKML() const
Definition: Geometry_impl.hpp:244
f2c::types::Geometry::get
T * get()
Definition: Geometry_impl.hpp:71
f2c::types::Geometry::exportToWkt
std::string exportToWkt() const
Definition: Geometry_impl.hpp:227
f2c::types::Geometry::importFromWkt
void importFromWkt(const std::string &text)
Definition: Geometry_impl.hpp:260
f2c::types::Geometry::exportToGML
std::string exportToGML() const
Definition: Geometry_impl.hpp:236
f2c::types::EmptyDestructor
Definition: Geometry.h:23
f2c::types::Geometry::Geometry
Geometry()
Definition: Geometry_impl.hpp:21

fields2cover
Author(s):
autogenerated on Fri Apr 25 2025 02:18:31