1 /******************************************************************************
2  *
3  * Project:  Raster Matrix Format
4  * Purpose:  Read/write raster files used in GIS "Integratsia"
5  *           (also known as "Panorama" GIS).
6  * Author:   Andrey Kiselev, dron@ak4719.spb.edu
7  *
8  ******************************************************************************
9  * Copyright (c) 2005, Andrey Kiselev <dron@ak4719.spb.edu>
10  * Copyright (c) 2007-2012, Even Rouault <even dot rouault at spatialys.com>
11  *
12  * Permission is hereby granted, free of charge, to any person obtaining a
13  * copy of this software and associated documentation files (the "Software"),
14  * to deal in the Software without restriction, including without limitation
15  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
16  * and/or sell copies of the Software, and to permit persons to whom the
17  * Software is furnished to do so, subject to the following conditions:
18  *
19  * The above copyright notice and this permission notice shall be included
20  * in all copies or substantial portions of the Software.
21  *
22  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
23  * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
24  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
25  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
26  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
27  * DEALINGS IN THE SOFTWARE.
28  ****************************************************************************/
29 #include <algorithm>
30 
31 #include "cpl_string.h"
32 #include "gdal_frmts.h"
33 #include "ogr_spatialref.h"
34 
35 #include "rmfdataset.h"
36 
37 #include "cpl_safemaths.hpp"
38 
39 CPL_CVSID("$Id: rmfdataset.cpp 86933038c3926cd4dc3ff37c431b317abb69e602 2021-03-27 23:20:49 +0100 Even Rouault $")
40 
41 constexpr int RMF_DEFAULT_BLOCKXSIZE = 256;
42 constexpr int RMF_DEFAULT_BLOCKYSIZE = 256;
43 
44 static const char RMF_SigRSW[] = { 'R', 'S', 'W', '\0' };
45 static const char RMF_SigRSW_BE[] = { '\0', 'W', 'S', 'R' };
46 static const char RMF_SigMTW[] = { 'M', 'T', 'W', '\0' };
47 
48 static const char RMF_UnitsEmpty[] = "";
49 static const char RMF_UnitsM[] = "m";
50 static const char RMF_UnitsCM[] = "cm";
51 static const char RMF_UnitsDM[] = "dm";
52 static const char RMF_UnitsMM[] = "mm";
53 
54 constexpr double RMF_DEFAULT_SCALE = 10000.0;
55 constexpr double RMF_DEFAULT_RESOLUTION = 100.0;
56 
57 /* -------------------------------------------------------------------- */
58 /*  Note: Due to the fact that in the early versions of RMF             */
59 /*  format the field of the iEPSGCode was marked as a 'reserved',       */
60 /*  in the header on its place in many cases garbage values were written.*/
61 /*  Most of them can be weeded out by the minimum EPSG code value.      */
62 /*                                                                      */
63 /*  see: Surveying and Positioning Guidance Note Number 7, part 1       */
64 /*       Using the EPSG Geodetic Parameter Dataset p. 22                */
65 /*       http://www.epsg.org/Portals/0/373-07-1.pdf                     */
66 /* -------------------------------------------------------------------- */
67 constexpr GInt32 RMF_EPSG_MIN_CODE = 1024;
68 
RMFUnitTypeToStr(GUInt32 iElevationUnit)69 static char* RMFUnitTypeToStr( GUInt32 iElevationUnit )
70 {
71     switch ( iElevationUnit )
72     {
73         case 0:
74             return CPLStrdup( RMF_UnitsM );
75         case 1:
76             return CPLStrdup( RMF_UnitsDM );
77         case 2:
78             return CPLStrdup( RMF_UnitsCM );
79         case 3:
80             return CPLStrdup( RMF_UnitsMM );
81         default:
82             return CPLStrdup( RMF_UnitsEmpty );
83     }
84 }
85 
RMFStrToUnitType(const char * pszUnit,int * pbSuccess=nullptr)86 static GUInt32 RMFStrToUnitType( const char* pszUnit, int *pbSuccess = nullptr )
87 {
88     if( pbSuccess != nullptr )
89     {
90         *pbSuccess = TRUE;
91     }
92     if( EQUAL(pszUnit, RMF_UnitsM) )
93         return 0;
94     else if( EQUAL(pszUnit, RMF_UnitsDM) )
95         return 1;
96     else if( EQUAL(pszUnit, RMF_UnitsCM) )
97         return 2;
98     else if( EQUAL(pszUnit, RMF_UnitsMM) )
99         return 3;
100     else
101     {
102         //There is no 'invalid unit' in RMF format. So meter is default...
103         if( pbSuccess != nullptr )
104         {
105             *pbSuccess = FALSE;
106         }
107         return 0;
108     }
109 }
110 
111 /************************************************************************/
112 /* ==================================================================== */
113 /*                            RMFRasterBand                             */
114 /* ==================================================================== */
115 /************************************************************************/
116 
117 /************************************************************************/
118 /*                           RMFRasterBand()                            */
119 /************************************************************************/
120 
RMFRasterBand(RMFDataset * poDSIn,int nBandIn,GDALDataType eType)121 RMFRasterBand::RMFRasterBand( RMFDataset *poDSIn, int nBandIn,
122                               GDALDataType eType ) :
123     nBytesPerPixel(poDSIn->sHeader.nBitDepth / 8),
124     nLastTileWidth(poDSIn->GetRasterXSize() % poDSIn->sHeader.nTileWidth),
125     nLastTileHeight(poDSIn->GetRasterYSize() % poDSIn->sHeader.nTileHeight),
126     nDataSize(GDALGetDataTypeSizeBytes( eType ))
127 {
128     poDS = poDSIn;
129     nBand = nBandIn;
130 
131     eDataType = eType;
132     nBlockXSize = poDSIn->sHeader.nTileWidth;
133     nBlockYSize = poDSIn->sHeader.nTileHeight;
134     nBlockSize = nBlockXSize * nBlockYSize;
135     nBlockBytes = nBlockSize * nDataSize;
136 
137 #ifdef DEBUG
138     CPLDebug( "RMF",
139               "Band %d: tile width is %d, tile height is %d, "
140               " last tile width %u, last tile height %u, "
141               "bytes per pixel is %d, data type size is %d",
142               nBand, nBlockXSize, nBlockYSize,
143               nLastTileWidth, nLastTileHeight,
144               nBytesPerPixel, nDataSize );
145 #endif
146 }
147 
148 /************************************************************************/
149 /*                           ~RMFRasterBand()                           */
150 /************************************************************************/
151 
~RMFRasterBand()152 RMFRasterBand::~RMFRasterBand() {}
153 
154 /************************************************************************/
155 /*                             IReadBlock()                             */
156 /************************************************************************/
157 
IReadBlock(int nBlockXOff,int nBlockYOff,void * pImage)158 CPLErr RMFRasterBand::IReadBlock( int nBlockXOff, int nBlockYOff,
159                                   void * pImage )
160 {
161     RMFDataset  *poGDS = reinterpret_cast<RMFDataset *>( poDS );
162 
163     CPLAssert( poGDS != nullptr
164                && nBlockXOff >= 0
165                && nBlockYOff >= 0
166                && pImage != nullptr );
167 
168     memset(pImage, 0, nBlockBytes);
169 
170     GUInt32 nRawXSize = nBlockXSize;
171     GUInt32 nRawYSize = nBlockYSize;
172 
173     if( nLastTileWidth && (GUInt32)nBlockXOff == poGDS->nXTiles - 1 )
174         nRawXSize = nLastTileWidth;
175 
176     if( nLastTileHeight && (GUInt32)nBlockYOff == poGDS->nYTiles - 1 )
177         nRawYSize = nLastTileHeight;
178 
179     GUInt32 nRawBytes = nRawXSize * nRawYSize * poGDS->sHeader.nBitDepth / 8;
180 
181     //Direct read optimization
182     if(poGDS->nBands == 1 && poGDS->sHeader.nBitDepth >= 8 &&
183        nRawXSize == static_cast<GUInt32>(nBlockXSize) &&
184        nRawYSize == static_cast<GUInt32>(nBlockYSize))
185     {
186         bool bNullTile = false;
187         if(CE_None != poGDS->ReadTile(nBlockXOff, nBlockYOff,
188                                       reinterpret_cast<GByte*>(pImage),
189                                       nRawBytes, nRawXSize, nRawYSize,
190                                       bNullTile))
191         {
192             CPLError(CE_Failure, CPLE_FileIO,
193                      "Failed to read tile xOff %d yOff %d",
194                      nBlockXOff, nBlockYOff);
195             return CE_Failure;
196         }
197         if(bNullTile)
198         {
199             const int nChunkSize = std::max(1, GDALGetDataTypeSizeBytes(eDataType));
200             const GPtrDiff_t nWords = static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
201             GDALCopyWords64(&poGDS->sHeader.dfNoData, GDT_Float64, 0,
202                             pImage, eDataType, nChunkSize, nWords);
203         }
204         return CE_None;
205     }
206 #ifdef DEBUG
207     CPLDebug("RMF", "IReadBlock nBand %d, RawSize [%d, %d], Bits %d",
208              nBand, nRawXSize, nRawYSize, (int)poGDS->sHeader.nBitDepth);
209 #endif //DEBUG
210     if(poGDS->pabyCurrentTile == nullptr ||
211        poGDS->nCurrentTileXOff != nBlockXOff ||
212        poGDS->nCurrentTileYOff != nBlockYOff ||
213        poGDS->nCurrentTileBytes != nRawBytes)
214     {
215         if(poGDS->pabyCurrentTile == nullptr)
216         {
217             GUInt32 nMaxTileBytes = poGDS->sHeader.nTileWidth *
218                                     poGDS->sHeader.nTileHeight *
219                                     poGDS->sHeader.nBitDepth / 8;
220             poGDS->pabyCurrentTile =
221                reinterpret_cast<GByte*>(VSIMalloc(std::max(1U, nMaxTileBytes)));
222             if(!poGDS->pabyCurrentTile)
223             {
224                 CPLError(CE_Failure, CPLE_OutOfMemory,
225                          "Can't allocate tile block of size %lu.\n%s",
226                          static_cast<unsigned long>(nMaxTileBytes),
227                          VSIStrerror(errno));
228                 poGDS->nCurrentTileBytes = 0;
229                 return CE_Failure;
230             }
231         }
232 
233         poGDS->nCurrentTileXOff = nBlockXOff;
234         poGDS->nCurrentTileYOff = nBlockYOff;
235         poGDS->nCurrentTileBytes = nRawBytes;
236 
237         if(CE_None != poGDS->ReadTile(nBlockXOff, nBlockYOff,
238                                       poGDS->pabyCurrentTile, nRawBytes,
239                                       nRawXSize, nRawYSize,
240                                       poGDS->bCurrentTileIsNull))
241         {
242             CPLError(CE_Failure, CPLE_FileIO,
243                      "Failed to read tile xOff %d yOff %d",
244                      nBlockXOff, nBlockYOff);
245             poGDS->nCurrentTileBytes = 0;
246             return CE_Failure;
247         }
248     }
249 
250 /* -------------------------------------------------------------------- */
251 /*  Deinterleave pixels from input buffer.                              */
252 /* -------------------------------------------------------------------- */
253 
254     if(poGDS->bCurrentTileIsNull)
255     {
256         const int nChunkSize = std::max(1, GDALGetDataTypeSizeBytes(eDataType));
257         const GPtrDiff_t nWords = static_cast<GPtrDiff_t>(nBlockXSize) * nBlockYSize;
258         GDALCopyWords64(&poGDS->sHeader.dfNoData, GDT_Float64, 0,
259                         pImage, eDataType, nChunkSize, nWords);
260         return CE_None;
261     }
262     else if((poGDS->eRMFType == RMFT_RSW &&
263              (poGDS->sHeader.nBitDepth == 8 ||
264               poGDS->sHeader.nBitDepth == 24 ||
265               poGDS->sHeader.nBitDepth == 32)) ||
266             (poGDS->eRMFType == RMFT_MTW))
267     {
268         size_t  nTilePixelSize = poGDS->sHeader.nBitDepth / 8;
269         size_t  nTileLineSize = nTilePixelSize * nRawXSize;
270         size_t  nBlockLineSize = nDataSize * nBlockXSize;
271         int     iDstBand = (poGDS->nBands - nBand);
272         for(GUInt32 iLine = 0; iLine != nRawYSize; ++iLine)
273         {
274             GByte* pabySrc;
275             GByte* pabyDst;
276             pabySrc = poGDS->pabyCurrentTile +
277                       iLine * nTileLineSize +
278                       iDstBand * nDataSize;
279             pabyDst = reinterpret_cast<GByte*>(pImage) +
280                       iLine * nBlockLineSize;
281             GDALCopyWords(pabySrc, eDataType, static_cast<int>(nTilePixelSize),
282                           pabyDst, eDataType, static_cast<int>(nDataSize),
283                           nRawXSize);
284         }
285         return CE_None;
286     }
287     else if(poGDS->eRMFType == RMFT_RSW &&
288             poGDS->sHeader.nBitDepth == 16 &&
289             poGDS->nBands == 3)
290     {
291         size_t  nTilePixelBits = poGDS->sHeader.nBitDepth;
292         size_t  nTileLineSize = nTilePixelBits * nRawXSize / 8;
293         size_t  nBlockLineSize = nDataSize * nBlockXSize;
294 
295         for(GUInt32 iLine = 0; iLine != nRawYSize; ++iLine)
296         {
297             GUInt16* pabySrc;
298             GByte* pabyDst;
299             pabySrc = reinterpret_cast<GUInt16*>(poGDS->pabyCurrentTile +
300                                                  iLine * nTileLineSize);
301             pabyDst = reinterpret_cast<GByte*>(pImage) +
302                       iLine * nBlockLineSize;
303 
304             for( GUInt32 i = 0; i < nRawXSize; i++ )
305             {
306                 switch( nBand )
307                 {
308                     case 1:
309                         pabyDst[i] =
310                             static_cast<GByte>((pabySrc[i] & 0x7c00) >> 7);
311                         break;
312                     case 2:
313                         pabyDst[i] =
314                             static_cast<GByte>((pabySrc[i] & 0x03e0) >> 2);
315                         break;
316                     case 3:
317                         pabyDst[i] =
318                             static_cast<GByte>((pabySrc[i] & 0x1F) << 3);
319                         break;
320                     default:
321                         break;
322                 }
323             }
324         }
325         return CE_None;
326     }
327     else if(poGDS->eRMFType == RMFT_RSW &&
328             poGDS->nBands == 1 &&
329             poGDS->sHeader.nBitDepth == 4)
330     {
331         if( poGDS->nCurrentTileBytes != (nBlockSize+1) / 2 )
332         {
333             CPLError(CE_Failure, CPLE_AppDefined,
334                      "Tile has %d bytes, %d were expected",
335                      poGDS->nCurrentTileBytes, (nBlockSize+1) / 2 );
336             return CE_Failure;
337         }
338 
339         size_t  nTilePixelBits = poGDS->sHeader.nBitDepth;
340         size_t  nTileLineSize = nTilePixelBits * nRawXSize / 8;
341         size_t  nBlockLineSize = nDataSize * nBlockXSize;
342 
343         for(GUInt32 iLine = 0; iLine != nRawYSize; ++iLine)
344         {
345             GByte* pabySrc;
346             GByte* pabyDst;
347             pabySrc = poGDS->pabyCurrentTile +
348                       iLine * nTileLineSize;
349             pabyDst = reinterpret_cast<GByte*>(pImage) +
350                       iLine * nBlockLineSize;
351             for( GUInt32 i = 0; i < nRawXSize; ++i )
352             {
353                 if( i & 0x01 )
354                     pabyDst[i] = (*pabySrc++ & 0xF0) >> 4;
355                 else
356                     pabyDst[i] = *pabySrc & 0x0F;
357             }
358         }
359         return CE_None;
360     }
361     else if(poGDS->eRMFType == RMFT_RSW &&
362             poGDS->nBands == 1 &&
363             poGDS->sHeader.nBitDepth == 1)
364     {
365         if( poGDS->nCurrentTileBytes != (nBlockSize+7) / 8 )
366         {
367             CPLError(CE_Failure, CPLE_AppDefined,
368                      "Tile has %d bytes, %d were expected",
369                      poGDS->nCurrentTileBytes, (nBlockSize+7) / 8 );
370                 return CE_Failure;
371         }
372 
373         size_t  nTilePixelBits = poGDS->sHeader.nBitDepth;
374         size_t  nTileLineSize = nTilePixelBits * nRawXSize / 8;
375         size_t  nBlockLineSize = nDataSize * nBlockXSize;
376 
377         for(GUInt32 iLine = 0; iLine != nRawYSize; ++iLine)
378         {
379             GByte* pabySrc;
380             GByte* pabyDst;
381             pabySrc = poGDS->pabyCurrentTile +
382                       iLine * nTileLineSize;
383             pabyDst = reinterpret_cast<GByte*>(pImage) +
384                       iLine * nBlockLineSize;
385 
386             for(GUInt32 i = 0; i < nRawXSize; ++i)
387             {
388                 switch(i & 0x7)
389                 {
390                     case 0:
391                         pabyDst[i] = (*pabySrc & 0x80) >> 7;
392                         break;
393                     case 1:
394                         pabyDst[i] = (*pabySrc & 0x40) >> 6;
395                         break;
396                     case 2:
397                         pabyDst[i] = (*pabySrc & 0x20) >> 5;
398                         break;
399                     case 3:
400                         pabyDst[i] = (*pabySrc & 0x10) >> 4;
401                         break;
402                     case 4:
403                         pabyDst[i] = (*pabySrc & 0x08) >> 3;
404                         break;
405                     case 5:
406                         pabyDst[i] = (*pabySrc & 0x04) >> 2;
407                         break;
408                     case 6:
409                         pabyDst[i] = (*pabySrc & 0x02) >> 1;
410                         break;
411                     case 7:
412                         pabyDst[i] = *pabySrc++ & 0x01;
413                         break;
414                     default:
415                         break;
416                 }
417             }
418         }
419         return CE_None;
420     }
421 
422     CPLError(CE_Failure, CPLE_AppDefined,
423              "Invalid block data type. BitDepth %d, nBands %d",
424              static_cast<int>(poGDS->sHeader.nBitDepth), poGDS->nBands);
425 
426     return CE_Failure;
427 }
428 
429 /************************************************************************/
430 /*                            IWriteBlock()                             */
431 /************************************************************************/
432 
IWriteBlock(int nBlockXOff,int nBlockYOff,void * pImage)433 CPLErr RMFRasterBand::IWriteBlock( int nBlockXOff, int nBlockYOff,
434                                    void * pImage )
435 {
436     CPLAssert(poDS != nullptr
437               && nBlockXOff >= 0
438               && nBlockYOff >= 0
439               && pImage != nullptr);
440 
441     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>(poDS);
442 
443     //First drop current tile read by IReadBlock
444     poGDS->nCurrentTileBytes = 0;
445 
446     GUInt32 nRawXSize = nBlockXSize;
447     GUInt32 nRawYSize = nBlockYSize;
448 
449     if(nLastTileWidth && static_cast<GUInt32>(nBlockXOff) == poGDS->nXTiles - 1)
450        nRawXSize = nLastTileWidth;
451 
452     if(nLastTileHeight && static_cast<GUInt32>(nBlockYOff) == poGDS->nYTiles - 1)
453        nRawYSize = nLastTileHeight;
454 
455     size_t  nTilePixelSize = nDataSize * poGDS->nBands;
456     size_t  nTileLineSize = nTilePixelSize * nRawXSize;
457     size_t  nTileSize = nTileLineSize * nRawYSize;
458     size_t  nBlockLineSize = nDataSize * nBlockXSize;
459 
460 #ifdef DEBUG
461     CPLDebug("RMF", "IWriteBlock BlockSize [%d, %d], RawSize [%d, %d], size %d, nBand %d",
462              nBlockXSize, nBlockYSize, nRawXSize, nRawYSize,
463              static_cast<int>(nTileSize), nBand);
464 #endif // DEBUG
465 
466     if(poGDS->nBands == 1 &&
467        nRawXSize == static_cast<GUInt32>(nBlockXSize) &&
468        nRawYSize == static_cast<GUInt32>(nBlockYSize))
469     {//Immediate write
470         return poGDS->WriteTile(nBlockXOff, nBlockYOff,
471                                 reinterpret_cast<GByte*>(pImage),
472                                 nRawXSize * nRawYSize * nDataSize,
473                                 nRawXSize, nRawYSize);
474     }
475     else
476     {//Try to construct full tile in memory and write later
477         const GUInt32 nTile = nBlockYOff * poGDS->nXTiles + nBlockXOff;
478 
479         // Find tile
480         auto    poTile(poGDS->oUnfinishedTiles.find(nTile));
481         if(poTile == poGDS->oUnfinishedTiles.end())
482         {
483             RMFTileData  oTile;
484             oTile.oData.resize(nTileSize);
485             // If not found, but exist on disk than read it
486             if(poGDS->paiTiles[2 * nTile + 1])
487             {
488                 CPLErr eRes;
489                 bool bNullTile = false;
490                 eRes = poGDS->ReadTile(nBlockXOff, nBlockYOff, oTile.oData.data(),
491                                        nTileSize, nRawXSize, nRawYSize, bNullTile);
492                 if(eRes != CE_None)
493                 {
494                     CPLError(CE_Failure, CPLE_FileIO,
495                             "Can't read block with offset [%d, %d]",
496                             nBlockXOff, nBlockYOff);
497                     return eRes;
498                 }
499             }
500             poTile = poGDS->oUnfinishedTiles.insert(poGDS->oUnfinishedTiles.end(),
501                                                 std::make_pair(nTile, oTile));
502         }
503 
504         GByte*  pabyTileData = poTile->second.oData.data();
505 
506         // Copy new data to a tile
507         int iDstBand = (poGDS->nBands - nBand);
508         for(GUInt32 iLine = 0; iLine != nRawYSize; ++iLine)
509         {
510             const GByte* pabySrc;
511             GByte* pabyDst;
512             pabySrc = reinterpret_cast<const GByte*>(pImage) +
513                       iLine * nBlockLineSize;
514             pabyDst = pabyTileData +
515                       iLine * nTileLineSize +
516                       iDstBand * nDataSize;
517             GDALCopyWords(pabySrc, eDataType, static_cast<int>(nDataSize),
518                           pabyDst, eDataType, static_cast<int>(nTilePixelSize),
519                           nRawXSize);
520         }
521         ++poTile->second.nBandsWritten;
522 
523         // Write to disk if tile is finished
524         if(poTile->second.nBandsWritten == poGDS->nBands)
525         {
526             poGDS->WriteTile(nBlockXOff, nBlockYOff,
527                              pabyTileData, nTileSize,
528                              nRawXSize, nRawYSize);
529             poGDS->oUnfinishedTiles.erase(poTile);
530         }
531 #ifdef DEBUG
532         CPLDebug("RMF", "poGDS->oUnfinishedTiles.size() %d",
533                  static_cast<int>(poGDS->oUnfinishedTiles.size()));
534 #endif //DEBUG
535     }
536 
537     return CE_None;
538 }
539 
540 /************************************************************************/
541 /*                          GetNoDataValue()                            */
542 /************************************************************************/
543 
GetNoDataValue(int * pbSuccess)544 double RMFRasterBand::GetNoDataValue( int *pbSuccess )
545 
546 {
547     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
548 
549     if( pbSuccess )
550         *pbSuccess = TRUE;
551 
552     return poGDS->sHeader.dfNoData;
553 }
554 
SetNoDataValue(double dfNoData)555 CPLErr RMFRasterBand::SetNoDataValue( double dfNoData )
556 {
557     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
558 
559     poGDS->sHeader.dfNoData = dfNoData;
560     poGDS->bHeaderDirty = true;
561 
562     return CE_None;
563 }
564 
565 /************************************************************************/
566 /*                            GetUnitType()                             */
567 /************************************************************************/
568 
GetUnitType()569 const char *RMFRasterBand::GetUnitType()
570 
571 {
572     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
573 
574     return (const char *)poGDS->pszUnitType;
575 }
576 
577 /************************************************************************/
578 /*                            SetUnitType()                             */
579 /************************************************************************/
580 
SetUnitType(const char * pszNewValue)581 CPLErr RMFRasterBand::SetUnitType( const char *pszNewValue )
582 
583 {
584     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
585     int         bSuccess = FALSE;
586     int         iNewUnit = RMFStrToUnitType(pszNewValue, &bSuccess);
587 
588     if( bSuccess )
589     {
590         CPLFree(poGDS->pszUnitType);
591         poGDS->pszUnitType = CPLStrdup( pszNewValue );
592         poGDS->sHeader.iElevationUnit = iNewUnit;
593         poGDS->bHeaderDirty = true;
594         return CE_None;
595     }
596     else
597     {
598         CPLError( CE_Warning, CPLE_NotSupported,
599                   "RMF driver does not support '%s' elevation units. "
600                   "Possible values are: m, dm, cm, mm.",
601                   pszNewValue );
602         return CE_Failure;
603     }
604 }
605 
606 /************************************************************************/
607 /*                           GetColorTable()                            */
608 /************************************************************************/
609 
GetColorTable()610 GDALColorTable *RMFRasterBand::GetColorTable()
611 {
612     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
613 
614     return poGDS->poColorTable;
615 }
616 
617 /************************************************************************/
618 /*                           SetColorTable()                            */
619 /************************************************************************/
620 
SetColorTable(GDALColorTable * poColorTable)621 CPLErr RMFRasterBand::SetColorTable( GDALColorTable *poColorTable )
622 {
623     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
624 
625     if( poColorTable )
626     {
627         if( poGDS->eRMFType == RMFT_RSW && poGDS->nBands == 1 )
628         {
629             if( !poGDS->pabyColorTable )
630                 return CE_Failure;
631 
632             GDALColorEntry  oEntry;
633             for( GUInt32 i = 0; i < poGDS->nColorTableSize; i++ )
634             {
635                 poColorTable->GetColorEntryAsRGB( i, &oEntry );
636                 poGDS->pabyColorTable[i * 4] = (GByte) oEntry.c1;     // Red
637                 poGDS->pabyColorTable[i * 4 + 1] = (GByte) oEntry.c2; // Green
638                 poGDS->pabyColorTable[i * 4 + 2] = (GByte) oEntry.c3; // Blue
639                 poGDS->pabyColorTable[i * 4 + 3] = 0;
640             }
641 
642             poGDS->bHeaderDirty = true;
643         }
644         return CE_None;
645     }
646 
647     return CE_Failure;
648 }
649 
GetOverviewCount()650 int RMFRasterBand::GetOverviewCount()
651 {
652     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
653     if( poGDS->poOvrDatasets.empty() )
654         return GDALRasterBand::GetOverviewCount();
655     else
656         return static_cast<int>( poGDS->poOvrDatasets.size() );
657 }
658 
GetOverview(int i)659 GDALRasterBand* RMFRasterBand::GetOverview(int i)
660 {
661     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
662     size_t      n = static_cast<size_t>( i );
663     if( poGDS->poOvrDatasets.empty() )
664         return GDALRasterBand::GetOverview(i);
665     else
666         return poGDS->poOvrDatasets[n]->GetRasterBand(nBand);
667 }
668 
IRasterIO(GDALRWFlag eRWFlag,int nXOff,int nYOff,int nXSize,int nYSize,void * pData,int nBufXSize,int nBufYSize,GDALDataType eType,GSpacing nPixelSpace,GSpacing nLineSpace,GDALRasterIOExtraArg * psExtraArg)669 CPLErr RMFRasterBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
670                                int nXSize, int nYSize, void* pData,
671                                int nBufXSize, int nBufYSize, GDALDataType eType,
672                                GSpacing nPixelSpace, GSpacing nLineSpace,
673                                GDALRasterIOExtraArg* psExtraArg)
674 {
675     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>(poDS);
676 
677     if(eRWFlag == GF_Read &&
678        poGDS->poCompressData != nullptr &&
679        poGDS->poCompressData->oThreadPool.GetThreadCount() > 0)
680     {
681         poGDS->poCompressData->oThreadPool.WaitCompletion();
682     }
683 
684     return GDALRasterBand::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
685                                      pData, nBufXSize, nBufYSize, eType,
686                                      nPixelSpace, nLineSpace, psExtraArg);
687 }
688 
689 /************************************************************************/
690 /*                       GetColorInterpretation()                       */
691 /************************************************************************/
692 
GetColorInterpretation()693 GDALColorInterp RMFRasterBand::GetColorInterpretation()
694 {
695     RMFDataset *poGDS = reinterpret_cast<RMFDataset *>( poDS );
696 
697     if( poGDS->nBands == 3 )
698     {
699         if( nBand == 1 )
700             return GCI_RedBand;
701         else if( nBand == 2 )
702             return GCI_GreenBand;
703         else if( nBand == 3 )
704             return GCI_BlueBand;
705 
706         return GCI_Undefined;
707     }
708 
709     if( poGDS->eRMFType == RMFT_RSW )
710         return GCI_PaletteIndex;
711 
712     return GCI_Undefined;
713 }
714 
715 /************************************************************************/
716 /* ==================================================================== */
717 /*                              RMFDataset                              */
718 /* ==================================================================== */
719 /************************************************************************/
720 
721 /************************************************************************/
722 /*                           RMFDataset()                               */
723 /************************************************************************/
724 
RMFDataset()725 RMFDataset::RMFDataset() :
726     eRMFType(RMFT_RSW),
727     nXTiles(0),
728     nYTiles(0),
729     paiTiles(nullptr),
730     pabyDecompressBuffer(nullptr),
731     pabyCurrentTile(nullptr),
732     bCurrentTileIsNull(false),
733     nCurrentTileXOff(-1),
734     nCurrentTileYOff(-1),
735     nCurrentTileBytes(0),
736     nColorTableSize(0),
737     pabyColorTable(nullptr),
738     poColorTable(nullptr),
739     pszProjection(CPLStrdup( "" )),
740     pszUnitType(CPLStrdup( RMF_UnitsEmpty )),
741     bBigEndian(false),
742     bHeaderDirty(false),
743     fp(nullptr),
744     Decompress(nullptr),
745     Compress(nullptr),
746     nHeaderOffset(0),
747     poParentDS(nullptr)
748 {
749     nBands = 0;
750     adfGeoTransform[0] = 0.0;
751     adfGeoTransform[1] = 1.0;
752     adfGeoTransform[2] = 0.0;
753     adfGeoTransform[3] = 0.0;
754     adfGeoTransform[4] = 0.0;
755     adfGeoTransform[5] = 1.0;
756     memset( &sHeader, 0, sizeof(sHeader) );
757     memset( &sExtHeader, 0, sizeof(sExtHeader) );
758 }
759 
760 /************************************************************************/
761 /*                            ~RMFDataset()                             */
762 /************************************************************************/
763 
~RMFDataset()764 RMFDataset::~RMFDataset()
765 {
766     RMFDataset::FlushCache();
767     for( size_t n = 0; n != poOvrDatasets.size(); ++n )
768     {
769         poOvrDatasets[n]->RMFDataset::FlushCache();
770     }
771 
772     VSIFree( paiTiles );
773     VSIFree( pabyDecompressBuffer );
774     VSIFree( pabyCurrentTile );
775     CPLFree( pszProjection );
776     CPLFree( pszUnitType );
777     CPLFree( pabyColorTable );
778     if( poColorTable != nullptr )
779         delete poColorTable;
780 
781     for( size_t n = 0; n != poOvrDatasets.size(); ++n )
782     {
783         GDALClose( poOvrDatasets[n] );
784     }
785 
786     if( fp != nullptr && poParentDS == nullptr )
787     {
788         VSIFCloseL( fp );
789     }
790 }
791 
792 /************************************************************************/
793 /*                          GetGeoTransform()                           */
794 /************************************************************************/
795 
GetGeoTransform(double * padfTransform)796 CPLErr RMFDataset::GetGeoTransform( double * padfTransform )
797 {
798     memcpy( padfTransform, adfGeoTransform, sizeof(adfGeoTransform[0]) * 6 );
799 
800     if( sHeader.iGeorefFlag )
801         return CE_None;
802 
803     return CE_Failure;
804 }
805 
806 /************************************************************************/
807 /*                          SetGeoTransform()                           */
808 /************************************************************************/
809 
SetGeoTransform(double * padfTransform)810 CPLErr RMFDataset::SetGeoTransform( double * padfTransform )
811 {
812     memcpy( adfGeoTransform, padfTransform, sizeof(double) * 6 );
813     sHeader.dfPixelSize = adfGeoTransform[1];
814     if( sHeader.dfPixelSize != 0.0 )
815         sHeader.dfResolution = sHeader.dfScale / sHeader.dfPixelSize;
816     sHeader.dfLLX = adfGeoTransform[0];
817     sHeader.dfLLY = adfGeoTransform[3] - nRasterYSize * sHeader.dfPixelSize;
818     sHeader.iGeorefFlag = 1;
819 
820     bHeaderDirty = true;
821 
822     return CE_None;
823 }
824 
825 /************************************************************************/
826 /*                          GetProjectionRef()                          */
827 /************************************************************************/
828 
_GetProjectionRef()829 const char *RMFDataset::_GetProjectionRef()
830 {
831     if( pszProjection )
832         return pszProjection;
833 
834     return "";
835 }
836 
837 /************************************************************************/
838 /*                           SetProjection()                            */
839 /************************************************************************/
840 
_SetProjection(const char * pszNewProjection)841 CPLErr RMFDataset::_SetProjection( const char * pszNewProjection )
842 
843 {
844     CPLFree( pszProjection );
845     pszProjection = CPLStrdup( (pszNewProjection) ? pszNewProjection : "" );
846 
847     bHeaderDirty = true;
848 
849     return CE_None;
850 }
851 
852 /************************************************************************/
853 /*                           WriteHeader()                              */
854 /************************************************************************/
855 
WriteHeader()856 CPLErr RMFDataset::WriteHeader()
857 {
858 /* -------------------------------------------------------------------- */
859 /*  Setup projection.                                                   */
860 /* -------------------------------------------------------------------- */
861     if( pszProjection && !EQUAL( pszProjection, "" ) )
862     {
863         OGRSpatialReference oSRS;
864         if( oSRS.importFromWkt( pszProjection ) == OGRERR_NONE )
865         {
866             long iProjection = 0;
867             long iDatum = 0;
868             long iEllips = 0;
869             long iZone = 0;
870             double adfPrjParams[7] = {};
871 
872             oSRS.exportToPanorama( &iProjection, &iDatum, &iEllips, &iZone,
873                                    adfPrjParams );
874             sHeader.iProjection = static_cast<int>(iProjection);
875             sHeader.dfStdP1 = adfPrjParams[0];
876             sHeader.dfStdP2 = adfPrjParams[1];
877             sHeader.dfCenterLat = adfPrjParams[2];
878             sHeader.dfCenterLong = adfPrjParams[3];
879             if( oSRS.GetAuthorityName(nullptr) != nullptr &&
880                 oSRS.GetAuthorityCode(nullptr) != nullptr &&
881                 EQUAL(oSRS.GetAuthorityName(nullptr), "EPSG") )
882             {
883                 sHeader.iEPSGCode = atoi(oSRS.GetAuthorityCode(nullptr));
884             }
885 
886             sExtHeader.nEllipsoid = static_cast<int>(iEllips);
887             sExtHeader.nDatum = static_cast<int>(iDatum);
888             sExtHeader.nZone = static_cast<int>(iZone);
889         }
890     }
891 
892 #define RMF_WRITE_LONG( ptr, value, offset )            \
893 do {                                                    \
894     GInt32  iLong = CPL_LSBWORD32( value );             \
895     memcpy( (ptr) + (offset), &iLong, 4 );              \
896 } while( false );
897 
898 #define RMF_WRITE_ULONG( ptr,value, offset )            \
899 do {                                                    \
900     GUInt32 iULong = CPL_LSBWORD32( value );            \
901     memcpy( (ptr) + (offset), &iULong, 4 );             \
902 } while( false );
903 
904 #define RMF_WRITE_DOUBLE( ptr,value, offset )           \
905 do {                                                    \
906     double  dfDouble = (value);                         \
907     CPL_LSBPTR64( &dfDouble );                          \
908     memcpy( (ptr) + (offset), &dfDouble, 8 );           \
909 } while( false );
910 
911     vsi_l_offset    iCurrentFileSize( GetLastOffset() );
912     sHeader.nFileSize0 = GetRMFOffset( iCurrentFileSize, &iCurrentFileSize );
913     sHeader.nSize = sHeader.nFileSize0 - GetRMFOffset( nHeaderOffset, nullptr );
914 /* -------------------------------------------------------------------- */
915 /*  Write out the main header.                                          */
916 /* -------------------------------------------------------------------- */
917     {
918         GByte abyHeader[RMF_HEADER_SIZE] = {};
919 
920         memcpy( abyHeader, sHeader.bySignature, RMF_SIGNATURE_SIZE );
921         RMF_WRITE_ULONG( abyHeader, sHeader.iVersion, 4 );
922         RMF_WRITE_ULONG( abyHeader, sHeader.nSize, 8 );
923         RMF_WRITE_ULONG( abyHeader, sHeader.nOvrOffset, 12 );
924         RMF_WRITE_ULONG( abyHeader, sHeader.iUserID, 16 );
925         memcpy( abyHeader + 20, sHeader.byName, RMF_NAME_SIZE );
926         RMF_WRITE_ULONG( abyHeader, sHeader.nBitDepth, 52 );
927         RMF_WRITE_ULONG( abyHeader, sHeader.nHeight, 56 );
928         RMF_WRITE_ULONG( abyHeader, sHeader.nWidth, 60 );
929         RMF_WRITE_ULONG( abyHeader, sHeader.nXTiles, 64 );
930         RMF_WRITE_ULONG( abyHeader, sHeader.nYTiles, 68 );
931         RMF_WRITE_ULONG( abyHeader, sHeader.nTileHeight, 72 );
932         RMF_WRITE_ULONG( abyHeader, sHeader.nTileWidth, 76 );
933         RMF_WRITE_ULONG( abyHeader, sHeader.nLastTileHeight, 80 );
934         RMF_WRITE_ULONG( abyHeader, sHeader.nLastTileWidth, 84 );
935         RMF_WRITE_ULONG( abyHeader, sHeader.nROIOffset, 88 );
936         RMF_WRITE_ULONG( abyHeader, sHeader.nROISize, 92 );
937         RMF_WRITE_ULONG( abyHeader, sHeader.nClrTblOffset, 96 );
938         RMF_WRITE_ULONG( abyHeader, sHeader.nClrTblSize, 100 );
939         RMF_WRITE_ULONG( abyHeader, sHeader.nTileTblOffset, 104 );
940         RMF_WRITE_ULONG( abyHeader, sHeader.nTileTblSize, 108 );
941         RMF_WRITE_LONG( abyHeader, sHeader.iMapType, 124 );
942         RMF_WRITE_LONG( abyHeader, sHeader.iProjection, 128 );
943         RMF_WRITE_LONG( abyHeader, sHeader.iEPSGCode, 132 );
944         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfScale, 136 );
945         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfResolution, 144 );
946         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfPixelSize, 152 );
947         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfLLY, 160 );
948         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfLLX, 168 );
949         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfStdP1, 176 );
950         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfStdP2, 184 );
951         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfCenterLong, 192 );
952         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfCenterLat, 200 );
953         *(abyHeader + 208) = sHeader.iCompression;
954         *(abyHeader + 209) = sHeader.iMaskType;
955         *(abyHeader + 210) = sHeader.iMaskStep;
956         *(abyHeader + 211) = sHeader.iFrameFlag;
957         RMF_WRITE_ULONG( abyHeader, sHeader.nFlagsTblOffset, 212 );
958         RMF_WRITE_ULONG( abyHeader, sHeader.nFlagsTblSize, 216 );
959         RMF_WRITE_ULONG( abyHeader, sHeader.nFileSize0, 220 );
960         RMF_WRITE_ULONG( abyHeader, sHeader.nFileSize1, 224 );
961         *(abyHeader + 228) = sHeader.iUnknown;
962         *(abyHeader + 244) = sHeader.iGeorefFlag;
963         *(abyHeader + 245) = sHeader.iInverse;
964         *(abyHeader + 246) = sHeader.iJpegQuality;
965         memcpy( abyHeader + 248, sHeader.abyInvisibleColors,
966                 sizeof(sHeader.abyInvisibleColors) );
967         RMF_WRITE_DOUBLE( abyHeader, sHeader.adfElevMinMax[0], 280 );
968         RMF_WRITE_DOUBLE( abyHeader, sHeader.adfElevMinMax[1], 288 );
969         RMF_WRITE_DOUBLE( abyHeader, sHeader.dfNoData, 296 );
970         RMF_WRITE_ULONG( abyHeader, sHeader.iElevationUnit, 304 );
971         *(abyHeader + 308) = sHeader.iElevationType;
972         RMF_WRITE_ULONG( abyHeader, sHeader.nExtHdrOffset, 312 );
973         RMF_WRITE_ULONG( abyHeader, sHeader.nExtHdrSize, 316 );
974 
975         VSIFSeekL( fp, nHeaderOffset, SEEK_SET );
976         VSIFWriteL( abyHeader, 1, sizeof(abyHeader), fp );
977     }
978 
979 /* -------------------------------------------------------------------- */
980 /*  Write out the extended header.                                      */
981 /* -------------------------------------------------------------------- */
982 
983     if( sHeader.nExtHdrOffset && sHeader.nExtHdrSize )
984     {
985         GByte *pabyExtHeader = reinterpret_cast<GByte *>(
986             CPLCalloc( sHeader.nExtHdrSize, 1 ) );
987 
988         RMF_WRITE_LONG( pabyExtHeader, sExtHeader.nEllipsoid, 24 );
989         RMF_WRITE_LONG( pabyExtHeader, sExtHeader.nVertDatum, 28 );
990         RMF_WRITE_LONG( pabyExtHeader, sExtHeader.nDatum, 32 );
991         RMF_WRITE_LONG( pabyExtHeader, sExtHeader.nZone, 36 );
992 
993         VSIFSeekL( fp, GetFileOffset( sHeader.nExtHdrOffset ), SEEK_SET );
994         VSIFWriteL( pabyExtHeader, 1, sHeader.nExtHdrSize, fp );
995 
996         CPLFree( pabyExtHeader );
997     }
998 
999 #undef RMF_WRITE_DOUBLE
1000 #undef RMF_WRITE_ULONG
1001 #undef RMF_WRITE_LONG
1002 
1003 /* -------------------------------------------------------------------- */
1004 /*  Write out the color table.                                          */
1005 /* -------------------------------------------------------------------- */
1006 
1007     if( sHeader.nClrTblOffset && sHeader.nClrTblSize )
1008     {
1009         VSIFSeekL( fp, GetFileOffset( sHeader.nClrTblOffset ), SEEK_SET );
1010         VSIFWriteL( pabyColorTable, 1, sHeader.nClrTblSize, fp );
1011     }
1012 
1013 /* -------------------------------------------------------------------- */
1014 /*  Write out the block table, swap if needed.                          */
1015 /* -------------------------------------------------------------------- */
1016 
1017     VSIFSeekL( fp, GetFileOffset( sHeader.nTileTblOffset ), SEEK_SET );
1018 
1019 #ifdef CPL_MSB
1020     GUInt32 *paiTilesSwapped = reinterpret_cast<GUInt32 *>(
1021         CPLMalloc( sHeader.nTileTblSize ) );
1022     if( !paiTilesSwapped )
1023         return CE_Failure;
1024 
1025     memcpy( paiTilesSwapped, paiTiles, sHeader.nTileTblSize );
1026     for( GUInt32 i = 0; i < sHeader.nTileTblSize / sizeof(GUInt32); i++ )
1027         CPL_SWAP32PTR( paiTilesSwapped + i );
1028     VSIFWriteL( paiTilesSwapped, 1, sHeader.nTileTblSize, fp );
1029 
1030     CPLFree( paiTilesSwapped );
1031 #else
1032     VSIFWriteL( paiTiles, 1, sHeader.nTileTblSize, fp );
1033 #endif
1034 
1035     bHeaderDirty = false;
1036 
1037     return CE_None;
1038 }
1039 
1040 /************************************************************************/
1041 /*                             FlushCache()                             */
1042 /************************************************************************/
1043 
FlushCache()1044 void RMFDataset::FlushCache()
1045 
1046 {
1047     GDALDataset::FlushCache();
1048 
1049     if(poCompressData != nullptr &&
1050        poCompressData->oThreadPool.GetThreadCount() > 0)
1051     {
1052         poCompressData->oThreadPool.WaitCompletion();
1053     }
1054 
1055     if( !bHeaderDirty )
1056         return;
1057 
1058     if( eRMFType == RMFT_MTW )
1059     {
1060         GDALRasterBand *poBand = GetRasterBand(1);
1061 
1062         if( poBand )
1063         {
1064             poBand->ComputeRasterMinMax( FALSE, sHeader.adfElevMinMax );
1065             bHeaderDirty = true;
1066         }
1067     }
1068     WriteHeader();
1069 }
1070 
1071 /************************************************************************/
1072 /*                              Identify()                              */
1073 /************************************************************************/
1074 
Identify(GDALOpenInfo * poOpenInfo)1075 int RMFDataset::Identify( GDALOpenInfo *poOpenInfo )
1076 
1077 {
1078     if( poOpenInfo->pabyHeader == nullptr)
1079         return FALSE;
1080 
1081     if( memcmp(poOpenInfo->pabyHeader, RMF_SigRSW, sizeof(RMF_SigRSW)) != 0
1082         && memcmp(poOpenInfo->pabyHeader, RMF_SigRSW_BE,
1083                   sizeof(RMF_SigRSW_BE)) != 0
1084         && memcmp(poOpenInfo->pabyHeader, RMF_SigMTW, sizeof(RMF_SigMTW)) != 0 )
1085         return FALSE;
1086 
1087     return TRUE;
1088 }
1089 
1090 /************************************************************************/
1091 /*                                Open()                                */
1092 /************************************************************************/
1093 
Open(GDALOpenInfo * poOpenInfo)1094 GDALDataset *RMFDataset::Open( GDALOpenInfo * poOpenInfo )
1095 {
1096     auto poDS = Open( poOpenInfo, nullptr, 0 );
1097     if( poDS == nullptr )
1098     {
1099         return nullptr;
1100     }
1101 
1102     RMFDataset* poCurrentLayer = poDS;
1103     RMFDataset* poParent = poCurrentLayer;
1104     const int   nMaxPossibleOvCount = 64;
1105 
1106     for( int iOv = 0; iOv < nMaxPossibleOvCount && poCurrentLayer != nullptr; ++iOv )
1107     {
1108         poCurrentLayer = poCurrentLayer->OpenOverview( poParent, poOpenInfo );
1109         if( poCurrentLayer == nullptr )
1110             break;
1111         poParent->poOvrDatasets.push_back( poCurrentLayer );
1112     }
1113 
1114     return poDS;
1115 }
1116 
Open(GDALOpenInfo * poOpenInfo,RMFDataset * poParentDS,vsi_l_offset nNextHeaderOffset)1117 RMFDataset *RMFDataset::Open(GDALOpenInfo * poOpenInfo,
1118                               RMFDataset* poParentDS,
1119                               vsi_l_offset nNextHeaderOffset )
1120 {
1121     if( !Identify(poOpenInfo) ||
1122         (poParentDS == nullptr && poOpenInfo->fpL == nullptr) )
1123         return nullptr;
1124 
1125 /* -------------------------------------------------------------------- */
1126 /*  Create a corresponding GDALDataset.                                 */
1127 /* -------------------------------------------------------------------- */
1128     RMFDataset *poDS = new RMFDataset();
1129 
1130     if( poParentDS == nullptr )
1131     {
1132         poDS->fp = poOpenInfo->fpL;
1133         poOpenInfo->fpL = nullptr;
1134         poDS->nHeaderOffset = 0;
1135         poDS->poParentDS = nullptr;
1136     }
1137     else
1138     {
1139         poDS->fp = poParentDS->fp;
1140         poDS->poParentDS = poParentDS;
1141         poDS->nHeaderOffset = nNextHeaderOffset;
1142     }
1143     poDS->eAccess = poOpenInfo->eAccess;
1144 
1145 #define RMF_READ_SHORT(ptr, value, offset)                              \
1146 do {                                                                    \
1147     memcpy(&(value), (GInt16*)((ptr) + (offset)), sizeof(GInt16));      \
1148     if( poDS->bBigEndian )                                              \
1149     {                                                                   \
1150         CPL_MSBPTR16(&(value));                                         \
1151     }                                                                   \
1152     else                                                                \
1153     {                                                                   \
1154         CPL_LSBPTR16(&(value));                                         \
1155     }                                                                   \
1156 } while( false );
1157 
1158 #define RMF_READ_ULONG(ptr, value, offset)                              \
1159 do {                                                                    \
1160     memcpy(&(value), (GUInt32*)((ptr) + (offset)), sizeof(GUInt32));    \
1161     if( poDS->bBigEndian )                                              \
1162     {                                                                   \
1163         CPL_MSBPTR32(&(value));                                         \
1164     }                                                                   \
1165     else                                                                \
1166     {                                                                   \
1167         CPL_LSBPTR32(&(value));                                         \
1168     }                                                                   \
1169 } while( false );
1170 
1171 #define RMF_READ_LONG(ptr, value, offset) RMF_READ_ULONG(ptr, value, offset)
1172 
1173 #define RMF_READ_DOUBLE(ptr, value, offset)                             \
1174 do {                                                                    \
1175     memcpy(&(value), (double*)((ptr) + (offset)), sizeof(double));      \
1176     if( poDS->bBigEndian )                                              \
1177     {                                                                   \
1178         CPL_MSBPTR64(&(value));                                         \
1179     }                                                                   \
1180     else                                                                \
1181     {                                                                   \
1182         CPL_LSBPTR64(&(value));                                         \
1183     }                                                                   \
1184 } while( false );
1185 
1186 /* -------------------------------------------------------------------- */
1187 /*  Read the main header.                                               */
1188 /* -------------------------------------------------------------------- */
1189 
1190     {
1191         GByte abyHeader[RMF_HEADER_SIZE] = {};
1192 
1193         VSIFSeekL( poDS->fp, nNextHeaderOffset, SEEK_SET );
1194         if( VSIFReadL( abyHeader, 1, sizeof(abyHeader),
1195                        poDS->fp ) != sizeof(abyHeader) )
1196         {
1197             delete poDS;
1198             return nullptr;
1199         }
1200 
1201         if( memcmp(abyHeader, RMF_SigMTW, sizeof(RMF_SigMTW)) == 0 )
1202         {
1203             poDS->eRMFType = RMFT_MTW;
1204         }
1205         else if( memcmp(abyHeader, RMF_SigRSW_BE, sizeof(RMF_SigRSW_BE)) == 0 )
1206         {
1207             poDS->eRMFType = RMFT_RSW;
1208             poDS->bBigEndian = true;
1209         }
1210         else
1211         {
1212             poDS->eRMFType = RMFT_RSW;
1213         }
1214 
1215         memcpy( poDS->sHeader.bySignature, abyHeader, RMF_SIGNATURE_SIZE );
1216         RMF_READ_ULONG( abyHeader, poDS->sHeader.iVersion, 4 );
1217         RMF_READ_ULONG( abyHeader, poDS->sHeader.nSize, 8 );
1218         RMF_READ_ULONG( abyHeader, poDS->sHeader.nOvrOffset, 12 );
1219         RMF_READ_ULONG( abyHeader, poDS->sHeader.iUserID, 16 );
1220         memcpy( poDS->sHeader.byName, abyHeader + 20,
1221                 sizeof(poDS->sHeader.byName) );
1222         poDS->sHeader.byName[sizeof(poDS->sHeader.byName) - 1] = '\0';
1223         RMF_READ_ULONG( abyHeader, poDS->sHeader.nBitDepth, 52 );
1224         RMF_READ_ULONG( abyHeader, poDS->sHeader.nHeight, 56 );
1225         RMF_READ_ULONG( abyHeader, poDS->sHeader.nWidth, 60 );
1226         RMF_READ_ULONG( abyHeader, poDS->sHeader.nXTiles, 64 );
1227         RMF_READ_ULONG( abyHeader, poDS->sHeader.nYTiles, 68 );
1228         RMF_READ_ULONG( abyHeader, poDS->sHeader.nTileHeight, 72 );
1229         RMF_READ_ULONG( abyHeader, poDS->sHeader.nTileWidth, 76 );
1230         RMF_READ_ULONG( abyHeader, poDS->sHeader.nLastTileHeight, 80 );
1231         RMF_READ_ULONG( abyHeader, poDS->sHeader.nLastTileWidth, 84 );
1232         RMF_READ_ULONG( abyHeader, poDS->sHeader.nROIOffset, 88 );
1233         RMF_READ_ULONG( abyHeader, poDS->sHeader.nROISize, 92 );
1234         RMF_READ_ULONG( abyHeader, poDS->sHeader.nClrTblOffset, 96 );
1235         RMF_READ_ULONG( abyHeader, poDS->sHeader.nClrTblSize, 100 );
1236         RMF_READ_ULONG( abyHeader, poDS->sHeader.nTileTblOffset, 104 );
1237         RMF_READ_ULONG( abyHeader, poDS->sHeader.nTileTblSize, 108 );
1238         RMF_READ_LONG( abyHeader, poDS->sHeader.iMapType, 124 );
1239         RMF_READ_LONG( abyHeader, poDS->sHeader.iProjection, 128 );
1240         RMF_READ_LONG( abyHeader, poDS->sHeader.iEPSGCode, 132 );
1241         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfScale, 136 );
1242         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfResolution, 144 );
1243         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfPixelSize, 152 );
1244         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfLLY, 160 );
1245         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfLLX, 168 );
1246         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfStdP1, 176 );
1247         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfStdP2, 184 );
1248         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfCenterLong, 192 );
1249         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.dfCenterLat, 200 );
1250         poDS->sHeader.iCompression = *(abyHeader + 208);
1251         poDS->sHeader.iMaskType = *(abyHeader + 209);
1252         poDS->sHeader.iMaskStep = *(abyHeader + 210);
1253         poDS->sHeader.iFrameFlag = *(abyHeader + 211);
1254         RMF_READ_ULONG( abyHeader, poDS->sHeader.nFlagsTblOffset, 212 );
1255         RMF_READ_ULONG( abyHeader, poDS->sHeader.nFlagsTblSize, 216 );
1256         RMF_READ_ULONG( abyHeader, poDS->sHeader.nFileSize0, 220 );
1257         RMF_READ_ULONG( abyHeader, poDS->sHeader.nFileSize1, 224 );
1258         poDS->sHeader.iUnknown = *(abyHeader + 228);
1259         poDS->sHeader.iGeorefFlag = *(abyHeader + 244);
1260         poDS->sHeader.iInverse = *(abyHeader + 245);
1261         poDS->sHeader.iJpegQuality = *(abyHeader + 246);
1262         memcpy( poDS->sHeader.abyInvisibleColors,
1263                 abyHeader + 248, sizeof(poDS->sHeader.abyInvisibleColors) );
1264         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.adfElevMinMax[0], 280 );
1265         RMF_READ_DOUBLE( abyHeader, poDS->sHeader.adfElevMinMax[1], 288 );
1266         RMF_READ_DOUBLE(abyHeader, poDS->sHeader.dfNoData, 296);
1267 
1268         RMF_READ_ULONG( abyHeader, poDS->sHeader.iElevationUnit, 304 );
1269         poDS->sHeader.iElevationType = *(abyHeader + 308);
1270         RMF_READ_ULONG( abyHeader, poDS->sHeader.nExtHdrOffset, 312 );
1271         RMF_READ_ULONG( abyHeader, poDS->sHeader.nExtHdrSize, 316 );
1272     }
1273 
1274     if(poDS->sHeader.nTileTblSize % (sizeof(GUInt32)*2))
1275     {
1276         CPLError( CE_Warning, CPLE_IllegalArg,
1277                   "Invalid tile table size." );
1278         delete poDS;
1279         return nullptr;
1280     }
1281 
1282     bool bInvalidTileSize;
1283     try
1284     {
1285         GUInt64 nMaxTileBits =
1286             (CPLSM(static_cast<GUInt64>(2)) *
1287              CPLSM(static_cast<GUInt64>(poDS->sHeader.nTileWidth)) *
1288              CPLSM(static_cast<GUInt64>(poDS->sHeader.nTileHeight)) *
1289              CPLSM(static_cast<GUInt64>(poDS->sHeader.nBitDepth))).v();
1290         bInvalidTileSize = (
1291             nMaxTileBits > static_cast<GUInt64>(std::numeric_limits<GUInt32>::max()));
1292     }
1293     catch( ... )
1294     {
1295         bInvalidTileSize = true;
1296     }
1297     if( bInvalidTileSize )
1298     {
1299         CPLError(CE_Warning, CPLE_IllegalArg,
1300                  "Invalid tile size. Width %lu, height %lu, bit depth %lu.",
1301                  static_cast<unsigned long>(poDS->sHeader.nTileWidth),
1302                  static_cast<unsigned long>(poDS->sHeader.nTileHeight),
1303                  static_cast<unsigned long>(poDS->sHeader.nBitDepth));
1304         delete poDS;
1305         return nullptr;
1306     }
1307 
1308     if(poDS->sHeader.nLastTileWidth > poDS->sHeader.nTileWidth ||
1309        poDS->sHeader.nLastTileHeight > poDS->sHeader.nTileHeight)
1310     {
1311         CPLError(CE_Warning, CPLE_IllegalArg,
1312                  "Invalid last tile size %lu x %lu. "
1313                  "It can't be greater than %lu x %lu.",
1314                  static_cast<unsigned long>(poDS->sHeader.nLastTileWidth),
1315                  static_cast<unsigned long>(poDS->sHeader.nLastTileHeight),
1316                  static_cast<unsigned long>(poDS->sHeader.nTileWidth),
1317                  static_cast<unsigned long>(poDS->sHeader.nTileHeight));
1318         delete poDS;
1319         return nullptr;
1320     }
1321 
1322     if( poParentDS != nullptr )
1323     {
1324         if( 0 != memcmp( poDS->sHeader.bySignature,
1325                          poParentDS->sHeader.bySignature,
1326                          RMF_SIGNATURE_SIZE ) )
1327         {
1328             CPLError( CE_Warning, CPLE_IllegalArg,
1329                       "Invalid subheader signature." );
1330             delete poDS;
1331             return nullptr;
1332         }
1333     }
1334 
1335 /* -------------------------------------------------------------------- */
1336 /*  Read the extended header.                                           */
1337 /* -------------------------------------------------------------------- */
1338 
1339     if( poDS->sHeader.nExtHdrOffset && poDS->sHeader.nExtHdrSize )
1340     {
1341         if( poDS->sHeader.nExtHdrSize > 1000000 )
1342         {
1343             delete poDS;
1344             return nullptr;
1345         }
1346         GByte *pabyExtHeader = reinterpret_cast<GByte *>(
1347             CPLCalloc( poDS->sHeader.nExtHdrSize, 1 ) );
1348         if( pabyExtHeader == nullptr )
1349         {
1350             delete poDS;
1351             return nullptr;
1352         }
1353 
1354         VSIFSeekL( poDS->fp, poDS->GetFileOffset( poDS->sHeader.nExtHdrOffset ),
1355                    SEEK_SET );
1356         VSIFReadL( pabyExtHeader, 1, poDS->sHeader.nExtHdrSize, poDS->fp );
1357 
1358         if( poDS->sHeader.nExtHdrSize >= 36 + 4 )
1359         {
1360             RMF_READ_LONG( pabyExtHeader, poDS->sExtHeader.nEllipsoid, 24 );
1361             RMF_READ_LONG( pabyExtHeader, poDS->sExtHeader.nVertDatum, 28 );
1362             RMF_READ_LONG( pabyExtHeader, poDS->sExtHeader.nDatum, 32 );
1363             RMF_READ_LONG( pabyExtHeader, poDS->sExtHeader.nZone, 36 );
1364         }
1365 
1366         CPLFree( pabyExtHeader );
1367     }
1368 
1369 #undef RMF_READ_DOUBLE
1370 #undef RMF_READ_LONG
1371 #undef RMF_READ_ULONG
1372 
1373     CPLDebug( "RMF", "Version %d", poDS->sHeader.iVersion );
1374 
1375 #ifdef DEBUG
1376 
1377     CPLDebug( "RMF", "%s image has width %d, height %d, bit depth %d, "
1378               "compression scheme %d, %s, nodata %f",
1379               (poDS->eRMFType == RMFT_MTW) ? "MTW" : "RSW",
1380               poDS->sHeader.nWidth, poDS->sHeader.nHeight,
1381               poDS->sHeader.nBitDepth, poDS->sHeader.iCompression,
1382               poDS->bBigEndian ? "big endian" : "little endian",
1383               poDS->sHeader.dfNoData );
1384     CPLDebug( "RMF", "Size %d, offset to overview %#lx, user ID %d, "
1385               "ROI offset %#lx, ROI size %d",
1386               poDS->sHeader.nSize,
1387               static_cast<unsigned long>( poDS->sHeader.nOvrOffset ),
1388               poDS->sHeader.iUserID,
1389               static_cast<unsigned long>( poDS->sHeader.nROIOffset ),
1390               poDS->sHeader.nROISize );
1391     CPLDebug( "RMF", "Map type %d, projection %d, scale %f, resolution %f, ",
1392               poDS->sHeader.iMapType, poDS->sHeader.iProjection,
1393               poDS->sHeader.dfScale, poDS->sHeader.dfResolution );
1394     CPLDebug( "RMF", "EPSG %d ", (int)poDS->sHeader.iEPSGCode );
1395     CPLDebug( "RMF", "Georeferencing: pixel size %f, LLX %f, LLY %f",
1396               poDS->sHeader.dfPixelSize,
1397               poDS->sHeader.dfLLX, poDS->sHeader.dfLLY );
1398     if( poDS->sHeader.nROIOffset && poDS->sHeader.nROISize )
1399     {
1400         GInt32 nValue = 0;
1401 
1402         CPLDebug( "RMF", "ROI coordinates:" );
1403         /* coverity[tainted_data] */
1404         for( GUInt32 i = 0; i < poDS->sHeader.nROISize; i += sizeof(nValue) )
1405         {
1406             if( VSIFSeekL( poDS->fp,
1407                            poDS->GetFileOffset( poDS->sHeader.nROIOffset + i ),
1408                            SEEK_SET ) != 0 ||
1409                 VSIFReadL( &nValue, 1, sizeof(nValue),
1410                            poDS->fp ) != sizeof(nValue) )
1411             {
1412                 CPLDebug("RMF", "Cannot read ROI at index %u", i);
1413                 break;
1414                 //delete poDS;
1415                 //return nullptr;
1416             }
1417 
1418             CPLDebug( "RMF", "%d", nValue );
1419         }
1420     }
1421 #endif
1422     if( poDS->sHeader.nWidth >= INT_MAX ||
1423         poDS->sHeader.nHeight >= INT_MAX ||
1424         !GDALCheckDatasetDimensions(poDS->sHeader.nWidth, poDS->sHeader.nHeight) )
1425     {
1426         delete poDS;
1427         return nullptr;
1428     }
1429 
1430 /* -------------------------------------------------------------------- */
1431 /*  Read array of blocks offsets/sizes.                                 */
1432 /* -------------------------------------------------------------------- */
1433 
1434     // To avoid useless excessive memory allocation
1435     if( poDS->sHeader.nTileTblSize > 1000000 )
1436     {
1437         VSIFSeekL( poDS->fp, 0, SEEK_END );
1438         vsi_l_offset nFileSize = VSIFTellL( poDS->fp );
1439         if( nFileSize < poDS->sHeader.nTileTblSize )
1440         {
1441             delete poDS;
1442             return nullptr;
1443         }
1444     }
1445 
1446     if( VSIFSeekL( poDS->fp,
1447                    poDS->GetFileOffset( poDS->sHeader.nTileTblOffset ),
1448                    SEEK_SET ) < 0 )
1449     {
1450         delete poDS;
1451         return nullptr;
1452     }
1453 
1454     poDS->paiTiles = reinterpret_cast<GUInt32 *>(
1455         VSIMalloc( poDS->sHeader.nTileTblSize ) );
1456     if( !poDS->paiTiles )
1457     {
1458         delete poDS;
1459         return nullptr;
1460     }
1461 
1462     if( VSIFReadL( poDS->paiTiles, 1, poDS->sHeader.nTileTblSize,
1463                    poDS->fp ) < poDS->sHeader.nTileTblSize )
1464     {
1465         CPLDebug( "RMF", "Can't read tiles offsets/sizes table." );
1466         delete poDS;
1467         return nullptr;
1468     }
1469 
1470 #ifdef CPL_MSB
1471     if( !poDS->bBigEndian )
1472     {
1473         for( GUInt32 i = 0;
1474              i < poDS->sHeader.nTileTblSize / sizeof(GUInt32);
1475              i++ )
1476             CPL_SWAP32PTR( poDS->paiTiles + i );
1477     }
1478 #else
1479     if( poDS->bBigEndian )
1480     {
1481         for( GUInt32 i = 0;
1482              i < poDS->sHeader.nTileTblSize / sizeof(GUInt32);
1483              i++ )
1484             CPL_SWAP32PTR( poDS->paiTiles + i );
1485     }
1486 #endif
1487 
1488 #ifdef DEBUG
1489     CPLDebug( "RMF", "List of block offsets/sizes:" );
1490 
1491     for( GUInt32 i = 0;
1492          i < poDS->sHeader.nTileTblSize / sizeof(GUInt32);
1493          i += 2 )
1494     {
1495         CPLDebug( "RMF", "    %u / %u",
1496                   poDS->paiTiles[i], poDS->paiTiles[i + 1] );
1497     }
1498 #endif
1499 
1500 /* -------------------------------------------------------------------- */
1501 /*  Set up essential image parameters.                                  */
1502 /* -------------------------------------------------------------------- */
1503     GDALDataType eType = GDT_Byte;
1504 
1505     poDS->nRasterXSize = poDS->sHeader.nWidth;
1506     poDS->nRasterYSize = poDS->sHeader.nHeight;
1507 
1508     if( poDS->eRMFType == RMFT_RSW )
1509     {
1510         switch( poDS->sHeader.nBitDepth )
1511         {
1512             case 32:
1513             case 24:
1514             case 16:
1515                 poDS->nBands = 3;
1516                 break;
1517             case 1:
1518             case 4:
1519             case 8:
1520                 if( poParentDS != nullptr && poParentDS->poColorTable != nullptr )
1521                 {
1522                     poDS->poColorTable = poParentDS->poColorTable->Clone();
1523                 }
1524                 else
1525                 {
1526                     // Allocate memory for colour table and read it
1527                     poDS->nColorTableSize = 1 << poDS->sHeader.nBitDepth;
1528                     GUInt32 nExpectedColorTableBytes = poDS->nColorTableSize * 4;
1529                     if(nExpectedColorTableBytes > poDS->sHeader.nClrTblSize )
1530                     {
1531                         // We could probably test for strict equality in
1532                         // the above test ???
1533                         CPLDebug( "RMF",
1534                                   "Wrong color table size. "
1535                                   "Expected %u, got %u.",
1536                                   nExpectedColorTableBytes,
1537                                   poDS->sHeader.nClrTblSize );
1538                         delete poDS;
1539                         return nullptr;
1540                     }
1541                     poDS->pabyColorTable = reinterpret_cast<GByte *>(
1542                         VSIMalloc( nExpectedColorTableBytes ) );
1543                     if( poDS->pabyColorTable == nullptr )
1544                     {
1545                         CPLDebug( "RMF", "Can't allocate color table." );
1546                         delete poDS;
1547                         return nullptr;
1548                     }
1549                     if( VSIFSeekL( poDS->fp,
1550                                    poDS->GetFileOffset( poDS->sHeader.nClrTblOffset ),
1551                                    SEEK_SET ) < 0 )
1552                     {
1553                         CPLDebug( "RMF",
1554                                   "Can't seek to color table location." );
1555                         delete poDS;
1556                         return nullptr;
1557                     }
1558                     if( VSIFReadL( poDS->pabyColorTable, 1,
1559                                    nExpectedColorTableBytes, poDS->fp )
1560                         < nExpectedColorTableBytes )
1561                     {
1562                         CPLDebug( "RMF", "Can't read color table." );
1563                         delete poDS;
1564                         return nullptr;
1565                     }
1566 
1567                     poDS->poColorTable = new GDALColorTable();
1568                     for( GUInt32 i = 0; i < poDS->nColorTableSize; i++ )
1569                     {
1570                         const GDALColorEntry oEntry = {
1571                             poDS->pabyColorTable[i * 4],     // Red
1572                             poDS->pabyColorTable[i * 4 + 1], // Green
1573                             poDS->pabyColorTable[i * 4 + 2], // Blue
1574                             255                              // Alpha
1575                         };
1576 
1577                         poDS->poColorTable->SetColorEntry( i, &oEntry );
1578                     }
1579                 }
1580                 poDS->nBands = 1;
1581                 break;
1582             default:
1583                 CPLError(CE_Warning, CPLE_IllegalArg,
1584                          "Invalid RSW bit depth %lu.",
1585                          static_cast<unsigned long>(poDS->sHeader.nBitDepth));
1586                 delete poDS;
1587                 return nullptr;
1588         }
1589         eType = GDT_Byte;
1590     }
1591     else
1592     {
1593         poDS->nBands = 1;
1594         if( poDS->sHeader.nBitDepth == 8 )
1595         {
1596             eType = GDT_Byte;
1597         }
1598         else if( poDS->sHeader.nBitDepth == 16 )
1599         {
1600             eType = GDT_Int16;
1601         }
1602         else if( poDS->sHeader.nBitDepth == 32 )
1603         {
1604             eType = GDT_Int32;
1605         }
1606         else if( poDS->sHeader.nBitDepth == 64 )
1607         {
1608             eType = GDT_Float64;
1609         }
1610         else
1611         {
1612             CPLError(CE_Warning, CPLE_IllegalArg,
1613                      "Invalid MTW bit depth %lu.",
1614                      static_cast<unsigned long>(poDS->sHeader.nBitDepth));
1615             delete poDS;
1616             return nullptr;
1617         }
1618     }
1619 
1620     if( poDS->sHeader.nTileWidth == 0 || poDS->sHeader.nTileWidth > INT_MAX ||
1621         poDS->sHeader.nTileHeight == 0 || poDS->sHeader.nTileHeight > INT_MAX )
1622     {
1623         CPLDebug("RMF", "Invalid tile dimension : %u x %u",
1624                  poDS->sHeader.nTileWidth, poDS->sHeader.nTileHeight);
1625         delete poDS;
1626         return nullptr;
1627     }
1628 
1629     const int nDataSize = GDALGetDataTypeSizeBytes( eType );
1630     const int nBlockXSize = static_cast<int>(poDS->sHeader.nTileWidth);
1631     const int nBlockYSize = static_cast<int>(poDS->sHeader.nTileHeight);
1632     if( nDataSize == 0 ||
1633         nBlockXSize > INT_MAX / nBlockYSize ||
1634         nBlockYSize > INT_MAX / nDataSize ||
1635         nBlockXSize > INT_MAX / (nBlockYSize * nDataSize) )
1636     {
1637         CPLDebug ("RMF", "Too big raster / tile dimension");
1638         delete poDS;
1639         return nullptr;
1640     }
1641 
1642     poDS->nXTiles = DIV_ROUND_UP( poDS->nRasterXSize, nBlockXSize );
1643     poDS->nYTiles = DIV_ROUND_UP( poDS->nRasterYSize, nBlockYSize );
1644 
1645 #ifdef DEBUG
1646     CPLDebug( "RMF", "Image is %d tiles wide, %d tiles long",
1647               poDS->nXTiles, poDS->nYTiles );
1648 #endif
1649 
1650 /* -------------------------------------------------------------------- */
1651 /*  Choose compression scheme.                                          */
1652 /* -------------------------------------------------------------------- */
1653     if(CE_None != poDS->SetupCompression(eType, poOpenInfo->pszFilename))
1654     {
1655         delete poDS;
1656         return nullptr;
1657     }
1658 
1659     if(poOpenInfo->eAccess == GA_Update)
1660     {
1661         if(poParentDS == nullptr)
1662         {
1663             if(CE_None != poDS->InitCompressorData(poOpenInfo->papszOpenOptions))
1664             {
1665                 delete poDS;
1666                 return nullptr;
1667             }
1668         }
1669         else
1670         {
1671             poDS->poCompressData = poParentDS->poCompressData;
1672         }
1673     }
1674 /* -------------------------------------------------------------------- */
1675 /*  Create band information objects.                                    */
1676 /* -------------------------------------------------------------------- */
1677     for( int iBand = 1; iBand <= poDS->nBands; iBand++ )
1678         poDS->SetBand( iBand, new RMFRasterBand( poDS, iBand, eType ) );
1679 
1680     poDS->SetupNBits();
1681 
1682     if(poDS->nBands > 1)
1683     {
1684         poDS->SetMetadataItem("INTERLEAVE", "PIXEL", "IMAGE_STRUCTURE");
1685     }
1686 /* -------------------------------------------------------------------- */
1687 /*  Set up projection.                                                  */
1688 /*                                                                      */
1689 /*  XXX: If projection value is not specified, but image still have     */
1690 /*  georeferencing information, assume Gauss-Kruger projection.         */
1691 /* -------------------------------------------------------------------- */
1692     if(poDS->sHeader.iEPSGCode > RMF_EPSG_MIN_CODE ||
1693        poDS->sHeader.iProjection > 0 ||
1694        (poDS->sHeader.dfPixelSize != 0.0 &&
1695         poDS->sHeader.dfLLX != 0.0 &&
1696         poDS->sHeader.dfLLY != 0.0))
1697     {
1698         OGRSpatialReference oSRS;
1699         GInt32 nProj =
1700             (poDS->sHeader.iProjection) ? poDS->sHeader.iProjection : 1;
1701         double padfPrjParams[8] = {
1702             poDS->sHeader.dfStdP1,
1703             poDS->sHeader.dfStdP2,
1704             poDS->sHeader.dfCenterLat,
1705             poDS->sHeader.dfCenterLong,
1706             1.0,
1707             0.0,
1708             0.0,
1709             0.0
1710           };
1711 
1712         // XXX: Compute zone number for Gauss-Kruger (Transverse Mercator)
1713         // projection if it is not specified.
1714         if( nProj == 1L && poDS->sHeader.dfCenterLong == 0.0 )
1715         {
1716             if( poDS->sExtHeader.nZone == 0 )
1717             {
1718                 double centerXCoord = poDS->sHeader.dfLLX +
1719                     (poDS->nRasterXSize * poDS->sHeader.dfPixelSize / 2.0);
1720                 padfPrjParams[7] =
1721                     floor((centerXCoord - 500000.0 ) / 1000000.0);
1722             }
1723             else
1724             {
1725                 padfPrjParams[7] = poDS->sExtHeader.nZone;
1726             }
1727         }
1728 
1729         OGRErr  res = OGRERR_FAILURE;
1730         if(nProj >= 0 &&
1731            (poDS->sExtHeader.nDatum >= 0 || poDS->sExtHeader.nEllipsoid >= 0))
1732         {
1733             res = oSRS.importFromPanorama( nProj, poDS->sExtHeader.nDatum,
1734                                  poDS->sExtHeader.nEllipsoid, padfPrjParams );
1735         }
1736 
1737         if(poDS->sHeader.iEPSGCode > RMF_EPSG_MIN_CODE &&
1738            (OGRERR_NONE != res || oSRS.IsLocal()))
1739         {
1740             res = oSRS.importFromEPSG(poDS->sHeader.iEPSGCode);
1741         }
1742 
1743         const char* pszSetVertCS =
1744             CSLFetchNameValueDef(poOpenInfo->papszOpenOptions,
1745                                 "RMF_SET_VERTCS",
1746                                  CPLGetConfigOption("RMF_SET_VERTCS", "NO"));
1747         if(CPLTestBool(pszSetVertCS) && res == OGRERR_NONE &&
1748            poDS->sExtHeader.nVertDatum > 0)
1749         {
1750             oSRS.importVertCSFromPanorama(poDS->sExtHeader.nVertDatum);
1751         }
1752 
1753         if( poDS->pszProjection )
1754             CPLFree( poDS->pszProjection );
1755         oSRS.exportToWkt( &poDS->pszProjection );
1756     }
1757 
1758 /* -------------------------------------------------------------------- */
1759 /*  Set up georeferencing.                                              */
1760 /* -------------------------------------------------------------------- */
1761     if( (poDS->eRMFType == RMFT_RSW && poDS->sHeader.iGeorefFlag) ||
1762         (poDS->eRMFType == RMFT_MTW && poDS->sHeader.dfPixelSize != 0.0) )
1763     {
1764         poDS->adfGeoTransform[0] = poDS->sHeader.dfLLX;
1765         poDS->adfGeoTransform[3] = poDS->sHeader.dfLLY
1766             + poDS->nRasterYSize * poDS->sHeader.dfPixelSize;
1767         poDS->adfGeoTransform[1] = poDS->sHeader.dfPixelSize;
1768         poDS->adfGeoTransform[5] = - poDS->sHeader.dfPixelSize;
1769         poDS->adfGeoTransform[2] = 0.0;
1770         poDS->adfGeoTransform[4] = 0.0;
1771     }
1772 
1773 /* -------------------------------------------------------------------- */
1774 /*  Set units.                                                          */
1775 /* -------------------------------------------------------------------- */
1776 
1777     if( poDS->eRMFType == RMFT_MTW )
1778     {
1779         CPLFree(poDS->pszUnitType);
1780         poDS->pszUnitType = RMFUnitTypeToStr(poDS->sHeader.iElevationUnit);
1781     }
1782 
1783 /* -------------------------------------------------------------------- */
1784 /*  Report some other dataset related information.                      */
1785 /* -------------------------------------------------------------------- */
1786 
1787     if( poDS->eRMFType == RMFT_MTW )
1788     {
1789         char szTemp[256] = {};
1790 
1791         snprintf(szTemp, sizeof(szTemp), "%g", poDS->sHeader.adfElevMinMax[0]);
1792         poDS->SetMetadataItem( "ELEVATION_MINIMUM", szTemp );
1793 
1794         snprintf(szTemp, sizeof(szTemp), "%g", poDS->sHeader.adfElevMinMax[1]);
1795         poDS->SetMetadataItem( "ELEVATION_MAXIMUM", szTemp );
1796 
1797         poDS->SetMetadataItem( "ELEVATION_UNITS", poDS->pszUnitType );
1798 
1799         snprintf(szTemp, sizeof(szTemp), "%d", poDS->sHeader.iElevationType);
1800         poDS->SetMetadataItem( "ELEVATION_TYPE", szTemp );
1801     }
1802 
1803 /* -------------------------------------------------------------------- */
1804 /*      Check for overviews.                                            */
1805 /* -------------------------------------------------------------------- */
1806     if( nNextHeaderOffset == 0 && poParentDS == nullptr )
1807     {
1808         poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
1809     }
1810 
1811     return poDS;
1812 }
1813 
1814 /************************************************************************/
1815 /*                               Create()                               */
1816 /************************************************************************/
Create(const char * pszFilename,int nXSize,int nYSize,int nBands,GDALDataType eType,char ** papszParamList)1817 GDALDataset *RMFDataset::Create( const char * pszFilename,
1818                                  int nXSize, int nYSize, int nBands,
1819                                  GDALDataType eType, char **papszParamList )
1820 {
1821     return Create( pszFilename, nXSize, nYSize, nBands,
1822                    eType, papszParamList, nullptr, 1.0 );
1823 }
1824 
Create(const char * pszFilename,int nXSize,int nYSize,int nBands,GDALDataType eType,char ** papszParamList,RMFDataset * poParentDS,double dfOvFactor)1825 GDALDataset *RMFDataset::Create( const char * pszFilename,
1826                                  int nXSize, int nYSize, int nBands,
1827                                  GDALDataType eType, char **papszParamList,
1828                                  RMFDataset* poParentDS, double dfOvFactor )
1829 
1830 {
1831     if( nBands != 1 && nBands != 3 )
1832     {
1833         CPLError( CE_Failure, CPLE_NotSupported,
1834                   "RMF driver doesn't support %d bands. Must be 1 or 3.",
1835                   nBands );
1836 
1837         return nullptr;
1838     }
1839 
1840     if( nBands == 1
1841         && eType != GDT_Byte
1842         && eType != GDT_Int16
1843         && eType != GDT_Int32
1844         && eType != GDT_Float64 )
1845     {
1846          CPLError(
1847              CE_Failure, CPLE_AppDefined,
1848              "Attempt to create RMF dataset with an illegal data type (%s), "
1849              "only Byte, Int16, Int32 and Float64 types supported "
1850              "by the format for single-band images.",
1851              GDALGetDataTypeName(eType) );
1852 
1853         return nullptr;
1854     }
1855 
1856     if( nBands == 3 && eType != GDT_Byte )
1857     {
1858          CPLError(
1859              CE_Failure, CPLE_AppDefined,
1860              "Attempt to create RMF dataset with an illegal data type (%s), "
1861              "only Byte type supported by the format for three-band images.",
1862              GDALGetDataTypeName(eType) );
1863 
1864         return nullptr;
1865     }
1866 
1867 /* -------------------------------------------------------------------- */
1868 /*  Create the dataset.                                                 */
1869 /* -------------------------------------------------------------------- */
1870     RMFDataset *poDS = new RMFDataset();
1871 
1872     GUInt32 nBlockXSize =
1873         ( nXSize < RMF_DEFAULT_BLOCKXSIZE ) ? nXSize : RMF_DEFAULT_BLOCKXSIZE;
1874     GUInt32 nBlockYSize =
1875         ( nYSize < RMF_DEFAULT_BLOCKYSIZE ) ? nYSize : RMF_DEFAULT_BLOCKYSIZE;
1876     double dfScale;
1877     double dfResolution;
1878     double dfPixelSize;
1879     if( poParentDS == nullptr )
1880     {
1881         poDS->fp = VSIFOpenL( pszFilename, "w+b" );
1882         if( poDS->fp == nullptr )
1883         {
1884             CPLError( CE_Failure, CPLE_OpenFailed, "Unable to create file %s.",
1885                       pszFilename );
1886             delete poDS;
1887             return nullptr;
1888         }
1889 
1890         dfScale = RMF_DEFAULT_SCALE;
1891         dfResolution = RMF_DEFAULT_RESOLUTION;
1892         dfPixelSize = 1;
1893 
1894         if( CPLFetchBool( papszParamList, "MTW", false) )
1895             poDS->eRMFType = RMFT_MTW;
1896         else
1897             poDS->eRMFType = RMFT_RSW;
1898 
1899         GUInt32 iVersion = RMF_VERSION;
1900         const char *pszRMFHUGE = CSLFetchNameValue(papszParamList, "RMFHUGE");
1901 
1902         if( pszRMFHUGE == nullptr )
1903             pszRMFHUGE = "NO";// Keep old behavior by default
1904 
1905         if( EQUAL(pszRMFHUGE,"NO") )
1906         {
1907             iVersion = RMF_VERSION;
1908         }
1909         else if( EQUAL(pszRMFHUGE,"YES") )
1910         {
1911             iVersion = RMF_VERSION_HUGE;
1912         }
1913         else if( EQUAL(pszRMFHUGE,"IF_SAFER") )
1914         {
1915             const double dfImageSize =
1916                 static_cast<double>(nXSize) *
1917                 static_cast<double>(nYSize) *
1918                 static_cast<double>(nBands) *
1919                 static_cast<double>(GDALGetDataTypeSizeBytes(eType));
1920             if( dfImageSize > 3.0*1024.0*1024.0*1024.0 )
1921             {
1922                 iVersion = RMF_VERSION_HUGE;
1923             }
1924             else
1925             {
1926                 iVersion = RMF_VERSION;
1927             }
1928         }
1929 
1930         const char *pszValue = CSLFetchNameValue(papszParamList,"BLOCKXSIZE");
1931         if( pszValue != nullptr )
1932             nBlockXSize = atoi( pszValue );
1933         if( static_cast<int>(nBlockXSize) <= 0 )
1934             nBlockXSize = RMF_DEFAULT_BLOCKXSIZE;
1935 
1936         pszValue = CSLFetchNameValue(papszParamList,"BLOCKYSIZE");
1937         if( pszValue != nullptr )
1938             nBlockYSize = atoi( pszValue );
1939         if( static_cast<int>(nBlockYSize) <= 0 )
1940             nBlockYSize = RMF_DEFAULT_BLOCKXSIZE;
1941 
1942         if( poDS->eRMFType == RMFT_MTW )
1943             memcpy( poDS->sHeader.bySignature, RMF_SigMTW, RMF_SIGNATURE_SIZE );
1944         else
1945             memcpy( poDS->sHeader.bySignature, RMF_SigRSW, RMF_SIGNATURE_SIZE );
1946         poDS->sHeader.iVersion = iVersion;
1947         poDS->sHeader.nOvrOffset = 0x00;
1948     }
1949     else
1950     {
1951         poDS->fp = poParentDS->fp;
1952         memcpy( poDS->sHeader.bySignature, poParentDS->sHeader.bySignature,
1953                 RMF_SIGNATURE_SIZE );
1954         poDS->sHeader.iVersion = poParentDS->sHeader.iVersion;
1955         poDS->eRMFType = poParentDS->eRMFType;
1956         nBlockXSize = poParentDS->sHeader.nTileWidth;
1957         nBlockYSize = poParentDS->sHeader.nTileHeight;
1958         dfScale = poParentDS->sHeader.dfScale;
1959         dfResolution = poParentDS->sHeader.dfResolution/dfOvFactor;
1960         dfPixelSize = poParentDS->sHeader.dfPixelSize*dfOvFactor;
1961 
1962         poDS->nHeaderOffset = poParentDS->GetLastOffset();
1963         poParentDS->sHeader.nOvrOffset =
1964                poDS->GetRMFOffset( poDS->nHeaderOffset, &poDS->nHeaderOffset );
1965         poParentDS->bHeaderDirty = true;
1966         VSIFSeekL( poDS->fp, poDS->nHeaderOffset, SEEK_SET );
1967         poDS->poParentDS = poParentDS;
1968         CPLDebug( "RMF",
1969                   "Create overview subfile at " CPL_FRMT_GUIB
1970                   " with size %dx%d, parent overview offset %d",
1971                   poDS->nHeaderOffset, nXSize, nYSize,
1972                   poParentDS->sHeader.nOvrOffset );
1973     }
1974 /* -------------------------------------------------------------------- */
1975 /*  Fill the RMFHeader                                                  */
1976 /* -------------------------------------------------------------------- */
1977     CPLDebug( "RMF", "Version %d", poDS->sHeader.iVersion );
1978 
1979     poDS->sHeader.iUserID = 0x00;
1980     memset( poDS->sHeader.byName, 0, sizeof(poDS->sHeader.byName) );
1981     poDS->sHeader.nBitDepth = GDALGetDataTypeSizeBits( eType ) * nBands;
1982     poDS->sHeader.nHeight = nYSize;
1983     poDS->sHeader.nWidth = nXSize;
1984     poDS->sHeader.nTileWidth = nBlockXSize;
1985     poDS->sHeader.nTileHeight = nBlockYSize;
1986 
1987     poDS->nXTiles = poDS->sHeader.nXTiles =
1988         ( nXSize + poDS->sHeader.nTileWidth - 1 ) / poDS->sHeader.nTileWidth;
1989     poDS->nYTiles = poDS->sHeader.nYTiles =
1990         ( nYSize + poDS->sHeader.nTileHeight - 1 ) / poDS->sHeader.nTileHeight;
1991     poDS->sHeader.nLastTileHeight = nYSize % poDS->sHeader.nTileHeight;
1992     if( !poDS->sHeader.nLastTileHeight )
1993         poDS->sHeader.nLastTileHeight = poDS->sHeader.nTileHeight;
1994     poDS->sHeader.nLastTileWidth = nXSize % poDS->sHeader.nTileWidth;
1995     if( !poDS->sHeader.nLastTileWidth )
1996         poDS->sHeader.nLastTileWidth = poDS->sHeader.nTileWidth;
1997 
1998     poDS->sHeader.nROIOffset = 0x00;
1999     poDS->sHeader.nROISize = 0x00;
2000 
2001     vsi_l_offset nCurPtr = poDS->nHeaderOffset + RMF_HEADER_SIZE;
2002 
2003     // Extended header
2004     poDS->sHeader.nExtHdrOffset = poDS->GetRMFOffset( nCurPtr, &nCurPtr );
2005     poDS->sHeader.nExtHdrSize = RMF_EXT_HEADER_SIZE;
2006     nCurPtr += poDS->sHeader.nExtHdrSize;
2007 
2008     // Color table
2009     if( poDS->eRMFType == RMFT_RSW && nBands == 1 )
2010     {
2011         if( poDS->sHeader.nBitDepth > 8 )
2012         {
2013             CPLError( CE_Failure, CPLE_AppDefined,
2014                       "Cannot create color table of RSW with nBitDepth = %d. "
2015                       "Retry with MTW ?",
2016                       poDS->sHeader.nBitDepth );
2017             delete poDS;
2018             return nullptr;
2019         }
2020 
2021         poDS->sHeader.nClrTblOffset = poDS->GetRMFOffset( nCurPtr, &nCurPtr );
2022         poDS->nColorTableSize = 1 << poDS->sHeader.nBitDepth;
2023         poDS->sHeader.nClrTblSize = poDS->nColorTableSize * 4;
2024         poDS->pabyColorTable = reinterpret_cast<GByte *>(
2025             VSI_MALLOC_VERBOSE( poDS->sHeader.nClrTblSize ) );
2026         if( poDS->pabyColorTable == nullptr )
2027         {
2028             delete poDS;
2029             return nullptr;
2030         }
2031         for( GUInt32 i = 0; i < poDS->nColorTableSize; i++ )
2032         {
2033             poDS->pabyColorTable[i * 4] =
2034                 poDS->pabyColorTable[i * 4 + 1] =
2035                 poDS->pabyColorTable[i * 4 + 2] = (GByte) i;
2036                 poDS->pabyColorTable[i * 4 + 3] = 0;
2037         }
2038         nCurPtr += poDS->sHeader.nClrTblSize;
2039     }
2040     else
2041     {
2042         poDS->sHeader.nClrTblOffset = 0x00;
2043         poDS->sHeader.nClrTblSize = 0x00;
2044     }
2045 
2046     // Blocks table
2047     poDS->sHeader.nTileTblOffset = poDS->GetRMFOffset( nCurPtr, &nCurPtr );
2048     poDS->sHeader.nTileTblSize =
2049         poDS->sHeader.nXTiles * poDS->sHeader.nYTiles * 4 * 2;
2050     poDS->paiTiles = reinterpret_cast<GUInt32 *>(
2051         CPLCalloc( poDS->sHeader.nTileTblSize, 1 ) );
2052     // nCurPtr += poDS->sHeader.nTileTblSize;
2053     const GUInt32 nTileSize =
2054         poDS->sHeader.nTileWidth * poDS->sHeader.nTileHeight
2055         * GDALGetDataTypeSizeBytes( eType );
2056     poDS->sHeader.nSize =
2057         poDS->paiTiles[poDS->sHeader.nTileTblSize / 4 - 2] + nTileSize;
2058 
2059     // Elevation units
2060     poDS->sHeader.iElevationUnit = RMFStrToUnitType(poDS->pszUnitType);
2061 
2062     poDS->sHeader.iMapType = -1;
2063     poDS->sHeader.iProjection = -1;
2064     poDS->sHeader.iEPSGCode = -1;
2065     poDS->sHeader.dfScale = dfScale;
2066     poDS->sHeader.dfResolution = dfResolution;
2067     poDS->sHeader.dfPixelSize = dfPixelSize;
2068     poDS->sHeader.iMaskType = 0;
2069     poDS->sHeader.iMaskStep = 0;
2070     poDS->sHeader.iFrameFlag = 0;
2071     poDS->sHeader.nFlagsTblOffset = 0x00;
2072     poDS->sHeader.nFlagsTblSize = 0x00;
2073     poDS->sHeader.nFileSize0 = 0x00;
2074     poDS->sHeader.nFileSize1 = 0x00;
2075     poDS->sHeader.iUnknown = 0;
2076     poDS->sHeader.iGeorefFlag = 0;
2077     poDS->sHeader.iInverse = 0;
2078     poDS->sHeader.iJpegQuality = 0;
2079     memset( poDS->sHeader.abyInvisibleColors, 0,
2080             sizeof(poDS->sHeader.abyInvisibleColors) );
2081     poDS->sHeader.iElevationType = 0;
2082 
2083     poDS->nRasterXSize = nXSize;
2084     poDS->nRasterYSize = nYSize;
2085     poDS->eAccess = GA_Update;
2086     poDS->nBands = nBands;
2087 
2088     if(poParentDS == nullptr)
2089     {
2090         poDS->sHeader.adfElevMinMax[0] = 0.0;
2091         poDS->sHeader.adfElevMinMax[1] = 0.0;
2092         poDS->sHeader.dfNoData = 0.0;
2093         poDS->sHeader.iCompression = GetCompressionType(
2094                                         CSLFetchNameValue(papszParamList,
2095                                                           "COMPRESS"));
2096         if(CE_None != poDS->InitCompressorData(papszParamList))
2097         {
2098             delete poDS;
2099             return nullptr;
2100         }
2101 
2102         if(poDS->sHeader.iCompression == RMF_COMPRESSION_JPEG)
2103         {
2104             const char* pszJpegQuality = CSLFetchNameValue(papszParamList,
2105                                                            "JPEG_QUALITY");
2106             if(pszJpegQuality == nullptr)
2107             {
2108                 poDS->sHeader.iJpegQuality = 75;
2109             }
2110             else
2111             {
2112                 int iJpegQuality = atoi(pszJpegQuality);
2113                 if( iJpegQuality < 10 || iJpegQuality > 100)
2114                 {
2115                     CPLError(CE_Failure, CPLE_IllegalArg,
2116                              "JPEG_QUALITY=%s is not a legal value in the range 10-100.\n"
2117                              "Defaulting to 75",
2118                              pszJpegQuality);
2119                     iJpegQuality = 75;
2120                 }
2121                 poDS->sHeader.iJpegQuality = static_cast<GByte>(iJpegQuality);
2122             }
2123         }
2124 
2125         if(CE_None != poDS->SetupCompression(eType, pszFilename))
2126         {
2127             delete poDS;
2128             return nullptr;
2129         }
2130     }
2131     else
2132     {
2133         poDS->sHeader.adfElevMinMax[0] = poParentDS->sHeader.adfElevMinMax[0];
2134         poDS->sHeader.adfElevMinMax[1] = poParentDS->sHeader.adfElevMinMax[1];
2135         poDS->sHeader.dfNoData = poParentDS->sHeader.dfNoData;
2136         poDS->sHeader.iCompression = poParentDS->sHeader.iCompression;
2137         poDS->sHeader.iJpegQuality = poParentDS->sHeader.iJpegQuality;
2138         poDS->Decompress = poParentDS->Decompress;
2139         poDS->Compress = poParentDS->Compress;
2140         poDS->poCompressData = poParentDS->poCompressData;
2141     }
2142 
2143     if(nBands > 1)
2144     {
2145         poDS->SetMetadataItem("INTERLEAVE", "PIXEL", "IMAGE_STRUCTURE");
2146     }
2147 
2148     poDS->WriteHeader();
2149 
2150 /* -------------------------------------------------------------------- */
2151 /*      Create band information objects.                                */
2152 /* -------------------------------------------------------------------- */
2153     for( int iBand = 1; iBand <= poDS->nBands; iBand++ )
2154         poDS->SetBand( iBand, new RMFRasterBand( poDS, iBand, eType ) );
2155 
2156     poDS->SetupNBits();
2157 
2158     return reinterpret_cast<GDALDataset *>( poDS );
2159 }
2160 
2161 //GIS Panorama 11 was introduced new format for huge files (greater than 3 Gb)
GetFileOffset(GUInt32 iRMFOffset) const2162 vsi_l_offset RMFDataset::GetFileOffset( GUInt32 iRMFOffset ) const
2163 {
2164     if( sHeader.iVersion >= RMF_VERSION_HUGE )
2165     {
2166         return ((vsi_l_offset)iRMFOffset) * RMF_HUGE_OFFSET_FACTOR;
2167     }
2168 
2169     return (vsi_l_offset)iRMFOffset;
2170 }
2171 
GetRMFOffset(vsi_l_offset nFileOffset,vsi_l_offset * pnNewFileOffset) const2172 GUInt32 RMFDataset::GetRMFOffset( vsi_l_offset nFileOffset,
2173                                   vsi_l_offset* pnNewFileOffset ) const
2174 {
2175     if( sHeader.iVersion >= RMF_VERSION_HUGE )
2176     {
2177         //Round offset to next RMF_HUGE_OFFSET_FACTOR
2178         const GUInt32 iRMFOffset = static_cast<GUInt32>(
2179             (nFileOffset + (RMF_HUGE_OFFSET_FACTOR-1) ) /
2180             RMF_HUGE_OFFSET_FACTOR );
2181         if( pnNewFileOffset != nullptr )
2182         {
2183             *pnNewFileOffset = GetFileOffset( iRMFOffset );
2184         }
2185         return iRMFOffset;
2186     }
2187 
2188     if( pnNewFileOffset != nullptr )
2189     {
2190         *pnNewFileOffset = nFileOffset;
2191     }
2192     return static_cast<GUInt32>(nFileOffset);
2193 }
2194 
OpenOverview(RMFDataset * poParent,GDALOpenInfo * poOpenInfo)2195 RMFDataset* RMFDataset::OpenOverview(RMFDataset* poParent, GDALOpenInfo* poOpenInfo)
2196 {
2197     if( sHeader.nOvrOffset == 0 )
2198     {
2199         return nullptr;
2200     }
2201 
2202     if( poParent == nullptr )
2203     {
2204         return nullptr;
2205     }
2206 
2207     vsi_l_offset nSubOffset = GetFileOffset(sHeader.nOvrOffset);
2208 
2209     CPLDebug( "RMF",
2210               "Try to open overview subfile at " CPL_FRMT_GUIB " for '%s'",
2211               nSubOffset, poOpenInfo->pszFilename );
2212 
2213     if( !poParent->poOvrDatasets.empty() )
2214     {
2215         if( poParent->GetFileOffset( poParent->sHeader.nOvrOffset ) ==
2216             nSubOffset )
2217         {
2218             CPLError( CE_Warning, CPLE_IllegalArg,
2219                       "Recursive subdataset list is detected. "
2220                       "Overview open failed." );
2221             return nullptr;
2222         }
2223 
2224         for( size_t n = 0; n != poParent->poOvrDatasets.size() - 1; ++n )
2225         {
2226             RMFDataset* poOvr( poParent->poOvrDatasets[n] );
2227 
2228             if( poOvr == nullptr )
2229                 continue;
2230             if( poOvr->GetFileOffset( poOvr->sHeader.nOvrOffset ) ==
2231                 nSubOffset )
2232             {
2233                 CPLError( CE_Warning, CPLE_IllegalArg,
2234                           "Recursive subdataset list is detected. "
2235                           "Overview open failed." );
2236                 return nullptr;
2237             }
2238         }
2239     }
2240 
2241     size_t nHeaderSize( RMF_HEADER_SIZE );
2242     GByte * pabyNewHeader;
2243     pabyNewHeader = static_cast<GByte *>( CPLRealloc(poOpenInfo->pabyHeader,
2244                                                      nHeaderSize + 1) );
2245     if( pabyNewHeader == nullptr )
2246     {
2247         CPLError( CE_Warning, CPLE_OutOfMemory,
2248                   "Can't allocate buffer for overview header" );
2249         return nullptr;
2250     }
2251 
2252     poOpenInfo->pabyHeader = pabyNewHeader;
2253     memset( poOpenInfo->pabyHeader, 0, nHeaderSize + 1 );
2254     VSIFSeekL( fp, nSubOffset, SEEK_SET );
2255     poOpenInfo->nHeaderBytes = static_cast<int>( VSIFReadL( poOpenInfo->pabyHeader,
2256                                                  1, nHeaderSize, fp ) );
2257 
2258     RMFDataset* poSub = (RMFDataset*)Open( poOpenInfo, poParent, nSubOffset );
2259 
2260     if( poSub == nullptr )
2261     {
2262         return nullptr;
2263     }
2264 
2265     return poSub;
2266 }
2267 
IBuildOverviews(const char * pszResampling,int nOverviews,int * panOverviewList,int nBandsIn,int * panBandList,GDALProgressFunc pfnProgress,void * pProgressData)2268 CPLErr RMFDataset::IBuildOverviews( const char* pszResampling,
2269                                     int nOverviews, int* panOverviewList,
2270                                     int nBandsIn, int* panBandList,
2271                                     GDALProgressFunc pfnProgress,
2272                                     void* pProgressData )
2273 {
2274     bool bUseGenericHandling = false;
2275 
2276     if( GetAccess() != GA_Update )
2277     {
2278         CPLDebug( "RMF",
2279                   "File open for read-only accessing, "
2280                   "creating overviews externally." );
2281 
2282         bUseGenericHandling = true;
2283     }
2284 
2285     if( bUseGenericHandling )
2286     {
2287         if( !poOvrDatasets.empty() )
2288         {
2289             CPLError(
2290                 CE_Failure, CPLE_NotSupported,
2291                 "Cannot add external overviews when there are already "
2292                 "internal overviews" );
2293             return CE_Failure;
2294         }
2295 
2296         return GDALDataset::IBuildOverviews(
2297             pszResampling, nOverviews, panOverviewList,
2298             nBandsIn, panBandList, pfnProgress, pProgressData );
2299     }
2300 
2301     if( nBandsIn != GetRasterCount() )
2302     {
2303         CPLError( CE_Failure, CPLE_NotSupported,
2304                   "Generation of overviews in RMF is only "
2305                   "supported when operating on all bands.  "
2306                   "Operation failed." );
2307         return CE_Failure;
2308     }
2309 
2310     if( nOverviews == 0 )
2311     {
2312         if( poOvrDatasets.empty() )
2313         {
2314             return GDALDataset::IBuildOverviews(
2315                 pszResampling, nOverviews, panOverviewList,
2316                 nBandsIn, panBandList, pfnProgress, pProgressData );
2317         }
2318         return CleanOverviews();
2319     }
2320 
2321     // First destroy old overviews
2322     if( CE_None != CleanOverviews() )
2323     {
2324         return CE_Failure;
2325     }
2326 
2327     CPLDebug( "RMF", "Build overviews on dataset %d x %d size",
2328               GetRasterXSize(), GetRasterYSize() );
2329 
2330     GDALDataType    eMainType = GetRasterBand(1)->GetRasterDataType();
2331     RMFDataset*     poParent = this;
2332     double          prevOvLevel = 1.0;
2333     for( int n = 0; n != nOverviews; ++n )
2334     {
2335         int nOvLevel = panOverviewList[n];
2336         const int nOXSize = (GetRasterXSize() + nOvLevel - 1) / nOvLevel;
2337         const int nOYSize = (GetRasterYSize() + nOvLevel - 1) / nOvLevel;
2338         CPLDebug( "RMF", "\tCreate overview #%d size %d x %d",
2339                   nOvLevel, nOXSize, nOYSize );
2340 
2341         RMFDataset*    poOvrDataset;
2342         poOvrDataset = static_cast<RMFDataset*>(
2343                     RMFDataset::Create( nullptr, nOXSize, nOYSize,
2344                                         GetRasterCount(), eMainType,
2345                                         nullptr, poParent, nOvLevel / prevOvLevel ) );
2346 
2347         if( poOvrDataset == nullptr )
2348         {
2349             CPLError( CE_Failure, CPLE_AppDefined,
2350                       "Can't create overview dataset #%d size %d x %d",
2351                       nOvLevel, nOXSize, nOYSize );
2352             return CE_Failure;
2353         }
2354 
2355         prevOvLevel = nOvLevel;
2356         poParent = poOvrDataset;
2357         poOvrDatasets.push_back( poOvrDataset );
2358     }
2359 
2360     GDALRasterBand ***papapoOverviewBands =
2361         static_cast<GDALRasterBand ***>(CPLCalloc(sizeof(void*),nBandsIn));
2362     GDALRasterBand **papoBandList =
2363         static_cast<GDALRasterBand **>(CPLCalloc(sizeof(void*),nBandsIn));
2364 
2365     for( int iBand = 0; iBand < nBandsIn; ++iBand )
2366     {
2367         GDALRasterBand* poBand = GetRasterBand( panBandList[iBand] );
2368 
2369         papoBandList[iBand] = poBand;
2370         papapoOverviewBands[iBand] =
2371             static_cast<GDALRasterBand **>( CPLCalloc(
2372                 sizeof(void*), poBand->GetOverviewCount()) );
2373 
2374         for( int i = 0; i < nOverviews; ++i )
2375         {
2376             papapoOverviewBands[iBand][i] = poBand->GetOverview( i );
2377         }
2378     }
2379 #ifdef DEBUG
2380     for( int iBand = 0; iBand < nBandsIn; ++iBand )
2381     {
2382         CPLDebug( "RMF",
2383                   "Try to create overview for #%d size %d x %d",
2384                   iBand + 1,
2385                   papoBandList[iBand]->GetXSize(),
2386                   papoBandList[iBand]->GetYSize() );
2387         for( int i = 0; i < nOverviews; ++i )
2388         {
2389             CPLDebug( "RMF",
2390                       "\t%d x %d",
2391                       papapoOverviewBands[iBand][i]->GetXSize(),
2392                       papapoOverviewBands[iBand][i]->GetYSize() );
2393         }
2394     }
2395 #endif //DEBUG
2396     CPLErr  res;
2397     res = GDALRegenerateOverviewsMultiBand( nBandsIn, papoBandList,
2398                                             nOverviews, papapoOverviewBands,
2399                                             pszResampling, pfnProgress,
2400                                             pProgressData );
2401 
2402     for( int iBand = 0; iBand < nBandsIn; ++iBand )
2403     {
2404         CPLFree(papapoOverviewBands[iBand]);
2405     }
2406 
2407     CPLFree(papapoOverviewBands);
2408     CPLFree(papoBandList);
2409 
2410     return res;
2411 }
2412 
IRasterIO(GDALRWFlag eRWFlag,int nXOff,int nYOff,int nXSize,int nYSize,void * pData,int nBufXSize,int nBufYSize,GDALDataType eBufType,int nBandCount,int * panBandMap,GSpacing nPixelSpace,GSpacing nLineSpace,GSpacing nBandSpace,GDALRasterIOExtraArg * psExtraArg)2413 CPLErr RMFDataset::IRasterIO(GDALRWFlag eRWFlag,
2414                              int nXOff, int nYOff, int nXSize, int nYSize,
2415                              void * pData, int nBufXSize, int nBufYSize,
2416                              GDALDataType eBufType, int nBandCount,
2417                              int *panBandMap, GSpacing nPixelSpace,
2418                              GSpacing nLineSpace, GSpacing nBandSpace,
2419                              GDALRasterIOExtraArg* psExtraArg)
2420 {
2421 #ifdef DEBUG
2422     CPLDebug("RMF", "Dataset %p, %s %d %d %d %d, %d %d",
2423              this, (eRWFlag == GF_Read ? "Read" : "Write"),
2424              nXOff, nYOff, nXSize, nYSize, nBufXSize, nBufYSize);
2425 #endif //DEBUG
2426     if(eRWFlag == GF_Read &&
2427        poCompressData != nullptr &&
2428        poCompressData->oThreadPool.GetThreadCount() > 0)
2429     {
2430         poCompressData->oThreadPool.WaitCompletion();
2431     }
2432 
2433     return GDALDataset::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
2434                                   pData, nBufXSize, nBufYSize, eBufType,
2435                                   nBandCount, panBandMap,
2436                                   nPixelSpace, nLineSpace, nBandSpace,
2437                                   psExtraArg);
2438 }
2439 
GetLastOffset() const2440 vsi_l_offset RMFDataset::GetLastOffset() const
2441 {
2442     vsi_l_offset    nLastTileOff = 0;
2443     GUInt32         nTiles( sHeader.nTileTblSize/sizeof(GUInt32) );
2444 
2445     for( GUInt32 n = 0; n < nTiles; n += 2 )
2446     {
2447         vsi_l_offset nTileOffset = GetFileOffset( paiTiles[n] );
2448         GUInt32      nTileBytes = paiTiles[n + 1];
2449         nLastTileOff = std::max( nLastTileOff, nTileOffset + nTileBytes );
2450     }
2451 
2452     nLastTileOff = std::max( nLastTileOff,
2453                              GetFileOffset( sHeader.nROIOffset ) +
2454                              sHeader.nROISize );
2455     nLastTileOff = std::max( nLastTileOff,
2456                              GetFileOffset( sHeader.nClrTblOffset ) +
2457                              sHeader.nClrTblSize );
2458     nLastTileOff = std::max( nLastTileOff,
2459                              GetFileOffset( sHeader.nTileTblOffset ) +
2460                              sHeader.nTileTblSize );
2461     nLastTileOff = std::max( nLastTileOff,
2462                              GetFileOffset( sHeader.nFlagsTblOffset ) +
2463                              sHeader.nFlagsTblSize );
2464     nLastTileOff = std::max( nLastTileOff,
2465                              GetFileOffset( sHeader.nExtHdrOffset ) +
2466                              sHeader.nExtHdrSize );
2467     return nLastTileOff;
2468 }
2469 
CleanOverviews()2470 CPLErr RMFDataset::CleanOverviews()
2471 {
2472     if( sHeader.nOvrOffset == 0 )
2473     {
2474         return CE_None;
2475     }
2476 
2477     if( GetAccess() != GA_Update )
2478     {
2479         CPLError( CE_Failure, CPLE_NotSupported,
2480                   "File open for read-only accessing, "
2481                   "overviews cleanup failed." );
2482         return CE_Failure;
2483     }
2484 
2485     if( poParentDS != nullptr )
2486     {
2487         CPLError( CE_Failure, CPLE_NotSupported,
2488                   "Overviews cleanup for non-root dataset is not possible." );
2489         return CE_Failure;
2490     }
2491 
2492     for( size_t n = 0; n != poOvrDatasets.size(); ++n )
2493     {
2494         GDALClose( poOvrDatasets[n] );
2495     }
2496     poOvrDatasets.clear();
2497 
2498     vsi_l_offset    nLastTileOff = GetLastOffset();
2499 
2500     if( 0 != VSIFSeekL( fp, 0, SEEK_END ) )
2501     {
2502         CPLError( CE_Failure, CPLE_FileIO,
2503                   "Failed to seek to end of file, "
2504                   "overviews cleanup failed." );
2505     }
2506 
2507     vsi_l_offset    nFileSize = VSIFTellL( fp );
2508     if( nFileSize < nLastTileOff )
2509     {
2510         CPLError( CE_Failure, CPLE_FileIO,
2511                   "Invalid file offset, "
2512                   "overviews cleanup failed." );
2513         return CE_Failure;
2514     }
2515 
2516     CPLDebug( "RMF", "Truncate to " CPL_FRMT_GUIB, nLastTileOff );
2517     CPLDebug( "RMF", "File size:  " CPL_FRMT_GUIB, nFileSize );
2518 
2519     if( 0 != VSIFTruncateL( fp, nLastTileOff ) )
2520     {
2521         CPLError( CE_Failure, CPLE_FileIO,
2522                   "Failed to truncate file, "
2523                   "overviews cleanup failed." );
2524         return CE_Failure;
2525     }
2526 
2527     sHeader.nOvrOffset = 0;
2528     bHeaderDirty = true;
2529 
2530     return CE_None;
2531 }
2532 
2533 /************************************************************************/
2534 /*                         GetCompressionType()                         */
2535 /************************************************************************/
2536 
GetCompressionType(const char * pszCompressName)2537 GByte RMFDataset::GetCompressionType(const char* pszCompressName)
2538 {
2539     if(pszCompressName == nullptr ||
2540        EQUAL(pszCompressName, "NONE"))
2541     {
2542         return RMF_COMPRESSION_NONE;
2543     }
2544     else if(EQUAL(pszCompressName, "LZW"))
2545     {
2546         return RMF_COMPRESSION_LZW;
2547     }
2548     else if(EQUAL(pszCompressName, "JPEG"))
2549     {
2550         return RMF_COMPRESSION_JPEG;
2551     }
2552     else if(EQUAL(pszCompressName, "RMF_DEM"))
2553     {
2554         return RMF_COMPRESSION_DEM;
2555     }
2556 
2557     CPLError(CE_Failure, CPLE_AppDefined,
2558              "RMF: Unknown compression scheme <%s>.\n"
2559              "Defaults to NONE compression.",
2560               pszCompressName);
2561     return RMF_COMPRESSION_NONE;
2562 }
2563 
2564 /************************************************************************/
2565 /*                        SetupCompression()                            */
2566 /************************************************************************/
2567 
SetupCompression(GDALDataType eType,const char * pszFilename)2568 int RMFDataset::SetupCompression(GDALDataType eType, const char* pszFilename)
2569 {
2570 /* -------------------------------------------------------------------- */
2571 /*  XXX: The DEM compression method seems to be only applicable         */
2572 /*  to Int32 data.                                                      */
2573 /* -------------------------------------------------------------------- */
2574     if( sHeader.iCompression == RMF_COMPRESSION_NONE )
2575     {
2576         Decompress = nullptr;
2577         Compress = nullptr;
2578     }
2579     else if( sHeader.iCompression == RMF_COMPRESSION_LZW )
2580     {
2581         Decompress = &LZWDecompress;
2582         Compress = &LZWCompress;
2583         SetMetadataItem("COMPRESSION", "LZW", "IMAGE_STRUCTURE");
2584     }
2585     else if(sHeader.iCompression == RMF_COMPRESSION_JPEG)
2586     {
2587         if(eType != GDT_Byte || nBands != RMF_JPEG_BAND_COUNT ||
2588            sHeader.nBitDepth != 24)
2589         {
2590             CPLError(CE_Failure, CPLE_AppDefined,
2591                     "RMF support only 24 bpp JPEG compressed files.");
2592             return CE_Failure;
2593         }
2594 #ifdef HAVE_LIBJPEG
2595         CPLString   oBuf;
2596         oBuf.Printf("%d", (int)sHeader.iJpegQuality);
2597         Decompress = &JPEGDecompress;
2598         Compress = &JPEGCompress;
2599         SetMetadataItem("JPEG_QUALITY", oBuf.c_str(), "IMAGE_STRUCTURE");
2600         SetMetadataItem("COMPRESSION", "JPEG", "IMAGE_STRUCTURE");
2601 #else //HAVE_LIBJPEG
2602          CPLError(CE_Failure, CPLE_AppDefined,
2603              "JPEG codec is needed to open <%s>.\n"
2604              "Please rebuild GDAL with libjpeg support.",
2605              pszFilename);
2606         return CE_Failure;
2607 #endif //HAVE_LIBJPEG
2608     }
2609     else if( sHeader.iCompression == RMF_COMPRESSION_DEM
2610              && eType == GDT_Int32 && nBands == RMF_DEM_BAND_COUNT )
2611     {
2612         Decompress = &DEMDecompress;
2613         Compress = &DEMCompress;
2614         SetMetadataItem("COMPRESSION", "RMF_DEM", "IMAGE_STRUCTURE");
2615     }
2616     else
2617     {
2618          CPLError(CE_Failure, CPLE_AppDefined,
2619              "Unknown compression #%d at file <%s>.",
2620              (int)sHeader.iCompression, pszFilename);
2621         return CE_Failure;
2622     }
2623 
2624     return CE_None;
2625 }
2626 
WriteTileJobFunc(void * pData)2627 void RMFDataset::WriteTileJobFunc(void* pData)
2628 {
2629     RMFCompressionJob* psJob = static_cast<RMFCompressionJob*>(pData);
2630     RMFDataset*        poDS = psJob->poDS;
2631 
2632     GByte*  pabyTileData;
2633     size_t  nTileSize;
2634 
2635     if(poDS->Compress)
2636     {
2637         // RMF doesn't store compressed tiles with size greater than 80% of
2638         // uncompressed size
2639         GUInt32 nMaxCompressedTileSize =
2640                     static_cast<GUInt32>((psJob->nUncompressedBytes*8)/10);
2641         size_t nCompressedBytes =
2642                 poDS->Compress(psJob->pabyUncompressedData,
2643                                static_cast<GUInt32>(psJob->nUncompressedBytes),
2644                                psJob->pabyCompressedData,
2645                                nMaxCompressedTileSize,
2646                                psJob->nXSize, psJob->nYSize, poDS);
2647         if(nCompressedBytes == 0)
2648         {
2649             pabyTileData = psJob->pabyUncompressedData;
2650             nTileSize = psJob->nUncompressedBytes;
2651         }
2652         else
2653         {
2654             pabyTileData = psJob->pabyCompressedData;
2655             nTileSize = nCompressedBytes;
2656         }
2657     }
2658     else
2659     {
2660         pabyTileData = psJob->pabyUncompressedData;
2661         nTileSize = psJob->nUncompressedBytes;
2662     }
2663 
2664     {
2665         CPLMutexHolder oHolder(poDS->poCompressData->hWriteTileMutex);
2666         psJob->eResult =
2667                 poDS->WriteRawTile(psJob->nBlockXOff, psJob->nBlockYOff,
2668                                    pabyTileData, nTileSize);
2669     }
2670     if(poDS->poCompressData->oThreadPool.GetThreadCount() > 0)
2671     {
2672         CPLMutexHolder  oHolder(poDS->poCompressData->hReadyJobMutex);
2673         poDS->poCompressData->asReadyJobs.push_back(psJob);
2674     }
2675 }
2676 
InitCompressorData(char ** papszParamList)2677 CPLErr RMFDataset::InitCompressorData(char **papszParamList)
2678 {
2679     const char* pszNumThreads = CSLFetchNameValue(papszParamList, "NUM_THREADS");
2680     if(pszNumThreads == nullptr)
2681         pszNumThreads = CPLGetConfigOption("GDAL_NUM_THREADS", nullptr);
2682 
2683     int nThreads = 0;
2684     if(pszNumThreads != nullptr)
2685     {
2686         nThreads = EQUAL(pszNumThreads, "ALL_CPUS") ?
2687                             CPLGetNumCPUs() :
2688                             atoi(pszNumThreads);
2689     }
2690 
2691     if(nThreads < 0)
2692     {
2693         nThreads = 0;
2694     }
2695     if( nThreads > 1024 )
2696     {
2697         nThreads = 1024;
2698     }
2699 
2700     poCompressData = std::make_shared<RMFCompressData>();
2701     if(nThreads > 0)
2702     {
2703         if(!poCompressData->oThreadPool.Setup(nThreads, nullptr, nullptr))
2704         {
2705             CPLError(CE_Failure, CPLE_AppDefined,
2706                      "Can't setup %d compressor threads", nThreads);
2707             return CE_Failure;
2708         }
2709     }
2710 
2711     poCompressData->asJobs.resize(nThreads + 1);
2712 
2713     size_t nMaxTileBytes = sHeader.nTileWidth * sHeader.nTileHeight *
2714                            sHeader.nBitDepth / 8;
2715     size_t nCompressBufferSize =
2716                            2 * nMaxTileBytes * poCompressData->asJobs.size();
2717     poCompressData->pabyBuffers =
2718                     reinterpret_cast<GByte*>(VSIMalloc(nCompressBufferSize));
2719 
2720     CPLDebug("RMF", "Setup %d compressor threads and allocate %lu bytes buffer",
2721              nThreads, (long int)nCompressBufferSize);
2722     if(poCompressData->pabyBuffers == nullptr)
2723     {
2724         CPLError( CE_Failure, CPLE_OutOfMemory,
2725                 "Can't allocate compress buffer of size %lu.",
2726                 (unsigned long) nCompressBufferSize);
2727         return CE_Failure;
2728     }
2729 
2730     for(size_t i = 0; i != poCompressData->asJobs.size(); ++i)
2731     {
2732         RMFCompressionJob& sJob(poCompressData->asJobs[i]);
2733         sJob.pabyCompressedData = poCompressData->pabyBuffers + 2*i*nMaxTileBytes;
2734         sJob.pabyUncompressedData = sJob.pabyCompressedData + nMaxTileBytes;
2735         poCompressData->asReadyJobs.push_back(&sJob);
2736     }
2737 
2738     if(nThreads > 0)
2739     {
2740         poCompressData->hReadyJobMutex = CPLCreateMutex();
2741         CPLReleaseMutex(poCompressData->hReadyJobMutex);
2742         poCompressData->hWriteTileMutex = CPLCreateMutex();
2743         CPLReleaseMutex(poCompressData->hWriteTileMutex);
2744     }
2745 
2746     return CE_None;
2747 }
2748 
WriteTile(int nBlockXOff,int nBlockYOff,GByte * pabyData,size_t nBytes,GUInt32 nRawXSize,GUInt32 nRawYSize)2749 CPLErr RMFDataset::WriteTile(int nBlockXOff, int nBlockYOff,
2750                              GByte* pabyData, size_t nBytes,
2751                              GUInt32 nRawXSize, GUInt32 nRawYSize)
2752 {
2753     RMFCompressionJob* poJob = nullptr;
2754     if(poCompressData == nullptr)
2755     {
2756         CPLError(CE_Failure, CPLE_AppDefined, "RMF: Compress data is null");
2757         return CE_Failure;
2758     }
2759 
2760     if(poCompressData->oThreadPool.GetThreadCount() > 0)
2761     {
2762         size_t             nJobs(poCompressData->asJobs.size());
2763 
2764         poCompressData->oThreadPool.WaitCompletion(
2765                                 static_cast<int>(nJobs - 1) );
2766 
2767         CPLMutexHolder oHolder(poCompressData->hReadyJobMutex);
2768         CPLAssert(!poCompressData->asReadyJobs.empty());
2769         poJob = poCompressData->asReadyJobs.front();
2770         poCompressData->asReadyJobs.pop_front();
2771     }
2772     else
2773     {
2774         poJob = poCompressData->asReadyJobs.front();
2775     }
2776 
2777     if(poJob->eResult != CE_None)
2778     {
2779         //One of the previous jobs is not done.
2780         //Detailed debug message is already emitted from WriteRawTile
2781         return poJob->eResult;
2782     }
2783     poJob->poDS = this;
2784     poJob->eResult = CE_Failure;
2785     poJob->nBlockXOff = nBlockXOff;
2786     poJob->nBlockYOff = nBlockYOff;
2787     poJob->nUncompressedBytes = nBytes;
2788     poJob->nXSize = nRawXSize;
2789     poJob->nYSize = nRawYSize;
2790 
2791     memcpy(poJob->pabyUncompressedData, pabyData, nBytes);
2792 
2793     if(poCompressData->oThreadPool.GetThreadCount() > 0)
2794     {
2795         if(!poCompressData->oThreadPool.SubmitJob(WriteTileJobFunc,
2796                                                   poJob))
2797         {
2798             CPLError(CE_Failure, CPLE_NotSupported,
2799                      "Can't submit job to thread pool.");
2800             return CE_Failure;
2801         }
2802     }
2803     else
2804     {
2805         WriteTileJobFunc(poJob);
2806         if(poJob->eResult != CE_None)
2807         {
2808             return poJob->eResult;
2809         }
2810     }
2811 
2812     return CE_None;
2813 }
2814 
WriteRawTile(int nBlockXOff,int nBlockYOff,GByte * pabyData,size_t nTileBytes)2815 CPLErr RMFDataset::WriteRawTile(int nBlockXOff, int nBlockYOff,
2816                                 GByte* pabyData, size_t nTileBytes)
2817 {
2818     CPLAssert(nBlockXOff >= 0
2819               && nBlockYOff >= 0
2820               && pabyData != nullptr
2821               && nTileBytes > 0);
2822 
2823     const GUInt32 nTile = nBlockYOff * nXTiles + nBlockXOff;
2824 
2825     vsi_l_offset nTileOffset = GetFileOffset( paiTiles[2 * nTile] );
2826     size_t       nTileSize = static_cast<size_t>(paiTiles[2 * nTile + 1]);
2827 
2828     if(nTileOffset && nTileSize <= nTileBytes)
2829     {
2830         if( VSIFSeekL( fp, nTileOffset, SEEK_SET ) < 0 )
2831         {
2832             CPLError( CE_Failure, CPLE_FileIO,
2833                 "Can't seek to offset %ld in output file to write data.\n%s",
2834                       static_cast<long>( nTileOffset ),
2835                       VSIStrerror( errno ) );
2836             return CE_Failure;
2837         }
2838     }
2839     else
2840     {
2841         if( VSIFSeekL( fp, 0, SEEK_END ) < 0 )
2842         {
2843             CPLError( CE_Failure, CPLE_FileIO,
2844                 "Can't seek to offset %ld in output file to write data.\n%s",
2845                       static_cast<long>( nTileOffset ),
2846                       VSIStrerror( errno ) );
2847             return CE_Failure;
2848         }
2849         nTileOffset = VSIFTellL( fp );
2850         vsi_l_offset nNewTileOffset = 0;
2851         paiTiles[2 * nTile] = GetRMFOffset( nTileOffset, &nNewTileOffset );
2852 
2853         if( nTileOffset != nNewTileOffset )
2854         {
2855             if( VSIFSeekL( fp, nNewTileOffset, SEEK_SET ) < 0 )
2856             {
2857                 CPLError( CE_Failure, CPLE_FileIO,
2858                           "Can't seek to offset %ld in output file to "
2859                           "write data.\n%s",
2860                           static_cast<long>( nNewTileOffset ),
2861                           VSIStrerror( errno ) );
2862                 return CE_Failure;
2863             }
2864         }
2865         bHeaderDirty = true;
2866     }
2867 
2868 #ifdef CPL_MSB
2869     // Compressed tiles are already with proper byte order
2870     if(eRMFType == RMFT_MTW &&
2871        sHeader.iCompression == RMF_COMPRESSION_NONE)
2872     {
2873         //Byte swap can be done in place
2874         if( sHeader.nBitDepth == 16 )
2875         {
2876             for( size_t i = 0; i < nTileBytes; i += 2 )
2877                 CPL_SWAP16PTR( pabyData + i );
2878         }
2879         else if( sHeader.nBitDepth == 32 )
2880         {
2881             for( size_t i = 0; i < nTileBytes; i += 4 )
2882                 CPL_SWAP32PTR( pabyData + i );
2883         }
2884         else if( sHeader.nBitDepth == 64 )
2885         {
2886             for( size_t i = 0; i < nTileBytes; i += 8 )
2887                 CPL_SWAPDOUBLE( pabyData + i );
2888         }
2889     }
2890 #endif
2891 
2892     bool bOk = (VSIFWriteL(pabyData, 1, nTileBytes, fp) == nTileBytes);
2893 
2894     if(!bOk)
2895     {
2896         CPLError(CE_Failure, CPLE_FileIO,
2897                  "Can't write tile with X offset %d and Y offset %d.\n%s",
2898                  nBlockXOff, nBlockYOff, VSIStrerror(errno));
2899         return CE_Failure;
2900     }
2901 
2902     paiTiles[2 * nTile + 1] = static_cast<GUInt32>(nTileBytes);
2903     bHeaderDirty = true;
2904 
2905     return CE_None;
2906 }
2907 
ReadTile(int nBlockXOff,int nBlockYOff,GByte * pabyData,size_t nRawBytes,GUInt32 nRawXSize,GUInt32 nRawYSize,bool & bNullTile)2908 CPLErr RMFDataset::ReadTile(int nBlockXOff, int nBlockYOff,
2909                             GByte* pabyData, size_t nRawBytes,
2910                             GUInt32 nRawXSize, GUInt32 nRawYSize,
2911                             bool& bNullTile)
2912 {
2913     bNullTile = false;
2914 
2915     const GUInt32 nTile = nBlockYOff * nXTiles + nBlockXOff;
2916     if(2 * nTile + 1 >= sHeader.nTileTblSize / sizeof(GUInt32))
2917     {
2918         return CE_Failure;
2919     }
2920     vsi_l_offset nTileOffset = GetFileOffset(paiTiles[2 * nTile]);
2921     GUInt32 nTileBytes = paiTiles[2 * nTile + 1];
2922     // RMF doesn't store compressed tiles with size greater than 80% of
2923     // uncompressed size. But just in case, select twice as many.
2924     GUInt32 nMaxTileBytes = 2 * sHeader.nTileWidth * sHeader.nTileHeight *
2925                             sHeader.nBitDepth / 8;
2926 
2927     if(nTileBytes >= nMaxTileBytes)
2928     {
2929         CPLError(CE_Failure, CPLE_AppDefined,
2930                  "Invalid tile size %lu at offset %ld. Must be less than %lu",
2931                  static_cast<unsigned long>(nTileBytes),
2932                  static_cast<long>(nTileOffset),
2933                  static_cast<unsigned long>(nMaxTileBytes));
2934         return CE_Failure;
2935     }
2936 
2937 
2938     if(nTileOffset == 0)
2939     {
2940         bNullTile = true;
2941         return CE_None;
2942     }
2943 
2944 #ifdef DEBUG
2945     CPLDebug("RMF", "Read RawSize [%d, %d], nTileBytes %d, nRawBytes %d",
2946              nRawXSize, nRawYSize,
2947              static_cast<int>(nTileBytes),
2948              static_cast<int>(nRawBytes));
2949 #endif // DEBUG
2950 
2951     if(VSIFSeekL(fp, nTileOffset, SEEK_SET) < 0)
2952     {
2953         // XXX: We will not report error here, because file just may be
2954         // in update state and data for this block will be available later
2955         if(eAccess == GA_Update)
2956             return CE_None;
2957 
2958         CPLError(CE_Failure, CPLE_FileIO,
2959                  "Can't seek to offset %ld in input file to read data.\n%s",
2960                  static_cast<long>(nTileOffset),
2961                  VSIStrerror( errno ));
2962         return CE_Failure;
2963     }
2964 
2965     if(Decompress == nullptr ||
2966        nTileBytes == nRawBytes)
2967     {
2968         if(nTileBytes != nRawBytes)
2969         {
2970             CPLError(CE_Failure, CPLE_AppDefined,
2971                      "RMF: Invalid tile size %lu, expected %lu",
2972                      static_cast<unsigned long>(nTileBytes),
2973                      static_cast<unsigned long>(nRawBytes));
2974             return CE_Failure;
2975         }
2976 
2977         if(VSIFReadL(pabyData, 1, nRawBytes, fp) < nRawBytes)
2978         {
2979             CPLError( CE_Failure, CPLE_FileIO,
2980                       "RMF: Can't read at offset %lu from input file.\n%s",
2981                       static_cast<unsigned long>(nTileOffset),
2982                       VSIStrerror(errno));
2983             return CE_Failure;
2984         }
2985 
2986 #ifdef CPL_MSB
2987     if(eRMFType == RMFT_MTW)
2988     {
2989         if(sHeader.nBitDepth == 16)
2990         {
2991             for(GUInt32 i = 0; i < nRawBytes; i += 2)
2992                 CPL_SWAP16PTR(pabyData + i);
2993         }
2994         else if(sHeader.nBitDepth == 32)
2995         {
2996             for(GUInt32 i = 0; i < nRawBytes; i += 4)
2997                 CPL_SWAP32PTR(pabyData + i);
2998         }
2999         else if(sHeader.nBitDepth == 64)
3000         {
3001             for( GUInt32 i = 0; i < nRawBytes; i += 8)
3002                 CPL_SWAPDOUBLE(pabyData + i);
3003         }
3004     }
3005 #endif
3006         return CE_None;
3007     }
3008 
3009     if(pabyDecompressBuffer == nullptr)
3010     {
3011         pabyDecompressBuffer =
3012            reinterpret_cast<GByte*>(VSIMalloc(std::max(1U, nMaxTileBytes)));
3013         if(!pabyDecompressBuffer)
3014         {
3015             CPLError(CE_Failure, CPLE_OutOfMemory,
3016                      "Can't allocate decompress buffer of size %lu.\n%s",
3017                      static_cast<unsigned long>(nMaxTileBytes),
3018                      VSIStrerror(errno));
3019             return CE_Failure;
3020         }
3021     }
3022 
3023     if(VSIFReadL(pabyDecompressBuffer, 1, nTileBytes, fp) < nTileBytes)
3024     {
3025         CPLError( CE_Failure, CPLE_FileIO,
3026                   "RMF: Can't read at offset %lu from input file.\n%s",
3027                   static_cast<unsigned long>(nTileOffset), VSIStrerror(errno));
3028         return CE_Failure;
3029     }
3030 
3031     size_t  nDecompressedSize =
3032                     Decompress(pabyDecompressBuffer, nTileBytes,
3033                                pabyData, static_cast<GUInt32>(nRawBytes),
3034                                nRawXSize, nRawYSize);
3035 
3036     if(nDecompressedSize != (size_t)nRawBytes)
3037     {
3038         CPLError( CE_Failure, CPLE_FileIO,
3039                   "Can't decompress tile xOff %d yOff %d. "
3040                   "Raw tile size is %lu but decompressed is %lu. "
3041                   "Compressed tile size is %lu",
3042                   nBlockXOff, nBlockYOff,
3043                   static_cast<unsigned long>(nRawBytes),
3044                   static_cast<unsigned long>(nDecompressedSize),
3045                   static_cast<unsigned long>(nTileBytes));
3046         return CE_Failure;
3047     }
3048     // We don't need to swap bytes here,
3049     // because decompressed data is in proper byte order
3050     return CE_None;
3051 }
3052 
SetupNBits()3053 void RMFDataset::SetupNBits()
3054 {
3055     int nBitDepth = 0;
3056     if(sHeader.nBitDepth < 8 && nBands == 1)
3057     {
3058         nBitDepth = static_cast<int>(sHeader.nBitDepth);
3059     }
3060     else
3061     if(sHeader.nBitDepth == 16 && nBands == 3 && eRMFType == RMFT_RSW )
3062     {
3063         nBitDepth = 5;
3064     }
3065 
3066     if(nBitDepth > 0)
3067     {
3068         char szNBits[32] = {};
3069         snprintf(szNBits, sizeof(szNBits), "%d", nBitDepth);
3070         for(int iBand = 1; iBand <= nBands; iBand++)
3071         {
3072             GetRasterBand(iBand)->SetMetadataItem("NBITS", szNBits,
3073                                                   "IMAGE_STRUCTURE");
3074         }
3075     }
3076 }
3077 
3078 /************************************************************************/
3079 /*                        GDALRegister_RMF()                            */
3080 /************************************************************************/
3081 
GDALRegister_RMF()3082 void GDALRegister_RMF()
3083 
3084 {
3085     if( GDALGetDriverByName( "RMF" ) != nullptr )
3086         return;
3087 
3088     GDALDriver *poDriver = new GDALDriver();
3089 
3090     poDriver->SetDescription( "RMF" );
3091     poDriver->SetMetadataItem( GDAL_DCAP_RASTER, "YES" );
3092     poDriver->SetMetadataItem( GDAL_DMD_LONGNAME, "Raster Matrix Format" );
3093     poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC, "drivers/raster/rmf.html" );
3094     poDriver->SetMetadataItem( GDAL_DMD_EXTENSION, "rsw" );
3095     poDriver->SetMetadataItem( GDAL_DMD_CREATIONDATATYPES,
3096                                "Byte Int16 Int32 Float64" );
3097     poDriver->SetMetadataItem( GDAL_DMD_CREATIONOPTIONLIST,
3098 "<CreationOptionList>"
3099 "   <Option name='MTW' type='boolean' description='Create MTW DEM matrix'/>"
3100 "   <Option name='BLOCKXSIZE' type='int' description='Tile Width'/>"
3101 "   <Option name='BLOCKYSIZE' type='int' description='Tile Height'/>"
3102 "   <Option name='RMFHUGE' type='string-select' description='Creation of huge RMF file (Supported by GIS Panorama since v11)'>"
3103 "     <Value>NO</Value>"
3104 "     <Value>YES</Value>"
3105 "     <Value>IF_SAFER</Value>"
3106 "   </Option>"
3107 "   <Option name='COMPRESS' type='string-select' default='NONE'>"
3108 "     <Value>NONE</Value>"
3109 "     <Value>LZW</Value>"
3110 "     <Value>JPEG</Value>"
3111 "     <Value>RMF_DEM</Value>"
3112 "   </Option>"
3113 "   <Option name='JPEG_QUALITY' type='int' description='JPEG quality 1-100' default='75'/>"
3114 "   <Option name='NUM_THREADS' type='string' description='Number of worker threads for compression. Can be set to ALL_CPUS' default='1'/>"
3115 "</CreationOptionList>" );
3116     poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
3117 
3118     poDriver->pfnIdentify = RMFDataset::Identify;
3119     poDriver->pfnOpen = RMFDataset::Open;
3120     poDriver->pfnCreate = RMFDataset::Create;
3121     poDriver->SetMetadataItem( GDAL_DMD_OPENOPTIONLIST,
3122         "<OpenOptionList>"
3123         "  <Option name='RMF_SET_VERTCS' type='string' description='Layers spatial reference will include vertical coordinate system description if exist' default='NO'/>"
3124         "</OpenOptionList>");
3125 
3126     GetGDALDriverManager()->RegisterDriver( poDriver );
3127 }
3128 
3129 /************************************************************************/
3130 /*                            RMFCompressData                           */
3131 /************************************************************************/
3132 
RMFCompressData()3133 RMFCompressData::RMFCompressData():
3134     pabyBuffers(nullptr)
3135 {
3136 }
3137 
~RMFCompressData()3138 RMFCompressData::~RMFCompressData()
3139 {
3140     if(pabyBuffers != nullptr)
3141     {
3142         VSIFree(pabyBuffers);
3143     }
3144 
3145     if(hWriteTileMutex != nullptr)
3146     {
3147         CPLDestroyMutex(hWriteTileMutex);
3148     }
3149 
3150     if(hReadyJobMutex != nullptr)
3151     {
3152         CPLDestroyMutex(hReadyJobMutex);
3153     }
3154 }
3155