Line data Source code
1 : /******************************************************************************
2 : * Project: ogr linear referencing utility
3 : * Purpose: main source file
4 : * Author: Dmitry Baryshnikov (aka Bishop), polimax@mail.ru
5 : *
6 : ******************************************************************************
7 : * Copyright (C) 2014 NextGIS
8 : *
9 : * Permission is hereby granted, free of charge, to any person obtaining a
10 : * copy of this software and associated documentation files (the "Software"),
11 : * to deal in the Software without restriction, including without limitation
12 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
13 : * and/or sell copies of the Software, and to permit persons to whom the
14 : * Software is furnished to do so, subject to the following conditions:
15 : *
16 : * The above copyright notice and this permission notice shall be included
17 : * in all copies or substantial portions of the Software.
18 : *
19 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
20 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
21 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
22 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
23 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
24 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
25 : * DEALINGS IN THE SOFTWARE.
26 : ****************************************************************************/
27 :
28 : #include "cpl_port.h"
29 :
30 : #include "commonutils.h"
31 : #include "cpl_conv.h"
32 : #include "cpl_error.h"
33 : #include "cpl_string.h"
34 : #include "gdal_version.h"
35 : #include "gdal.h"
36 : #include "gdal_alg.h"
37 : #include "ogr_api.h"
38 : #include "ogr_p.h"
39 : #include "ogrsf_frmts.h"
40 :
41 : #include <limits>
42 : #include <map>
43 : #include <set>
44 : #include <vector>
45 :
46 : #define FIELD_START "beg"
47 : #define FIELD_FINISH "end"
48 : #define FIELD_SCALE_FACTOR "scale"
49 : constexpr double DELTA = 0.00000001; // - delta
50 : #ifdef HAVE_GEOS
51 : constexpr double TOLERANCE_DEGREE = 0.00008983153;
52 : constexpr double TOLERANCE_METER = 10.0;
53 : #endif
54 :
55 : enum operation
56 : {
57 : op_unknown = 0,
58 : op_create,
59 : op_get_pos,
60 : op_get_coord,
61 : op_get_subline
62 : };
63 :
64 : typedef struct _curve_data
65 : {
66 : OGRLineString *pPart;
67 : double dfBeg, dfEnd, dfFactor;
68 :
69 70 : bool IsInside(const double &dfDist) const
70 : {
71 70 : return (dfDist + DELTA >= dfBeg) && (dfDist - DELTA <= dfEnd);
72 : }
73 : } CURVE_DATA;
74 :
75 : /************************************************************************/
76 : /* Usage() */
77 : /************************************************************************/
78 : static void Usage(bool bIsError, const char *pszAdditionalMsg = nullptr,
79 : bool bShort = true) CPL_NO_RETURN;
80 :
81 0 : static void Usage(bool bIsError, const char *pszAdditionalMsg, bool bShort)
82 : {
83 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
84 :
85 0 : fprintf(
86 : bIsError ? stderr : stdout,
87 : "Usage: ogrlineref [--help] [--help-general] [-progress] [-quiet]\n"
88 : " [-f <format_name>] [-dsco <NAME>=<VALUE>]... [-lco "
89 : "<NAME>=<VALUE>]...\n"
90 : " [-create]\n"
91 : " [-l <src_line_datasource_name>] [-ln <layer_name>] "
92 : "[-lf "
93 : "<field_name>]\n"
94 : " [-p <src_repers_datasource_name>] [-pn <layer_name>] "
95 : "[-pm <pos_field_name>] [-pf <field_name>]\n"
96 : " [-r <src_parts_datasource_name>] [-rn <layer_name>]\n"
97 : " [-o <dst_datasource_name>] [-on <layer_name>] [-of "
98 : "<field_name>] [-s <step>]\n"
99 : " [-get_pos] [-x <long>] [-y <lat>]\n"
100 : " [-get_coord] [-m <position>] \n"
101 : " [-get_subline] [-mb <position>] [-me <position>]\n");
102 :
103 0 : if (bShort)
104 : {
105 0 : fprintf(bIsError ? stderr : stdout,
106 : "\nNote: ogrlineref --long-usage for full help.\n");
107 0 : if (pszAdditionalMsg)
108 0 : fprintf(stderr, "\nFAILURE: %s\n", pszAdditionalMsg);
109 0 : exit(1);
110 : }
111 :
112 0 : fprintf(
113 : bIsError ? stderr : stdout,
114 : "\n -f format_name: output file format name, possible values are:\n");
115 :
116 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
117 : {
118 0 : GDALDriver *poDriver = poR->GetDriver(iDriver);
119 :
120 0 : if (CPLTestBool(CSLFetchNameValueDef(poDriver->GetMetadata(),
121 : GDAL_DCAP_CREATE, "FALSE")))
122 0 : fprintf(bIsError ? stderr : stdout, " -f \"%s\"\n",
123 0 : poDriver->GetDescription());
124 : }
125 :
126 0 : fprintf(bIsError ? stderr : stdout,
127 : " -progress: Display progress on terminal. Only works if input "
128 : "layers have the \n"
129 : " \"fast feature count\" "
130 : "capability\n"
131 : " -dsco NAME=VALUE: Dataset creation option (format specific)\n"
132 : " -lco NAME=VALUE: Layer creation option (format specific)\n"
133 : " -l src_line_datasource_name: Datasource of line path name\n"
134 : " -ln layer_name: Layer name in datasource (optional)\n"
135 : " -lf field_name: Field name for unique paths in layer (optional)\n"
136 : " -p src_repers_datasource_name: Datasource of repers name\n"
137 : " -pn layer_name: Layer name in datasource (optional)\n"
138 : " -pm pos_field_name: Line position field name\n"
139 : " -pf field_name: Field name for correspondence repers of separate "
140 : "paths in layer (optional)\n"
141 : " -r src_parts_datasource_name: Parts datasource name\n"
142 : " -rn layer_name: Layer name in datasource (optional)\n"
143 : " -o dst_datasource_name: Parts datasource name\n"
144 : " -on layer_name: Layer name in datasource (optional)\n"
145 : " -of field_name: Field name for correspondence parts of separate "
146 : "paths in layer (optional)\n"
147 : " -s step: part size in m\n");
148 :
149 0 : if (pszAdditionalMsg)
150 0 : fprintf(stderr, "\nFAILURE: %s\n", pszAdditionalMsg);
151 :
152 0 : exit(bIsError ? 1 : 0);
153 : }
154 :
155 : /************************************************************************/
156 : /* SetupTargetLayer() */
157 : /************************************************************************/
158 :
159 3 : static OGRLayer *SetupTargetLayer(OGRLayer *poSrcLayer, GDALDataset *poDstDS,
160 : char **papszLCO, const char *pszNewLayerName,
161 : const char *pszOutputSepFieldName = nullptr)
162 : {
163 : const CPLString szLayerName =
164 3 : pszNewLayerName == nullptr ? CPLGetBasename(poDstDS->GetDescription())
165 9 : : pszNewLayerName;
166 :
167 : /* -------------------------------------------------------------------- */
168 : /* Get other info. */
169 : /* -------------------------------------------------------------------- */
170 3 : OGRFeatureDefn *poSrcFDefn = poSrcLayer->GetLayerDefn();
171 :
172 : /* -------------------------------------------------------------------- */
173 : /* Find requested geometry fields. */
174 : /* -------------------------------------------------------------------- */
175 :
176 3 : OGRSpatialReference *poOutputSRS = poSrcLayer->GetSpatialRef();
177 :
178 : /* -------------------------------------------------------------------- */
179 : /* Find the layer. */
180 : /* -------------------------------------------------------------------- */
181 :
182 : // GetLayerByName() can instantiate layers that would have been
183 : // 'hidden' otherwise, for example, non-spatial tables in a
184 : // PostGIS-enabled database, so this apparently useless command is
185 : // not useless... (#4012)
186 3 : CPLPushErrorHandler(CPLQuietErrorHandler);
187 3 : OGRLayer *poDstLayer = poDstDS->GetLayerByName(szLayerName);
188 3 : CPLPopErrorHandler();
189 3 : CPLErrorReset();
190 :
191 3 : if (poDstLayer != nullptr)
192 : {
193 0 : const int nLayerCount = poDstDS->GetLayerCount();
194 0 : int iLayer = -1; // Used after for.
195 0 : for (iLayer = 0; iLayer < nLayerCount; iLayer++)
196 : {
197 0 : OGRLayer *poLayer = poDstDS->GetLayer(iLayer);
198 0 : if (poLayer == poDstLayer)
199 0 : break;
200 : }
201 :
202 0 : if (iLayer == nLayerCount)
203 : // Should not happen with an ideal driver.
204 0 : poDstLayer = nullptr;
205 : }
206 :
207 : /* -------------------------------------------------------------------- */
208 : /* If the layer does not exist, then create it. */
209 : /* -------------------------------------------------------------------- */
210 3 : if (poDstLayer == nullptr)
211 : {
212 3 : if (!poDstDS->TestCapability(ODsCCreateLayer))
213 : {
214 0 : fprintf(stderr,
215 : "Layer %s not found, and "
216 : "CreateLayer not supported by driver.\n",
217 : szLayerName.c_str());
218 0 : return nullptr;
219 : }
220 :
221 3 : OGRwkbGeometryType eGType = wkbLineString;
222 :
223 3 : CPLErrorReset();
224 :
225 3 : if (poDstDS->TestCapability(ODsCCreateGeomFieldAfterCreateLayer))
226 : {
227 0 : eGType = wkbNone;
228 : }
229 :
230 3 : poDstLayer = poDstDS->CreateLayer(
231 : szLayerName, poOutputSRS, static_cast<OGRwkbGeometryType>(eGType),
232 : papszLCO);
233 :
234 3 : if (poDstLayer == nullptr)
235 0 : return nullptr;
236 :
237 3 : if (poDstDS->TestCapability(ODsCCreateGeomFieldAfterCreateLayer))
238 : {
239 0 : OGRGeomFieldDefn oGFldDefn(poSrcFDefn->GetGeomFieldDefn(0));
240 0 : if (poOutputSRS != nullptr)
241 0 : oGFldDefn.SetSpatialRef(poOutputSRS);
242 0 : oGFldDefn.SetType(wkbLineString);
243 0 : poDstLayer->CreateGeomField(&oGFldDefn);
244 : }
245 : }
246 :
247 : /* -------------------------------------------------------------------- */
248 : /* Otherwise we will append to it, if append was requested. */
249 : /* -------------------------------------------------------------------- */
250 : else
251 : {
252 0 : fprintf(stderr, "FAILED: Layer %s already exists.\n",
253 : szLayerName.c_str());
254 0 : return nullptr;
255 : }
256 :
257 : // Create beg, end, scale factor fields.
258 6 : OGRFieldDefn oFieldDefn_Beg(FIELD_START, OFTReal);
259 3 : if (poDstLayer->CreateField(&oFieldDefn_Beg) != OGRERR_NONE)
260 : {
261 0 : CPLError(CE_Failure, CPLE_AppDefined, "Create %s field failed!",
262 : oFieldDefn_Beg.GetNameRef());
263 0 : return nullptr;
264 : }
265 :
266 6 : OGRFieldDefn oFieldDefn_End(FIELD_FINISH, OFTReal);
267 3 : if (poDstLayer->CreateField(&oFieldDefn_End) != OGRERR_NONE)
268 : {
269 0 : CPLError(CE_Failure, CPLE_AppDefined, "Create %s field failed!",
270 : oFieldDefn_End.GetNameRef());
271 0 : return nullptr;
272 : }
273 :
274 6 : OGRFieldDefn oFieldDefn_SF(FIELD_SCALE_FACTOR, OFTReal);
275 3 : if (poDstLayer->CreateField(&oFieldDefn_SF) != OGRERR_NONE)
276 : {
277 0 : CPLError(CE_Failure, CPLE_AppDefined, "Create %s field failed!",
278 : oFieldDefn_SF.GetNameRef());
279 0 : return nullptr;
280 : }
281 :
282 3 : if (pszOutputSepFieldName != nullptr)
283 : {
284 0 : OGRFieldDefn oSepField(pszOutputSepFieldName, OFTString);
285 0 : oSepField.SetWidth(254);
286 0 : if (poDstLayer->CreateField(&oSepField) != OGRERR_NONE)
287 : {
288 0 : CPLError(CE_Failure, CPLE_AppDefined, "Create %s field failed!",
289 : oSepField.GetNameRef());
290 0 : return nullptr;
291 : }
292 : }
293 :
294 : // Now that we've created a field, GetLayerDefn() won't return NULL.
295 3 : OGRFeatureDefn *poDstFDefn = poDstLayer->GetLayerDefn();
296 :
297 : // Sanity check: if it fails, the driver is buggy.
298 3 : if (poDstFDefn != nullptr && poDstFDefn->GetFieldCount() != 3)
299 : {
300 1 : CPLError(CE_Warning, CPLE_AppDefined,
301 : "The output driver has claimed to have added the %s field, "
302 : "but it did not!",
303 : oFieldDefn_Beg.GetNameRef());
304 : }
305 :
306 3 : return poDstLayer;
307 : }
308 :
309 : /* -------------------------------------------------------------------- */
310 : /* CheckDestDataSourceNameConsistency() */
311 : /* -------------------------------------------------------------------- */
312 :
313 3 : static void CheckDestDataSourceNameConsistency(const char *pszDestFilename,
314 : const char *pszDriverName)
315 : {
316 3 : char *pszDestExtension = CPLStrdup(CPLGetExtension(pszDestFilename));
317 :
318 : // TODO: Would be good to have driver metadata like for GDAL drivers.
319 : static const char *apszExtensions[][2] = {{"shp", "ESRI Shapefile"},
320 : {"dbf", "ESRI Shapefile"},
321 : {"sqlite", "SQLite"},
322 : {"db", "SQLite"},
323 : {"mif", "MapInfo File"},
324 : {"tab", "MapInfo File"},
325 : {"s57", "S57"},
326 : {"bna", "BNA"},
327 : {"csv", "CSV"},
328 : {"gml", "GML"},
329 : {"kml", "KML"},
330 : {"kmz", "LIBKML"},
331 : {"json", "GeoJSON"},
332 : {"geojson", "GeoJSON"},
333 : {"dxf", "DXF"},
334 : {"gdb", "FileGDB"},
335 : {"pix", "PCIDSK"},
336 : {"sql", "PGDump"},
337 : {"gtm", "GPSTrackMaker"},
338 : {"gmt", "GMT"},
339 : {"pdf", "PDF"},
340 : {nullptr, nullptr}};
341 : static const char *apszBeginName[][2] = {{"PG:", "PG"},
342 : {"MySQL:", "MySQL"},
343 : {"CouchDB:", "CouchDB"},
344 : {"GFT:", "GFT"},
345 : {"MSSQL:", "MSSQLSpatial"},
346 : {"ODBC:", "ODBC"},
347 : {"OCI:", "OCI"},
348 : {"SDE:", "SDE"},
349 : {"WFS:", "WFS"},
350 : {nullptr, nullptr}};
351 :
352 66 : for (int i = 0; apszExtensions[i][0] != nullptr; i++)
353 : {
354 63 : if (EQUAL(pszDestExtension, apszExtensions[i][0]) &&
355 3 : !EQUAL(pszDriverName, apszExtensions[i][1]))
356 : {
357 0 : fprintf(stderr,
358 : "Warning: The target file has a '%s' extension, "
359 : "which is normally used by the %s driver,\n"
360 : "but the requested output driver is %s. "
361 : "Is it really what you want ?\n",
362 : pszDestExtension, apszExtensions[i][1], pszDriverName);
363 0 : break;
364 : }
365 : }
366 :
367 30 : for (int i = 0; apszBeginName[i][0] != nullptr; i++)
368 : {
369 27 : if (EQUALN(pszDestFilename, apszBeginName[i][0],
370 0 : strlen(apszBeginName[i][0])) &&
371 0 : !EQUAL(pszDriverName, apszBeginName[i][1]))
372 : {
373 0 : fprintf(stderr,
374 : "Warning: The target file has a name which is normally "
375 : "recognized by the %s driver,\n"
376 : "but the requested output driver is %s. "
377 : "Is it really what you want ?\n",
378 : apszBeginName[i][1], pszDriverName);
379 0 : break;
380 : }
381 : }
382 :
383 3 : CPLFree(pszDestExtension);
384 3 : }
385 :
386 : //------------------------------------------------------------------------
387 : // AddFeature
388 : //------------------------------------------------------------------------
389 :
390 44 : static OGRErr AddFeature(OGRLayer *const poOutLayer, OGRLineString *pPart,
391 : double dfFrom, double dfTo, double dfScaleFactor,
392 : bool bQuiet,
393 : const char *pszOutputSepFieldName = nullptr,
394 : const char *pszOutputSepFieldValue = nullptr)
395 : {
396 : OGRFeature *poFeature =
397 44 : OGRFeature::CreateFeature(poOutLayer->GetLayerDefn());
398 :
399 44 : poFeature->SetField(FIELD_START, dfFrom);
400 44 : poFeature->SetField(FIELD_FINISH, dfTo);
401 44 : poFeature->SetField(FIELD_SCALE_FACTOR, dfScaleFactor);
402 :
403 44 : if (pszOutputSepFieldName != nullptr)
404 : {
405 0 : poFeature->SetField(pszOutputSepFieldName, pszOutputSepFieldValue);
406 : }
407 :
408 44 : poFeature->SetGeometryDirectly(pPart);
409 :
410 44 : if (poOutLayer->CreateFeature(poFeature) != OGRERR_NONE)
411 : {
412 0 : if (!bQuiet)
413 0 : printf("Failed to create feature in shapefile.\n");
414 0 : return OGRERR_FAILURE;
415 : }
416 :
417 44 : OGRFeature::DestroyFeature(poFeature);
418 :
419 44 : return OGRERR_NONE;
420 : }
421 :
422 : //------------------------------------------------------------------------
423 : // CreateSubline
424 : //------------------------------------------------------------------------
425 1 : static OGRErr CreateSubline(OGRLayer *const poPkLayer, double dfPosBeg,
426 : double dfPosEnd, OGRLayer *const poOutLayer,
427 : CPL_UNUSED int bDisplayProgress, bool bQuiet)
428 : {
429 : // Get step
430 1 : poPkLayer->ResetReading();
431 1 : OGRFeature *pFeature = poPkLayer->GetNextFeature();
432 1 : if (nullptr != pFeature)
433 : {
434 : // FIXME: Clang Static Analyzer rightly found that the following
435 : // code is dead
436 : // dfBeg = pFeature->GetFieldAsDouble(FIELD_START);
437 : // dfEnd = pFeature->GetFieldAsDouble(FIELD_FINISH);
438 1 : OGRFeature::DestroyFeature(pFeature);
439 : }
440 : else
441 : {
442 0 : fprintf(stderr, "Get step for positions %f - %f failed\n", dfPosBeg,
443 : dfPosEnd);
444 0 : return OGRERR_FAILURE;
445 : }
446 : // Get second part.
447 1 : double dfBeg = 0.0;
448 1 : double dfEnd = 0.0;
449 1 : pFeature = poPkLayer->GetNextFeature();
450 1 : if (nullptr != pFeature)
451 : {
452 1 : dfBeg = pFeature->GetFieldAsDouble(FIELD_START);
453 1 : dfEnd = pFeature->GetFieldAsDouble(FIELD_FINISH);
454 1 : OGRFeature::DestroyFeature(pFeature);
455 : }
456 : else
457 : {
458 0 : fprintf(stderr, "Get step for positions %f - %f failed\n", dfPosBeg,
459 : dfPosEnd);
460 0 : return OGRERR_FAILURE;
461 : }
462 1 : const double dfStep = dfEnd - dfBeg;
463 :
464 : // Round input to step
465 1 : const double dfPosBegLow = floor(dfPosBeg / dfStep) * dfStep;
466 1 : const double dfPosEndHigh = ceil(dfPosEnd / dfStep) * dfStep;
467 :
468 2 : CPLString szAttributeFilter;
469 : szAttributeFilter.Printf("%s >= %f AND %s <= %f", FIELD_START, dfPosBegLow,
470 1 : FIELD_FINISH, dfPosEndHigh);
471 : // TODO: ExecuteSQL should be faster.
472 1 : poPkLayer->SetAttributeFilter(szAttributeFilter);
473 1 : poPkLayer->ResetReading();
474 :
475 2 : std::map<double, OGRFeature *> moParts;
476 :
477 6 : while ((pFeature = poPkLayer->GetNextFeature()) != nullptr)
478 : {
479 5 : double dfStart = pFeature->GetFieldAsDouble(FIELD_START);
480 5 : moParts[dfStart] = pFeature;
481 : }
482 :
483 1 : if (moParts.empty())
484 : {
485 0 : fprintf(stderr, "Get parts for positions %f - %f failed\n", dfPosBeg,
486 : dfPosEnd);
487 0 : return OGRERR_FAILURE;
488 : }
489 :
490 2 : OGRLineString SubLine;
491 1 : if (moParts.size() == 1)
492 : {
493 0 : std::map<double, OGRFeature *>::iterator IT = moParts.begin();
494 0 : const double dfStart = IT->first;
495 0 : double dfPosBegCorr = dfPosBeg - dfStart;
496 0 : const double dfSF = IT->second->GetFieldAsDouble(FIELD_SCALE_FACTOR);
497 0 : dfPosBegCorr *= dfSF;
498 :
499 0 : const double dfPosEndCorr = (dfPosEnd - dfStart) * dfSF;
500 :
501 0 : OGRLineString *pLine = IT->second->GetGeometryRef()->toLineString();
502 :
503 : OGRLineString *pSubLine =
504 0 : pLine->getSubLine(dfPosBegCorr, dfPosEndCorr, FALSE);
505 :
506 0 : OGRFeature::DestroyFeature(IT->second);
507 : // Store.
508 0 : return AddFeature(poOutLayer, pSubLine, dfPosBeg, dfPosEnd, 1.0,
509 0 : bQuiet);
510 : }
511 : else
512 : {
513 1 : int nCounter = static_cast<int>(moParts.size());
514 1 : std::map<double, OGRFeature *>::iterator IT = moParts.begin();
515 1 : OGRLineString *pOutLine = new OGRLineString();
516 : // Get first part.
517 1 : const double dfStart = IT->first;
518 1 : double dfPosBegCorr = dfPosBeg - dfStart;
519 1 : double dfSF = IT->second->GetFieldAsDouble(FIELD_SCALE_FACTOR);
520 1 : dfPosBegCorr *= dfSF;
521 :
522 1 : OGRLineString *pLine = IT->second->GetGeometryRef()->toLineString();
523 :
524 : OGRLineString *pSubLine =
525 1 : pLine->getSubLine(dfPosBegCorr, pLine->get_Length(), FALSE);
526 :
527 1 : pOutLine->addSubLineString(pSubLine);
528 1 : delete pSubLine;
529 1 : OGRFeature::DestroyFeature(IT->second);
530 :
531 1 : ++IT;
532 1 : nCounter--;
533 :
534 4 : while (nCounter > 1)
535 : {
536 3 : pLine = IT->second->GetGeometryRef()->toLineString();
537 3 : pOutLine->addSubLineString(pLine);
538 3 : OGRFeature::DestroyFeature(IT->second);
539 3 : ++IT;
540 3 : nCounter--;
541 : }
542 :
543 : // Get last part
544 1 : double dfPosEndCorr = dfPosEnd - IT->first;
545 1 : dfSF = IT->second->GetFieldAsDouble(FIELD_SCALE_FACTOR);
546 1 : dfPosEndCorr *= dfSF;
547 :
548 1 : pLine = IT->second->GetGeometryRef()->toLineString();
549 :
550 1 : pSubLine = pLine->getSubLine(0, dfPosEndCorr, FALSE);
551 :
552 1 : pOutLine->addSubLineString(pSubLine);
553 1 : delete pSubLine;
554 :
555 1 : OGRFeature::DestroyFeature(IT->second);
556 : // Store
557 1 : return AddFeature(poOutLayer, pOutLine, dfPosBeg, dfPosEnd, 1.0,
558 1 : bQuiet);
559 : }
560 : }
561 :
562 : //------------------------------------------------------------------------
563 : // Project
564 : //------------------------------------------------------------------------
565 : #ifdef HAVE_GEOS
566 66 : static double Project(OGRLineString *pLine, OGRPoint *pPoint)
567 : {
568 66 : if (nullptr == pLine || nullptr == pPoint)
569 0 : return -1;
570 132 : OGRPoint TestPoint;
571 66 : pLine->StartPoint(&TestPoint);
572 66 : if (TestPoint.Equals(pPoint))
573 2 : return 0.0;
574 64 : pLine->EndPoint(&TestPoint);
575 64 : if (TestPoint.Equals(pPoint))
576 2 : return pLine->get_Length();
577 :
578 62 : return pLine->Project(pPoint);
579 : }
580 : #endif
581 :
582 : //------------------------------------------------------------------------
583 : // CreatePartsFromLineString
584 : //------------------------------------------------------------------------
585 : #ifdef HAVE_GEOS
586 2 : static OGRErr CreatePartsFromLineString(
587 : OGRLineString *pPathGeom, OGRLayer *const poPkLayer, int nMValField,
588 : double dfStep, OGRLayer *const poOutLayer, int bDisplayProgress,
589 : bool bQuiet, const char *pszOutputSepFieldName = nullptr,
590 : const char *pszOutputSepFieldValue = nullptr)
591 : {
592 : // Check repers/milestones/reference points type
593 2 : OGRwkbGeometryType eGeomType = poPkLayer->GetGeomType();
594 2 : if (wkbFlatten(eGeomType) != wkbPoint)
595 : {
596 0 : fprintf(stderr, "Unsupported geometry type %s for path\n",
597 : OGRGeometryTypeToName(eGeomType));
598 0 : return OGRERR_FAILURE;
599 : }
600 :
601 2 : double dfTolerance = 1.0;
602 2 : const OGRSpatialReference *pSpaRef = pPathGeom->getSpatialReference();
603 2 : if (pSpaRef->IsGeographic())
604 : {
605 2 : dfTolerance = TOLERANCE_DEGREE;
606 : }
607 : else
608 : {
609 0 : dfTolerance = TOLERANCE_METER;
610 : }
611 :
612 : // Create sorted list of repers.
613 4 : std::map<double, OGRPoint *> moRepers;
614 2 : poPkLayer->ResetReading();
615 2 : OGRFeature *pReperFeature = nullptr;
616 16 : while ((pReperFeature = poPkLayer->GetNextFeature()) != nullptr)
617 : {
618 14 : const double dfReperPos = pReperFeature->GetFieldAsDouble(nMValField);
619 14 : OGRGeometry *pGeom = pReperFeature->GetGeometryRef();
620 14 : if (nullptr != pGeom)
621 : {
622 14 : OGRPoint *pPt = pGeom->clone()->toPoint();
623 14 : if (!bQuiet)
624 : {
625 14 : if (moRepers.find(dfReperPos) != moRepers.end())
626 : {
627 0 : CPLError(
628 : CE_Warning, CPLE_AppDefined,
629 : "The distance %f is already present in repers file!",
630 : dfReperPos);
631 : }
632 : }
633 : // Check if reper is inside the path
634 14 : const double dfTestDistance = Project(pPathGeom, pPt);
635 14 : if (dfTestDistance < 0)
636 : {
637 0 : if (!bQuiet)
638 : {
639 0 : CPLError(CE_Warning, CPLE_AppDefined,
640 : "The distance %f is out of path!", dfReperPos);
641 : }
642 0 : delete pPt;
643 : }
644 : else
645 : {
646 14 : const double dfDist = pPathGeom->Distance(pPt);
647 14 : if (dfDist < dfTolerance)
648 6 : moRepers[dfReperPos] = pPt;
649 : else
650 8 : delete pPt;
651 : }
652 : }
653 14 : OGRFeature::DestroyFeature(pReperFeature);
654 : }
655 :
656 2 : if (moRepers.size() < 2)
657 : {
658 0 : fprintf(stderr, "Not enough repers to proceed\n");
659 0 : return OGRERR_FAILURE;
660 : }
661 :
662 : // Check direction.
663 2 : if (!bQuiet)
664 : {
665 2 : fprintf(stdout, "Check path direction\n");
666 : }
667 :
668 : // Get distance along path from pt1 and pt2.
669 : // If pt1 distance > pt2 distance, reverse path
670 2 : std::map<double, OGRPoint *>::const_iterator IT;
671 2 : IT = moRepers.begin();
672 2 : double dfPosition = IT->first;
673 2 : double dfBeginPosition = IT->first;
674 2 : OGRPoint *pt1 = IT->second;
675 2 : ++IT;
676 2 : OGRPoint *pt2 = IT->second;
677 :
678 2 : double dfDistance1 = Project(pPathGeom, pt1);
679 2 : double dfDistance2 = Project(pPathGeom, pt2);
680 :
681 2 : if (dfDistance1 > dfDistance2)
682 : {
683 0 : if (!bQuiet)
684 : {
685 0 : fprintf(stderr,
686 : "Warning: The path is opposite the repers direction. "
687 : "Let's reverse path\n");
688 : }
689 0 : pPathGeom->reversePoints();
690 :
691 0 : dfDistance1 = Project(pPathGeom, pt1);
692 0 : dfDistance2 = Project(pPathGeom, pt2);
693 : }
694 :
695 2 : OGRLineString *pPart = nullptr;
696 :
697 2 : std::vector<CURVE_DATA> astSubLines;
698 :
699 2 : if (!bQuiet)
700 : {
701 2 : fprintf(stdout, "Create parts\n");
702 : }
703 :
704 : // Get first part
705 : // If first point is not at the beginning of the path
706 : // The first part should be from the beginning of the path
707 : // to the first point. length == part.getLength
708 2 : OGRPoint *pPtBeg = nullptr;
709 2 : OGRPoint *pPtEnd = nullptr;
710 2 : double dfPtBegPosition = 0.0;
711 2 : double dfPtEndPosition = 0.0;
712 :
713 2 : if (dfDistance1 > DELTA)
714 : {
715 2 : pPart = pPathGeom->getSubLine(0, dfDistance1, FALSE);
716 2 : if (nullptr != pPart)
717 : {
718 2 : double dfLen = pPart->get_Length();
719 2 : if (pSpaRef->IsGeographic())
720 : {
721 : // convert to UTM/WGS84
722 4 : OGRPoint pt;
723 2 : pPart->Value(dfLen / 2, &pt);
724 : const int nZoneEnv =
725 2 : static_cast<int>(30 + (pt.getX() + 3.0) / 6.0 + 0.5);
726 : int nEPSG;
727 2 : if (pt.getY() > 0)
728 : {
729 2 : nEPSG = 32600 + nZoneEnv;
730 : }
731 : else
732 : {
733 0 : nEPSG = 32700 + nZoneEnv;
734 : }
735 4 : OGRSpatialReference SpatRef;
736 2 : SpatRef.importFromEPSG(nEPSG);
737 2 : OGRGeometry *pTransformPart = pPart->clone();
738 2 : if (pTransformPart->transformTo(&SpatRef) == OGRERR_NONE)
739 : {
740 : OGRLineString *pTransformPartLS =
741 2 : pTransformPart->toLineString();
742 2 : dfLen = pTransformPartLS->get_Length();
743 : }
744 :
745 2 : CURVE_DATA data = {pPart, dfPosition - dfLen, dfPosition,
746 2 : pPart->get_Length() / dfLen};
747 2 : astSubLines.push_back(data);
748 :
749 2 : pPtBeg = new OGRPoint();
750 2 : pPart->getPoint(0, pPtBeg);
751 2 : dfPtBegPosition = dfPosition - dfLen;
752 :
753 : // AddFeature(poOutLayer, pPart, dfPosition - dfLen, dfPosition,
754 : // pPart->get_Length() / dfLen, bQuiet);
755 2 : delete pTransformPart;
756 : }
757 : else
758 : {
759 0 : CURVE_DATA data = {pPart, dfPosition - dfLen, dfPosition, 1.0};
760 0 : astSubLines.push_back(data);
761 : // AddFeature(poOutLayer, pPart, dfPosition - dfLen, dfPosition,
762 : // 1.0, bQuiet);
763 0 : pPtBeg = new OGRPoint();
764 0 : pPart->getPoint(0, pPtBeg);
765 0 : dfPtBegPosition = dfPosition - dfLen;
766 : }
767 : }
768 : }
769 :
770 2 : if (dfDistance2 - dfDistance1 > DELTA)
771 : {
772 2 : pPart = pPathGeom->getSubLine(dfDistance1, dfDistance2, FALSE);
773 2 : if (nullptr != pPart)
774 : {
775 4 : CURVE_DATA data = {pPart, dfPosition, IT->first,
776 2 : pPart->get_Length() / (IT->first - dfPosition)};
777 2 : astSubLines.push_back(data);
778 : // AddFeature(poOutLayer, pPart, dfPosition, IT->first,
779 : // pPart->get_Length() / (IT->first - dfPosition),
780 : // bQuiet);
781 : }
782 : }
783 :
784 2 : GDALProgressFunc pfnProgress = nullptr;
785 2 : void *pProgressArg = nullptr;
786 :
787 2 : double dfFactor = 1.0 / moRepers.size();
788 2 : if (bDisplayProgress)
789 : {
790 0 : pfnProgress = GDALScaledProgress;
791 : pProgressArg =
792 0 : GDALCreateScaledProgress(0.0, 1.0, GDALTermProgress, nullptr);
793 : }
794 :
795 2 : int nCount = 2;
796 2 : dfDistance1 = dfDistance2;
797 2 : dfPosition = IT->first;
798 2 : ++IT; // Get third point
799 :
800 2 : double dfEndPosition = 0.0;
801 4 : while (IT != moRepers.end())
802 : {
803 2 : if (bDisplayProgress)
804 : {
805 0 : pfnProgress(nCount * dfFactor, "", pProgressArg);
806 0 : nCount++;
807 : }
808 :
809 2 : dfEndPosition = IT->first;
810 :
811 2 : dfDistance2 = Project(pPathGeom, IT->second);
812 :
813 2 : if (dfDistance2 - dfDistance1 > DELTA)
814 : {
815 2 : pPart = pPathGeom->getSubLine(dfDistance1, dfDistance2, FALSE);
816 2 : if (nullptr != pPart)
817 : {
818 4 : CURVE_DATA data = {pPart, dfPosition, IT->first,
819 6 : pPart->get_Length() /
820 2 : (IT->first - dfPosition)};
821 2 : astSubLines.push_back(data);
822 : // AddFeature(poOutLayer, pPart, dfPosition, IT->first,
823 : // pPart->get_Length() / (IT->first - dfPosition),
824 : // bQuiet);
825 2 : dfDistance1 = dfDistance2;
826 2 : dfPosition = IT->first;
827 : }
828 : }
829 :
830 2 : ++IT;
831 : }
832 :
833 : // Get last part
834 2 : if (pPathGeom->get_Length() - dfDistance1 > DELTA)
835 : {
836 : pPart =
837 2 : pPathGeom->getSubLine(dfDistance1, pPathGeom->get_Length(), FALSE);
838 2 : if (nullptr != pPart)
839 : {
840 2 : double dfLen = pPart->get_Length();
841 2 : if (pSpaRef->IsGeographic())
842 : {
843 : // convert to UTM/WGS84
844 4 : OGRPoint pt;
845 2 : pPart->Value(dfLen / 2, &pt);
846 : const int nZoneEnv =
847 2 : static_cast<int>(30 + (pt.getX() + 3.0) / 6.0 + 0.5);
848 : int nEPSG;
849 2 : if (pt.getY() > 0)
850 : {
851 2 : nEPSG = 32600 + nZoneEnv;
852 : }
853 : else
854 : {
855 0 : nEPSG = 32700 + nZoneEnv;
856 : }
857 4 : OGRSpatialReference SpatRef;
858 2 : SpatRef.importFromEPSG(nEPSG);
859 2 : OGRGeometry *pTransformPart = pPart->clone();
860 2 : if (pTransformPart->transformTo(&SpatRef) == OGRERR_NONE)
861 : {
862 : OGRLineString *pTransformPartLS =
863 2 : pTransformPart->toLineString();
864 2 : dfLen = pTransformPartLS->get_Length();
865 : }
866 2 : CURVE_DATA data = {pPart, dfPosition, dfPosition + dfLen,
867 2 : pPart->get_Length() / dfLen};
868 2 : astSubLines.push_back(data);
869 : // AddFeature(poOutLayer, pPart, dfPosition, dfPosition + dfLen,
870 : // pPart->get_Length() / dfLen, bQuiet);
871 :
872 2 : pPtEnd = new OGRPoint();
873 2 : pPart->getPoint(pPart->getNumPoints() - 1, pPtEnd);
874 2 : dfPtEndPosition = dfPosition + dfLen;
875 :
876 2 : delete pTransformPart;
877 : }
878 : else
879 : {
880 0 : CURVE_DATA data = {pPart, dfPosition, dfPosition + dfLen, 1.0};
881 0 : astSubLines.push_back(data);
882 : // AddFeature(poOutLayer, pPart, dfPosition - dfLen, dfPosition,
883 : // 1.0, bQuiet);
884 0 : pPtEnd = new OGRPoint();
885 0 : pPart->getPoint(pPart->getNumPoints() - 1, pPtEnd);
886 0 : dfPtEndPosition = dfPosition + dfLen;
887 : }
888 : }
889 : }
890 :
891 : // Create pickets
892 2 : if (!bQuiet)
893 : {
894 2 : fprintf(stdout, "\nCreate pickets\n");
895 : }
896 :
897 2 : const double dfRoundBeg = pPtBeg != nullptr
898 2 : ? ceil(dfPtBegPosition / dfStep) * dfStep
899 0 : : ceil(dfBeginPosition / dfStep) * dfStep;
900 :
901 2 : if (pPtEnd != nullptr)
902 2 : dfEndPosition = dfPtEndPosition;
903 :
904 2 : dfFactor = dfStep / (dfEndPosition - dfRoundBeg);
905 2 : nCount = 0;
906 8 : for (std::map<double, OGRPoint *>::iterator oIter = moRepers.begin();
907 14 : oIter != moRepers.end(); ++oIter)
908 : {
909 6 : delete oIter->second;
910 : }
911 :
912 2 : moRepers.clear();
913 :
914 2 : if (pPtBeg != nullptr)
915 2 : moRepers[dfPtBegPosition] = pPtBeg;
916 2 : if (pPtEnd != nullptr)
917 2 : moRepers[dfPtEndPosition] = pPtEnd;
918 :
919 43 : for (double dfDist = dfRoundBeg; dfDist <= dfEndPosition; dfDist += dfStep)
920 : {
921 41 : if (bDisplayProgress)
922 : {
923 0 : pfnProgress(nCount * dfFactor, "", pProgressArg);
924 0 : nCount++;
925 : }
926 :
927 70 : for (int j = 0; j < static_cast<int>(astSubLines.size()); j++)
928 : {
929 70 : if (astSubLines[j].IsInside(dfDist))
930 : {
931 : const double dfRealDist =
932 41 : (dfDist - astSubLines[j].dfBeg) * astSubLines[j].dfFactor;
933 41 : OGRPoint *pReperPoint = new OGRPoint();
934 41 : astSubLines[j].pPart->Value(dfRealDist, pReperPoint);
935 :
936 41 : moRepers[dfDist] = pReperPoint;
937 41 : break;
938 : }
939 : }
940 : }
941 :
942 10 : for (int i = 0; i < static_cast<int>(astSubLines.size()); i++)
943 : {
944 8 : delete astSubLines[i].pPart;
945 : }
946 2 : astSubLines.clear();
947 :
948 2 : if (!bQuiet)
949 : {
950 2 : fprintf(stdout, "\nCreate sublines\n");
951 : }
952 :
953 2 : IT = moRepers.begin();
954 2 : dfFactor = 1.0 / moRepers.size();
955 2 : nCount = 0;
956 2 : dfDistance1 = 0;
957 2 : dfPosition = IT->first;
958 :
959 47 : while (IT != moRepers.end())
960 : {
961 45 : if (bDisplayProgress)
962 : {
963 0 : pfnProgress(nCount * dfFactor, "", pProgressArg);
964 0 : nCount++;
965 : }
966 :
967 45 : dfDistance2 = Project(pPathGeom, IT->second);
968 :
969 45 : if (dfDistance2 - dfDistance1 > DELTA)
970 : {
971 43 : pPart = pPathGeom->getSubLine(dfDistance1, dfDistance2, FALSE);
972 43 : if (nullptr != pPart)
973 : {
974 43 : AddFeature(poOutLayer, pPart, dfPosition, IT->first,
975 43 : pPart->get_Length() / (IT->first - dfPosition),
976 : bQuiet, pszOutputSepFieldName,
977 : pszOutputSepFieldValue);
978 43 : dfDistance1 = dfDistance2;
979 43 : dfPosition = IT->first;
980 : }
981 : }
982 :
983 45 : ++IT;
984 : }
985 :
986 47 : for (std::map<double, OGRPoint *>::iterator oIter = moRepers.begin();
987 92 : oIter != moRepers.end(); ++oIter)
988 : {
989 45 : delete oIter->second;
990 : }
991 :
992 2 : if (!bQuiet)
993 : {
994 2 : fprintf(stdout, "\nSuccess!\n\n");
995 : }
996 :
997 2 : if (nullptr != pProgressArg)
998 : {
999 0 : GDALDestroyScaledProgress(pProgressArg);
1000 : }
1001 :
1002 2 : return OGRERR_NONE;
1003 : }
1004 : #endif
1005 :
1006 : //------------------------------------------------------------------------
1007 : // CreateParts
1008 : //------------------------------------------------------------------------
1009 : #ifdef HAVE_GEOS
1010 2 : static OGRErr CreateParts(OGRLayer *const poLnLayer, OGRLayer *const poPkLayer,
1011 : int nMValField, double dfStep,
1012 : OGRLayer *const poOutLayer, int bDisplayProgress,
1013 : bool bQuiet,
1014 : const char *pszOutputSepFieldName = nullptr,
1015 : const char *pszOutputSepFieldValue = nullptr)
1016 : {
1017 2 : OGRErr eRetCode = OGRERR_FAILURE;
1018 :
1019 : // check path and get first line
1020 2 : OGRwkbGeometryType eGeomType = poLnLayer->GetGeomType();
1021 2 : if (wkbFlatten(eGeomType) != wkbLineString &&
1022 0 : wkbFlatten(eGeomType) != wkbMultiLineString)
1023 : {
1024 0 : fprintf(stderr, "Unsupported geometry type %s for path\n",
1025 : OGRGeometryTypeToName(eGeomType));
1026 0 : return eRetCode;
1027 : }
1028 :
1029 2 : poLnLayer->ResetReading();
1030 : // Get first geometry
1031 : // TODO: Attribute filter for path geometry.
1032 2 : OGRFeature *pPathFeature = poLnLayer->GetNextFeature();
1033 2 : if (nullptr != pPathFeature)
1034 : {
1035 2 : OGRGeometry *pGeom = pPathFeature->GetGeometryRef();
1036 :
1037 4 : if (pGeom != nullptr &&
1038 2 : wkbFlatten(pGeom->getGeometryType()) == wkbMultiLineString)
1039 : {
1040 0 : if (!bQuiet)
1041 : {
1042 0 : fprintf(stdout,
1043 : "\nThe geometry " CPL_FRMT_GIB
1044 : " is wkbMultiLineString type\n",
1045 : pPathFeature->GetFID());
1046 : }
1047 :
1048 0 : OGRGeometryCollection *pGeomColl = pGeom->toGeometryCollection();
1049 0 : for (int i = 0; i < pGeomColl->getNumGeometries(); ++i)
1050 : {
1051 : OGRLineString *pPath =
1052 0 : pGeomColl->getGeometryRef(i)->clone()->toLineString();
1053 0 : pPath->assignSpatialReference(pGeomColl->getSpatialReference());
1054 0 : eRetCode = CreatePartsFromLineString(
1055 : pPath, poPkLayer, nMValField, dfStep, poOutLayer,
1056 : bDisplayProgress, bQuiet, pszOutputSepFieldName,
1057 : pszOutputSepFieldValue);
1058 :
1059 0 : if (eRetCode != OGRERR_NONE)
1060 : {
1061 0 : OGRFeature::DestroyFeature(pPathFeature);
1062 0 : return eRetCode;
1063 : }
1064 : }
1065 : }
1066 4 : else if (pGeom != nullptr &&
1067 2 : wkbFlatten(pGeom->getGeometryType()) == wkbLineString)
1068 : {
1069 2 : OGRLineString *pGeomClone = pGeom->clone()->toLineString();
1070 2 : eRetCode = CreatePartsFromLineString(
1071 : pGeomClone, poPkLayer, nMValField, dfStep, poOutLayer,
1072 : bDisplayProgress, bQuiet, pszOutputSepFieldName,
1073 : pszOutputSepFieldValue);
1074 2 : delete pGeomClone;
1075 : }
1076 :
1077 2 : OGRFeature::DestroyFeature(pPathFeature);
1078 : }
1079 :
1080 : // Should never reach
1081 :
1082 2 : return eRetCode;
1083 : }
1084 : #endif
1085 :
1086 : //------------------------------------------------------------------------
1087 : // CreatePartsMultiple
1088 : //------------------------------------------------------------------------
1089 : #ifdef HAVE_GEOS
1090 0 : static OGRErr CreatePartsMultiple(
1091 : OGRLayer *const poLnLayer, const char *pszLineSepFieldName,
1092 : OGRLayer *const poPkLayer, const char *pszPicketsSepFieldName,
1093 : int nMValField, double dfStep, OGRLayer *const poOutLayer,
1094 : const char *pszOutputSepFieldName, int bDisplayProgress, bool bQuiet)
1095 : {
1096 : // Read all separate field values into array
1097 0 : OGRFeatureDefn *pDefn = poLnLayer->GetLayerDefn();
1098 0 : const int nLineSepFieldInd = pDefn->GetFieldIndex(pszLineSepFieldName);
1099 0 : if (nLineSepFieldInd == -1)
1100 : {
1101 0 : fprintf(stderr, "The field %s not found\n", pszLineSepFieldName);
1102 0 : return OGRERR_FAILURE;
1103 : }
1104 :
1105 0 : poLnLayer->ResetReading();
1106 :
1107 0 : std::set<CPLString> asIDs;
1108 0 : OGRFeature *pFeature = nullptr;
1109 0 : while ((pFeature = poLnLayer->GetNextFeature()) != nullptr)
1110 : {
1111 0 : CPLString sID = pFeature->GetFieldAsString(nLineSepFieldInd);
1112 0 : asIDs.insert(sID);
1113 :
1114 0 : OGRFeature::DestroyFeature(pFeature);
1115 : }
1116 :
1117 0 : for (std::set<CPLString>::const_iterator it = asIDs.begin();
1118 0 : it != asIDs.end(); ++it)
1119 : {
1120 : // Create select clause
1121 0 : CPLString sLineWhere;
1122 0 : sLineWhere.Printf("%s = '%s'", pszLineSepFieldName, it->c_str());
1123 0 : poLnLayer->SetAttributeFilter(sLineWhere);
1124 :
1125 0 : CPLString sPkWhere;
1126 0 : sPkWhere.Printf("%s = '%s'", pszPicketsSepFieldName, it->c_str());
1127 0 : poPkLayer->SetAttributeFilter(sPkWhere);
1128 :
1129 0 : if (!bQuiet)
1130 : {
1131 0 : fprintf(stdout, "The %s %s\n", pszPicketsSepFieldName, it->c_str());
1132 : }
1133 :
1134 : // Don't check success as we want to try all paths
1135 0 : CreateParts(poLnLayer, poPkLayer, nMValField, dfStep, poOutLayer,
1136 0 : bDisplayProgress, bQuiet, pszOutputSepFieldName, *it);
1137 : }
1138 :
1139 0 : return OGRERR_NONE;
1140 : }
1141 : #endif
1142 :
1143 : //------------------------------------------------------------------------
1144 : // GetPosition
1145 : //------------------------------------------------------------------------
1146 : #ifdef HAVE_GEOS
1147 1 : static OGRErr GetPosition(OGRLayer *const poPkLayer, double dfX, double dfY,
1148 : int /* bDisplayProgress */, int bQuiet)
1149 : {
1150 : // Create point
1151 2 : OGRPoint pt;
1152 1 : pt.setX(dfX);
1153 1 : pt.setY(dfY);
1154 1 : pt.assignSpatialReference(poPkLayer->GetSpatialRef());
1155 :
1156 1 : poPkLayer->ResetReading();
1157 1 : OGRLineString *pCloserPart = nullptr;
1158 1 : double dfBeg = 0.0;
1159 1 : double dfScale = 0.0;
1160 1 : double dfMinDistance = std::numeric_limits<double>::max();
1161 1 : OGRFeature *pFeature = nullptr;
1162 10 : while ((pFeature = poPkLayer->GetNextFeature()) != nullptr)
1163 : {
1164 9 : OGRGeometry *pCurrentGeom = pFeature->GetGeometryRef();
1165 9 : if (pCurrentGeom != nullptr)
1166 : {
1167 9 : double dfCurrentDistance = pCurrentGeom->Distance(&pt);
1168 9 : if (dfCurrentDistance < dfMinDistance)
1169 : {
1170 5 : dfMinDistance = dfCurrentDistance;
1171 5 : if (pCloserPart != nullptr)
1172 4 : delete pCloserPart;
1173 5 : pCloserPart = pFeature->StealGeometry()->toLineString();
1174 5 : dfBeg = pFeature->GetFieldAsDouble(FIELD_START);
1175 5 : dfScale = pFeature->GetFieldAsDouble(FIELD_SCALE_FACTOR);
1176 : }
1177 : }
1178 9 : OGRFeature::DestroyFeature(pFeature);
1179 : }
1180 :
1181 1 : if (nullptr == pCloserPart)
1182 : {
1183 0 : fprintf(stderr, "Filed to find closest part\n");
1184 0 : return OGRERR_FAILURE;
1185 : }
1186 : // Now we have closest part
1187 : // Get real distance
1188 1 : const double dfRealDist = Project(pCloserPart, &pt);
1189 1 : delete pCloserPart;
1190 : // Compute reference distance
1191 1 : const double dfRefDist = dfBeg + dfRealDist / dfScale;
1192 1 : if (bQuiet)
1193 : {
1194 1 : fprintf(stdout, "%s", CPLSPrintf("%f\n", dfRefDist));
1195 : }
1196 : else
1197 : {
1198 0 : fprintf(
1199 : stdout, "%s",
1200 : CPLSPrintf("The position for coordinates lat:%f, long:%f is %f\n",
1201 : dfY, dfX, dfRefDist));
1202 : }
1203 :
1204 1 : return OGRERR_NONE;
1205 : }
1206 : #endif
1207 :
1208 : //------------------------------------------------------------------------
1209 : // GetCoordinates
1210 : //------------------------------------------------------------------------
1211 1 : static OGRErr GetCoordinates(OGRLayer *const poPkLayer, double dfPos,
1212 : /* CPL_UNUSED */ int /* bDisplayProgress */,
1213 : bool bQuiet)
1214 : {
1215 2 : CPLString szAttributeFilter;
1216 : szAttributeFilter.Printf("%s < %f AND %s > %f", FIELD_START, dfPos,
1217 1 : FIELD_FINISH, dfPos);
1218 : // TODO: ExecuteSQL should be faster.
1219 1 : poPkLayer->SetAttributeFilter(szAttributeFilter);
1220 1 : poPkLayer->ResetReading();
1221 :
1222 1 : bool bHaveCoords = false;
1223 2 : for (auto &pFeature : poPkLayer)
1224 : {
1225 1 : bHaveCoords = true;
1226 1 : const double dfStart = pFeature->GetFieldAsDouble(FIELD_START);
1227 1 : double dfPosCorr = dfPos - dfStart;
1228 1 : const double dfSF = pFeature->GetFieldAsDouble(FIELD_SCALE_FACTOR);
1229 1 : dfPosCorr *= dfSF;
1230 1 : OGRLineString *pLine = pFeature->GetGeometryRef()->toLineString();
1231 :
1232 2 : OGRPoint pt;
1233 1 : pLine->Value(dfPosCorr, &pt);
1234 :
1235 1 : if (bQuiet)
1236 : {
1237 1 : fprintf(stdout, "%s",
1238 : CPLSPrintf("%f,%f,%f\n", pt.getX(), pt.getY(), pt.getZ()));
1239 : }
1240 : else
1241 : {
1242 0 : fprintf(stdout, "%s",
1243 : CPLSPrintf("The position for distance %f is lat:%f, "
1244 : "long:%f, height:%f\n",
1245 : dfPos, pt.getY(), pt.getX(), pt.getZ()));
1246 : }
1247 : }
1248 :
1249 1 : if (bHaveCoords)
1250 : {
1251 1 : return OGRERR_NONE;
1252 : }
1253 : else
1254 : {
1255 0 : fprintf(stderr, "Get coordinates for position %f failed\n", dfPos);
1256 0 : return OGRERR_FAILURE;
1257 : }
1258 : }
1259 :
1260 : /************************************************************************/
1261 : /* main() */
1262 : /************************************************************************/
1263 :
1264 : #define CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(nExtraArg) \
1265 : do \
1266 : { \
1267 : if (iArg + nExtraArg >= nArgc) \
1268 : Usage(true, CPLSPrintf("%s option requires %d argument(s)", \
1269 : papszArgv[iArg], nExtraArg)); \
1270 : } while (false)
1271 :
1272 6 : MAIN_START(nArgc, papszArgv)
1273 :
1274 : {
1275 6 : OGRErr eErr = OGRERR_NONE;
1276 6 : bool bQuiet = false;
1277 6 : const char *pszFormat = "ESRI Shapefile";
1278 :
1279 6 : const char *pszOutputDataSource = nullptr;
1280 6 : const char *pszLineDataSource = nullptr;
1281 6 : const char *pszPicketsDataSource = nullptr;
1282 6 : const char *pszPartsDataSource = nullptr;
1283 6 : char *pszOutputLayerName = nullptr;
1284 6 : const char *pszLineLayerName = nullptr;
1285 : #ifdef HAVE_GEOS
1286 6 : const char *pszPicketsLayerName = nullptr;
1287 6 : const char *pszPicketsMField = nullptr;
1288 : #endif
1289 6 : const char *pszPartsLayerName = nullptr;
1290 :
1291 : #ifdef HAVE_GEOS
1292 6 : const char *pszLineSepFieldName = nullptr;
1293 6 : const char *pszPicketsSepFieldName = nullptr;
1294 6 : const char *pszOutputSepFieldName = "uniq_uid";
1295 : #endif
1296 :
1297 6 : char **papszDSCO = nullptr;
1298 6 : char **papszLCO = nullptr;
1299 :
1300 6 : operation stOper = op_unknown;
1301 : #ifdef HAVE_GEOS
1302 6 : double dfX = -100000000.0;
1303 6 : double dfY = -100000000.0;
1304 : #endif
1305 6 : double dfPos = -100000000.0;
1306 :
1307 6 : int bDisplayProgress = FALSE;
1308 :
1309 6 : double dfPosBeg = -100000000.0;
1310 6 : double dfPosEnd = -100000000.0;
1311 : #ifdef HAVE_GEOS
1312 6 : double dfStep = -100000000.0;
1313 : #endif
1314 :
1315 : // Check strict compilation and runtime library version as we use C++ API.
1316 6 : if (!GDAL_CHECK_VERSION(papszArgv[0]))
1317 0 : exit(1);
1318 :
1319 6 : EarlySetConfigOptions(nArgc, papszArgv);
1320 :
1321 6 : OGRRegisterAll();
1322 :
1323 : /* -------------------------------------------------------------------- */
1324 : /* Processing command line arguments. */
1325 : /* -------------------------------------------------------------------- */
1326 6 : nArgc = OGRGeneralCmdLineProcessor(nArgc, &papszArgv, 0);
1327 :
1328 6 : if (nArgc < 1)
1329 0 : exit(-nArgc);
1330 :
1331 33 : for (int iArg = 1; iArg < nArgc; iArg++)
1332 : {
1333 28 : if (EQUAL(papszArgv[iArg], "--utility_version"))
1334 : {
1335 1 : printf("%s was compiled against GDAL %s and "
1336 : "is running against GDAL %s\n",
1337 : papszArgv[0], GDAL_RELEASE_NAME,
1338 : GDALVersionInfo("RELEASE_NAME"));
1339 1 : CSLDestroy(papszArgv);
1340 1 : return 0;
1341 : }
1342 27 : else if (EQUAL(papszArgv[iArg], "--help"))
1343 : {
1344 0 : Usage(false);
1345 : }
1346 27 : else if (EQUAL(papszArgv[iArg], "--long-usage"))
1347 : {
1348 0 : Usage(false, nullptr, false);
1349 : }
1350 :
1351 27 : else if (EQUAL(papszArgv[iArg], "-q") ||
1352 27 : EQUAL(papszArgv[iArg], "-quiet"))
1353 : {
1354 2 : bQuiet = true;
1355 : }
1356 25 : else if ((EQUAL(papszArgv[iArg], "-f") ||
1357 24 : EQUAL(papszArgv[iArg], "-of")))
1358 : {
1359 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1360 : // bFormatExplicitlySet = TRUE;
1361 : // coverity[tainted_data]
1362 1 : pszFormat = papszArgv[++iArg];
1363 : }
1364 24 : else if (EQUAL(papszArgv[iArg], "-dsco"))
1365 : {
1366 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1367 : // coverity[tainted_data]
1368 0 : papszDSCO = CSLAddString(papszDSCO, papszArgv[++iArg]);
1369 : }
1370 24 : else if (EQUAL(papszArgv[iArg], "-lco"))
1371 : {
1372 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1373 : // coverity[tainted_data]
1374 0 : papszLCO = CSLAddString(papszLCO, papszArgv[++iArg]);
1375 : }
1376 24 : else if (EQUAL(papszArgv[iArg], "-create"))
1377 : {
1378 2 : stOper = op_create;
1379 : }
1380 22 : else if (EQUAL(papszArgv[iArg], "-get_pos"))
1381 : {
1382 1 : stOper = op_get_pos;
1383 : }
1384 21 : else if (EQUAL(papszArgv[iArg], "-get_coord"))
1385 : {
1386 1 : stOper = op_get_coord;
1387 : }
1388 20 : else if (EQUAL(papszArgv[iArg], "-get_subline"))
1389 : {
1390 1 : stOper = op_get_subline;
1391 : }
1392 19 : else if (EQUAL(papszArgv[iArg], "-l"))
1393 : {
1394 2 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1395 : // coverity[tainted_data]
1396 2 : pszLineDataSource = papszArgv[++iArg];
1397 : }
1398 17 : else if (EQUAL(papszArgv[iArg], "-ln"))
1399 : {
1400 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1401 : // coverity[tainted_data]
1402 0 : pszLineLayerName = papszArgv[++iArg];
1403 : }
1404 17 : else if (EQUAL(papszArgv[iArg], "-lf"))
1405 : {
1406 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1407 : #ifdef HAVE_GEOS
1408 : // coverity[tainted_data]
1409 0 : pszLineSepFieldName = papszArgv[++iArg];
1410 : #else
1411 : fprintf(stderr,
1412 : "GEOS support not enabled or incompatible version.\n");
1413 : exit(1);
1414 : #endif
1415 : }
1416 17 : else if (EQUAL(papszArgv[iArg], "-p"))
1417 : {
1418 2 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1419 : // coverity[tainted_data]
1420 2 : pszPicketsDataSource = papszArgv[++iArg];
1421 : }
1422 15 : else if (EQUAL(papszArgv[iArg], "-pn"))
1423 : {
1424 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1425 : #ifdef HAVE_GEOS
1426 : // coverity[tainted_data]
1427 0 : pszPicketsLayerName = papszArgv[++iArg];
1428 : #else
1429 : fprintf(stderr,
1430 : "GEOS support not enabled or incompatible version.\n");
1431 : exit(1);
1432 : #endif
1433 : }
1434 15 : else if (EQUAL(papszArgv[iArg], "-pm"))
1435 : {
1436 2 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1437 : #ifdef HAVE_GEOS
1438 : // coverity[tainted_data]
1439 2 : pszPicketsMField = papszArgv[++iArg];
1440 : #else
1441 : fprintf(stderr,
1442 : "GEOS support not enabled or incompatible version.\n");
1443 : exit(1);
1444 : #endif
1445 : }
1446 13 : else if (EQUAL(papszArgv[iArg], "-pf"))
1447 : {
1448 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1449 : #ifdef HAVE_GEOS
1450 : // coverity[tainted_data]
1451 0 : pszPicketsSepFieldName = papszArgv[++iArg];
1452 : #else
1453 : fprintf(stderr,
1454 : "GEOS support not enabled or incompatible version.\n");
1455 : exit(1);
1456 : #endif
1457 : }
1458 13 : else if (EQUAL(papszArgv[iArg], "-r"))
1459 : {
1460 3 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1461 : // coverity[tainted_data]
1462 3 : pszPartsDataSource = papszArgv[++iArg];
1463 : }
1464 10 : else if (EQUAL(papszArgv[iArg], "-rn"))
1465 : {
1466 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1467 : // coverity[tainted_data]
1468 0 : pszPartsLayerName = papszArgv[++iArg];
1469 : }
1470 10 : else if (EQUAL(papszArgv[iArg], "-o"))
1471 : {
1472 3 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1473 : // coverity[tainted_data]
1474 3 : pszOutputDataSource = papszArgv[++iArg];
1475 : }
1476 7 : else if (EQUAL(papszArgv[iArg], "-on"))
1477 : {
1478 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1479 : // coverity[tainted_data]
1480 0 : pszOutputLayerName = CPLStrdup(papszArgv[++iArg]);
1481 : }
1482 7 : else if (EQUAL(papszArgv[iArg], "-of"))
1483 : {
1484 0 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1485 : #ifdef HAVE_GEOS
1486 : // coverity[tainted_data]
1487 0 : pszOutputSepFieldName = papszArgv[++iArg];
1488 : #else
1489 : fprintf(stderr,
1490 : "GEOS support not enabled or incompatible version.\n");
1491 : exit(1);
1492 : #endif
1493 : }
1494 7 : else if (EQUAL(papszArgv[iArg], "-x"))
1495 : {
1496 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1497 : #ifdef HAVE_GEOS
1498 : // coverity[tainted_data]
1499 1 : dfX = CPLAtofM(papszArgv[++iArg]);
1500 : #else
1501 : fprintf(stderr,
1502 : "GEOS support not enabled or incompatible version.\n");
1503 : exit(1);
1504 : #endif
1505 : }
1506 6 : else if (EQUAL(papszArgv[iArg], "-y"))
1507 : {
1508 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1509 : #ifdef HAVE_GEOS
1510 : // coverity[tainted_data]
1511 1 : dfY = CPLAtofM(papszArgv[++iArg]);
1512 : #else
1513 : fprintf(stderr,
1514 : "GEOS support not enabled or incompatible version.\n");
1515 : exit(1);
1516 : #endif
1517 : }
1518 5 : else if (EQUAL(papszArgv[iArg], "-m"))
1519 : {
1520 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1521 : // coverity[tainted_data]
1522 1 : dfPos = CPLAtofM(papszArgv[++iArg]);
1523 : }
1524 4 : else if (EQUAL(papszArgv[iArg], "-mb"))
1525 : {
1526 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1527 : // coverity[tainted_data]
1528 1 : dfPosBeg = CPLAtofM(papszArgv[++iArg]);
1529 : }
1530 3 : else if (EQUAL(papszArgv[iArg], "-me"))
1531 : {
1532 1 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1533 : // coverity[tainted_data]
1534 1 : dfPosEnd = CPLAtofM(papszArgv[++iArg]);
1535 : }
1536 2 : else if (EQUAL(papszArgv[iArg], "-s"))
1537 : {
1538 2 : CHECK_HAS_ENOUGH_ADDITIONAL_ARGS(1);
1539 : #ifdef HAVE_GEOS
1540 : // coverity[tainted_data]
1541 2 : dfStep = CPLAtofM(papszArgv[++iArg]);
1542 : #else
1543 : fprintf(stderr,
1544 : "GEOS support not enabled or incompatible version.\n");
1545 : exit(1);
1546 : #endif
1547 : }
1548 0 : else if (EQUAL(papszArgv[iArg], "-progress"))
1549 : {
1550 0 : bDisplayProgress = TRUE;
1551 : }
1552 0 : else if (papszArgv[iArg][0] == '-')
1553 : {
1554 0 : Usage(true,
1555 0 : CPLSPrintf("Unknown option name '%s'", papszArgv[iArg]));
1556 : }
1557 : }
1558 :
1559 5 : if (stOper == op_create)
1560 : {
1561 : #ifdef HAVE_GEOS
1562 2 : if (pszOutputDataSource == nullptr)
1563 0 : Usage(true, "no output datasource provided");
1564 2 : else if (pszLineDataSource == nullptr)
1565 0 : Usage(true, "no path datasource provided");
1566 2 : else if (pszPicketsDataSource == nullptr)
1567 0 : Usage(true, "no repers datasource provided");
1568 2 : else if (pszPicketsMField == nullptr)
1569 0 : Usage(true, "no position field provided");
1570 2 : else if (dfStep == -100000000.0)
1571 0 : Usage(true, "no step provided");
1572 :
1573 : /* --------------------------------------------------------------------
1574 : */
1575 : /* Open data source. */
1576 : /* --------------------------------------------------------------------
1577 : */
1578 :
1579 : GDALDataset *poLnDS =
1580 2 : GDALDataset::FromHandle(OGROpen(pszLineDataSource, FALSE, nullptr));
1581 :
1582 : /* --------------------------------------------------------------------
1583 : */
1584 : /* Report failure */
1585 : /* --------------------------------------------------------------------
1586 : */
1587 2 : if (poLnDS == nullptr)
1588 : {
1589 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1590 :
1591 0 : fprintf(stderr,
1592 : "FAILURE:\n"
1593 : "Unable to open path datasource `%s' with "
1594 : "the following drivers.\n",
1595 : pszLineDataSource);
1596 :
1597 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1598 : {
1599 0 : fprintf(stderr, " -> %s\n",
1600 0 : poR->GetDriver(iDriver)->GetDescription());
1601 : }
1602 :
1603 0 : exit(1);
1604 : }
1605 :
1606 2 : GDALDataset *poPkDS = GDALDataset::FromHandle(
1607 2 : OGROpen(pszPicketsDataSource, FALSE, nullptr));
1608 : /* --------------------------------------------------------------------
1609 : */
1610 : /* Report failure */
1611 : /* --------------------------------------------------------------------
1612 : */
1613 2 : if (poPkDS == nullptr)
1614 : {
1615 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1616 :
1617 0 : fprintf(stderr,
1618 : "FAILURE:\n"
1619 : "Unable to open repers datasource `%s' "
1620 : "with the following drivers.\n",
1621 : pszPicketsDataSource);
1622 :
1623 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1624 : {
1625 0 : fprintf(stderr, " -> %s\n",
1626 0 : poR->GetDriver(iDriver)->GetDescription());
1627 : }
1628 :
1629 0 : exit(1);
1630 : }
1631 :
1632 : /* --------------------------------------------------------------------
1633 : */
1634 : /* Find the output driver. */
1635 : /* --------------------------------------------------------------------
1636 : */
1637 :
1638 2 : if (!bQuiet)
1639 2 : CheckDestDataSourceNameConsistency(pszOutputDataSource, pszFormat);
1640 :
1641 2 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1642 :
1643 2 : GDALDriver *poDriver = poR->GetDriverByName(pszFormat);
1644 2 : if (poDriver == nullptr)
1645 : {
1646 0 : fprintf(stderr, "Unable to find driver `%s'.\n", pszFormat);
1647 0 : fprintf(stderr, "The following drivers are available:\n");
1648 :
1649 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1650 : {
1651 0 : fprintf(stderr, " -> `%s'\n",
1652 0 : poR->GetDriver(iDriver)->GetDescription());
1653 : }
1654 0 : exit(1);
1655 : }
1656 :
1657 2 : if (!CPLTestBool(CSLFetchNameValueDef(poDriver->GetMetadata(),
1658 : GDAL_DCAP_CREATE, "FALSE")))
1659 : {
1660 0 : fprintf(stderr,
1661 : "%s driver does not support data source creation.\n",
1662 : pszFormat);
1663 0 : exit(1);
1664 : }
1665 :
1666 : /* --------------------------------------------------------------------
1667 : */
1668 : /* Create the output data source. */
1669 : /* --------------------------------------------------------------------
1670 : */
1671 2 : GDALDataset *poODS = poDriver->Create(pszOutputDataSource, 0, 0, 0,
1672 : GDT_Unknown, papszDSCO);
1673 2 : if (poODS == nullptr)
1674 : {
1675 0 : fprintf(stderr, "%s driver failed to create %s\n", pszFormat,
1676 : pszOutputDataSource);
1677 0 : exit(1);
1678 : }
1679 :
1680 : OGRLayer *poLnLayer = pszLineLayerName == nullptr
1681 2 : ? poLnDS->GetLayer(0)
1682 0 : : poLnDS->GetLayerByName(pszLineLayerName);
1683 :
1684 2 : if (poLnLayer == nullptr)
1685 : {
1686 0 : fprintf(stderr, "Get path layer failed.\n");
1687 0 : exit(1);
1688 : }
1689 :
1690 : OGRLayer *poPkLayer = pszPicketsLayerName == nullptr
1691 2 : ? poPkDS->GetLayer(0)
1692 0 : : poPkDS->GetLayerByName(pszPicketsLayerName);
1693 :
1694 2 : if (poPkLayer == nullptr)
1695 : {
1696 0 : fprintf(stderr, "Get repers layer failed.\n");
1697 0 : exit(1);
1698 : }
1699 :
1700 2 : OGRFeatureDefn *poPkFDefn = poPkLayer->GetLayerDefn();
1701 2 : int nMValField = poPkFDefn->GetFieldIndex(pszPicketsMField);
1702 :
1703 2 : OGRLayer *poOutLayer = nullptr;
1704 2 : if (pszLineSepFieldName != nullptr && pszPicketsSepFieldName != nullptr)
1705 : {
1706 : poOutLayer =
1707 0 : SetupTargetLayer(poLnLayer, poODS, papszLCO, pszOutputLayerName,
1708 : pszOutputSepFieldName);
1709 0 : if (poOutLayer == nullptr)
1710 : {
1711 0 : fprintf(stderr, "Create output layer failed.\n");
1712 0 : exit(1);
1713 : }
1714 :
1715 : // Do the work
1716 0 : eErr = CreatePartsMultiple(
1717 : poLnLayer, pszLineSepFieldName, poPkLayer,
1718 : pszPicketsSepFieldName, nMValField, dfStep, poOutLayer,
1719 : pszOutputSepFieldName, bDisplayProgress, bQuiet);
1720 : }
1721 : else
1722 : {
1723 2 : poOutLayer = SetupTargetLayer(poLnLayer, poODS, papszLCO,
1724 : pszOutputLayerName);
1725 2 : if (poOutLayer == nullptr)
1726 : {
1727 0 : fprintf(stderr, "Create output layer failed.\n");
1728 0 : exit(1);
1729 : }
1730 :
1731 : // Do the work
1732 2 : eErr = CreateParts(poLnLayer, poPkLayer, nMValField, dfStep,
1733 : poOutLayer, bDisplayProgress, bQuiet);
1734 : }
1735 :
1736 2 : GDALClose(poLnDS);
1737 2 : GDALClose(poPkDS);
1738 2 : if (GDALClose(poODS) != CE_None)
1739 0 : eErr = CE_Failure;
1740 :
1741 2 : if (nullptr != pszOutputLayerName)
1742 0 : CPLFree(pszOutputLayerName);
1743 : #else // HAVE_GEOS
1744 : fprintf(stderr, "GEOS support not enabled or incompatible version.\n");
1745 : exit(1);
1746 : #endif // HAVE_GEOS
1747 : }
1748 3 : else if (stOper == op_get_pos)
1749 : {
1750 : #ifdef HAVE_GEOS
1751 1 : if (pszPartsDataSource == nullptr)
1752 0 : Usage(true, "no parts datasource provided");
1753 1 : else if (dfX == -100000000.0 || dfY == -100000000.0)
1754 0 : Usage(true, "no coordinates provided");
1755 :
1756 1 : GDALDataset *poPartsDS = GDALDataset::FromHandle(
1757 1 : OGROpen(pszPartsDataSource, FALSE, nullptr));
1758 : /* --------------------------------------------------------------------
1759 : */
1760 : /* Report failure */
1761 : /* --------------------------------------------------------------------
1762 : */
1763 1 : if (poPartsDS == nullptr)
1764 : {
1765 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1766 :
1767 0 : fprintf(stderr,
1768 : "FAILURE:\n"
1769 : "Unable to open parts datasource `%s' with "
1770 : "the following drivers.\n",
1771 : pszPicketsDataSource);
1772 :
1773 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1774 : {
1775 0 : fprintf(stderr, " -> %s\n",
1776 0 : poR->GetDriver(iDriver)->GetDescription());
1777 : }
1778 :
1779 0 : exit(1);
1780 : }
1781 :
1782 : OGRLayer *poPartsLayer =
1783 : pszPartsLayerName == nullptr
1784 1 : ? poPartsDS->GetLayer(0)
1785 0 : : poPartsDS->GetLayerByName(pszPartsLayerName);
1786 :
1787 1 : if (poPartsLayer == nullptr)
1788 : {
1789 0 : fprintf(stderr, "Get parts layer failed.\n");
1790 0 : exit(1);
1791 : }
1792 :
1793 : // Do the work
1794 1 : eErr = GetPosition(poPartsLayer, dfX, dfY, bDisplayProgress, bQuiet);
1795 :
1796 1 : GDALClose(poPartsDS);
1797 :
1798 : #else // HAVE_GEOS
1799 : fprintf(stderr, "GEOS support not enabled or incompatible version.\n");
1800 : exit(1);
1801 : #endif // HAVE_GEOS
1802 : }
1803 2 : else if (stOper == op_get_coord)
1804 : {
1805 1 : if (pszPartsDataSource == nullptr)
1806 0 : Usage(true, "no parts datasource provided");
1807 1 : else if (dfPos == -100000000.0)
1808 0 : Usage(true, "no position provided");
1809 :
1810 1 : GDALDataset *poPartsDS = GDALDataset::FromHandle(
1811 1 : OGROpen(pszPartsDataSource, FALSE, nullptr));
1812 : /* --------------------------------------------------------------------
1813 : */
1814 : /* Report failure */
1815 : /* --------------------------------------------------------------------
1816 : */
1817 1 : if (poPartsDS == nullptr)
1818 : {
1819 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1820 :
1821 0 : fprintf(stderr,
1822 : "FAILURE:\n"
1823 : "Unable to open parts datasource `%s' with "
1824 : "the following drivers.\n",
1825 : pszPicketsDataSource);
1826 :
1827 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1828 : {
1829 0 : fprintf(stderr, " -> %s\n",
1830 0 : poR->GetDriver(iDriver)->GetDescription());
1831 : }
1832 :
1833 0 : exit(1);
1834 : }
1835 :
1836 : OGRLayer *poPartsLayer =
1837 : pszPartsLayerName == nullptr
1838 1 : ? poPartsDS->GetLayer(0)
1839 0 : : poPartsDS->GetLayerByName(pszPartsLayerName);
1840 :
1841 1 : if (poPartsLayer == nullptr)
1842 : {
1843 0 : fprintf(stderr, "Get parts layer failed.\n");
1844 0 : exit(1);
1845 : }
1846 : // Do the work
1847 1 : eErr = GetCoordinates(poPartsLayer, dfPos, bDisplayProgress, bQuiet);
1848 :
1849 1 : GDALClose(poPartsDS);
1850 : }
1851 1 : else if (stOper == op_get_subline)
1852 : {
1853 1 : if (pszOutputDataSource == nullptr)
1854 0 : Usage(true, "no output datasource provided");
1855 1 : else if (pszPartsDataSource == nullptr)
1856 0 : Usage(true, "no parts datasource provided");
1857 1 : else if (dfPosBeg == -100000000.0)
1858 0 : Usage(true, "no begin position provided");
1859 1 : else if (dfPosEnd == -100000000.0)
1860 0 : Usage(true, "no end position provided");
1861 :
1862 : // Open data source.
1863 1 : GDALDataset *poPartsDS = GDALDataset::FromHandle(
1864 1 : OGROpen(pszPartsDataSource, FALSE, nullptr));
1865 :
1866 : // Report failure.
1867 1 : if (poPartsDS == nullptr)
1868 : {
1869 0 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1870 :
1871 0 : fprintf(stderr,
1872 : "FAILURE:\n"
1873 : "Unable to open parts datasource `%s' with "
1874 : "the following drivers.\n",
1875 : pszLineDataSource);
1876 :
1877 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1878 : {
1879 0 : fprintf(stderr, " -> %s\n",
1880 0 : poR->GetDriver(iDriver)->GetDescription());
1881 : }
1882 :
1883 0 : exit(1);
1884 : }
1885 :
1886 : // Find the output driver.
1887 :
1888 1 : if (!bQuiet)
1889 1 : CheckDestDataSourceNameConsistency(pszOutputDataSource, pszFormat);
1890 :
1891 1 : OGRSFDriverRegistrar *poR = OGRSFDriverRegistrar::GetRegistrar();
1892 :
1893 1 : GDALDriver *poDriver = poR->GetDriverByName(pszFormat);
1894 1 : if (poDriver == nullptr)
1895 : {
1896 0 : fprintf(stderr, "Unable to find driver `%s'.\n", pszFormat);
1897 0 : fprintf(stderr, "The following drivers are available:\n");
1898 :
1899 0 : for (int iDriver = 0; iDriver < poR->GetDriverCount(); iDriver++)
1900 : {
1901 0 : fprintf(stderr, " -> `%s'\n",
1902 0 : poR->GetDriver(iDriver)->GetDescription());
1903 : }
1904 0 : exit(1);
1905 : }
1906 :
1907 1 : if (!CPLTestBool(CSLFetchNameValueDef(poDriver->GetMetadata(),
1908 : GDAL_DCAP_CREATE, "FALSE")))
1909 : {
1910 0 : fprintf(stderr,
1911 : "%s driver does not support data source creation.\n",
1912 : pszFormat);
1913 0 : exit(1);
1914 : }
1915 :
1916 : // Create the output data source.
1917 1 : GDALDataset *poODS = poDriver->Create(pszOutputDataSource, 0, 0, 0,
1918 : GDT_Unknown, papszDSCO);
1919 1 : if (poODS == nullptr)
1920 : {
1921 0 : fprintf(stderr, "%s driver failed to create %s\n", pszFormat,
1922 : pszOutputDataSource);
1923 0 : exit(1);
1924 : }
1925 :
1926 : OGRLayer *poPartsLayer =
1927 : pszLineLayerName == nullptr
1928 1 : ? poPartsDS->GetLayer(0)
1929 0 : : poPartsDS->GetLayerByName(pszLineLayerName);
1930 :
1931 1 : if (poPartsLayer == nullptr)
1932 : {
1933 0 : fprintf(stderr, "Get parts layer failed.\n");
1934 0 : exit(1);
1935 : }
1936 :
1937 : OGRLayer *poOutLayer =
1938 1 : SetupTargetLayer(poPartsLayer, poODS, papszLCO, pszOutputLayerName);
1939 1 : if (poOutLayer == nullptr)
1940 : {
1941 0 : fprintf(stderr, "Create output layer failed.\n");
1942 0 : exit(1);
1943 : }
1944 :
1945 : // Do the work.
1946 1 : eErr = CreateSubline(poPartsLayer, dfPosBeg, dfPosEnd, poOutLayer,
1947 : bDisplayProgress, bQuiet);
1948 :
1949 1 : GDALClose(poPartsDS);
1950 1 : if (GDALClose(poODS) != CE_None)
1951 0 : eErr = CE_Failure;
1952 :
1953 1 : if (nullptr != pszOutputLayerName)
1954 0 : CPLFree(pszOutputLayerName);
1955 : }
1956 : else
1957 : {
1958 0 : Usage(true, "no operation provided");
1959 : }
1960 :
1961 5 : CSLDestroy(papszArgv);
1962 5 : CSLDestroy(papszDSCO);
1963 5 : CSLDestroy(papszLCO);
1964 :
1965 5 : OGRCleanupAll();
1966 :
1967 : #ifdef DBMALLOC
1968 : malloc_dump(1);
1969 : #endif
1970 :
1971 5 : return eErr == OGRERR_NONE ? 0 : 1;
1972 : }
1973 :
1974 0 : MAIN_END
|
