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