Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: OpenGIS Simple Features Reference Implementation
4 : * Purpose: Implements OGRGmtLayer class.
5 : * Author: Frank Warmerdam, warmerdam@pobox.com
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2007, Frank Warmerdam <warmerdam@pobox.com>
9 : *
10 : * SPDX-License-Identifier: MIT
11 : ****************************************************************************/
12 :
13 : #include "ogr_gmt.h"
14 : #include "cpl_conv.h"
15 : #include "ogr_p.h"
16 :
17 : #include <algorithm>
18 :
19 : /************************************************************************/
20 : /* OGRGmtLayer() */
21 : /************************************************************************/
22 :
23 107 : OGRGmtLayer::OGRGmtLayer(GDALDataset *poDS, const char *pszFilename,
24 : VSILFILE *fp, const OGRSpatialReference *poSRS,
25 107 : int bUpdateIn)
26 : : m_poDS(poDS), poFeatureDefn(nullptr), iNextFID(0),
27 107 : bUpdate(CPL_TO_BOOL(bUpdateIn)),
28 : // Assume header complete in readonly mode.
29 107 : bHeaderComplete(CPL_TO_BOOL(!bUpdate)), bRegionComplete(false),
30 : nRegionOffset(0),
31 107 : m_fp(fp ? fp : VSIFOpenL(pszFilename, (bUpdateIn ? "r+" : "r"))),
32 214 : papszKeyedValues(nullptr), bValidFile(false)
33 : {
34 107 : if (m_fp == nullptr)
35 35 : return;
36 :
37 : /* -------------------------------------------------------------------- */
38 : /* Create the feature definition */
39 : /* -------------------------------------------------------------------- */
40 72 : poFeatureDefn = new OGRFeatureDefn(CPLGetBasenameSafe(pszFilename).c_str());
41 72 : SetDescription(poFeatureDefn->GetName());
42 72 : poFeatureDefn->Reference();
43 :
44 : /* -------------------------------------------------------------------- */
45 : /* Read the header. */
46 : /* -------------------------------------------------------------------- */
47 72 : if (!STARTS_WITH(pszFilename, "/vsistdout"))
48 : {
49 142 : CPLString osFieldNames;
50 142 : CPLString osFieldTypes;
51 142 : CPLString osGeometryType;
52 142 : CPLString osRegion;
53 142 : CPLString osWKT;
54 142 : CPLString osProj4;
55 142 : CPLString osEPSG;
56 71 : vsi_l_offset nStartOfLine = 0;
57 :
58 71 : VSIFSeekL(m_fp, 0, SEEK_SET);
59 :
60 258 : while (ReadLine() && osLine[0] == '#')
61 : {
62 207 : if (strstr(osLine, "FEATURE_DATA"))
63 : {
64 20 : bHeaderComplete = true;
65 20 : ReadLine();
66 20 : break;
67 : }
68 :
69 187 : if (STARTS_WITH_CI(osLine, "# REGION_STUB "))
70 34 : nRegionOffset = nStartOfLine;
71 :
72 390 : for (int iKey = 0; papszKeyedValues != nullptr &&
73 356 : papszKeyedValues[iKey] != nullptr;
74 : iKey++)
75 : {
76 203 : if (papszKeyedValues[iKey][0] == 'N')
77 20 : osFieldNames = papszKeyedValues[iKey] + 1;
78 203 : if (papszKeyedValues[iKey][0] == 'T')
79 20 : osFieldTypes = papszKeyedValues[iKey] + 1;
80 203 : if (papszKeyedValues[iKey][0] == 'G')
81 52 : osGeometryType = papszKeyedValues[iKey] + 1;
82 203 : if (papszKeyedValues[iKey][0] == 'R')
83 35 : osRegion = papszKeyedValues[iKey] + 1;
84 203 : if (papszKeyedValues[iKey][0] == 'J' &&
85 6 : papszKeyedValues[iKey][1] != 0 &&
86 6 : papszKeyedValues[iKey][2] != 0)
87 : {
88 12 : std::string osArg = papszKeyedValues[iKey] + 2;
89 10 : if (osArg[0] == '"' && osArg.size() >= 2 &&
90 4 : osArg.back() == '"')
91 : {
92 4 : osArg = osArg.substr(1, osArg.length() - 2);
93 4 : char *pszArg = CPLUnescapeString(
94 : osArg.c_str(), nullptr, CPLES_BackslashQuotable);
95 4 : osArg = pszArg;
96 4 : CPLFree(pszArg);
97 : }
98 :
99 6 : if (papszKeyedValues[iKey][1] == 'e')
100 2 : osEPSG = std::move(osArg);
101 6 : if (papszKeyedValues[iKey][1] == 'p')
102 2 : osProj4 = std::move(osArg);
103 6 : if (papszKeyedValues[iKey][1] == 'w')
104 2 : osWKT = std::move(osArg);
105 : }
106 : }
107 :
108 187 : nStartOfLine = VSIFTellL(m_fp);
109 : }
110 :
111 : /* --------------------------------------------------------------------
112 : */
113 : /* Handle coordinate system. */
114 : /* --------------------------------------------------------------------
115 : */
116 71 : if (osWKT.length())
117 : {
118 2 : m_poSRS = new OGRSpatialReference();
119 2 : m_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
120 2 : if (m_poSRS->importFromWkt(osWKT.c_str()) != OGRERR_NONE)
121 : {
122 0 : delete m_poSRS;
123 0 : m_poSRS = nullptr;
124 : }
125 : }
126 69 : else if (osEPSG.length())
127 : {
128 0 : m_poSRS = new OGRSpatialReference();
129 0 : m_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
130 0 : if (m_poSRS->importFromEPSG(atoi(osEPSG)) != OGRERR_NONE)
131 : {
132 0 : delete m_poSRS;
133 0 : m_poSRS = nullptr;
134 : }
135 : }
136 69 : else if (osProj4.length())
137 : {
138 0 : m_poSRS = new OGRSpatialReference();
139 0 : m_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
140 0 : if (m_poSRS->importFromProj4(osProj4) != OGRERR_NONE)
141 : {
142 0 : delete m_poSRS;
143 0 : m_poSRS = nullptr;
144 : }
145 : }
146 :
147 71 : if (osGeometryType == "POINT")
148 8 : poFeatureDefn->SetGeomType(wkbPoint);
149 63 : else if (osGeometryType == "MULTIPOINT")
150 8 : poFeatureDefn->SetGeomType(wkbMultiPoint);
151 55 : else if (osGeometryType == "LINESTRING")
152 8 : poFeatureDefn->SetGeomType(wkbLineString);
153 47 : else if (osGeometryType == "MULTILINESTRING")
154 9 : poFeatureDefn->SetGeomType(wkbMultiLineString);
155 38 : else if (osGeometryType == "POLYGON")
156 10 : poFeatureDefn->SetGeomType(wkbPolygon);
157 28 : else if (osGeometryType == "MULTIPOLYGON")
158 9 : poFeatureDefn->SetGeomType(wkbMultiPolygon);
159 :
160 : /* --------------------------------------------------------------------
161 : */
162 : /* Process a region line. */
163 : /* --------------------------------------------------------------------
164 : */
165 71 : if (osRegion.length() > 0)
166 : {
167 : char **papszTokens =
168 35 : CSLTokenizeStringComplex(osRegion.c_str(), "/", FALSE, FALSE);
169 :
170 35 : if (CSLCount(papszTokens) == 4)
171 : {
172 35 : sRegion.MinX = CPLAtofM(papszTokens[0]);
173 35 : sRegion.MaxX = CPLAtofM(papszTokens[1]);
174 35 : sRegion.MinY = CPLAtofM(papszTokens[2]);
175 35 : sRegion.MaxY = CPLAtofM(papszTokens[3]);
176 : }
177 :
178 35 : bRegionComplete = true;
179 :
180 35 : CSLDestroy(papszTokens);
181 : }
182 :
183 : /* --------------------------------------------------------------------
184 : */
185 : /* Process fields. */
186 : /* --------------------------------------------------------------------
187 : */
188 71 : if (osFieldNames.length() || osFieldTypes.length())
189 : {
190 : char **papszFN =
191 20 : CSLTokenizeStringComplex(osFieldNames, "|", TRUE, TRUE);
192 : char **papszFT =
193 20 : CSLTokenizeStringComplex(osFieldTypes, "|", TRUE, TRUE);
194 20 : const int nFNCount = CSLCount(papszFN);
195 20 : const int nFTCount = CSLCount(papszFT);
196 20 : const int nFieldCount = std::max(nFNCount, nFTCount);
197 :
198 110 : for (int iField = 0; iField < nFieldCount; iField++)
199 : {
200 180 : OGRFieldDefn oField("", OFTString);
201 :
202 90 : if (iField < nFNCount)
203 90 : oField.SetName(papszFN[iField]);
204 : else
205 0 : oField.SetName(CPLString().Printf("Field_%d", iField + 1));
206 :
207 90 : if (iField < nFTCount)
208 : {
209 90 : if (EQUAL(papszFT[iField], "integer"))
210 19 : oField.SetType(OFTInteger);
211 71 : else if (EQUAL(papszFT[iField], "double"))
212 19 : oField.SetType(OFTReal);
213 52 : else if (EQUAL(papszFT[iField], "datetime"))
214 17 : oField.SetType(OFTDateTime);
215 : }
216 :
217 90 : poFeatureDefn->AddFieldDefn(&oField);
218 : }
219 :
220 20 : CSLDestroy(papszFN);
221 20 : CSLDestroy(papszFT);
222 : }
223 : }
224 : else
225 : {
226 1 : if (poSRS)
227 : {
228 1 : m_poSRS = poSRS->Clone();
229 1 : m_poSRS->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
230 : }
231 : }
232 :
233 72 : poFeatureDefn->GetGeomFieldDefn(0)->SetSpatialRef(m_poSRS);
234 :
235 72 : bValidFile = true;
236 : }
237 :
238 : /************************************************************************/
239 : /* ~OGRGmtLayer() */
240 : /************************************************************************/
241 :
242 214 : OGRGmtLayer::~OGRGmtLayer()
243 :
244 : {
245 107 : if (m_nFeaturesRead > 0 && poFeatureDefn != nullptr)
246 : {
247 10 : CPLDebug("Gmt", "%d features read on layer '%s'.",
248 5 : static_cast<int>(m_nFeaturesRead), poFeatureDefn->GetName());
249 : }
250 :
251 : /* -------------------------------------------------------------------- */
252 : /* Write out the region bounds if we know where they go, and we */
253 : /* are in update mode. */
254 : /* -------------------------------------------------------------------- */
255 107 : if (nRegionOffset != 0 && bUpdate)
256 : {
257 34 : VSIFSeekL(m_fp, nRegionOffset, SEEK_SET);
258 34 : VSIFPrintfL(m_fp, "# @R%.12g/%.12g/%.12g/%.12g", sRegion.MinX,
259 : sRegion.MaxX, sRegion.MinY, sRegion.MaxY);
260 : }
261 :
262 : /* -------------------------------------------------------------------- */
263 : /* Clean up. */
264 : /* -------------------------------------------------------------------- */
265 107 : CSLDestroy(papszKeyedValues);
266 :
267 107 : if (poFeatureDefn)
268 72 : poFeatureDefn->Release();
269 :
270 107 : if (m_poSRS)
271 3 : m_poSRS->Release();
272 :
273 107 : if (m_fp != nullptr)
274 72 : VSIFCloseL(m_fp);
275 214 : }
276 :
277 : /************************************************************************/
278 : /* ReadLine() */
279 : /* */
280 : /* Read a line into osLine. If it is a comment line with @ */
281 : /* keyed values, parse out the keyed values into */
282 : /* papszKeyedValues. */
283 : /************************************************************************/
284 :
285 1244 : bool OGRGmtLayer::ReadLine()
286 :
287 : {
288 : /* -------------------------------------------------------------------- */
289 : /* Clear last line. */
290 : /* -------------------------------------------------------------------- */
291 1244 : osLine.erase();
292 1244 : if (papszKeyedValues)
293 : {
294 286 : CSLDestroy(papszKeyedValues);
295 286 : papszKeyedValues = nullptr;
296 : }
297 :
298 : /* -------------------------------------------------------------------- */
299 : /* Read newline. */
300 : /* -------------------------------------------------------------------- */
301 1244 : const char *pszLine = CPLReadLineL(m_fp);
302 1244 : if (pszLine == nullptr)
303 56 : return false; // end of file.
304 :
305 1188 : osLine = pszLine;
306 :
307 : /* -------------------------------------------------------------------- */
308 : /* If this is a comment line with keyed values, parse them. */
309 : /* -------------------------------------------------------------------- */
310 :
311 1188 : if (osLine[0] != '#' || osLine.find_first_of('@') == std::string::npos)
312 898 : return true;
313 :
314 290 : CPLStringList aosKeyedValues;
315 3163 : for (size_t i = 0; i < osLine.length(); i++)
316 : {
317 2873 : if (osLine[i] == '@' && i + 2 <= osLine.size())
318 : {
319 341 : bool bInQuotes = false;
320 :
321 341 : size_t iValEnd = i + 2; // Used after for.
322 5660 : for (; iValEnd < osLine.length(); iValEnd++)
323 : {
324 9610 : if (!bInQuotes &&
325 4204 : isspace(static_cast<unsigned char>(osLine[iValEnd])))
326 87 : break;
327 :
328 6511 : if (bInQuotes && iValEnd < osLine.length() - 1 &&
329 1192 : osLine[iValEnd] == '\\')
330 : {
331 88 : iValEnd++;
332 : }
333 5231 : else if (osLine[iValEnd] == '"')
334 26 : bInQuotes = !bInQuotes;
335 : }
336 :
337 682 : const CPLString osValue = osLine.substr(i + 2, iValEnd - i - 2);
338 :
339 : // Unecape contents
340 : char *pszUEValue =
341 341 : CPLUnescapeString(osValue, nullptr, CPLES_BackslashQuotable);
342 :
343 341 : CPLString osKeyValue = osLine.substr(i + 1, 1);
344 341 : osKeyValue += pszUEValue;
345 341 : CPLFree(pszUEValue);
346 341 : aosKeyedValues.AddString(osKeyValue);
347 :
348 341 : i = iValEnd;
349 : }
350 : }
351 290 : papszKeyedValues = aosKeyedValues.StealList();
352 :
353 290 : return true;
354 : }
355 :
356 : /************************************************************************/
357 : /* ResetReading() */
358 : /************************************************************************/
359 :
360 20 : void OGRGmtLayer::ResetReading()
361 :
362 : {
363 20 : if (iNextFID == 0)
364 18 : return;
365 :
366 2 : iNextFID = 0;
367 2 : VSIFSeekL(m_fp, 0, SEEK_SET);
368 2 : ReadLine();
369 : }
370 :
371 : /************************************************************************/
372 : /* ScanAheadForHole() */
373 : /* */
374 : /* Scan ahead to see if the next geometry is a hole. If so */
375 : /* return true, otherwise seek back to where we were and return */
376 : /* false. */
377 : /************************************************************************/
378 :
379 30 : bool OGRGmtLayer::ScanAheadForHole()
380 :
381 : {
382 60 : CPLString osSavedLine = osLine;
383 30 : const vsi_l_offset nSavedLocation = VSIFTellL(m_fp);
384 :
385 87 : while (ReadLine() && osLine[0] == '#')
386 : {
387 58 : if (papszKeyedValues != nullptr && papszKeyedValues[0][0] == 'H')
388 1 : return true;
389 : }
390 :
391 29 : VSIFSeekL(m_fp, nSavedLocation, SEEK_SET);
392 29 : osLine = std::move(osSavedLine);
393 :
394 : // We do not actually restore papszKeyedValues, but we
395 : // assume it does not matter since this method is only called
396 : // when processing the '>' line.
397 :
398 29 : return false;
399 : }
400 :
401 : /************************************************************************/
402 : /* NextIsFeature() */
403 : /* */
404 : /* Returns true if the next line is a feature attribute line. */
405 : /* This generally indicates the end of a multilinestring or */
406 : /* multipolygon feature. */
407 : /************************************************************************/
408 :
409 5 : bool OGRGmtLayer::NextIsFeature()
410 :
411 : {
412 5 : CPLString osSavedLine = osLine;
413 5 : const vsi_l_offset nSavedLocation = VSIFTellL(m_fp);
414 5 : bool bReturn = false;
415 :
416 5 : ReadLine();
417 :
418 5 : if (osLine[0] == '#' && strstr(osLine, "@D") != nullptr)
419 2 : bReturn = true;
420 :
421 5 : VSIFSeekL(m_fp, nSavedLocation, SEEK_SET);
422 5 : osLine = std::move(osSavedLine);
423 :
424 : // We do not actually restore papszKeyedValues, but we
425 : // assume it does not matter since this method is only called
426 : // when processing the '>' line.
427 :
428 10 : return bReturn;
429 : }
430 :
431 : /************************************************************************/
432 : /* GetNextRawFeature() */
433 : /************************************************************************/
434 :
435 55 : OGRFeature *OGRGmtLayer::GetNextRawFeature()
436 :
437 : {
438 : #if 0
439 : bool bMultiVertex =
440 : poFeatureDefn->GetGeomType() != wkbPoint
441 : && poFeatureDefn->GetGeomType() != wkbUnknown;
442 : #endif
443 110 : CPLString osFieldData;
444 55 : OGRGeometry *poGeom = nullptr;
445 :
446 : /* -------------------------------------------------------------------- */
447 : /* Read lines associated with this feature. */
448 : /* -------------------------------------------------------------------- */
449 871 : for (; true; ReadLine())
450 : {
451 926 : if (osLine.length() == 0)
452 25 : break;
453 :
454 901 : if (osLine[0] == '>')
455 : {
456 67 : OGRwkbGeometryType eType = wkbUnknown;
457 67 : if (poGeom)
458 33 : eType = wkbFlatten(poGeom->getGeometryType());
459 67 : if (eType == wkbMultiPolygon)
460 : {
461 3 : OGRMultiPolygon *poMP = poGeom->toMultiPolygon();
462 3 : if (ScanAheadForHole())
463 : {
464 : // Add a hole to the current polygon.
465 1 : poMP->getGeometryRef(poMP->getNumGeometries() - 1)
466 1 : ->addRingDirectly(new OGRLinearRing());
467 : }
468 2 : else if (!NextIsFeature())
469 : {
470 1 : OGRPolygon *poPoly = new OGRPolygon();
471 :
472 1 : poPoly->addRingDirectly(new OGRLinearRing());
473 :
474 1 : poMP->addGeometryDirectly(poPoly);
475 : }
476 : else
477 1 : break; /* done geometry */
478 : }
479 64 : else if (eType == wkbPolygon)
480 : {
481 27 : if (ScanAheadForHole())
482 0 : poGeom->toPolygon()->addRingDirectly(new OGRLinearRing());
483 : else
484 27 : break; /* done geometry */
485 : }
486 37 : else if (eType == wkbMultiLineString && !NextIsFeature())
487 : {
488 4 : poGeom->toMultiLineString()->addGeometryDirectly(
489 2 : new OGRLineString());
490 : }
491 35 : else if (poGeom != nullptr)
492 : {
493 1 : break;
494 : }
495 34 : else if (poFeatureDefn->GetGeomType() == wkbUnknown)
496 : {
497 0 : poFeatureDefn->SetGeomType(wkbLineString);
498 : // bMultiVertex = true;
499 : }
500 : }
501 834 : else if (osLine[0] == '#')
502 : {
503 145 : for (int i = 0;
504 145 : papszKeyedValues != nullptr && papszKeyedValues[i] != nullptr;
505 : i++)
506 : {
507 72 : if (papszKeyedValues[i][0] == 'D')
508 34 : osFieldData = papszKeyedValues[i] + 1;
509 : }
510 : }
511 : else
512 : {
513 : // Parse point line.
514 761 : double dfX = 0.0;
515 761 : double dfY = 0.0;
516 761 : double dfZ = 0.0;
517 761 : const int nDim = CPLsscanf(osLine, "%lf %lf %lf", &dfX, &dfY, &dfZ);
518 :
519 761 : if (nDim >= 2)
520 : {
521 761 : if (poGeom == nullptr)
522 : {
523 35 : switch (poFeatureDefn->GetGeomType())
524 : {
525 0 : case wkbLineString:
526 0 : poGeom = new OGRLineString();
527 0 : break;
528 :
529 30 : case wkbPolygon:
530 : {
531 30 : OGRPolygon *poPoly = new OGRPolygon();
532 30 : poGeom = poPoly;
533 30 : poPoly->addRingDirectly(new OGRLinearRing());
534 30 : break;
535 : }
536 :
537 2 : case wkbMultiPolygon:
538 : {
539 2 : OGRPolygon *poPoly = new OGRPolygon();
540 2 : poPoly->addRingDirectly(new OGRLinearRing());
541 :
542 2 : OGRMultiPolygon *poMP = new OGRMultiPolygon();
543 2 : poGeom = poMP;
544 2 : poMP->addGeometryDirectly(poPoly);
545 : }
546 2 : break;
547 :
548 0 : case wkbMultiPoint:
549 0 : poGeom = new OGRMultiPoint();
550 0 : break;
551 :
552 2 : case wkbMultiLineString:
553 : {
554 : OGRMultiLineString *poMLS =
555 2 : new OGRMultiLineString();
556 2 : poGeom = poMLS;
557 2 : poMLS->addGeometryDirectly(new OGRLineString());
558 2 : break;
559 : }
560 :
561 1 : case wkbPoint:
562 : case wkbUnknown:
563 : default:
564 1 : poGeom = new OGRPoint();
565 1 : break;
566 : }
567 : }
568 :
569 761 : CPLAssert(poGeom != nullptr);
570 : // cppcheck-suppress nullPointerRedundantCheck
571 761 : switch (wkbFlatten(poGeom->getGeometryType()))
572 : {
573 1 : case wkbPoint:
574 : {
575 1 : OGRPoint *poPoint = poGeom->toPoint();
576 1 : poPoint->setX(dfX);
577 1 : poPoint->setY(dfY);
578 1 : if (nDim == 3)
579 1 : poPoint->setZ(dfZ);
580 1 : break;
581 : }
582 :
583 0 : case wkbLineString:
584 : {
585 0 : OGRLineString *poLS = poGeom->toLineString();
586 0 : if (nDim == 3)
587 0 : poLS->addPoint(dfX, dfY, dfZ);
588 : else
589 0 : poLS->addPoint(dfX, dfY);
590 0 : break;
591 : }
592 :
593 752 : case wkbPolygon:
594 : case wkbMultiPolygon:
595 : {
596 752 : OGRPolygon *poPoly = nullptr;
597 :
598 752 : if (wkbFlatten(poGeom->getGeometryType()) ==
599 : wkbMultiPolygon)
600 : {
601 17 : OGRMultiPolygon *poMP = poGeom->toMultiPolygon();
602 17 : poPoly = poMP->getGeometryRef(
603 17 : poMP->getNumGeometries() - 1);
604 : }
605 : else
606 735 : poPoly = poGeom->toPolygon();
607 :
608 752 : OGRLinearRing *poRing = nullptr;
609 752 : if (poPoly->getNumInteriorRings() == 0)
610 748 : poRing = poPoly->getExteriorRing();
611 : else
612 4 : poRing = poPoly->getInteriorRing(
613 4 : poPoly->getNumInteriorRings() - 1);
614 :
615 752 : if (nDim == 3)
616 0 : poRing->addPoint(dfX, dfY, dfZ);
617 : else
618 752 : poRing->addPoint(dfX, dfY);
619 : }
620 752 : break;
621 :
622 8 : case wkbMultiLineString:
623 : {
624 8 : OGRMultiLineString *poML = poGeom->toMultiLineString();
625 : OGRLineString *poLine =
626 8 : poML->getGeometryRef(poML->getNumGeometries() - 1);
627 :
628 8 : if (nDim == 3)
629 0 : poLine->addPoint(dfX, dfY, dfZ);
630 : else
631 8 : poLine->addPoint(dfX, dfY);
632 : }
633 8 : break;
634 :
635 0 : default:
636 0 : CPLAssert(false);
637 : }
638 : }
639 : }
640 :
641 872 : if (poGeom && wkbFlatten(poGeom->getGeometryType()) == wkbPoint)
642 : {
643 1 : ReadLine();
644 1 : break;
645 : }
646 871 : }
647 :
648 55 : if (poGeom == nullptr)
649 20 : return nullptr;
650 :
651 : /* -------------------------------------------------------------------- */
652 : /* Create feature. */
653 : /* -------------------------------------------------------------------- */
654 35 : OGRFeature *poFeature = new OGRFeature(poFeatureDefn);
655 35 : poGeom->assignSpatialReference(m_poSRS);
656 35 : poFeature->SetGeometryDirectly(poGeom);
657 35 : poFeature->SetFID(iNextFID++);
658 :
659 : /* -------------------------------------------------------------------- */
660 : /* Process field values. */
661 : /* -------------------------------------------------------------------- */
662 35 : char **papszFD = CSLTokenizeStringComplex(osFieldData, "|", TRUE, TRUE);
663 :
664 133 : for (int iField = 0; papszFD != nullptr && papszFD[iField] != nullptr;
665 : iField++)
666 : {
667 98 : if (iField >= poFeatureDefn->GetFieldCount())
668 0 : break;
669 :
670 98 : poFeature->SetField(iField, papszFD[iField]);
671 : }
672 :
673 35 : CSLDestroy(papszFD);
674 :
675 35 : m_nFeaturesRead++;
676 :
677 35 : return poFeature;
678 : }
679 :
680 : /************************************************************************/
681 : /* CompleteHeader() */
682 : /* */
683 : /* Finish writing out the header with field definitions and the */
684 : /* layer geometry type. */
685 : /************************************************************************/
686 :
687 19 : OGRErr OGRGmtLayer::CompleteHeader(OGRGeometry *poThisGeom)
688 :
689 : {
690 : /* -------------------------------------------------------------------- */
691 : /* If we do not already have a geometry type, try to work one */
692 : /* out and write it now. */
693 : /* -------------------------------------------------------------------- */
694 19 : if (poFeatureDefn->GetGeomType() == wkbUnknown && poThisGeom != nullptr)
695 : {
696 2 : poFeatureDefn->SetGeomType(wkbFlatten(poThisGeom->getGeometryType()));
697 :
698 2 : const char *pszGeom = nullptr;
699 2 : switch (wkbFlatten(poFeatureDefn->GetGeomType()))
700 : {
701 0 : case wkbPoint:
702 0 : pszGeom = " @GPOINT";
703 0 : break;
704 0 : case wkbLineString:
705 0 : pszGeom = " @GLINESTRING";
706 0 : break;
707 1 : case wkbPolygon:
708 1 : pszGeom = " @GPOLYGON";
709 1 : break;
710 0 : case wkbMultiPoint:
711 0 : pszGeom = " @GMULTIPOINT";
712 0 : break;
713 0 : case wkbMultiLineString:
714 0 : pszGeom = " @GMULTILINESTRING";
715 0 : break;
716 1 : case wkbMultiPolygon:
717 1 : pszGeom = " @GMULTIPOLYGON";
718 1 : break;
719 0 : default:
720 0 : pszGeom = "";
721 0 : break;
722 : }
723 :
724 2 : VSIFPrintfL(m_fp, "#%s\n", pszGeom);
725 : }
726 :
727 : /* -------------------------------------------------------------------- */
728 : /* Prepare and write the field names and types. */
729 : /* -------------------------------------------------------------------- */
730 38 : CPLString osFieldNames;
731 19 : CPLString osFieldTypes;
732 :
733 106 : for (int iField = 0; iField < poFeatureDefn->GetFieldCount(); iField++)
734 : {
735 87 : if (iField > 0)
736 : {
737 68 : osFieldNames += "|";
738 68 : osFieldTypes += "|";
739 : }
740 :
741 87 : osFieldNames += poFeatureDefn->GetFieldDefn(iField)->GetNameRef();
742 87 : switch (poFeatureDefn->GetFieldDefn(iField)->GetType())
743 : {
744 18 : case OFTInteger:
745 18 : osFieldTypes += "integer";
746 18 : break;
747 :
748 19 : case OFTReal:
749 19 : osFieldTypes += "double";
750 19 : break;
751 :
752 16 : case OFTDateTime:
753 16 : osFieldTypes += "datetime";
754 16 : break;
755 :
756 34 : default:
757 34 : osFieldTypes += "string";
758 34 : break;
759 : }
760 : }
761 :
762 19 : if (poFeatureDefn->GetFieldCount() > 0)
763 : {
764 19 : VSIFPrintfL(m_fp, "# @N%s\n", osFieldNames.c_str());
765 19 : VSIFPrintfL(m_fp, "# @T%s\n", osFieldTypes.c_str());
766 : }
767 :
768 : /* -------------------------------------------------------------------- */
769 : /* Mark the end of the header, and start of feature data. */
770 : /* -------------------------------------------------------------------- */
771 19 : VSIFPrintfL(m_fp, "# FEATURE_DATA\n");
772 :
773 19 : bHeaderComplete = true;
774 19 : bRegionComplete = true; // no feature written, so we know them all!
775 :
776 38 : return OGRERR_NONE;
777 : }
778 :
779 : /************************************************************************/
780 : /* ICreateFeature() */
781 : /************************************************************************/
782 :
783 86 : OGRErr OGRGmtLayer::ICreateFeature(OGRFeature *poFeature)
784 :
785 : {
786 86 : if (!bUpdate)
787 : {
788 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
789 : "Cannot create features on read-only dataset.");
790 0 : return OGRERR_FAILURE;
791 : }
792 :
793 : /* -------------------------------------------------------------------- */
794 : /* Do we need to write the header describing the fields? */
795 : /* -------------------------------------------------------------------- */
796 86 : if (!bHeaderComplete)
797 : {
798 19 : OGRErr eErr = CompleteHeader(poFeature->GetGeometryRef());
799 :
800 19 : if (eErr != OGRERR_NONE)
801 0 : return eErr;
802 : }
803 :
804 : /* -------------------------------------------------------------------- */
805 : /* Write out the feature */
806 : /* -------------------------------------------------------------------- */
807 86 : OGRGeometry *poGeom = poFeature->GetGeometryRef();
808 :
809 86 : if (poGeom == nullptr)
810 : {
811 33 : CPLError(CE_Failure, CPLE_AppDefined,
812 : "Features without geometry not supported by GMT writer.");
813 33 : return OGRERR_FAILURE;
814 : }
815 :
816 53 : if (poFeatureDefn->GetGeomType() == wkbUnknown)
817 4 : poFeatureDefn->SetGeomType(wkbFlatten(poGeom->getGeometryType()));
818 :
819 : // Do we need a vertex collection marker grouping vertices.
820 53 : if (poFeatureDefn->GetGeomType() != wkbPoint)
821 47 : VSIFPrintfL(m_fp, ">\n");
822 :
823 : /* -------------------------------------------------------------------- */
824 : /* Write feature properties() */
825 : /* -------------------------------------------------------------------- */
826 53 : if (poFeatureDefn->GetFieldCount() > 0)
827 : {
828 106 : CPLString osFieldData;
829 :
830 270 : for (int iField = 0; iField < poFeatureDefn->GetFieldCount(); iField++)
831 : {
832 : OGRFieldType eFType =
833 217 : poFeatureDefn->GetFieldDefn(iField)->GetType();
834 217 : const char *pszRawValue = poFeature->GetFieldAsString(iField);
835 :
836 217 : if (iField > 0)
837 164 : osFieldData += "|";
838 :
839 : // We do not want prefix spaces for numeric values.
840 217 : if (eFType == OFTInteger || eFType == OFTReal)
841 104 : while (*pszRawValue == ' ')
842 0 : pszRawValue++;
843 :
844 217 : if (strchr(pszRawValue, ' ') || strchr(pszRawValue, '|') ||
845 186 : strchr(pszRawValue, '\t') || strchr(pszRawValue, '\n'))
846 : {
847 31 : osFieldData += "\"";
848 :
849 : char *pszEscapedVal =
850 31 : CPLEscapeString(pszRawValue, -1, CPLES_BackslashQuotable);
851 31 : osFieldData += pszEscapedVal;
852 31 : CPLFree(pszEscapedVal);
853 :
854 31 : osFieldData += "\"";
855 : }
856 : else
857 186 : osFieldData += pszRawValue;
858 : }
859 :
860 53 : VSIFPrintfL(m_fp, "# @D%s\n", osFieldData.c_str());
861 : }
862 :
863 : /* -------------------------------------------------------------------- */
864 : /* Write Geometry */
865 : /* -------------------------------------------------------------------- */
866 53 : return WriteGeometry(OGRGeometry::ToHandle(poGeom), true);
867 : }
868 :
869 : /************************************************************************/
870 : /* WriteGeometry() */
871 : /* */
872 : /* Write a geometry to the file. If bHaveAngle is true it */
873 : /* means the angle bracket preceding the point stream has */
874 : /* already been written out. */
875 : /* */
876 : /* We use the C API for geometry access because of its */
877 : /* simplified access to vertices and children geometries. */
878 : /************************************************************************/
879 :
880 116 : OGRErr OGRGmtLayer::WriteGeometry(OGRGeometryH hGeom, bool bHaveAngle)
881 :
882 : {
883 : /* -------------------------------------------------------------------- */
884 : /* This is a geometry with sub-geometries. */
885 : /* -------------------------------------------------------------------- */
886 116 : if (OGR_G_GetGeometryCount(hGeom) > 0)
887 : {
888 53 : OGRErr eErr = OGRERR_NONE;
889 :
890 53 : for (int iGeom = 0;
891 116 : iGeom < OGR_G_GetGeometryCount(hGeom) && eErr == OGRERR_NONE;
892 : iGeom++)
893 : {
894 : // We need to emit polygon @P and @H items while we still
895 : // know this is a polygon and which is the outer and inner
896 : // ring.
897 63 : if (wkbFlatten(OGR_G_GetGeometryType(hGeom)) == wkbPolygon)
898 : {
899 36 : if (!bHaveAngle)
900 : {
901 6 : VSIFPrintfL(m_fp, ">\n");
902 6 : bHaveAngle = true;
903 : }
904 36 : if (iGeom == 0)
905 35 : VSIFPrintfL(m_fp, "# @P\n");
906 : else
907 1 : VSIFPrintfL(m_fp, "# @H\n");
908 : }
909 :
910 : eErr =
911 63 : WriteGeometry(OGR_G_GetGeometryRef(hGeom, iGeom), bHaveAngle);
912 63 : bHaveAngle = false;
913 : }
914 53 : return eErr;
915 : }
916 :
917 : /* -------------------------------------------------------------------- */
918 : /* If this is not a point we need to have an angle bracket to */
919 : /* mark the vertex list. */
920 : /* -------------------------------------------------------------------- */
921 63 : if (wkbFlatten(OGR_G_GetGeometryType(hGeom)) != wkbPoint && !bHaveAngle)
922 4 : VSIFPrintfL(m_fp, ">\n");
923 :
924 : /* -------------------------------------------------------------------- */
925 : /* Dump vertices. */
926 : /* -------------------------------------------------------------------- */
927 63 : const int nPointCount = OGR_G_GetPointCount(hGeom);
928 63 : const int nDim = OGR_G_GetCoordinateDimension(hGeom);
929 : // For testing only. Ticket #6453
930 : const bool bUseTab =
931 63 : CPLTestBool(CPLGetConfigOption("GMT_USE_TAB", "FALSE"));
932 :
933 671 : for (int iPoint = 0; iPoint < nPointCount; iPoint++)
934 : {
935 608 : const double dfX = OGR_G_GetX(hGeom, iPoint);
936 608 : const double dfY = OGR_G_GetY(hGeom, iPoint);
937 608 : const double dfZ = OGR_G_GetZ(hGeom, iPoint);
938 :
939 608 : sRegion.Merge(dfX, dfY);
940 : char szLine[128];
941 608 : OGRMakeWktCoordinate(szLine, dfX, dfY, dfZ, nDim);
942 608 : if (bUseTab)
943 : {
944 60 : for (char *szPtr = szLine; *szPtr != '\0'; ++szPtr)
945 : {
946 47 : if (*szPtr == ' ')
947 13 : *szPtr = '\t';
948 : }
949 : }
950 608 : if (VSIFPrintfL(m_fp, "%s\n", szLine) < 1)
951 : {
952 0 : CPLError(CE_Failure, CPLE_FileIO, "Gmt write failure: %s",
953 0 : VSIStrerror(errno));
954 0 : return OGRERR_FAILURE;
955 : }
956 : }
957 :
958 63 : return OGRERR_NONE;
959 : }
960 :
961 : /************************************************************************/
962 : /* IGetExtent() */
963 : /* */
964 : /* Fetch extent of the data currently stored in the dataset. */
965 : /* The bForce flag has no effect on SHO files since that value */
966 : /* is always in the header. */
967 : /* */
968 : /* Returns OGRERR_NONE/OGRRERR_FAILURE. */
969 : /************************************************************************/
970 :
971 0 : OGRErr OGRGmtLayer::IGetExtent(int iGeomField, OGREnvelope *psExtent,
972 : bool bForce)
973 :
974 : {
975 0 : if (bRegionComplete && sRegion.IsInit())
976 : {
977 0 : *psExtent = sRegion;
978 0 : return OGRERR_NONE;
979 : }
980 :
981 0 : return OGRLayer::IGetExtent(iGeomField, psExtent, bForce);
982 : }
983 :
984 : /************************************************************************/
985 : /* TestCapability() */
986 : /************************************************************************/
987 :
988 72 : int OGRGmtLayer::TestCapability(const char *pszCap)
989 :
990 : {
991 72 : if (EQUAL(pszCap, OLCRandomRead))
992 0 : return FALSE;
993 :
994 72 : if (EQUAL(pszCap, OLCSequentialWrite))
995 16 : return TRUE;
996 :
997 56 : if (EQUAL(pszCap, OLCFastSpatialFilter))
998 0 : return FALSE;
999 :
1000 56 : if (EQUAL(pszCap, OLCFastGetExtent))
1001 0 : return bRegionComplete;
1002 :
1003 56 : if (EQUAL(pszCap, OLCCreateField))
1004 16 : return TRUE;
1005 :
1006 40 : if (EQUAL(pszCap, OLCZGeometries))
1007 0 : return TRUE;
1008 :
1009 40 : return FALSE;
1010 : }
1011 :
1012 : /************************************************************************/
1013 : /* CreateField() */
1014 : /************************************************************************/
1015 :
1016 87 : OGRErr OGRGmtLayer::CreateField(const OGRFieldDefn *poField, int bApproxOK)
1017 :
1018 : {
1019 87 : if (!bUpdate)
1020 : {
1021 0 : CPLError(CE_Failure, CPLE_NoWriteAccess,
1022 : "Cannot create fields on read-only dataset.");
1023 0 : return OGRERR_FAILURE;
1024 : }
1025 :
1026 87 : if (bHeaderComplete)
1027 : {
1028 0 : CPLError(CE_Failure, CPLE_AppDefined,
1029 : "Unable to create fields after features have been created.");
1030 0 : return OGRERR_FAILURE;
1031 : }
1032 :
1033 87 : switch (poField->GetType())
1034 : {
1035 71 : case OFTInteger:
1036 : case OFTReal:
1037 : case OFTString:
1038 : case OFTDateTime:
1039 71 : poFeatureDefn->AddFieldDefn(poField);
1040 71 : return OGRERR_NONE;
1041 : break;
1042 :
1043 : break;
1044 :
1045 16 : default:
1046 16 : if (!bApproxOK)
1047 : {
1048 0 : CPLError(CE_Failure, CPLE_AppDefined,
1049 : "Field %s is of unsupported type %s.",
1050 : poField->GetNameRef(),
1051 : poField->GetFieldTypeName(poField->GetType()));
1052 0 : return OGRERR_FAILURE;
1053 : }
1054 16 : else if (poField->GetType() == OFTDate ||
1055 0 : poField->GetType() == OFTTime)
1056 : {
1057 16 : OGRFieldDefn oModDef(poField);
1058 16 : oModDef.SetType(OFTDateTime);
1059 16 : poFeatureDefn->AddFieldDefn(poField);
1060 16 : return OGRERR_NONE;
1061 : }
1062 : else
1063 : {
1064 0 : OGRFieldDefn oModDef(poField);
1065 0 : oModDef.SetType(OFTString);
1066 0 : poFeatureDefn->AddFieldDefn(poField);
1067 0 : return OGRERR_NONE;
1068 : }
1069 : }
1070 : }
|