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