Changeset 11894

Timestamp:

08/17/07 15:43:44 (1 year ago)

Author:

rouault

Message:

* (bug 1735 : BNA driver)
* add two methods to OGREnvelope
* make the BNA parser a bit more tolerant again
* complete change in the handling of complex polygons
* add write support
* add HTML documentation

Files:

• trunk/gdal/ogr/ogr_core.h (modified) (1 diff)
• trunk/gdal/ogr/ogrsf_frmts/bna/GNUmakefile (modified) (1 diff)
• trunk/gdal/ogr/ogrsf_frmts/bna/drv_bna.html (added)
• trunk/gdal/ogr/ogrsf_frmts/bna/makefile.vc (modified) (1 diff)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogr_bna.h (modified) (5 diffs)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnadatasource.cpp (modified) (7 diffs)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnadriver.cpp (modified) (2 diffs)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnalayer.cpp (modified) (19 diffs)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnaparser.cpp (modified) (9 diffs)
• trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnaparser.h (modified) (1 diff)
• trunk/gdal/ogr/ogrsf_frmts/ogr_formats.html (modified) (1 diff)

trunk/gdal/ogr/ogr_core.h

@@ -86 +86 @@
-            MinY = MaxY = dfY;
-        }
+    }
+
+    int Intersects(OGREnvelope& other)
+    {
+        return MinX <= other.MaxX && MaxX >= other.MinX && 
+               MinY <= other.MaxY && MaxY >= other.MinY;
+    }
+
+    int Contains(OGREnvelope& other)
+    {
+        return MinX <= other.MinX && MinY <= other.MinY &&
+               MaxX >= other.MaxX && MaxY >= other.MaxY;
-    }
-};

trunk/gdal/ogr/ogrsf_frmts/bna/GNUmakefile

@@ -4 +4 @@
-
-OBJ     =       ogrbnadriver.o ogrbnadatasource.o ogrbnalayer.o ogrbnaparser.o
+
+ifeq ($(HAVE_GEOS),yes)
+CPPFLAGS        :=      -DHAVE_GEOS=1 $(GEOS_CFLAGS) $(CPPFLAGS)
+endif
-
-CPPFLAGS        :=      -I.. -I../.. $(GDAL_INCLUDE) $(CPPFLAGS)

trunk/gdal/ogr/ogrsf_frmts/bna/makefile.vc

@@ -3 +3 @@
-                ogrbnaparser.obj
-
-
+ $(GEOS_CFLAGS)
-
-GDAL_ROOT       =       ..\..\..

trunk/gdal/ogr/ogrsf_frmts/bna/ogr_bna.h

@@ -49 +49 @@
-class OGRBNALayer : public OGRLayer
-{
-
+
+
+
+
-
+    int                nIDs;
-    int                eof;
-    int                failed;
-    int                curLine;
-    int                nNextFID;
-               *
+*              
-    int                nFeatures;
-    int                partialIndexTable;
@@ -61 +65 @@
-
-    BNAFeatureType     bnaFeatureType;
-
+
-    OGRFeature *       BuildFeatureFromBNARecord (BNARecord* record, long fid);
-    void               FastParseUntil ( int interestFID);
+    void               WriteFeatureAttributes(FILE* fp, OGRFeature *poFeature);
+    
-  public:
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
-
-    void                ResetReading();
-    OGRFeature *        GetNextFeature();
+    
+    OGRErr              CreateFeature( OGRFeature *poFeature );
+    OGRErr              CreateField( OGRFieldDefn *poField, int bApproxOK );
-
-    OGRFeatureDefn *    GetLayerDefn() { return poFeatureDefn; }
-    
-    int                 GetFeatureCount( int );
-    
-    OGRFeature *        GetFeature( long nFID );
@@ -90 +104 @@
-class OGRBNADataSource : public OGRDataSource
-{
-                *
+*               
-
-       **
+**       
-    int                 nLayers;
-
-    int                 bUpdate;
+    
+    /*  Export related */
+    FILE                *fpOutput;
+    int                 bUseCRLF;
+    int                 bMultiLine;
+    int                 nbOutID;
+    int                 bEllipsesAsEllipses;
+    int                 nbPairPerLine;
+    int                 coordinatePrecision;
+    char*               pszCoordinateSeparator;
-    
-  public:
@@ -101 +125 @@
-                        ~OGRBNADataSource();
-
+    FILE                *GetOutputFP() { return fpOutput; }
+    int                 GetUseCRLF() { return bUseCRLF; }
+    int                 GetMultiLine() { return bMultiLine; }
+    int                 GetNbOutId() { return nbOutID; }
+    int                 GetEllipsesAsEllipses() { return bEllipsesAsEllipses; }
+    int                 GetNbPairPerLine() { return nbPairPerLine; }
+    int                 GetCoordinatePrecision() { return coordinatePrecision; }
+    const char*         GetCoordinateSeparator() { return pszCoordinateSeparator; }
+
-    int                 Open( const char * pszFilename,
-                              int bUpdate );
-    
-
+
+
+
+
-
-    int                 GetLayerCount() { return nLayers; }
-
+
+
+
+
+
+
-
-    int                 TestCapability( const char * );
@@ -121 +162 @@
-                ~OGRBNADriver();
-
-
-
-
+
+
+
+
+
-    
-};

trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnadatasource.cpp

@@ -43 +43 @@
-    papoLayers = NULL;
-    nLayers = 0;
+    
+    fpOutput = NULL;
-
-    pszName = NULL;
+    
+    pszCoordinateSeparator = NULL;
-
-    bUpdate = FALSE;
@@ -56 +60 @@
-
-{
+    if ( fpOutput != NULL )
+    {
+        VSIFClose( fpOutput);
+    }
+
-    for( int i = 0; i < nLayers; i++ )
-        delete papoLayers[i];
-    CPLFree( papoLayers );
+    
+    CPLFree( pszCoordinateSeparator );
-
-    CPLFree( pszName );
@@ -71 +82 @@
-{
-    if( EQUAL(pszCap,ODsCCreateLayer) )
-FALS
+TRU
-    else if( EQUAL(pszCap,ODsCDeleteLayer) )
-        return FALSE;
@@ -89 +100 @@
-    else
-        return papoLayers[iLayer];
+}
+
+/************************************************************************/
+/*                            CreateLayer()                             */
+/************************************************************************/
+
+OGRLayer * OGRBNADataSource::CreateLayer( const char * pszLayerName,
+                                          OGRSpatialReference *poSRS,
+                                          OGRwkbGeometryType eType,
+                                          char ** papszOptions )
+
+{
+    BNAFeatureType bnaFeatureType;
+    
+    switch(eType)
+    {
+        case wkbPolygon:
+        case wkbPolygon25D:
+        case wkbMultiPolygon:
+        case wkbMultiPolygon25D:
+            bnaFeatureType = BNA_POLYGON;
+            break;
+        
+        case wkbPoint:
+        case wkbPoint25D:
+            bnaFeatureType = BNA_POINT;
+            break;
+            
+        case wkbLineString:
+        case wkbLineString25D:
+            bnaFeatureType = BNA_POLYLINE;
+            break;
+            
+        default:
+            CPLError( CE_Failure, CPLE_NotSupported,
+                    "Geometry type of `%s' not supported in BNAs.\n",
+                    OGRGeometryTypeToName(eType) );
+            return NULL;
+    }
+    
+    nLayers++;
+    papoLayers = (OGRBNALayer **) CPLRealloc(papoLayers, nLayers * sizeof(OGRBNALayer*));
+    papoLayers[nLayers-1] = new OGRBNALayer( pszName, pszLayerName, bnaFeatureType, eType, TRUE, this );
+    
+    return papoLayers[nLayers-1];
-}
-
@@ -116 +172 @@
-    if( !EQUAL( CPLGetExtension(pszFilename), "bna" ) )
-        return FALSE;
-
- = VSIFOpen(pszFilename, "rb
-
+
+p = VSIFOpen(pszFilename, "r
+p
-    {
-        BNARecord* record;
-        int curLine = 0;
-
+
+
+
+
+
-        BNA_FreeRecord(record);
-
@@ -130 +190 @@
-
-            papoLayers = (OGRBNALayer **) CPLMalloc(nLayers * sizeof(OGRBNALayer*));
-
-
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-    }
-
@@ -142 +263 @@
-}
-
+
+/************************************************************************/
+/*                               Create()                               */
+/************************************************************************/
+
+int OGRBNADataSource::Create( const char *pszFilename, 
+                              char **papszOptions )
+{
+    if( fpOutput != NULL)
+    {
+        CPLAssert( FALSE );
+        return FALSE;
+    }
+
+/* -------------------------------------------------------------------- */
+/*     Do not override exiting file.                                    */
+/* -------------------------------------------------------------------- */
+    VSIStatBufL sStatBuf;
+
+    if( VSIStatL( pszFilename, &sStatBuf ) == 0 )
+        return FALSE;
+    
+/* -------------------------------------------------------------------- */
+/*      Create the output file.                                         */
+/* -------------------------------------------------------------------- */
+    pszName = CPLStrdup( pszFilename );
+
+    if( EQUAL(pszFilename,"stdout") )
+        fpOutput = stdout;
+    else
+        fpOutput = VSIFOpen( pszFilename, "w" );
+    if( fpOutput == NULL )
+    {
+        CPLError( CE_Failure, CPLE_OpenFailed, 
+                  "Failed to create BNA file %s.", 
+                  pszFilename );
+        return FALSE;
+    }
+    
+    /* EOL token */
+    const char *pszCRLFFormat = CSLFetchNameValue( papszOptions, "LINEFORMAT");
+
+    if( pszCRLFFormat == NULL )
+    {
+#ifdef WIN32
+        bUseCRLF = TRUE;
+#else
+        bUseCRLF = FALSE;
+#endif
+    }
+    else if( EQUAL(pszCRLFFormat,"CRLF") )
+        bUseCRLF = TRUE;
+    else if( EQUAL(pszCRLFFormat,"LF") )
+        bUseCRLF = FALSE;
+    else
+    {
+        CPLError( CE_Warning, CPLE_AppDefined, 
+                  "LINEFORMAT=%s not understood, use one of CRLF or LF.",
+                  pszCRLFFormat );
+#ifdef WIN32
+        bUseCRLF = TRUE;
+#else
+        bUseCRLF = FALSE;
+#endif
+    }
+
+    /* Multi line or single line format ? */
+    bMultiLine = CSLFetchBoolean( papszOptions, "MULTILINE", TRUE);
+    
+    /* Number of identifiers per record */
+    const char* pszNbOutID = CSLFetchNameValue ( papszOptions, "NB_IDS");
+    if (pszNbOutID == NULL)
+    {
+        nbOutID = NB_MIN_BNA_IDS;
+    }
+    else if (EQUAL(pszNbOutID, "NB_SOURCE_FIELDS"))
+    {
+        nbOutID = -1;
+    }
+    else
+    {
+        nbOutID = atoi(pszNbOutID);
+        if (nbOutID <= 0)
+        {
+            CPLError( CE_Warning, CPLE_AppDefined, 
+                  "NB_ID=%s not understood. Must be >=%d and <=%d or equal to NB_SOURCE_FIELDS",
+                  pszNbOutID, NB_MIN_BNA_IDS, NB_MAX_BNA_IDS );
+            nbOutID = NB_MIN_BNA_IDS;
+        }
+        if (nbOutID > NB_MAX_BNA_IDS)
+        {
+            CPLError( CE_Warning, CPLE_AppDefined, 
+                  "NB_ID=%s not understood. Must be >=%d and <=%d or equal to NB_SOURCE_FIELDS",
+                  pszNbOutID, NB_MIN_BNA_IDS, NB_MAX_BNA_IDS );
+            nbOutID = NB_MAX_BNA_IDS;
+        }
+    }
+    
+    /* Ellipses export as ellipses or polygons ? */
+    bEllipsesAsEllipses = CSLFetchBoolean( papszOptions, "ELLIPSES_AS_ELLIPSES", TRUE);
+    
+    /* Number of coordinate pairs per line */
+    const char* pszNbPairPerLine = CSLFetchNameValue( papszOptions, "NB_PAIRS_PER_LINE");
+    if (pszNbPairPerLine)
+    {
+        nbPairPerLine = atoi(pszNbPairPerLine);
+        if (nbPairPerLine <= 0)
+            nbPairPerLine = (bMultiLine == FALSE) ? 1000000000 : 1;
+        if (bMultiLine == FALSE)
+        {
+            CPLError( CE_Warning, CPLE_AppDefined, "NB_PAIR_PER_LINE option is ignored when MULTILINE=NO");
+        }
+    }
+    else
+    {
+        nbPairPerLine = (bMultiLine == FALSE) ? 1000000000 : 1;
+    }
+    
+    /* Coordinate precision */
+    const char* pszCoordinatePrecision = CSLFetchNameValue( papszOptions, "COORDINATE_PRECISION");
+    if (pszCoordinatePrecision)
+    {
+        coordinatePrecision = atoi(pszCoordinatePrecision);
+        if (coordinatePrecision <= 0)
+            coordinatePrecision = 0;
+        else if (coordinatePrecision >= 20)
+            coordinatePrecision = 20;
+    }
+    else
+    {
+        coordinatePrecision = 10;
+    }
+    
+    pszCoordinateSeparator = (char*)CSLFetchNameValue( papszOptions, "COORDINATE_SEPARATOR");
+    if (pszCoordinateSeparator == NULL)
+        pszCoordinateSeparator = CPLStrdup(",");
+    else
+        pszCoordinateSeparator = CPLStrdup(pszCoordinateSeparator);
+
+    return TRUE;
+}

trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnadriver.cpp

@@ -69 +69 @@
-
-/************************************************************************/
+/*                          CreateDataSource()                          */
+/************************************************************************/
+
+OGRDataSource *OGRBNADriver::CreateDataSource( const char * pszName,
+                                                 char **papszOptions )
+
+{
+    OGRBNADataSource   *poDS = new OGRBNADataSource();
+
+    if( !poDS->Create( pszName, papszOptions ) )
+    {
+        delete poDS;
+        poDS = NULL;
+    }
+
+    return poDS;
+}
+
+/************************************************************************/
+/*                          DeleteDataSource()                          */
+/************************************************************************/
+
+OGRErr OGRBNADriver::DeleteDataSource( const char *pszFilename )
+
+{
+    if( VSIUnlink( pszFilename ) == 0 )
+        return OGRERR_NONE;
+    else
+        return OGRERR_FAILURE;
+}
+
+/************************************************************************/
-/*                           TestCapability()                           */
-/************************************************************************/
@@ -76 +108 @@
-{
-    if( EQUAL(pszCap,ODrCCreateDataSource) )
-FALS
+TRU
-    else if( EQUAL(pszCap,ODrCDeleteDataSource) )
-FALS
+TRU
-    else
-        return FALSE;

trunk/gdal/ogr/ogrsf_frmts/bna/ogrbnalayer.cpp

@@ -1 +1 @@
-/******************************************************************************
- 10646 2007-01-18 02:38:10Z warmerdam $
+
- *
- * Project:  BNA Translator
@@ -47 +47 @@
-                          const char* layerName,
-                          BNAFeatureType bnaFeatureType,
-
+
+
+
+
-
-{
@@ -54 +57 @@
-    curLine = 0;
-    nNextFID = 0;
+    
+    this->bWriter = bWriter;
+    this->poDS = poDS;
+    this->nIDs = nIDs;
-
-    nFeatures = 0;
@@ -69 +76 @@
-    this->bnaFeatureType = bnaFeatureType;
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
-
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
-}
-
@@ -113 +127 @@
-    if (fpBNA)
-        VSIFClose( fpBNA );
+}
+
+/************************************************************************/
+/*                         SetFeatureIndexTable()                       */
+/************************************************************************/
+void  OGRBNALayer::SetFeatureIndexTable(int nFeatures, OffsetAndLine* offsetAndLineFeaturesTable, int partialIndexTable)
+{
+    this->nFeatures = nFeatures;
+    this->offsetAndLineFeaturesTable = offsetAndLineFeaturesTable;
+    this->partialIndexTable = partialIndexTable;
-}
-
@@ -155 +179 @@
-        if (ok == FALSE)
-        {
+            BNA_FreeRecord(record);
-            failed = TRUE;
-            return NULL;
@@ -162 +187 @@
-            /* end of file */
-            eof = TRUE;
-
+
-            /* and we have finally build the whole index table */
-            partialIndexTable = FALSE;
@@ -170 +195 @@
-        if (record->featureType == bnaFeatureType)
-        {
+            if (nNextFID >= nFeatures)
+            {
+                nFeatures++;
+                offsetAndLineFeaturesTable =
+                    (OffsetAndLine*)CPLRealloc(offsetAndLineFeaturesTable, nFeatures * sizeof(OffsetAndLine));
+                offsetAndLineFeaturesTable[nFeatures-1].offset = offset;
+                offsetAndLineFeaturesTable[nFeatures-1].line = line;
+            }
+
+            poFeature = BuildFeatureFromBNARecord(record, nNextFID++);
+
+            BNA_FreeRecord(record);
+
+            if(   (m_poFilterGeom == NULL
+                || FilterGeometry( poFeature->GetGeometryRef() ) )
+            && (m_poAttrQuery == NULL
+                || m_poAttrQuery->Evaluate( poFeature )) )
+            {
+                 return poFeature;
+            }
+
+            delete poFeature;
+        }
+        else
+        {
+            BNA_FreeRecord(record);
+        }
+    }
+}
+
+
+void OGRBNALayer::WriteFeatureAttributes(FILE* fp, OGRFeature *poFeature )
+{
+    int i;
+    OGRFieldDefn *poField;
+    int nbOutID = poDS->GetNbOutId();
+    if (nbOutID < 0)
+        nbOutID = poFeatureDefn->GetFieldCount();
+    for(i=0;i<nbOutID;i++)
+    { 
+        if (i < poFeatureDefn->GetFieldCount())
+        {
+            poField = poFeatureDefn->GetFieldDefn( i );
+            if( poFeature->IsFieldSet( i ) )
+            {
+                const char *pszRaw = poFeature->GetFieldAsString( i );
+                VSIFPrintf( fp, "\"%s\",", pszRaw);
+            }
+            else
+            {
+                VSIFPrintf( fp, "\"\",");
+            }
+        }
+        else
+        {
+            VSIFPrintf( fp, "\"\",");
+        }
+    }
+}
+
+/************************************************************************/
+/*                           CreateFeature()                            */
+/************************************************************************/
+
+OGRErr OGRBNALayer::CreateFeature( OGRFeature *poFeature )
+
+{
+    int i,j,k,n;
+    OGRGeometry     *poGeom = poFeature->GetGeometryRef();
+    char eol[3];
+    const char* partialEol = (poDS->GetMultiLine()) ? eol : poDS->GetCoordinateSeparator();
+    
+    if (poDS->GetUseCRLF())
+    {
+        eol[0] = 13;
+        eol[1] = 10;
+        eol[2] = 0;
+    }
+    else
+    {
+        eol[0] = 10;
+        eol[1] = 0;
+    }
+    
+    if ( ! bWriter )
+    {
+        return OGRERR_FAILURE;
+    }
+    
+    if( poFeature->GetFID() == OGRNullFID )
+        poFeature->SetFID( nFeatures++ );
+    
+    FILE* fp = poDS->GetOutputFP();
+    int nbPairPerLine = poDS->GetNbPairPerLine();
+    char formatCoordinates[32];
+    sprintf(formatCoordinates, "%%s%%.%df%s%%.%df",
+            poDS->GetCoordinatePrecision(), poDS->GetCoordinateSeparator(), poDS->GetCoordinatePrecision());
+
+    switch( poGeom->getGeometryType() )
+    {
+        case wkbPoint:
+        case wkbPoint25D:
+        {
+            OGRPoint* point = (OGRPoint*)poGeom;
+            WriteFeatureAttributes(fp, poFeature);
+            VSIFPrintf( fp, "1");
+            VSIFPrintf( fp, formatCoordinates, partialEol, point->getX(), point->getY());
+            VSIFPrintf( fp, "%s", eol);
-            break;
-        }
+            
+        case wkbPolygon:
+        case wkbPolygon25D:
+        {
+            OGRPolygon* polygon = (OGRPolygon*)poGeom;
+            OGRLinearRing* ring = polygon->getExteriorRing();
+            double firstX = ring->getX(0);
+            double firstY = ring->getY(0);
+            int nBNAPoints = ring->getNumPoints();
+            int is_ellipse = FALSE;
+            
+            /* This code tries to detect an ellipse in a polygon geometry */
+            /* This will only work presumably on ellipses already read from a BNA file */
+            /* Mostly a BNA to BNA feature... */
+            if (poDS->GetEllipsesAsEllipses() &&
+                polygon->getNumInteriorRings() == 0 &&
+                nBNAPoints == 361)
+            {
+                double oppositeX = ring->getX(180);
+                double oppositeY = ring->getY(180);
+                double quarterX = ring->getX(90);
+                double quarterY = ring->getY(90);
+                double antiquarterX = ring->getX(270);
+                double antiquarterY = ring->getY(270);
+                double center1X = 0.5*(firstX + oppositeX);
+                double center1Y = 0.5*(firstY + oppositeY);
+                double center2X = 0.5*(quarterX + antiquarterX);
+                double center2Y = 0.5*(quarterY + antiquarterY);
+                if (fabs(center1X - center2X) < 1e-5 && fabs(center1Y - center2Y) < 1e-5 &&
+                    fabs(oppositeY - firstY) < 1e-5 &&
+                    fabs(quarterX - antiquarterX) < 1e-5)
+                {
+                    double major_radius = fabs(firstX - center1X);
+                    double minor_radius = fabs(quarterY - center1Y);
+                    is_ellipse = TRUE;
+                    for(i=0;i<360;i++)
+                    {
+                        if (!(fabs(center1X + major_radius * cos(i * (M_PI / 180)) - ring->getX(i)) < 1e-5 &&
+                              fabs(center1Y + minor_radius * sin(i * (M_PI / 180)) - ring->getY(i)) < 1e-5))
+                        {
+                            is_ellipse = FALSE;
+                            break;
+                        }
+                    }
+                    if ( is_ellipse == TRUE )
+                    {
+                        WriteFeatureAttributes(fp, poFeature);
+                        VSIFPrintf( fp, "2");
+                        VSIFPrintf( fp, formatCoordinates, partialEol, center1X, center1Y);
+                        VSIFPrintf( fp, formatCoordinates,  partialEol, major_radius, minor_radius);
+                        VSIFPrintf( fp, "%s", eol);
+                    }
+                }
+            }
+
+            if ( is_ellipse == FALSE)
+            {
+                int nInteriorRings = polygon->getNumInteriorRings();
+                for(i=0;i<nInteriorRings;i++)
+                {
+                    nBNAPoints += polygon->getInteriorRing(i)->getNumPoints() + 1;
+                }
+                if (nBNAPoints <= 3)
+                {
+                    CPLError( CE_Failure, CPLE_AppDefined, "Invalid geometry" );
+                    return OGRERR_FAILURE;
+                }
+                WriteFeatureAttributes(fp, poFeature);
+                VSIFPrintf( fp, "%d", nBNAPoints);
+                n = ring->getNumPoints();
+                int nbPair = 0;
+                for(i=0;i<n;i++)
+                {
+                    VSIFPrintf( fp, formatCoordinates,
+                                ((nbPair % nbPairPerLine) == 0) ? partialEol : " ", ring->getX(i), ring->getY(i));
+                    nbPair++;
+                }
+                for(i=0;i<nInteriorRings;i++)
+                {
+                    ring = polygon->getInteriorRing(i);
+                    n = ring->getNumPoints();
+                    for(j=0;i<n;i++)
+                    {
+                        VSIFPrintf( fp, formatCoordinates,
+                                    ((nbPair % nbPairPerLine) == 0) ? partialEol : " ", ring->getX(j), ring->getY(j));
+                        nbPair++;
+                    }
+                    VSIFPrintf( fp, formatCoordinates,
+                                ((nbPair % nbPairPerLine) == 0) ? partialEol : " ", firstX, firstY);
+                    nbPair++;
+                }
+                VSIFPrintf( fp, "%s", eol);
+            }
+            break;
+        }
+
+        case wkbMultiPolygon:
+        case wkbMultiPolygon25D:
+        {
+            OGRMultiPolygon* multipolygon = (OGRMultiPolygon*)poGeom;
+            int N = multipolygon->getNumGeometries();
+            int nBNAPoints = 0;
+            double firstX = 0, firstY = 0; 
+            for(i=0;i<N;i++)
+            {
+                OGRPolygon* polygon = (OGRPolygon*)multipolygon->getGeometryRef(i);
+                OGRLinearRing* ring = polygon->getExteriorRing();
+                if (nBNAPoints)
+                    nBNAPoints ++;
+                else
+                {
+                    firstX = ring->getX(0);
+                    firstY = ring->getY(0);
+                }
+                nBNAPoints += ring->getNumPoints();
+                int nInteriorRings = polygon->getNumInteriorRings();
+                for(j=0;j<nInteriorRings;j++)
+                {
+                    nBNAPoints += polygon->getInteriorRing(j)->getNumPoints() + 1;
+                }
+            }
+            if (nBNAPoints <= 3)
+            {
+                CPLError( CE_Failure, CPLE_AppDefined, "Invalid geometry" );
+