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