Line data    Source code

       1             : /******************************************************************************
       2             :  * Project:  Selafin importer
       3             :  * Purpose:  Implementation of OGRSelafinLayer class.
       4             :  * Author:   François Hissel, francois.hissel@gmail.com
       5             :  *
       6             :  ******************************************************************************
       7             :  * Copyright (c) 2014,  François Hissel <francois.hissel@gmail.com>
       8             :  *
       9             :  * SPDX-License-Identifier: MIT
      10             :  ****************************************************************************/
      11             : 
      12             : #include <cstdlib>
      13             : #include "cpl_conv.h"
      14             : #include "cpl_string.h"
      15             : #include "ogr_p.h"
      16             : #include "io_selafin.h"
      17             : #include "ogr_selafin.h"
      18             : #include "cpl_error.h"
      19             : #include "cpl_quad_tree.h"
      20             : #include "cpl_vsi_virtual.h"
      21             : 
      22             : /************************************************************************/
      23             : /*                           Utilities functions                        */
      24             : /************************************************************************/
      25          74 : static void MoveOverwrite(VSILFILE *fpDest, VSILFILE *fpSource)
      26             : {
      27          74 :     VSIRewindL(fpSource);
      28          74 :     VSIRewindL(fpDest);
      29          74 :     VSIFTruncateL(fpDest, 0);
      30             :     char anBuf[0x10000];
      31         148 :     while (!fpSource->Eof() && !fpSource->Error())
      32             :     {
      33          74 :         size_t nSize = VSIFReadL(anBuf, 1, sizeof(anBuf), fpSource);
      34          74 :         size_t nLeft = nSize;
      35         148 :         while (nLeft > 0)
      36          74 :             nLeft -= VSIFWriteL(anBuf + nSize - nLeft, 1, nLeft, fpDest);
      37             :     }
      38          74 :     VSIFCloseL(fpSource);
      39          74 :     VSIFFlushL(fpDest);
      40          74 : }
      41             : 
      42             : /************************************************************************/
      43             : /*                            OGRSelafinLayer()                         */
      44             : /*       Note that no operation on OGRSelafinLayer is thread-safe       */
      45             : /************************************************************************/
      46             : 
      47          16 : OGRSelafinLayer::OGRSelafinLayer(GDALDataset *poDS, const char *pszLayerNameP,
      48             :                                  int bUpdateP,
      49             :                                  const OGRSpatialReference *poSpatialRefP,
      50             :                                  Selafin::Header *poHeaderP, int nStepNumberP,
      51          16 :                                  SelafinTypeDef eTypeP)
      52          16 :     : m_poDS(poDS), eType(eTypeP), bUpdate(CPL_TO_BOOL(bUpdateP)),
      53             :       nStepNumber(nStepNumberP), poHeader(poHeaderP),
      54             :       poFeatureDefn(
      55          16 :           new OGRFeatureDefn(CPLGetBasenameSafe(pszLayerNameP).c_str())),
      56          32 :       poSpatialRef(nullptr), nCurrentId(-1)
      57             : {
      58             : #ifdef DEBUG_VERBOSE
      59             :     CPLDebug("Selafin", "Opening layer %s", pszLayerNameP);
      60             : #endif
      61          16 :     SetDescription(poFeatureDefn->GetName());
      62          16 :     poFeatureDefn->Reference();
      63          16 :     if (eType == POINTS)
      64           8 :         poFeatureDefn->SetGeomType(wkbPoint);
      65             :     else
      66           8 :         poFeatureDefn->SetGeomType(wkbPolygon);
      67          16 :     if (poSpatialRefP)
      68             :     {
      69           8 :         poSpatialRef = poSpatialRefP->Clone();
      70           8 :         poSpatialRef->SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
      71             :     }
      72          26 :     for (int i = 0; i < poHeader->nVar; ++i)
      73             :     {
      74          20 :         OGRFieldDefn oFieldDefn(poHeader->papszVariables[i], OFTReal);
      75          10 :         poFeatureDefn->AddFieldDefn(&oFieldDefn);
      76             :     }
      77          16 : }
      78             : 
      79             : /************************************************************************/
      80             : /*                           ~OGRSelafinLayer()                         */
      81             : /************************************************************************/
      82          32 : OGRSelafinLayer::~OGRSelafinLayer()
      83             : {
      84             : #ifdef DEBUG_VERBOSE
      85             :     CPLDebug("Selafin", "Closing layer %s", GetName());
      86             : #endif
      87          16 :     poFeatureDefn->Release();
      88          16 :     if (poSpatialRef)
      89           8 :         poSpatialRef->Release();
      90             :     // poHeader->nRefCount--;
      91             :     // if (poHeader->nRefCount==0) delete poHeader;
      92          32 : }
      93             : 
      94             : /************************************************************************/
      95             : /*                           GetNextFeature()                           */
      96             : /************************************************************************/
      97           1 : OGRFeature *OGRSelafinLayer::GetNextFeature()
      98             : {
      99             :     // CPLDebug("Selafin","GetNextFeature(%li)",nCurrentId+1);
     100             :     while (true)
     101             :     {
     102           1 :         OGRFeature *poFeature = GetFeature(++nCurrentId);
     103           1 :         if (poFeature == nullptr)
     104           0 :             return nullptr;
     105           2 :         if ((m_poFilterGeom == nullptr ||
     106           2 :              FilterGeometry(poFeature->GetGeometryRef())) &&
     107           1 :             (m_poAttrQuery == nullptr || m_poAttrQuery->Evaluate(poFeature)))
     108           1 :             return poFeature;
     109           0 :         delete poFeature;
     110           0 :     }
     111             : }
     112             : 
     113             : /************************************************************************/
     114             : /*                            ResetReading()                            */
     115             : /************************************************************************/
     116           1 : void OGRSelafinLayer::ResetReading()
     117             : {
     118             :     // CPLDebug("Selafin","ResetReading()");
     119           1 :     nCurrentId = -1;
     120           1 : }
     121             : 
     122             : /************************************************************************/
     123             : /*                           SetNextByIndex()                           */
     124             : /************************************************************************/
     125           0 : OGRErr OGRSelafinLayer::SetNextByIndex(GIntBig nIndex)
     126             : {
     127             :     // CPLDebug("Selafin","SetNexByIndex(%li)",nIndex);
     128           0 :     if (nIndex < 0 || nIndex >= poHeader->nPoints)
     129           0 :         return OGRERR_FAILURE;
     130           0 :     nCurrentId = nIndex - 1;
     131           0 :     return OGRERR_NONE;
     132             : }
     133             : 
     134             : /************************************************************************/
     135             : /*                           TestCapability()                           */
     136             : /************************************************************************/
     137          14 : int OGRSelafinLayer::TestCapability(const char *pszCap)
     138             : {
     139             :     // CPLDebug("Selafin","TestCapability(%s)",pszCap);
     140          14 :     if (EQUAL(pszCap, OLCRandomRead))
     141           0 :         return TRUE;
     142          14 :     if (EQUAL(pszCap, OLCSequentialWrite))
     143           1 :         return bUpdate;
     144          13 :     if (EQUAL(pszCap, OLCRandomWrite))
     145           0 :         return bUpdate;
     146          13 :     if (EQUAL(pszCap, OLCFastSpatialFilter))
     147           0 :         return FALSE;
     148          13 :     if (EQUAL(pszCap, OLCFastFeatureCount))
     149           0 :         return TRUE;
     150          13 :     if (EQUAL(pszCap, OLCFastGetExtent))
     151           0 :         return TRUE;
     152          13 :     if (EQUAL(pszCap, OLCFastSetNextByIndex))
     153           0 :         return TRUE;
     154          13 :     if (EQUAL(pszCap, OLCCreateField))
     155           1 :         return bUpdate;
     156          12 :     if (EQUAL(pszCap, OLCCreateGeomField))
     157           0 :         return FALSE;
     158          12 :     if (EQUAL(pszCap, OLCDeleteField))
     159           0 :         return bUpdate;
     160          12 :     if (EQUAL(pszCap, OLCReorderFields))
     161           0 :         return bUpdate;
     162          12 :     if (EQUAL(pszCap, OLCAlterFieldDefn))
     163           0 :         return bUpdate;
     164          12 :     if (EQUAL(pszCap, OLCDeleteFeature))
     165           0 :         return bUpdate;
     166          12 :     if (EQUAL(pszCap, OLCStringsAsUTF8))
     167           0 :         return FALSE;
     168          12 :     if (EQUAL(pszCap, OLCTransactions))
     169           0 :         return FALSE;
     170          12 :     if (EQUAL(pszCap, OLCIgnoreFields))
     171           0 :         return FALSE;
     172          12 :     return FALSE;
     173             : }
     174             : 
     175             : /************************************************************************/
     176             : /*                            GetFeature()                              */
     177             : /************************************************************************/
     178          29 : OGRFeature *OGRSelafinLayer::GetFeature(GIntBig nFID)
     179             : {
     180          29 :     CPLDebug("Selafin", "GetFeature(" CPL_FRMT_GIB ")", nFID);
     181          29 :     if (nFID < 0)
     182           0 :         return nullptr;
     183          29 :     if (eType == POINTS)
     184             :     {
     185          27 :         if (nFID >= poHeader->nPoints)
     186           0 :             return nullptr;
     187          27 :         OGRFeature *poFeature = new OGRFeature(poFeatureDefn);
     188          27 :         poFeature->SetGeometryDirectly(new OGRPoint(
     189          27 :             poHeader->paadfCoords[0][nFID], poHeader->paadfCoords[1][nFID]));
     190          27 :         poFeature->SetFID(nFID);
     191          80 :         for (int i = 0; i < poHeader->nVar; ++i)
     192             :         {
     193          53 :             VSIFSeekL(poHeader->fp,
     194          53 :                       poHeader->getPosition(nStepNumber, (int)nFID, i),
     195             :                       SEEK_SET);
     196          53 :             double nData = 0.0;
     197          53 :             if (Selafin::read_float(poHeader->fp, nData) == 1)
     198          53 :                 poFeature->SetField(i, nData);
     199             :         }
     200          27 :         return poFeature;
     201             :     }
     202             :     else
     203             :     {
     204           2 :         if (nFID >= poHeader->nElements)
     205           0 :             return nullptr;
     206             :         double *anData =
     207           2 :             (double *)VSI_MALLOC2_VERBOSE(sizeof(double), poHeader->nVar);
     208           2 :         if (poHeader->nVar > 0 && anData == nullptr)
     209           0 :             return nullptr;
     210           4 :         for (int i = 0; i < poHeader->nVar; ++i)
     211           2 :             anData[i] = 0;
     212           2 :         OGRFeature *poFeature = new OGRFeature(poFeatureDefn);
     213           2 :         poFeature->SetFID(nFID);
     214           2 :         OGRPolygon *poPolygon = new OGRPolygon();
     215           2 :         OGRLinearRing *poLinearRing = new OGRLinearRing();
     216          10 :         for (int j = 0; j < poHeader->nPointsPerElement; ++j)
     217             :         {
     218           8 :             int nPointNum =
     219           8 :                 poHeader
     220           8 :                     ->panConnectivity[nFID * poHeader->nPointsPerElement + j] -
     221             :                 1;
     222           8 :             poLinearRing->addPoint(poHeader->paadfCoords[0][nPointNum],
     223           8 :                                    poHeader->paadfCoords[1][nPointNum]);
     224          16 :             for (int i = 0; i < poHeader->nVar; ++i)
     225             :             {
     226           8 :                 VSIFSeekL(poHeader->fp,
     227           8 :                           poHeader->getPosition(nStepNumber, nPointNum, i),
     228             :                           SEEK_SET);
     229           8 :                 double nData = 0.0;
     230           8 :                 if (Selafin::read_float(poHeader->fp, nData) == 1)
     231           8 :                     anData[i] += nData;
     232             :             }
     233             :         }
     234           2 :         poPolygon->addRingDirectly(poLinearRing);
     235           2 :         poPolygon->closeRings();
     236           2 :         poFeature->SetGeometryDirectly(poPolygon);
     237           2 :         if (poHeader->nPointsPerElement)
     238             :         {
     239           4 :             for (int i = 0; i < poHeader->nVar; ++i)
     240           2 :                 poFeature->SetField(i, anData[i] / poHeader->nPointsPerElement);
     241             :         }
     242           2 :         CPLFree(anData);
     243           2 :         return poFeature;
     244             :     }
     245             : }
     246             : 
     247             : /************************************************************************/
     248             : /*                           GetFeatureCount()                          */
     249             : /************************************************************************/
     250           4 : GIntBig OGRSelafinLayer::GetFeatureCount(int bForce)
     251             : {
     252             :     // CPLDebug("Selafin","GetFeatureCount(%i)",bForce);
     253           4 :     if (m_poFilterGeom == nullptr && m_poAttrQuery == nullptr)
     254           4 :         return (eType == POINTS) ? poHeader->nPoints : poHeader->nElements;
     255           0 :     if (!bForce)
     256           0 :         return -1;
     257           0 :     int i = 0;
     258           0 :     int nFeatureCount = 0;
     259           0 :     const int nMax = eType == POINTS ? poHeader->nPoints : poHeader->nElements;
     260           0 :     while (i < nMax)
     261             :     {
     262           0 :         OGRFeature *poFeature = GetFeature(i++);
     263           0 :         if ((m_poFilterGeom == nullptr ||
     264           0 :              FilterGeometry(poFeature->GetGeometryRef())) &&
     265           0 :             (m_poAttrQuery == nullptr || m_poAttrQuery->Evaluate(poFeature)))
     266           0 :             ++nFeatureCount;
     267           0 :         delete poFeature;
     268             :     }
     269           0 :     return nFeatureCount;
     270             : }
     271             : 
     272             : /************************************************************************/
     273             : /*                            IGetExtent()                              */
     274             : /************************************************************************/
     275           0 : OGRErr OGRSelafinLayer::IGetExtent(int /* iGeomField*/, OGREnvelope *psExtent,
     276             :                                    bool /* bForce*/)
     277             : {
     278             :     // CPLDebug("Selafin","GetExtent(%i)",bForce);
     279           0 :     if (poHeader->nPoints == 0)
     280           0 :         return OGRERR_NONE;
     281           0 :     CPLRectObj *poObj = poHeader->getBoundingBox();
     282           0 :     psExtent->MinX = poObj->minx;
     283           0 :     psExtent->MaxX = poObj->maxx;
     284           0 :     psExtent->MinY = poObj->miny;
     285           0 :     psExtent->MaxY = poObj->maxy;
     286           0 :     delete poObj;
     287           0 :     return OGRERR_NONE;
     288             : }
     289             : 
     290             : /************************************************************************/
     291             : /*                             ISetFeature()                             */
     292             : /************************************************************************/
     293          25 : OGRErr OGRSelafinLayer::ISetFeature(OGRFeature *poFeature)
     294             : {
     295          25 :     OGRGeometry *poGeom = poFeature->GetGeometryRef();
     296          25 :     if (poGeom == nullptr)
     297           0 :         return OGRERR_FAILURE;
     298          25 :     if (eType == POINTS)
     299             :     {
     300             :         // If it is a point layer, it is the "easy" case: we change the
     301             :         // coordinates and attributes of the feature and update the file
     302          25 :         if (poGeom->getGeometryType() != wkbPoint)
     303             :         {
     304           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     305             :                      "The new feature should be of the same Point geometry as "
     306             :                      "the existing ones in the layer.");
     307           0 :             return OGRERR_FAILURE;
     308             :         }
     309          25 :         OGRPoint *poPoint = poGeom->toPoint();
     310          25 :         GIntBig nFID = poFeature->GetFID();
     311          25 :         poHeader->paadfCoords[0][nFID] = poPoint->getX();
     312          25 :         poHeader->paadfCoords[1][nFID] = poPoint->getY();
     313          25 :         CPLDebug("Selafin", "SetFeature(" CPL_FRMT_GIB ",%f,%f)", nFID,
     314          25 :                  poHeader->paadfCoords[0][nFID],
     315          25 :                  poHeader->paadfCoords[1][nFID]);
     316          50 :         if (VSIFSeekL(
     317          25 :                 poHeader->fp,
     318          50 :                 88 + 16 + 40 * poHeader->nVar + 48 +
     319          25 :                     ((poHeader->panStartDate != nullptr) ? 32 : 0) + 24 +
     320          25 :                     (poHeader->nElements * poHeader->nPointsPerElement + 2) *
     321          25 :                         4 +
     322          25 :                     (poHeader->nPoints + 2) * 4 + 4 + nFID * 4,
     323          25 :                 SEEK_SET) != 0)
     324           0 :             return OGRERR_FAILURE;
     325          25 :         CPLDebug("Selafin", "Write_float(" CPL_FRMT_GUIB ",%f)",
     326          25 :                  VSIFTellL(poHeader->fp),
     327          25 :                  poHeader->paadfCoords[0][nFID] - poHeader->adfOrigin[0]);
     328          50 :         if (Selafin::write_float(poHeader->fp, poHeader->paadfCoords[0][nFID] -
     329          25 :                                                    poHeader->adfOrigin[0]) == 0)
     330           0 :             return OGRERR_FAILURE;
     331          50 :         if (VSIFSeekL(
     332          25 :                 poHeader->fp,
     333          50 :                 88 + 16 + 40 * poHeader->nVar + 48 +
     334          25 :                     ((poHeader->panStartDate != nullptr) ? 32 : 0) + 24 +
     335          25 :                     (poHeader->nElements * poHeader->nPointsPerElement + 2) *
     336          25 :                         4 +
     337          25 :                     (poHeader->nPoints + 2) * 4 + (poHeader->nPoints + 2) * 4 +
     338          25 :                     4 + nFID * 4,
     339          25 :                 SEEK_SET) != 0)
     340           0 :             return OGRERR_FAILURE;
     341          25 :         CPLDebug("Selafin", "Write_float(" CPL_FRMT_GUIB ",%f)",
     342          25 :                  VSIFTellL(poHeader->fp),
     343          25 :                  poHeader->paadfCoords[1][nFID] - poHeader->adfOrigin[1]);
     344          50 :         if (Selafin::write_float(poHeader->fp, poHeader->paadfCoords[1][nFID] -
     345          25 :                                                    poHeader->adfOrigin[1]) == 0)
     346           0 :             return OGRERR_FAILURE;
     347          75 :         for (int i = 0; i < poHeader->nVar; ++i)
     348             :         {
     349          50 :             double nData = poFeature->GetFieldAsDouble(i);
     350         100 :             if (VSIFSeekL(poHeader->fp,
     351          50 :                           poHeader->getPosition(nStepNumber, (int)nFID, i),
     352          50 :                           SEEK_SET) != 0)
     353           0 :                 return OGRERR_FAILURE;
     354          50 :             if (Selafin::write_float(poHeader->fp, nData) == 0)
     355           0 :                 return OGRERR_FAILURE;
     356             :         }
     357             :     }
     358             :     else
     359             :     {
     360             :         // Else, we have a layer of polygonal elements. Here we consider that
     361             :         // the vertices are moved when we change the geometry (which will also
     362             :         // lead to a modification in the corresponding point layer). The
     363             :         // attributes table can't be changed, because attributes are calculated
     364             :         // from those of the vertices. First we check that the new feature is a
     365             :         // polygon with the right number of vertices
     366           0 :         if (poGeom->getGeometryType() != wkbPolygon)
     367             :         {
     368           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     369             :                      "The new feature should be of the same Polygon geometry "
     370             :                      "as the existing ones in the layer.");
     371           0 :             return OGRERR_FAILURE;
     372             :         }
     373           0 :         OGRLinearRing *poLinearRing = poGeom->toPolygon()->getExteriorRing();
     374           0 :         if (poLinearRing->getNumPoints() != poHeader->nPointsPerElement + 1)
     375             :         {
     376           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     377             :                      "The new feature should have the same number of vertices "
     378             :                      "%d as the existing ones in the layer.",
     379           0 :                      poHeader->nPointsPerElement);
     380           0 :             return OGRERR_FAILURE;
     381             :         }
     382           0 :         CPLError(CE_Warning, CPLE_AppDefined,
     383             :                  "The attributes of elements layer in Selafin files can't be "
     384             :                  "updated.");
     385           0 :         CPLDebug("Selafin", "SetFeature(" CPL_FRMT_GIB ",%f,%f,%f,%f,%f,%f)",
     386             :                  poFeature->GetFID(), poLinearRing->getX(0),
     387             :                  poLinearRing->getY(0), poLinearRing->getX(1),
     388             :                  poLinearRing->getY(1), poLinearRing->getX(2),
     389             :                  poLinearRing->getY(
     390             :                      2));  //!< This is not safe as we can't be sure there are
     391             :                            //!< at least three vertices in the linear ring, but
     392             :                            //!< we can assume that for a debug mode
     393           0 :         int nFID = static_cast<int>(poFeature->GetFID());
     394             :         // Now we change the coordinates of points in the layer based on the
     395             :         // vertices of the new polygon. We don't look at the order of points and
     396             :         // we assume that it is the same as in the original layer.
     397           0 :         for (int i = 0; i < poHeader->nPointsPerElement; ++i)
     398             :         {
     399           0 :             int nPointId =
     400           0 :                 poHeader
     401           0 :                     ->panConnectivity[nFID * poHeader->nPointsPerElement + i] -
     402             :                 1;
     403           0 :             poHeader->paadfCoords[0][nPointId] = poLinearRing->getX(i);
     404           0 :             poHeader->paadfCoords[1][nPointId] = poLinearRing->getY(i);
     405           0 :             if (VSIFSeekL(
     406           0 :                     poHeader->fp,
     407           0 :                     88 + 16 + 40 * poHeader->nVar + 48 +
     408           0 :                         ((poHeader->panStartDate != nullptr) ? 32 : 0) + 24 +
     409           0 :                         (poHeader->nElements * poHeader->nPointsPerElement +
     410           0 :                          2) *
     411           0 :                             4 +
     412           0 :                         (poHeader->nPoints + 2) * 4 + 4 + nPointId * 4,
     413           0 :                     SEEK_SET) != 0)
     414           0 :                 return OGRERR_FAILURE;
     415           0 :             CPLDebug("Selafin", "Write_float(" CPL_FRMT_GUIB ",%f)",
     416           0 :                      VSIFTellL(poHeader->fp),
     417           0 :                      poHeader->paadfCoords[0][nPointId] -
     418           0 :                          poHeader->adfOrigin[0]);
     419           0 :             if (Selafin::write_float(poHeader->fp,
     420           0 :                                      poHeader->paadfCoords[0][nPointId] -
     421           0 :                                          poHeader->adfOrigin[0]) == 0)
     422           0 :                 return OGRERR_FAILURE;
     423           0 :             if (VSIFSeekL(
     424           0 :                     poHeader->fp,
     425           0 :                     88 + 16 + 40 * poHeader->nVar + 48 +
     426           0 :                         ((poHeader->panStartDate != nullptr) ? 32 : 0) + 24 +
     427           0 :                         (poHeader->nElements * poHeader->nPointsPerElement +
     428           0 :                          2) *
     429           0 :                             4 +
     430           0 :                         (poHeader->nPoints + 2) * 4 +
     431           0 :                         (poHeader->nPoints + 2) * 4 + 4 + nPointId * 4,
     432           0 :                     SEEK_SET) != 0)
     433           0 :                 return OGRERR_FAILURE;
     434           0 :             CPLDebug("Selafin", "Write_float(" CPL_FRMT_GUIB ",%f)",
     435           0 :                      VSIFTellL(poHeader->fp),
     436           0 :                      poHeader->paadfCoords[1][nPointId] -
     437           0 :                          poHeader->adfOrigin[1]);
     438           0 :             if (Selafin::write_float(poHeader->fp,
     439           0 :                                      poHeader->paadfCoords[1][nPointId] -
     440           0 :                                          poHeader->adfOrigin[1]) == 0)
     441           0 :                 return OGRERR_FAILURE;
     442             :         }
     443             :     }
     444          25 :     VSIFFlushL(poHeader->fp);
     445          25 :     poHeader->UpdateFileSize();
     446          25 :     return OGRERR_NONE;
     447             : }
     448             : 
     449             : /************************************************************************/
     450             : /*                           ICreateFeature()                            */
     451             : /************************************************************************/
     452          72 : OGRErr OGRSelafinLayer::ICreateFeature(OGRFeature *poFeature)
     453             : {
     454          72 :     OGRGeometry *poGeom = poFeature->GetGeometryRef();
     455          72 :     if (poGeom == nullptr)
     456           2 :         return OGRERR_FAILURE;
     457          70 :     if (VSIFSeekL(poHeader->fp, poHeader->getPosition(0), SEEK_SET) != 0)
     458           0 :         return OGRERR_FAILURE;
     459          70 :     if (eType == POINTS)
     460             :     {
     461             :         // If it is a point layer, it is the "easy" case: we add a new point
     462             :         // feature and update the file
     463          53 :         if (poGeom->getGeometryType() != wkbPoint)
     464             :         {
     465           1 :             CPLError(CE_Failure, CPLE_AppDefined,
     466             :                      "The new feature should be of the same Point geometry as "
     467             :                      "the existing ones in the layer.");
     468           1 :             return OGRERR_FAILURE;
     469             :         }
     470          52 :         OGRPoint *poPoint = poGeom->toPoint();
     471          52 :         poFeature->SetFID(poHeader->nPoints);
     472          52 :         CPLDebug("Selafin", "CreateFeature(%d,%f,%f)", poHeader->nPoints,
     473             :                  poPoint->getX(), poPoint->getY());
     474             :         // Change the header to add the new feature
     475          52 :         poHeader->addPoint(poPoint->getX(), poPoint->getY());
     476             :     }
     477             :     else
     478             :     {
     479             :         // This is the most difficult case. The user wants to add a polygon
     480             :         // element. First we check that it has the same number of vertices as
     481             :         // the other polygon elements in the file. If there is no other element,
     482             :         // then we define the number of vertices. Every vertex in the layer
     483             :         // should have a corresponding point in the corresponding point layer.
     484             :         // So if we add a polygon element, we also have to add points in the
     485             :         // corresponding layer. The function tries to add as few new points as
     486             :         // possible, reusing already existing points. This is generally what the
     487             :         // user will expect.
     488             : 
     489             :         // First we check that we have the required geometry
     490          17 :         if (poGeom->getGeometryType() != wkbPolygon)
     491             :         {
     492           0 :             CPLError(CE_Failure, CPLE_AppDefined,
     493             :                      "The new feature should be of the same Polygon geometry "
     494             :                      "as the existing ones in the layer.");
     495           0 :             return OGRERR_FAILURE;
     496             :         }
     497             : 
     498             :         // Now we check that we have the right number of vertices, or if this
     499             :         // number was not defined yet (0), we define it at once
     500          17 :         OGRLinearRing *poLinearRing = poGeom->toPolygon()->getExteriorRing();
     501          17 :         poFeature->SetFID(poHeader->nElements);
     502          17 :         CPLDebug("Selafin", "CreateFeature(" CPL_FRMT_GIB ",%f,%f,%f,%f,%f,%f)",
     503             :                  poFeature->GetFID(), poLinearRing->getX(0),
     504             :                  poLinearRing->getY(0), poLinearRing->getX(1),
     505             :                  poLinearRing->getY(1), poLinearRing->getX(2),
     506             :                  poLinearRing->getY(
     507             :                      2));  //!< This is not safe as we can't be sure there are
     508             :                            //!< at least three vertices in the linear ring, but
     509             :                            //!< we can assume that for a debug mode
     510          17 :         int nNum = poLinearRing->getNumPoints();
     511          17 :         if (poHeader->nPointsPerElement == 0)
     512             :         {
     513           1 :             if (nNum < 4)
     514             :             {
     515           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
     516             :                          "The new feature should have at least 3 vertices.");
     517           0 :                 return OGRERR_FAILURE;
     518             :             }
     519           1 :             poHeader->nPointsPerElement = nNum - 1;
     520           1 :             if (poHeader->nElements > 0)
     521             :             {
     522             :                 int *panConnectivity =
     523           0 :                     reinterpret_cast<int *>(VSI_REALLOC_VERBOSE(
     524             :                         poHeader->panConnectivity,
     525             :                         static_cast<size_t>(poHeader->nElements) *
     526             :                             poHeader->nPointsPerElement));
     527           0 :                 if (panConnectivity == nullptr)
     528             :                 {
     529           0 :                     VSIFree(poHeader->panConnectivity);
     530           0 :                     poHeader->panConnectivity = nullptr;
     531           0 :                     return OGRERR_FAILURE;
     532             :                 }
     533           0 :                 poHeader->panConnectivity = panConnectivity;
     534             :             }
     535             :         }
     536             :         else
     537             :         {
     538          16 :             if (poLinearRing->getNumPoints() != poHeader->nPointsPerElement + 1)
     539             :             {
     540           0 :                 CPLError(CE_Failure, CPLE_AppDefined,
     541             :                          "The new feature should have the same number of "
     542             :                          "vertices %d as the existing ones in the layer.",
     543           0 :                          poHeader->nPointsPerElement);
     544           0 :                 return OGRERR_FAILURE;
     545             :             }
     546             :         }
     547             : 
     548             :         // Now we look for vertices that are already referenced as points in the
     549             :         // file
     550          17 :         int *anMap = (int *)VSI_MALLOC2_VERBOSE(sizeof(int),
     551             :                                                 poHeader->nPointsPerElement);
     552          17 :         if (anMap == nullptr)
     553             :         {
     554           0 :             return OGRERR_FAILURE;
     555             :         }
     556          85 :         for (int i = 0; i < poHeader->nPointsPerElement; ++i)
     557          68 :             anMap[i] = -1;
     558          17 :         if (poHeader->nPoints > 0)
     559             :         {
     560          17 :             CPLRectObj *poBB = poHeader->getBoundingBox();
     561             :             double dfMaxDist =
     562          17 :                 (poBB->maxx - poBB->minx) / sqrt((double)(poHeader->nPoints)) /
     563          17 :                 1000.0;  //!< Heuristic approach to estimate a maximum distance
     564             :                          //!< such that two points are considered equal if they
     565             :                          //!< are closer from each other
     566          17 :             dfMaxDist *= dfMaxDist;
     567          17 :             delete poBB;
     568          85 :             for (int i = 0; i < poHeader->nPointsPerElement; ++i)
     569         204 :                 anMap[i] = poHeader->getClosestPoint(
     570          68 :                     poLinearRing->getX(i), poLinearRing->getY(i), dfMaxDist);
     571             :         }
     572             : 
     573             :         // We add new points if needed only
     574          85 :         for (int i = 0; i < poHeader->nPointsPerElement; ++i)
     575          68 :             if (anMap[i] == -1)
     576             :             {
     577           3 :                 poHeader->addPoint(poLinearRing->getX(i),
     578           3 :                                    poLinearRing->getY(i));
     579           3 :                 anMap[i] = poHeader->nPoints - 1;
     580             :             }
     581             : 
     582             :         // And we update the connectivity table to add the new element
     583          17 :         poHeader->nElements++;
     584          34 :         poHeader->panConnectivity = (int *)CPLRealloc(
     585          17 :             poHeader->panConnectivity,
     586          17 :             sizeof(int) * poHeader->nPointsPerElement * poHeader->nElements);
     587          85 :         for (int i = 0; i < poHeader->nPointsPerElement; ++i)
     588             :         {
     589          68 :             poHeader->panConnectivity[poHeader->nPointsPerElement *
     590          68 :                                           (poHeader->nElements - 1) +
     591          68 :                                       i] = anMap[i] + 1;
     592             :         }
     593          17 :         poHeader->setUpdated();
     594          17 :         CPLFree(anMap);
     595             :     }
     596             : 
     597             :     // Now comes the real insertion. Since values have to be inserted nearly
     598             :     // everywhere in the file and we don't want to store everything in memory to
     599             :     // overwrite it, we create a new copy of it where we write the new values
     600         138 :     const std::string osTempfile = CPLGenerateTempFilenameSafe(nullptr);
     601          69 :     VSILFILE *fpNew = VSIFOpenL(osTempfile.c_str(), "wb+");
     602          69 :     if (fpNew == nullptr)
     603             :     {
     604           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     605             :                  "Failed to open temporary file %s with write access, %s.",
     606           0 :                  osTempfile.c_str(), VSIStrerror(errno));
     607           0 :         return OGRERR_FAILURE;
     608             :     }
     609          69 :     if (Selafin::write_header(fpNew, poHeader) == 0)
     610             :     {
     611           0 :         VSIFCloseL(fpNew);
     612           0 :         VSIUnlink(osTempfile.c_str());
     613           0 :         return OGRERR_FAILURE;
     614             :     }
     615         139 :     for (int i = 0; i < poHeader->nSteps; ++i)
     616             :     {
     617          70 :         int nLen = 0;
     618          70 :         double dfDate = 0.0;
     619          70 :         if (Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     620          70 :             Selafin::read_float(poHeader->fp, dfDate) == 0 ||
     621          70 :             Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     622          70 :             Selafin::write_integer(fpNew, 4) == 0 ||
     623         210 :             Selafin::write_float(fpNew, dfDate) == 0 ||
     624          70 :             Selafin::write_integer(fpNew, 4) == 0)
     625             :         {
     626           0 :             VSIFCloseL(fpNew);
     627           0 :             VSIUnlink(osTempfile.c_str());
     628           0 :             return OGRERR_FAILURE;
     629             :         }
     630         140 :         for (int j = 0; j < poHeader->nVar; ++j)
     631             :         {
     632          70 :             double *padfValues = nullptr;
     633         140 :             if (Selafin::read_floatarray(poHeader->fp, &padfValues,
     634          70 :                                          poHeader->nFileSize) == -1)
     635             :             {
     636           0 :                 VSIFCloseL(fpNew);
     637           0 :                 VSIUnlink(osTempfile.c_str());
     638           0 :                 return OGRERR_FAILURE;
     639             :             }
     640         140 :             padfValues = (double *)CPLRealloc(
     641          70 :                 padfValues, sizeof(double) * poHeader->nPoints);
     642          70 :             if (padfValues == nullptr)
     643             :             {
     644           0 :                 VSIFCloseL(fpNew);
     645           0 :                 VSIUnlink(osTempfile.c_str());
     646           0 :                 return OGRERR_FAILURE;
     647             :             }
     648          70 :             if (eType == POINTS)
     649         106 :                 padfValues[poHeader->nPoints - 1] =
     650          53 :                     poFeature->GetFieldAsDouble(j);
     651             :             else
     652          17 :                 padfValues[poHeader->nPoints - 1] = 0;
     653         140 :             if (Selafin::write_floatarray(fpNew, padfValues,
     654          70 :                                           poHeader->nPoints) == 0)
     655             :             {
     656           0 :                 CPLFree(padfValues);
     657           0 :                 VSIFCloseL(fpNew);
     658           0 :                 VSIUnlink(osTempfile.c_str());
     659           0 :                 return OGRERR_FAILURE;
     660             :             }
     661          70 :             CPLFree(padfValues);
     662             :         }
     663             :     }
     664             : 
     665             :     // If everything went fine, we overwrite the new file with the content of
     666             :     // the old one. This way, even if something goes bad, we can still recover
     667             :     // the layer. The copy process is format-agnostic.
     668          69 :     MoveOverwrite(poHeader->fp, fpNew);
     669          69 :     VSIUnlink(osTempfile.c_str());
     670          69 :     poHeader->UpdateFileSize();
     671          69 :     return OGRERR_NONE;
     672             : }
     673             : 
     674             : /************************************************************************/
     675             : /*                           CreateField()                              */
     676             : /************************************************************************/
     677           9 : OGRErr OGRSelafinLayer::CreateField(const OGRFieldDefn *poField,
     678             :                                     CPL_UNUSED int bApproxOK)
     679             : {
     680           9 :     CPLDebug("Selafin", "CreateField(%s,%s)", poField->GetNameRef(),
     681             :              OGRFieldDefn::GetFieldTypeName(poField->GetType()));
     682             :     // Test if the field does not exist yet
     683           9 :     if (poFeatureDefn->GetFieldIndex(poField->GetNameRef()) != -1)
     684             :     {
     685             :         // Those two lines are copied from CSV driver, but I am not quite sure
     686             :         // what they actually do
     687           0 :         if (poFeatureDefn->GetGeomFieldIndex(poField->GetNameRef()) != -1)
     688           0 :             return OGRERR_NONE;
     689           0 :         if (poFeatureDefn->GetGeomFieldIndex(
     690           0 :                 CPLSPrintf("geom_%s", poField->GetNameRef())) != -1)
     691           0 :             return OGRERR_NONE;
     692           0 :         CPLError(CE_Failure, CPLE_AppDefined,
     693             :                  "Attempt to create field %s, but a field with this name "
     694             :                  "already exists.",
     695             :                  poField->GetNameRef());
     696           0 :         return OGRERR_FAILURE;
     697             :     }
     698             :     // Test if the field type is legal (only double precision values are
     699             :     // allowed)
     700           9 :     if (poField->GetType() != OFTReal)
     701             :     {
     702           5 :         CPLError(
     703             :             CE_Failure, CPLE_AppDefined,
     704             :             "Attempt to create field of type %s, but this is not supported for "
     705             :             "Selafin files (only double precision fields are allowed).",
     706             :             poField->GetFieldTypeName(poField->GetType()));
     707           5 :         return OGRERR_FAILURE;
     708             :     }
     709           4 :     if (VSIFSeekL(poHeader->fp, poHeader->getPosition(0), SEEK_SET) != 0)
     710           0 :         return OGRERR_FAILURE;
     711             :     // Change the header to add the new field
     712           4 :     poHeader->nVar++;
     713           4 :     poHeader->setUpdated();
     714           8 :     poHeader->papszVariables = (char **)CPLRealloc(
     715           4 :         poHeader->papszVariables, sizeof(char *) * poHeader->nVar);
     716           8 :     poHeader->papszVariables[poHeader->nVar - 1] =
     717           4 :         (char *)VSI_MALLOC2_VERBOSE(sizeof(char), 33);
     718           4 :     strncpy(poHeader->papszVariables[poHeader->nVar - 1], poField->GetNameRef(),
     719             :             32);
     720           4 :     poHeader->papszVariables[poHeader->nVar - 1][32] = 0;
     721           4 :     poFeatureDefn->AddFieldDefn(poField);
     722             : 
     723             :     // Now comes the real insertion. Since values have to be inserted nearly
     724             :     // everywhere in the file and we don't want to store everything in memory to
     725             :     // overwrite it, we create a new copy of it where we write the new values
     726           8 :     const std::string osTempfile = CPLGenerateTempFilenameSafe(nullptr);
     727           4 :     VSILFILE *fpNew = VSIFOpenL(osTempfile.c_str(), "wb+");
     728           4 :     if (fpNew == nullptr)
     729             :     {
     730           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     731             :                  "Failed to open temporary file %s with write access, %s.",
     732           0 :                  osTempfile.c_str(), VSIStrerror(errno));
     733           0 :         return OGRERR_FAILURE;
     734             :     }
     735           4 :     if (Selafin::write_header(fpNew, poHeader) == 0)
     736             :     {
     737           0 :         VSIFCloseL(fpNew);
     738           0 :         VSIUnlink(osTempfile.c_str());
     739           0 :         return OGRERR_FAILURE;
     740             :     }
     741           8 :     for (int i = 0; i < poHeader->nSteps; ++i)
     742             :     {
     743           4 :         int nLen = 0;
     744           4 :         double dfDate = 0.0;
     745           4 :         if (Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     746           4 :             Selafin::read_float(poHeader->fp, dfDate) == 0 ||
     747           4 :             Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     748           4 :             Selafin::write_integer(fpNew, 4) == 0 ||
     749          12 :             Selafin::write_float(fpNew, dfDate) == 0 ||
     750           4 :             Selafin::write_integer(fpNew, 4) == 0)
     751             :         {
     752           0 :             VSIFCloseL(fpNew);
     753           0 :             VSIUnlink(osTempfile.c_str());
     754           0 :             return OGRERR_FAILURE;
     755             :         }
     756           4 :         double *padfValues = nullptr;
     757           5 :         for (int j = 0; j < poHeader->nVar - 1; ++j)
     758             :         {
     759           2 :             if (Selafin::read_floatarray(poHeader->fp, &padfValues,
     760           1 :                                          poHeader->nFileSize) == -1)
     761             :             {
     762           0 :                 VSIFCloseL(fpNew);
     763           0 :                 VSIUnlink(osTempfile.c_str());
     764           0 :                 return OGRERR_FAILURE;
     765             :             }
     766           2 :             if (Selafin::write_floatarray(fpNew, padfValues,
     767           1 :                                           poHeader->nPoints) == 0)
     768             :             {
     769           0 :                 CPLFree(padfValues);
     770           0 :                 VSIFCloseL(fpNew);
     771           0 :                 VSIUnlink(osTempfile.c_str());
     772           0 :                 return OGRERR_FAILURE;
     773             :             }
     774           1 :             CPLFree(padfValues);
     775             :         }
     776           4 :         padfValues =
     777           4 :             (double *)VSI_MALLOC2_VERBOSE(sizeof(double), poHeader->nPoints);
     778          29 :         for (int k = 0; k < poHeader->nPoints; ++k)
     779          25 :             padfValues[k] = 0;
     780           4 :         if (Selafin::write_floatarray(fpNew, padfValues, poHeader->nPoints) ==
     781             :             0)
     782             :         {
     783           0 :             CPLFree(padfValues);
     784           0 :             VSIFCloseL(fpNew);
     785           0 :             VSIUnlink(osTempfile.c_str());
     786           0 :             return OGRERR_FAILURE;
     787             :         }
     788           4 :         CPLFree(padfValues);
     789             :     }
     790           4 :     MoveOverwrite(poHeader->fp, fpNew);
     791           4 :     VSIUnlink(osTempfile.c_str());
     792           4 :     poHeader->UpdateFileSize();
     793           4 :     return OGRERR_NONE;
     794             : }
     795             : 
     796             : /************************************************************************/
     797             : /*                           DeleteField()                              */
     798             : /************************************************************************/
     799           0 : OGRErr OGRSelafinLayer::DeleteField(int iField)
     800             : {
     801           0 :     CPLDebug("Selafin", "DeleteField(%i)", iField);
     802           0 :     if (VSIFSeekL(poHeader->fp, poHeader->getPosition(0), SEEK_SET) != 0)
     803           0 :         return OGRERR_FAILURE;
     804             :     // Change the header to remove the field
     805           0 :     poHeader->nVar--;
     806           0 :     poHeader->setUpdated();
     807           0 :     CPLFree(poHeader->papszVariables[iField]);
     808           0 :     for (int i = iField; i < poHeader->nVar; ++i)
     809           0 :         poHeader->papszVariables[i] = poHeader->papszVariables[i + 1];
     810           0 :     poHeader->papszVariables = (char **)CPLRealloc(
     811           0 :         poHeader->papszVariables, sizeof(char *) * poHeader->nVar);
     812           0 :     poFeatureDefn->DeleteFieldDefn(iField);
     813             : 
     814             :     // Now comes the real deletion. Since values have to be deleted nearly
     815             :     // everywhere in the file and we don't want to store everything in memory to
     816             :     // overwrite it, we create a new copy of it where we write the new values
     817           0 :     const std::string osTempfile = CPLGenerateTempFilenameSafe(nullptr);
     818           0 :     VSILFILE *fpNew = VSIFOpenL(osTempfile.c_str(), "wb+");
     819           0 :     if (fpNew == nullptr)
     820             :     {
     821           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     822             :                  "Failed to open temporary file %s with write access, %s.",
     823           0 :                  osTempfile.c_str(), VSIStrerror(errno));
     824           0 :         return OGRERR_FAILURE;
     825             :     }
     826           0 :     if (Selafin::write_header(fpNew, poHeader) == 0)
     827             :     {
     828           0 :         VSIFCloseL(fpNew);
     829           0 :         VSIUnlink(osTempfile.c_str());
     830           0 :         return OGRERR_FAILURE;
     831             :     }
     832           0 :     for (int i = 0; i < poHeader->nSteps; ++i)
     833             :     {
     834           0 :         int nLen = 0;
     835           0 :         double dfDate = 0.0;
     836           0 :         if (Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     837           0 :             Selafin::read_float(poHeader->fp, dfDate) == 0 ||
     838           0 :             Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     839           0 :             Selafin::write_integer(fpNew, 4) == 0 ||
     840           0 :             Selafin::write_float(fpNew, dfDate) == 0 ||
     841           0 :             Selafin::write_integer(fpNew, 4) == 0)
     842             :         {
     843           0 :             VSIFCloseL(fpNew);
     844           0 :             VSIUnlink(osTempfile.c_str());
     845           0 :             return OGRERR_FAILURE;
     846             :         }
     847           0 :         for (int j = 0; j < poHeader->nVar; ++j)
     848             :         {
     849           0 :             double *padfValues = nullptr;
     850           0 :             if (Selafin::read_floatarray(poHeader->fp, &padfValues,
     851           0 :                                          poHeader->nFileSize) == -1)
     852             :             {
     853           0 :                 VSIFCloseL(fpNew);
     854           0 :                 VSIUnlink(osTempfile.c_str());
     855           0 :                 return OGRERR_FAILURE;
     856             :             }
     857           0 :             if (j != iField)
     858             :             {
     859           0 :                 if (Selafin::write_floatarray(fpNew, padfValues,
     860           0 :                                               poHeader->nPoints) == 0)
     861             :                 {
     862           0 :                     CPLFree(padfValues);
     863           0 :                     VSIFCloseL(fpNew);
     864           0 :                     VSIUnlink(osTempfile.c_str());
     865           0 :                     return OGRERR_FAILURE;
     866             :                 }
     867             :             }
     868           0 :             CPLFree(padfValues);
     869             :         }
     870             :     }
     871           0 :     MoveOverwrite(poHeader->fp, fpNew);
     872           0 :     VSIUnlink(osTempfile.c_str());
     873           0 :     poHeader->UpdateFileSize();
     874           0 :     return OGRERR_NONE;
     875             : }
     876             : 
     877             : /************************************************************************/
     878             : /*                          ReorderFields()                             */
     879             : /************************************************************************/
     880           1 : OGRErr OGRSelafinLayer::ReorderFields(int *panMap)
     881             : {
     882           1 :     CPLDebug("Selafin", "ReorderFields()");
     883           1 :     if (VSIFSeekL(poHeader->fp, poHeader->getPosition(0), SEEK_SET) != 0)
     884           0 :         return OGRERR_FAILURE;
     885             :     // Change the header according to the map
     886             :     char **papszNew =
     887           1 :         (char **)VSI_MALLOC2_VERBOSE(sizeof(char *), poHeader->nVar);
     888           3 :     for (int i = 0; i < poHeader->nVar; ++i)
     889           2 :         papszNew[i] = poHeader->papszVariables[panMap[i]];
     890           1 :     CPLFree(poHeader->papszVariables);
     891           1 :     poHeader->papszVariables = papszNew;
     892           1 :     poFeatureDefn->ReorderFieldDefns(panMap);
     893             : 
     894             :     // Now comes the real change.
     895           2 :     const std::string osTempfile = CPLGenerateTempFilenameSafe(nullptr);
     896           1 :     VSILFILE *fpNew = VSIFOpenL(osTempfile.c_str(), "wb+");
     897           1 :     if (fpNew == nullptr)
     898             :     {
     899           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
     900             :                  "Failed to open temporary file %s with write access, %s.",
     901           0 :                  osTempfile.c_str(), VSIStrerror(errno));
     902           0 :         return OGRERR_FAILURE;
     903             :     }
     904           1 :     if (Selafin::write_header(fpNew, poHeader) == 0)
     905             :     {
     906           0 :         VSIFCloseL(fpNew);
     907           0 :         VSIUnlink(osTempfile.c_str());
     908           0 :         return OGRERR_FAILURE;
     909             :     }
     910           1 :     double *padfValues = nullptr;
     911           2 :     for (int i = 0; i < poHeader->nSteps; ++i)
     912             :     {
     913           1 :         int nLen = 0;
     914           1 :         double dfDate = 0.0;
     915           1 :         if (Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     916           1 :             Selafin::read_float(poHeader->fp, dfDate) == 0 ||
     917           1 :             Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
     918           1 :             Selafin::write_integer(fpNew, 4) == 0 ||
     919           3 :             Selafin::write_float(fpNew, dfDate) == 0 ||
     920           1 :             Selafin::write_integer(fpNew, 4) == 0)
     921             :         {
     922           0 :             VSIFCloseL(fpNew);
     923           0 :             VSIUnlink(osTempfile.c_str());
     924           0 :             return OGRERR_FAILURE;
     925             :         }
     926           3 :         for (int j = 0; j < poHeader->nVar; ++j)
     927             :         {
     928           2 :             if (VSIFSeekL(poHeader->fp, poHeader->getPosition(i, -1, panMap[j]),
     929           4 :                           SEEK_SET) != 0 ||
     930           2 :                 Selafin::read_floatarray(poHeader->fp, &padfValues,
     931           2 :                                          poHeader->nFileSize) == -1)
     932             :             {
     933           0 :                 VSIFCloseL(fpNew);
     934           0 :                 VSIUnlink(osTempfile.c_str());
     935           0 :                 return OGRERR_FAILURE;
     936             :             }
     937           4 :             if (Selafin::write_floatarray(fpNew, padfValues,
     938           2 :                                           poHeader->nPoints) == 0)
     939             :             {
     940           0 :                 CPLFree(padfValues);
     941           0 :                 VSIFCloseL(fpNew);
     942           0 :                 VSIUnlink(osTempfile.c_str());
     943           0 :                 return OGRERR_FAILURE;
     944             :             }
     945           2 :             CPLFree(padfValues);
     946             :         }
     947             :     }
     948           1 :     MoveOverwrite(poHeader->fp, fpNew);
     949           1 :     VSIUnlink(osTempfile.c_str());
     950           1 :     poHeader->UpdateFileSize();
     951           1 :     return OGRERR_NONE;
     952             : }
     953             : 
     954             : /************************************************************************/
     955             : /*                         AlterFieldDefn()                             */
     956             : /************************************************************************/
     957           1 : OGRErr OGRSelafinLayer::AlterFieldDefn(int iField, OGRFieldDefn *poNewFieldDefn,
     958             :                                        int /* nFlagsIn */)
     959             : {
     960           1 :     CPLDebug("Selafin", "AlterFieldDefn(%i,%s,%s)", iField,
     961             :              poNewFieldDefn->GetNameRef(),
     962             :              OGRFieldDefn::GetFieldTypeName(poNewFieldDefn->GetType()));
     963             :     // Test if the field type is legal (only double precision values are
     964             :     // allowed)
     965           1 :     if (poNewFieldDefn->GetType() != OFTReal)
     966             :     {
     967           0 :         CPLError(
     968             :             CE_Failure, CPLE_AppDefined,
     969             :             "Attempt to update field with type %s, but this is not supported "
     970             :             "for Selafin files (only double precision fields are allowed).",
     971             :             poNewFieldDefn->GetFieldTypeName(poNewFieldDefn->GetType()));
     972           0 :         return OGRERR_FAILURE;
     973             :     }
     974             :     // Since the field type can't change, only the field name is changed. We
     975             :     // change it in the header
     976           1 :     CPLFree(poHeader->papszVariables[iField]);
     977           2 :     poHeader->papszVariables[iField] =
     978           1 :         (char *)VSI_MALLOC2_VERBOSE(sizeof(char), 33);
     979           1 :     strncpy(poHeader->papszVariables[iField], poNewFieldDefn->GetNameRef(), 32);
     980           1 :     poHeader->papszVariables[iField][32] = 0;
     981             :     // And we update the file
     982           1 :     if (VSIFSeekL(poHeader->fp, 88 + 16 + 40 * iField, SEEK_SET) != 0)
     983           0 :         return OGRERR_FAILURE;
     984           1 :     if (Selafin::write_string(poHeader->fp, poHeader->papszVariables[iField],
     985           1 :                               32) == 0)
     986           0 :         return OGRERR_FAILURE;
     987           1 :     VSIFFlushL(poHeader->fp);
     988           1 :     poHeader->UpdateFileSize();
     989           1 :     return OGRERR_NONE;
     990             : }
     991             : 
     992             : /************************************************************************/
     993             : /*                          DeleteFeature()                             */
     994             : /************************************************************************/
     995           0 : OGRErr OGRSelafinLayer::DeleteFeature(GIntBig nFID)
     996             : {
     997           0 :     CPLDebug("Selafin", "DeleteFeature(" CPL_FRMT_GIB ")", nFID);
     998           0 :     if (VSIFSeekL(poHeader->fp, poHeader->getPosition(0), SEEK_SET) != 0)
     999           0 :         return OGRERR_FAILURE;
    1000             :     // Change the header to delete the feature
    1001           0 :     if (eType == POINTS)
    1002           0 :         poHeader->removePoint((int)nFID);
    1003             :     else
    1004             :     {
    1005             :         // For elements layer, we only delete the element and not the vertices
    1006           0 :         poHeader->nElements--;
    1007           0 :         for (int i = (int)nFID; i < poHeader->nElements; ++i)
    1008           0 :             for (int j = 0; j < poHeader->nPointsPerElement; ++j)
    1009           0 :                 poHeader->panConnectivity[poHeader->nPointsPerElement * i + j] =
    1010           0 :                     poHeader->panConnectivity[poHeader->nPointsPerElement *
    1011           0 :                                                   (i + 1) +
    1012           0 :                                               j];
    1013           0 :         poHeader->panConnectivity = (int *)CPLRealloc(
    1014           0 :             poHeader->panConnectivity,
    1015           0 :             sizeof(int) * poHeader->nPointsPerElement * poHeader->nElements);
    1016           0 :         poHeader->setUpdated();
    1017             :     }
    1018             : 
    1019             :     // Now we perform the deletion by creating a new temporary layer
    1020           0 :     const std::string osTempfile = CPLGenerateTempFilenameSafe(nullptr);
    1021           0 :     VSILFILE *fpNew = VSIFOpenL(osTempfile.c_str(), "wb+");
    1022           0 :     if (fpNew == nullptr)
    1023             :     {
    1024           0 :         CPLError(CE_Failure, CPLE_OpenFailed,
    1025             :                  "Failed to open temporary file %s with write access, %s.",
    1026           0 :                  osTempfile.c_str(), VSIStrerror(errno));
    1027           0 :         return OGRERR_FAILURE;
    1028             :     }
    1029           0 :     if (Selafin::write_header(fpNew, poHeader) == 0)
    1030             :     {
    1031           0 :         VSIFCloseL(fpNew);
    1032           0 :         VSIUnlink(osTempfile.c_str());
    1033           0 :         return OGRERR_FAILURE;
    1034             :     }
    1035           0 :     for (int i = 0; i < poHeader->nSteps; ++i)
    1036             :     {
    1037           0 :         int nLen = 0;
    1038           0 :         double dfDate = 0.0;
    1039           0 :         if (Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
    1040           0 :             Selafin::read_float(poHeader->fp, dfDate) == 0 ||
    1041           0 :             Selafin::read_integer(poHeader->fp, nLen, true) == 0 ||
    1042           0 :             Selafin::write_integer(fpNew, 4) == 0 ||
    1043           0 :             Selafin::write_float(fpNew, dfDate) == 0 ||
    1044           0 :             Selafin::write_integer(fpNew, 4) == 0)
    1045             :         {
    1046           0 :             VSIFCloseL(fpNew);
    1047           0 :             VSIUnlink(osTempfile.c_str());
    1048           0 :             return OGRERR_FAILURE;
    1049             :         }
    1050           0 :         for (int j = 0; j < poHeader->nVar; ++j)
    1051             :         {
    1052           0 :             double *padfValues = nullptr;
    1053           0 :             if (Selafin::read_floatarray(poHeader->fp, &padfValues,
    1054           0 :                                          poHeader->nFileSize) == -1)
    1055             :             {
    1056           0 :                 VSIFCloseL(fpNew);
    1057           0 :                 VSIUnlink(osTempfile.c_str());
    1058           0 :                 return OGRERR_FAILURE;
    1059             :             }
    1060           0 :             if (eType == POINTS)
    1061             :             {
    1062           0 :                 for (int k = (int)nFID; k <= poHeader->nPoints; ++k)
    1063           0 :                     padfValues[k - 1] = padfValues[k];
    1064             :             }
    1065           0 :             if (Selafin::write_floatarray(fpNew, padfValues,
    1066           0 :                                           poHeader->nPoints) == 0)
    1067             :             {
    1068           0 :                 CPLFree(padfValues);
    1069           0 :                 VSIFCloseL(fpNew);
    1070           0 :                 VSIUnlink(osTempfile.c_str());
    1071           0 :                 return OGRERR_FAILURE;
    1072             :             }
    1073           0 :             CPLFree(padfValues);
    1074             :         }
    1075             :     }
    1076             : 
    1077             :     // If everything went fine, we overwrite the new file with the content of
    1078             :     // the old one. This way, even if something goes bad, we can still recover
    1079             :     // the layer. The copy process is format-agnostic.
    1080           0 :     MoveOverwrite(poHeader->fp, fpNew);
    1081           0 :     VSIUnlink(osTempfile.c_str());
    1082           0 :     poHeader->UpdateFileSize();
    1083             : 
    1084           0 :     return OGRERR_NONE;
    1085             : }