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