Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: OpenGIS Simple Features Reference Implementation
4 : * Purpose: Implements OGRGeomCoordinatePrecision.
5 : * Author: Even Rouault <even dot rouault at spatialys.com>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2024, Even Rouault <even dot rouault at spatialys.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "ogr_core.h"
14 : #include "ogr_api.h"
15 : #include "ogr_spatialref.h"
16 : #include "ogr_geomcoordinateprecision.h"
17 :
18 : #include <algorithm>
19 : #include <cmath>
20 :
21 : /************************************************************************/
22 : /* OGRGeomCoordinatePrecisionCreate() */
23 : /************************************************************************/
24 :
25 : /** Creates a new instance of OGRGeomCoordinatePrecision.
26 : *
27 : * The default X,Y,Z,M resolutions are set to OGR_GEOM_COORD_PRECISION_UNKNOWN.
28 : *
29 : * @since GDAL 3.9
30 : */
31 21 : OGRGeomCoordinatePrecisionH OGRGeomCoordinatePrecisionCreate(void)
32 : {
33 : static_assert(OGR_GEOM_COORD_PRECISION_UNKNOWN ==
34 : OGRGeomCoordinatePrecision::UNKNOWN);
35 :
36 21 : return new OGRGeomCoordinatePrecision();
37 : }
38 :
39 : /************************************************************************/
40 : /* OGRGeomCoordinatePrecisionDestroy() */
41 : /************************************************************************/
42 :
43 : /** Destroy a OGRGeomCoordinatePrecision.
44 : *
45 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance or nullptr
46 : * @since GDAL 3.9
47 : */
48 21 : void OGRGeomCoordinatePrecisionDestroy(
49 : OGRGeomCoordinatePrecisionH hGeomCoordPrec)
50 : {
51 21 : delete hGeomCoordPrec;
52 21 : }
53 :
54 : /************************************************************************/
55 : /* OGRGeomCoordinatePrecisionGetXYResolution() */
56 : /************************************************************************/
57 :
58 : /** Get the X/Y resolution of a OGRGeomCoordinatePrecision
59 : *
60 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
61 : * @return the the X/Y resolution of a OGRGeomCoordinatePrecision or
62 : * OGR_GEOM_COORD_PRECISION_UNKNOWN
63 : * @since GDAL 3.9
64 : */
65 47 : double OGRGeomCoordinatePrecisionGetXYResolution(
66 : OGRGeomCoordinatePrecisionH hGeomCoordPrec)
67 : {
68 47 : VALIDATE_POINTER1(hGeomCoordPrec,
69 : "OGRGeomCoordinatePrecisionGetXYResolution", 0);
70 47 : return hGeomCoordPrec->dfXYResolution;
71 : }
72 :
73 : /************************************************************************/
74 : /* OGRGeomCoordinatePrecisionGetZResolution() */
75 : /************************************************************************/
76 :
77 : /** Get the Z resolution of a OGRGeomCoordinatePrecision
78 : *
79 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
80 : * @return the the Z resolution of a OGRGeomCoordinatePrecision or
81 : * OGR_GEOM_COORD_PRECISION_UNKNOWN
82 : * @since GDAL 3.9
83 : */
84 44 : double OGRGeomCoordinatePrecisionGetZResolution(
85 : OGRGeomCoordinatePrecisionH hGeomCoordPrec)
86 : {
87 44 : VALIDATE_POINTER1(hGeomCoordPrec,
88 : "OGRGeomCoordinatePrecisionGetZResolution", 0);
89 44 : return hGeomCoordPrec->dfZResolution;
90 : }
91 :
92 : /************************************************************************/
93 : /* OGRGeomCoordinatePrecisionGetMResolution() */
94 : /************************************************************************/
95 :
96 : /** Get the M resolution of a OGRGeomCoordinatePrecision
97 : *
98 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
99 : * @return the the M resolution of a OGRGeomCoordinatePrecision or
100 : * OGR_GEOM_COORD_PRECISION_UNKNOWN
101 : * @since GDAL 3.9
102 : */
103 29 : double OGRGeomCoordinatePrecisionGetMResolution(
104 : OGRGeomCoordinatePrecisionH hGeomCoordPrec)
105 : {
106 29 : VALIDATE_POINTER1(hGeomCoordPrec,
107 : "OGRGeomCoordinatePrecisionGetMResolution", 0);
108 29 : return hGeomCoordPrec->dfMResolution;
109 : }
110 :
111 : /************************************************************************/
112 : /* OGRGeomCoordinatePrecisionGetFormats() */
113 : /************************************************************************/
114 :
115 : /** Get the list of format names for coordinate precision format specific
116 : * options.
117 : *
118 : * An example of a supported value for pszFormatName is
119 : * "FileGeodatabase" for layers of the OpenFileGDB driver.
120 : *
121 : * The returned values may be used for the pszFormatName argument of
122 : * OGRGeomCoordinatePrecisionGetFormatSpecificOptions().
123 : *
124 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
125 : * @return a null-terminated list to free with CSLDestroy(), or nullptr.
126 : * @since GDAL 3.9
127 : */
128 : char **
129 6 : OGRGeomCoordinatePrecisionGetFormats(OGRGeomCoordinatePrecisionH hGeomCoordPrec)
130 : {
131 6 : VALIDATE_POINTER1(hGeomCoordPrec, "OGRGeomCoordinatePrecisionGetFormats",
132 : nullptr);
133 12 : CPLStringList aosFormats;
134 11 : for (const auto &kv : hGeomCoordPrec->oFormatSpecificOptions)
135 : {
136 5 : aosFormats.AddString(kv.first.c_str());
137 : }
138 6 : return aosFormats.StealList();
139 : }
140 :
141 : /************************************************************************/
142 : /* OGRGeomCoordinatePrecisionGetFormatSpecificOptions() */
143 : /************************************************************************/
144 :
145 : /** Get format specific coordinate precision options.
146 : *
147 : * An example of a supported value for pszFormatName is
148 : * "FileGeodatabase" for layers of the OpenFileGDB driver.
149 : *
150 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
151 : * @param pszFormatName A format name (one of those returned by
152 : * OGRGeomCoordinatePrecisionGetFormats())
153 : * @return a null-terminated list, or nullptr. The list must *not* be freed,
154 : * and is owned by hGeomCoordPrec
155 : * @since GDAL 3.9
156 : */
157 7 : CSLConstList OGRGeomCoordinatePrecisionGetFormatSpecificOptions(
158 : OGRGeomCoordinatePrecisionH hGeomCoordPrec, const char *pszFormatName)
159 : {
160 7 : VALIDATE_POINTER1(hGeomCoordPrec,
161 : "OGRGeomCoordinatePrecisionGetFormatSpecificOptions",
162 : nullptr);
163 : const auto oIter =
164 7 : hGeomCoordPrec->oFormatSpecificOptions.find(pszFormatName);
165 7 : if (oIter == hGeomCoordPrec->oFormatSpecificOptions.end())
166 : {
167 1 : return nullptr;
168 : }
169 6 : return oIter->second.List();
170 : }
171 :
172 : /************************************************************************/
173 : /* OGRGeomCoordinatePrecisionSetFormatSpecificOptions() */
174 : /************************************************************************/
175 :
176 : /** Set format specific coordinate precision options.
177 : *
178 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
179 : * @param pszFormatName A format name (must not be null)
180 : * @param papszOptions null-terminated list of options.
181 : * @since GDAL 3.9
182 : */
183 1 : void OGRGeomCoordinatePrecisionSetFormatSpecificOptions(
184 : OGRGeomCoordinatePrecisionH hGeomCoordPrec, const char *pszFormatName,
185 : CSLConstList papszOptions)
186 : {
187 1 : VALIDATE_POINTER0(hGeomCoordPrec,
188 : "OGRGeomCoordinatePrecisionSetFormatSpecificOptions");
189 1 : hGeomCoordPrec->oFormatSpecificOptions[pszFormatName] = papszOptions;
190 : }
191 :
192 : /************************************************************************/
193 : /* OGRGeomCoordinatePrecisionSet() */
194 : /************************************************************************/
195 :
196 : /**
197 : * \brief Set the resolution of the geometry coordinate components.
198 : *
199 : * For the X, Y and Z ordinates, the precision should be expressed in the units
200 : * of the CRS of the geometry. So typically degrees for geographic CRS, or
201 : * meters/feet/US-feet for projected CRS.
202 : * Users might use OGRGeomCoordinatePrecisionSetFromMeter() for an even more
203 : * convenient interface.
204 : *
205 : * For a projected CRS with meters as linear unit, 1e-3 corresponds to a
206 : * millimetric precision.
207 : * For a geographic CRS in 8.9e-9 corresponds to a millimetric precision
208 : * (for a Earth CRS)
209 : *
210 : * Resolution should be stricty positive, or set to
211 : * OGR_GEOM_COORD_PRECISION_UNKNOWN when unknown.
212 : *
213 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
214 : * @param dfXYResolution Resolution for for X and Y coordinates.
215 : * @param dfZResolution Resolution for for Z coordinates.
216 : * @param dfMResolution Resolution for for M coordinates.
217 : * @since GDAL 3.9
218 : */
219 :
220 21 : void OGRGeomCoordinatePrecisionSet(OGRGeomCoordinatePrecisionH hGeomCoordPrec,
221 : double dfXYResolution, double dfZResolution,
222 : double dfMResolution)
223 : {
224 21 : VALIDATE_POINTER0(hGeomCoordPrec, "OGRGeomCoordinatePrecisionSet");
225 21 : hGeomCoordPrec->dfXYResolution = dfXYResolution;
226 21 : hGeomCoordPrec->dfZResolution = dfZResolution;
227 21 : hGeomCoordPrec->dfMResolution = dfMResolution;
228 : }
229 :
230 : /************************************************************************/
231 : /* OGRGeomCoordinatePrecisionSetFromMeter() */
232 : /************************************************************************/
233 :
234 : /**
235 : * \brief Set the resolution of the geometry coordinate components.
236 : *
237 : * For the X, Y and Z ordinates, the precision should be expressed in meter,
238 : * e.g 1e-3 for millimetric precision.
239 : *
240 : * Resolution should be stricty positive, or set to
241 : * OGR_GEOM_COORD_PRECISION_UNKNOWN when unknown.
242 : *
243 : * @param hGeomCoordPrec OGRGeomCoordinatePrecision instance (must not be null)
244 : * @param hSRS Spatial reference system, used for metric to SRS unit conversion
245 : * (must not be null)
246 : * @param dfXYMeterResolution Resolution for for X and Y coordinates, in meter.
247 : * @param dfZMeterResolution Resolution for for Z coordinates, in meter.
248 : * @param dfMResolution Resolution for for M coordinates.
249 : * @since GDAL 3.9
250 : */
251 3 : void OGRGeomCoordinatePrecisionSetFromMeter(
252 : OGRGeomCoordinatePrecisionH hGeomCoordPrec, OGRSpatialReferenceH hSRS,
253 : double dfXYMeterResolution, double dfZMeterResolution, double dfMResolution)
254 : {
255 3 : VALIDATE_POINTER0(hGeomCoordPrec, "OGRGeomCoordinatePrecisionSet");
256 3 : VALIDATE_POINTER0(hSRS, "OGRGeomCoordinatePrecisionSet");
257 3 : return hGeomCoordPrec->SetFromMeter(OGRSpatialReference::FromHandle(hSRS),
258 : dfXYMeterResolution, dfZMeterResolution,
259 3 : dfMResolution);
260 : }
261 :
262 : /************************************************************************/
263 : /* GetConversionFactors() */
264 : /************************************************************************/
265 :
266 1860 : static void GetConversionFactors(const OGRSpatialReference *poSRS,
267 : double &dfXYFactor, double &dfZFactor)
268 : {
269 1860 : dfXYFactor = 1;
270 1860 : dfZFactor = 1;
271 :
272 1860 : if (poSRS)
273 : {
274 1216 : if (poSRS->IsGeographic())
275 : {
276 789 : dfXYFactor = poSRS->GetSemiMajor(nullptr) * M_PI / 180;
277 : }
278 : else
279 : {
280 427 : dfXYFactor = poSRS->GetLinearUnits(nullptr);
281 : }
282 :
283 1216 : if (poSRS->GetAxesCount() == 3)
284 : {
285 40 : poSRS->GetAxis(nullptr, 2, nullptr, &dfZFactor);
286 : }
287 : }
288 1860 : }
289 :
290 : /************************************************************************/
291 : /* OGRGeomCoordinatePrecision::SetFromMeter() */
292 : /************************************************************************/
293 :
294 : /**
295 : * \brief Set the resolution of the geometry coordinate components.
296 : *
297 : * For the X, Y and Z coordinates, the precision should be expressed in meter,
298 : * e.g 1e-3 for millimetric precision.
299 : *
300 : * Resolution should be stricty positive, or set to
301 : * OGRGeomCoordinatePrecision::UNKNOWN when unknown.
302 : *
303 : * @param poSRS Spatial reference system, used for metric to SRS unit conversion
304 : * (must not be null)
305 : * @param dfXYMeterResolution Resolution for for X and Y coordinates, in meter.
306 : * @param dfZMeterResolution Resolution for for Z coordinates, in meter.
307 : * @param dfMResolutionIn Resolution for for M coordinates.
308 : * @since GDAL 3.9
309 : */
310 12 : void OGRGeomCoordinatePrecision::SetFromMeter(const OGRSpatialReference *poSRS,
311 : double dfXYMeterResolution,
312 : double dfZMeterResolution,
313 : double dfMResolutionIn)
314 : {
315 12 : double dfXYFactor = 1;
316 12 : double dfZFactor = 1;
317 12 : GetConversionFactors(poSRS, dfXYFactor, dfZFactor);
318 :
319 12 : dfXYResolution = dfXYMeterResolution / dfXYFactor;
320 12 : dfZResolution = dfZMeterResolution / dfZFactor;
321 12 : dfMResolution = dfMResolutionIn;
322 12 : }
323 :
324 : /************************************************************************/
325 : /* OGRGeomCoordinatePrecision::ConvertToOtherSRS() */
326 : /************************************************************************/
327 :
328 : /**
329 : * \brief Return equivalent coordinate precision setting taking into account
330 : * a change of SRS.
331 : *
332 : * @param poSRSSrc Spatial reference system of the current instance
333 : * (if null, meter unit is assumed)
334 : * @param poSRSDst Spatial reference system of the returned instance
335 : * (if null, meter unit is assumed)
336 : * @return a new OGRGeomCoordinatePrecision instance, with a poSRSDst SRS.
337 : * @since GDAL 3.9
338 : */
339 924 : OGRGeomCoordinatePrecision OGRGeomCoordinatePrecision::ConvertToOtherSRS(
340 : const OGRSpatialReference *poSRSSrc,
341 : const OGRSpatialReference *poSRSDst) const
342 : {
343 924 : double dfXYFactorSrc = 1;
344 924 : double dfZFactorSrc = 1;
345 924 : GetConversionFactors(poSRSSrc, dfXYFactorSrc, dfZFactorSrc);
346 :
347 924 : double dfXYFactorDst = 1;
348 924 : double dfZFactorDst = 1;
349 924 : GetConversionFactors(poSRSDst, dfXYFactorDst, dfZFactorDst);
350 :
351 924 : OGRGeomCoordinatePrecision oNewPrec;
352 924 : oNewPrec.dfXYResolution = dfXYResolution * dfXYFactorSrc / dfXYFactorDst;
353 924 : oNewPrec.dfZResolution = dfZResolution * dfZFactorSrc / dfZFactorDst;
354 924 : oNewPrec.dfMResolution = dfMResolution;
355 :
356 : // Only preserve source forma specific options if no reprojection is
357 : // involved
358 1391 : if ((!poSRSSrc && !poSRSDst) ||
359 467 : (poSRSSrc && poSRSDst && poSRSSrc->IsSame(poSRSDst)))
360 : {
361 503 : oNewPrec.oFormatSpecificOptions = oFormatSpecificOptions;
362 : }
363 :
364 1848 : return oNewPrec;
365 : }
366 :
367 : /************************************************************************/
368 : /* OGRGeomCoordinatePrecision::ResolutionToPrecision() */
369 : /************************************************************************/
370 :
371 : /**
372 : * \brief Return the number of decimal digits after the decimal point to
373 : * get the specified resolution.
374 : *
375 : * @since GDAL 3.9
376 : */
377 :
378 : /* static */
379 57 : int OGRGeomCoordinatePrecision::ResolutionToPrecision(double dfResolution)
380 : {
381 : return static_cast<int>(
382 57 : std::ceil(std::log10(1. / std::min(1.0, dfResolution))));
383 : }
|
