1 /******************************************************************************
2  *
3  * Project:  ROI_PAC Raster Reader
4  * Purpose:  Implementation of the ROI_PAC raster reader
5  * Author:   Matthieu Volat (ISTerre), matthieu.volat@ujf-grenoble.fr
6  *
7  ******************************************************************************
8  * Copyright (c) 2014, Matthieu Volat <matthieu.volat@ujf-grenoble.fr>
9  *
10  * Permission is hereby granted, free of charge, to any person obtaining a
11  * copy of this software and associated documentation files (the "Software"),
12  * to deal in the Software without restriction, including without limitation
13  * the rights to use, copy, modify, merge, publish, distribute, sublicense,
14  * and/or sell copies of the Software, and to permit persons to whom the
15  * Software is furnished to do so, subject to the following conditions:
16  *
17  * The above copyright notice and this permission notice shall be included
18  * in all copies or substantial portions of the Software.
19  *
20  * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
21  * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
22  * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
23  * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
24  * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
25  * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
26  * DEALINGS IN THE SOFTWARE.
27  ****************************************************************************/
28 
29 #include "gdal_frmts.h"
30 #include "ogr_spatialref.h"
31 #include "rawdataset.h"
32 
33 CPL_CVSID("$Id: roipacdataset.cpp f6099e5ed704166bf5cc113a053dd1b2725cb391 2020-03-22 11:20:10 +0100 Kai Pastor $")
34 
35 /************************************************************************/
36 /* ==================================================================== */
37 /*                             ROIPACDataset                            */
38 /* ==================================================================== */
39 /************************************************************************/
40 
41 class ROIPACRasterBand;
42 
43 class ROIPACDataset final: public RawDataset
44 {
45     friend class ROIPACRasterBand;
46 
47     VSILFILE    *fpImage;
48     VSILFILE    *fpRsc;
49 
50     char        *pszRscFilename;
51 
52     double      adfGeoTransform[6];
53     bool        bValidGeoTransform;
54     char        *pszProjection;
55 
56     CPL_DISALLOW_COPY_ASSIGN(ROIPACDataset)
57 
58   public:
59     ROIPACDataset();
60     ~ROIPACDataset() override;
61 
62     static GDALDataset *Open( GDALOpenInfo *poOpenInfo );
63     static int          Identify( GDALOpenInfo *poOpenInfo );
64     static GDALDataset *Create( const char *pszFilename,
65                                 int nXSize, int nYSize, int nBands,
66                                 GDALDataType eType, char **papszOptions );
67 
68     void        FlushCache() override;
69     CPLErr              GetGeoTransform( double *padfTransform ) override;
70     CPLErr      SetGeoTransform( double *padfTransform ) override;
71     const char *_GetProjectionRef( void ) override;
72     CPLErr      _SetProjection( const char *pszNewProjection ) override;
GetSpatialRef() const73     const OGRSpatialReference* GetSpatialRef() const override {
74         return GetSpatialRefFromOldGetProjectionRef();
75     }
SetSpatialRef(const OGRSpatialReference * poSRS)76     CPLErr SetSpatialRef(const OGRSpatialReference* poSRS) override {
77         return OldSetProjectionFromSetSpatialRef(poSRS);
78     }
79 
80     char      **GetFileList() override;
81 };
82 
83 /************************************************************************/
84 /* ==================================================================== */
85 /*                           ROIPACRasterBand                           */
86 /* ==================================================================== */
87 /************************************************************************/
88 
89 class ROIPACRasterBand final: public RawRasterBand
90 {
91     CPL_DISALLOW_COPY_ASSIGN(ROIPACRasterBand)
92 
93     public:
94                 ROIPACRasterBand( GDALDataset *poDS, int nBand, VSILFILE *fpRaw,
95                                   vsi_l_offset nImgOffset, int nPixelOffset,
96                                   int nLineOffset,
97                                   GDALDataType eDataType, int bNativeOrder );
98 };
99 
100 /************************************************************************/
101 /*                           getRscFilename()                           */
102 /************************************************************************/
103 
getRscFilename(GDALOpenInfo * poOpenInfo)104 static CPLString getRscFilename( GDALOpenInfo *poOpenInfo )
105 {
106     char **papszSiblingFiles = poOpenInfo->GetSiblingFiles();
107     if ( papszSiblingFiles == nullptr )
108     {
109         const CPLString osRscFilename =
110             CPLFormFilename( nullptr, poOpenInfo->pszFilename,
111                              "rsc" );
112         VSIStatBufL psRscStatBuf;
113         if ( VSIStatL( osRscFilename, &psRscStatBuf ) != 0 )
114         {
115             return "";
116         }
117         return osRscFilename;
118     }
119 
120     /* ------------------------------------------------------------ */
121     /*      We need to tear apart the filename to form a .rsc       */
122     /*      filename.                                               */
123     /* ------------------------------------------------------------ */
124     const CPLString osPath = CPLGetPath( poOpenInfo->pszFilename );
125     const CPLString osName = CPLGetFilename( poOpenInfo->pszFilename );
126 
127     int iFile = CSLFindString( papszSiblingFiles,
128                                CPLFormFilename( nullptr, osName, "rsc" ) );
129     if( iFile >= 0 )
130     {
131         return CPLFormFilename( osPath,
132                                 papszSiblingFiles[iFile],
133                                 nullptr );
134     }
135 
136     return "";
137 }
138 
139 /************************************************************************/
140 /*                            ROIPACDataset()                           */
141 /************************************************************************/
142 
ROIPACDataset()143 ROIPACDataset::ROIPACDataset() :
144     fpImage(nullptr),
145     fpRsc(nullptr),
146     pszRscFilename(nullptr),
147     bValidGeoTransform(false),
148     pszProjection(nullptr)
149 {
150     adfGeoTransform[0] =  0.0;
151     adfGeoTransform[1] =  1.0;
152     adfGeoTransform[2] =  0.0;
153     adfGeoTransform[3] =  0.0;
154     adfGeoTransform[4] =  0.0;
155     adfGeoTransform[5] =  1.0;
156 }
157 
158 /************************************************************************/
159 /*                            ~ROIPACDataset()                          */
160 /************************************************************************/
161 
~ROIPACDataset()162 ROIPACDataset::~ROIPACDataset()
163 {
164     ROIPACDataset::FlushCache();
165     if ( fpRsc != nullptr && VSIFCloseL( fpRsc ) != 0 )
166     {
167         CPLError( CE_Failure, CPLE_FileIO, "I/O error" );
168     }
169     if ( fpImage != nullptr && VSIFCloseL( fpImage ) != 0 )
170     {
171         CPLError( CE_Failure, CPLE_FileIO, "I/O error" );
172     }
173     CPLFree( pszRscFilename );
174     CPLFree( pszProjection );
175 }
176 
177 /************************************************************************/
178 /*                                Open()                                */
179 /************************************************************************/
180 
Open(GDALOpenInfo * poOpenInfo)181 GDALDataset *ROIPACDataset::Open( GDALOpenInfo *poOpenInfo )
182 {
183 /* -------------------------------------------------------------------- */
184 /*      Confirm that the header is compatible with a ROIPAC dataset.    */
185 /* -------------------------------------------------------------------- */
186     if ( !Identify(poOpenInfo) || poOpenInfo->fpL == nullptr )
187     {
188         return nullptr;
189     }
190 
191 /* -------------------------------------------------------------------- */
192 /*      Open the .rsc file                                              */
193 /* -------------------------------------------------------------------- */
194     CPLString osRscFilename = getRscFilename( poOpenInfo );
195     if ( osRscFilename.empty() )
196     {
197         return nullptr;
198     }
199     VSILFILE *fpRsc = nullptr;
200     if ( poOpenInfo->eAccess == GA_Update )
201     {
202         fpRsc = VSIFOpenL( osRscFilename, "r+" );
203     }
204     else
205     {
206         fpRsc = VSIFOpenL( osRscFilename, "r" );
207     }
208     if ( fpRsc == nullptr )
209     {
210         return nullptr;
211     }
212 
213 /* -------------------------------------------------------------------- */
214 /*      Load the .rsc information.                                      */
215 /* -------------------------------------------------------------------- */
216     char **papszRsc = nullptr;
217     while ( true )
218     {
219         const char *pszLine = CPLReadLineL( fpRsc );
220         if (pszLine == nullptr)
221         {
222             break;
223         }
224 
225         char **papszTokens = CSLTokenizeString2( pszLine, " \t",
226                                                  CSLT_STRIPLEADSPACES
227                                                  | CSLT_STRIPENDSPACES
228                                                  | CSLT_PRESERVEQUOTES
229                                                  | CSLT_PRESERVEESCAPES );
230         if ( papszTokens == nullptr
231              || papszTokens[0] == nullptr || papszTokens[1] == nullptr )
232         {
233             CSLDestroy ( papszTokens );
234             break;
235         }
236         papszRsc = CSLSetNameValue( papszRsc,
237                                     papszTokens[0], papszTokens[1] );
238 
239         CSLDestroy ( papszTokens );
240     }
241 
242 /* -------------------------------------------------------------------- */
243 /*      Fetch required fields.                                          */
244 /* -------------------------------------------------------------------- */
245     if( CSLFetchNameValue( papszRsc, "WIDTH" ) == nullptr
246         || CSLFetchNameValue( papszRsc, "FILE_LENGTH" ) == nullptr )
247     {
248         CSLDestroy( papszRsc );
249         CPL_IGNORE_RET_VAL(VSIFCloseL( fpRsc ));
250         return nullptr;
251     }
252     const int nWidth = atoi( CSLFetchNameValue( papszRsc, "WIDTH" ) );
253     const int nFileLength =
254         atoi( CSLFetchNameValue( papszRsc, "FILE_LENGTH" ) );
255 
256     if (!GDALCheckDatasetDimensions(nWidth, nFileLength) )
257     {
258         CSLDestroy( papszRsc );
259         CPL_IGNORE_RET_VAL(VSIFCloseL( fpRsc ));
260         return nullptr;
261     }
262 
263 /* -------------------------------------------------------------------- */
264 /*      Create a corresponding GDALDataset.                             */
265 /* -------------------------------------------------------------------- */
266     ROIPACDataset *poDS = new ROIPACDataset();
267     poDS->nRasterXSize = nWidth;
268     poDS->nRasterYSize = nFileLength;
269     poDS->eAccess = poOpenInfo->eAccess;
270     poDS->fpRsc = fpRsc;
271     poDS->pszRscFilename = CPLStrdup( osRscFilename.c_str() );
272     poDS->fpImage = poOpenInfo->fpL;
273     poOpenInfo->fpL = nullptr;
274 
275 /* -------------------------------------------------------------------- */
276 /*      Create band information objects.                                */
277 /* -------------------------------------------------------------------- */
278     GDALDataType eDataType = GDT_Unknown;
279     int nBands = 0;
280     enum Interleave { UNKNOWN, LINE, PIXEL } eInterleave = UNKNOWN;
281     const char *pszExtension = CPLGetExtension(poOpenInfo->pszFilename);
282     if ( strcmp( pszExtension, "raw" ) == 0 )
283     {
284         /* ------------------------------------------------------------ */
285         /* TODO: ROI_PAC raw images are what would be GDT_CInt8 typed,  */
286         /* but since that type do not exist, we will have to implement  */
287         /* a specific case in the RasterBand to convert it to           */
288         /* GDT_CInt16 for example                                       */
289         /* ------------------------------------------------------------ */
290 #if 0
291         eDataType = GDT_CInt8;
292         nBands = 1;
293         eInterleave = PIXEL;
294 #else
295         CPLError( CE_Failure, CPLE_NotSupported,
296                   "Reading ROI_PAC raw files is not supported yet." );
297         delete poDS;
298         CSLDestroy( papszRsc );
299         return nullptr;
300 #endif
301     }
302     else if ( strcmp( pszExtension, "int" ) == 0
303               || strcmp( pszExtension, "slc" ) == 0 )
304     {
305         eDataType = GDT_CFloat32;
306         nBands = 1;
307         eInterleave = PIXEL;
308     }
309     else if ( strcmp( pszExtension, "amp" ) == 0 )
310     {
311         eDataType = GDT_Float32;
312         nBands = 2;
313         eInterleave = PIXEL;
314     }
315     else if ( strcmp( pszExtension, "cor" ) == 0
316               || strcmp( pszExtension, "hgt" ) == 0
317               || strcmp( pszExtension, "unw" ) == 0
318               || strcmp( pszExtension, "msk" ) == 0
319               || strcmp( pszExtension, "trans" ) == 0 )
320     {
321         eDataType = GDT_Float32;
322         nBands = 2;
323         eInterleave = LINE;
324     }
325     else if ( strcmp( pszExtension, "dem" ) == 0 )
326     {
327         eDataType = GDT_Int16;
328         nBands = 1;
329         eInterleave = PIXEL;
330     }
331     else if ( strcmp( pszExtension, "flg" ) == 0 )
332     {
333         eDataType = GDT_Byte;
334         nBands = 1;
335         eInterleave = PIXEL;
336     }
337     else { /* Eeek */
338         delete poDS;
339         CSLDestroy( papszRsc );
340         return nullptr;
341     }
342 
343     int nPixelOffset = 0;
344     int nLineOffset = 0;
345     vsi_l_offset nBandOffset = 0;
346     const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
347     bool bIntOverflow =  false;
348     if (eInterleave == LINE)
349     {
350         nPixelOffset = nDTSize;
351         if( nWidth > INT_MAX / (nPixelOffset * nBands) )
352             bIntOverflow = true;
353         else
354         {
355             nLineOffset = nPixelOffset * nWidth * nBands;
356             nBandOffset = static_cast<vsi_l_offset>(nDTSize) * nWidth;
357         }
358     }
359     else { /* PIXEL */
360         nPixelOffset = nDTSize * nBands;
361         if( nWidth > INT_MAX / nPixelOffset )
362             bIntOverflow = true;
363         else
364         {
365             nLineOffset = nPixelOffset * nWidth;
366             nBandOffset = nDTSize;
367 
368             if( nBands > 1 )
369             {
370                 // GDAL 2.0.[0-3] and 2.1.0  had a value of nLineOffset that was
371                 // equal to the theoretical nLineOffset multiplied by nBands.
372                 VSIFSeekL( poDS->fpImage, 0, SEEK_END );
373                 const GUIntBig nWrongFileSize =
374                     nDTSize *
375                     nWidth * (static_cast<GUIntBig>(nFileLength - 1) *
376                             nBands * nBands + nBands);
377                 if( VSIFTellL( poDS->fpImage ) == nWrongFileSize )
378                 {
379                     CPLError(CE_Warning, CPLE_AppDefined,
380                             "This file has been incorrectly generated by an older "
381                             "GDAL version whose line offset computation was "
382                             "erroneous.  Taking that into account, "
383                             "but the file should be re-encoded ideally.");
384                     nLineOffset = nLineOffset * nBands;
385                 }
386             }
387         }
388     }
389 
390     if (bIntOverflow)
391     {
392         CPLError( CE_Failure, CPLE_AppDefined,
393                   "Int overflow occurred." );
394         delete poDS;
395         CSLDestroy( papszRsc );
396         return nullptr;
397     }
398 
399 #ifdef CPL_LSB
400     const bool bNativeOrder = true;
401 #else
402     const bool bNativeOrder = false;
403 #endif
404     poDS->nBands = nBands;
405     for( int b = 0; b < nBands; b++ )
406     {
407         poDS->SetBand( b + 1,
408                        new ROIPACRasterBand( poDS, b + 1, poDS->fpImage,
409                                              nBandOffset * b,
410                                              nPixelOffset, nLineOffset,
411                                              eDataType, bNativeOrder ) );
412     }
413 
414 /* -------------------------------------------------------------------- */
415 /*      Interpret georeferencing, if present.                           */
416 /* -------------------------------------------------------------------- */
417     if ( CSLFetchNameValue( papszRsc, "X_FIRST" ) != nullptr
418          && CSLFetchNameValue( papszRsc, "X_STEP" ) != nullptr
419          && CSLFetchNameValue( papszRsc, "Y_FIRST" ) != nullptr
420          && CSLFetchNameValue( papszRsc, "Y_STEP" ) != nullptr )
421     {
422         poDS->adfGeoTransform[0] = CPLAtof( CSLFetchNameValue( papszRsc,
423                                                                "X_FIRST" ) );
424         poDS->adfGeoTransform[1] = CPLAtof( CSLFetchNameValue( papszRsc,
425                                                                "X_STEP" ) );
426         poDS->adfGeoTransform[2] = 0.0;
427         poDS->adfGeoTransform[3] = CPLAtof( CSLFetchNameValue( papszRsc,
428                                                                "Y_FIRST" ) );
429         poDS->adfGeoTransform[4] = 0.0;
430         poDS->adfGeoTransform[5] = CPLAtof( CSLFetchNameValue( papszRsc,
431                                                                "Y_STEP" ) );
432         poDS->bValidGeoTransform = true;
433     }
434     if ( CSLFetchNameValue( papszRsc, "PROJECTION" ) != nullptr )
435     {
436         /* ------------------------------------------------------------ */
437         /* In ROI_PAC, images are georeferenced either with lat/long or */
438         /* UTM projection. However, using UTM projection is dangerous   */
439         /* because there is no North/South field, or use of latitude    */
440         /* bands!                                                       */
441         /* ------------------------------------------------------------ */
442         OGRSpatialReference oSRS;
443         if ( strcmp( CSLFetchNameValue( papszRsc, "PROJECTION" ), "LL" ) == 0 )
444         {
445             if ( CSLFetchNameValue( papszRsc, "DATUM" ) != nullptr )
446             {
447                 oSRS.SetWellKnownGeogCS( CSLFetchNameValue( papszRsc,
448                                                             "DATUM" ) );
449             }
450             else {
451                 oSRS.SetWellKnownGeogCS( "WGS84" );
452             }
453         }
454         else if( STARTS_WITH(CSLFetchNameValue( papszRsc, "PROJECTION" ),
455                              "UTM") )
456         {
457             const char *pszZone = CSLFetchNameValue( papszRsc,
458                                                      "PROJECTION" ) + 3;
459             oSRS.SetUTM( atoi( pszZone ), TRUE ); /* FIXME: north/south? */
460             if ( CSLFetchNameValue( papszRsc, "DATUM" ) != nullptr )
461             {
462                 oSRS.SetWellKnownGeogCS( CSLFetchNameValue( papszRsc,
463                                                             "DATUM" ) );
464             }
465             else {
466                 oSRS.SetWellKnownGeogCS( "NAD27" );
467             }
468         }
469         oSRS.exportToWkt( &poDS->pszProjection );
470     }
471     if ( CSLFetchNameValue( papszRsc, "Z_OFFSET" ) != nullptr )
472     {
473         const double dfOffset =
474             strtod( CSLFetchNameValue( papszRsc, "Z_OFFSET" ), nullptr);
475         for( int b = 1; b <= nBands; b++ )
476         {
477             GDALRasterBand *poBand = poDS->GetRasterBand(b);
478             poBand->SetOffset( dfOffset );
479         }
480     }
481     if ( CSLFetchNameValue( papszRsc, "Z_SCALE" ) != nullptr )
482     {
483         const double dfScale =
484             strtod( CSLFetchNameValue( papszRsc, "Z_SCALE" ), nullptr);
485         for( int b = 1; b <= nBands; b++ )
486         {
487             GDALRasterBand *poBand = poDS->GetRasterBand(b);
488             poBand->SetScale( dfScale );
489         }
490     }
491 
492 /* -------------------------------------------------------------------- */
493 /*      Set all the other header metadata into the ROI_PAC domain       */
494 /* -------------------------------------------------------------------- */
495     for( int i = 0; papszRsc != nullptr && papszRsc[i] != nullptr; i++ )
496     {
497         char **papszTokens = CSLTokenizeString2( papszRsc[i],
498                                                  "=",
499                                                  CSLT_STRIPLEADSPACES
500                                                  | CSLT_STRIPENDSPACES);
501         if ( CSLCount(papszTokens) < 2
502               || strcmp( papszTokens[0], "WIDTH" ) == 0
503               || strcmp( papszTokens[0], "FILE_LENGTH" ) == 0
504               || strcmp( papszTokens[0], "X_FIRST" ) == 0
505               || strcmp( papszTokens[0], "X_STEP" ) == 0
506               || strcmp( papszTokens[0], "Y_FIRST" ) == 0
507               || strcmp( papszTokens[0], "Y_STEP" ) == 0
508               || strcmp( papszTokens[0], "PROJECTION" ) == 0
509               || strcmp( papszTokens[0], "DATUM" ) == 0
510               || strcmp( papszTokens[0], "Z_OFFSET" ) == 0
511               || strcmp( papszTokens[0], "Z_SCALE" ) == 0 )
512         {
513             CSLDestroy( papszTokens );
514             continue;
515         }
516         poDS->SetMetadataItem( papszTokens[0], papszTokens[1], "ROI_PAC" );
517         CSLDestroy( papszTokens );
518     }
519 
520 /* -------------------------------------------------------------------- */
521 /*      Free papszRsc                                                   */
522 /* -------------------------------------------------------------------- */
523     CSLDestroy( papszRsc );
524 
525 /* -------------------------------------------------------------------- */
526 /*      Initialize any PAM information.                                 */
527 /* -------------------------------------------------------------------- */
528     poDS->SetDescription( poOpenInfo->pszFilename );
529     poDS->TryLoadXML();
530 
531 /* -------------------------------------------------------------------- */
532 /*      Check for overviews.                                            */
533 /* -------------------------------------------------------------------- */
534     poDS->oOvManager.Initialize( poDS, poOpenInfo->pszFilename );
535 
536     return poDS;
537 }
538 
539 /************************************************************************/
540 /*                             Identify()                               */
541 /************************************************************************/
542 
Identify(GDALOpenInfo * poOpenInfo)543 int ROIPACDataset::Identify( GDALOpenInfo *poOpenInfo )
544 {
545 /* -------------------------------------------------------------------- */
546 /*      Check if:                                                       */
547 /*      * 1. The data file extension is known                           */
548 /* -------------------------------------------------------------------- */
549     const char *pszExtension = CPLGetExtension( poOpenInfo->pszFilename );
550     if ( strcmp( pszExtension, "raw" ) == 0 )
551     {
552         /* Since gdal do not read natively CInt8, more work is needed
553          * to read raw files */
554         return false;
555     }
556     const bool bExtensionIsValid =
557         strcmp( pszExtension, "int" ) == 0
558         || strcmp( pszExtension, "slc" ) == 0
559         || strcmp( pszExtension, "amp" ) == 0
560         || strcmp( pszExtension, "cor" ) == 0
561         || strcmp( pszExtension, "hgt" ) == 0
562         || strcmp( pszExtension, "unw" ) == 0
563         || strcmp( pszExtension, "msk" ) == 0
564         || strcmp( pszExtension, "trans" ) == 0
565         || strcmp( pszExtension, "dem" ) == 0
566         || strcmp( pszExtension, "flg" ) == 0;
567     if ( !bExtensionIsValid )
568     {
569         return false;
570     }
571 
572 /* -------------------------------------------------------------------- */
573 /*      * 2. there is a .rsc file                                      */
574 /* -------------------------------------------------------------------- */
575     CPLString osRscFilename = getRscFilename( poOpenInfo );
576     if ( osRscFilename.empty() )
577     {
578         return false;
579     }
580 
581     return true;
582 }
583 
584 /************************************************************************/
585 /*                              Create()                                */
586 /************************************************************************/
587 
Create(const char * pszFilename,int nXSize,int nYSize,int nBands,GDALDataType eType,char **)588 GDALDataset *ROIPACDataset::Create( const char *pszFilename,
589                                     int nXSize, int nYSize, int nBands,
590                                     GDALDataType eType,
591                                     char ** /* papszOptions */ )
592 {
593 /* -------------------------------------------------------------------- */
594 /*      Verify input options.                                           */
595 /* -------------------------------------------------------------------- */
596     const char *pszExtension = CPLGetExtension(pszFilename);
597     if ( strcmp( pszExtension, "int" ) == 0
598                 || strcmp( pszExtension, "slc" ) == 0 )
599     {
600         if ( nBands != 1 || eType != GDT_CFloat32 )
601         {
602             CPLError( CE_Failure, CPLE_AppDefined,
603                       "Attempt to create ROI_PAC %s dataset with an illegal "
604                       "number of bands (%d) and/or data type (%s).",
605                       pszExtension, nBands, GDALGetDataTypeName(eType) );
606             return nullptr;
607         }
608     }
609     else if ( strcmp( pszExtension, "amp" ) == 0 )
610     {
611         if ( nBands != 2 || eType != GDT_Float32 )
612         {
613             CPLError( CE_Failure, CPLE_AppDefined,
614                       "Attempt to create ROI_PAC %s dataset with an illegal "
615                       "number of bands (%d) and/or data type (%s).",
616                       pszExtension, nBands, GDALGetDataTypeName(eType) );
617             return nullptr;
618         }
619     }
620     else if ( strcmp( pszExtension, "cor" ) == 0
621                 || strcmp( pszExtension, "hgt" ) == 0
622                 || strcmp( pszExtension, "unw" ) == 0
623                 || strcmp( pszExtension, "msk" ) == 0
624                 || strcmp( pszExtension, "trans" ) == 0 )
625     {
626         if ( nBands != 2 || eType != GDT_Float32 )
627         {
628             CPLError( CE_Failure, CPLE_AppDefined,
629                       "Attempt to create ROI_PAC %s dataset with an illegal "
630                           "number of bands (%d) and/or data type (%s).",
631                       pszExtension, nBands, GDALGetDataTypeName(eType) );
632             return nullptr;
633         }
634     }
635     else if ( strcmp( pszExtension, "dem" ) == 0 )
636     {
637         if ( nBands != 1 || eType != GDT_Int16 )
638         {
639             CPLError( CE_Failure, CPLE_AppDefined,
640                       "Attempt to create ROI_PAC %s dataset with an illegal "
641                       "number of bands (%d) and/or data type (%s).",
642                       pszExtension, nBands, GDALGetDataTypeName(eType) );
643             return nullptr;
644         }
645     }
646     else if ( strcmp( pszExtension, "flg" ) == 0 )
647     {
648         if ( nBands != 1 || eType != GDT_Byte )
649         {
650             CPLError( CE_Failure, CPLE_AppDefined,
651                       "Attempt to create ROI_PAC %s dataset with an illegal "
652                       "number of bands (%d) and/or data type (%s).",
653                       pszExtension, nBands, GDALGetDataTypeName(eType) );
654             return nullptr;
655         }
656     }
657     else { /* Eeek */
658         CPLError( CE_Failure, CPLE_AppDefined,
659                   "Attempt to create ROI_PAC dataset with an unknown type (%s)",
660                   pszExtension );
661         return nullptr;
662     }
663 
664 /* -------------------------------------------------------------------- */
665 /*      Try to create the file.                                         */
666 /* -------------------------------------------------------------------- */
667     VSILFILE *fp = VSIFOpenL( pszFilename, "wb" );
668     if( fp == nullptr )
669     {
670         CPLError( CE_Failure, CPLE_OpenFailed,
671                   "Attempt to create file `%s' failed.",
672                   pszFilename );
673         return nullptr;
674     }
675 
676 /* -------------------------------------------------------------------- */
677 /*      Just write out a couple of bytes to establish the binary        */
678 /*      file, and then close it.                                        */
679 /* -------------------------------------------------------------------- */
680     CPL_IGNORE_RET_VAL(VSIFWriteL( "\0\0", 2, 1, fp ));
681     CPL_IGNORE_RET_VAL(VSIFCloseL( fp ));
682 
683 /* -------------------------------------------------------------------- */
684 /*      Open the RSC file.                                              */
685 /* -------------------------------------------------------------------- */
686     const char *pszRSCFilename = CPLFormFilename( nullptr, pszFilename, "rsc" );
687     fp = VSIFOpenL( pszRSCFilename, "wt" );
688     if( fp == nullptr )
689     {
690         CPLError( CE_Failure, CPLE_OpenFailed,
691                   "Attempt to create file `%s' failed.",
692                   pszRSCFilename );
693         return nullptr;
694     }
695 
696 /* -------------------------------------------------------------------- */
697 /*      Write out the header.                                           */
698 /* -------------------------------------------------------------------- */
699     CPL_IGNORE_RET_VAL(VSIFPrintfL( fp, "%-40s %d\n", "WIDTH", nXSize ));
700     CPL_IGNORE_RET_VAL(VSIFPrintfL( fp, "%-40s %d\n", "FILE_LENGTH", nYSize ));
701     CPL_IGNORE_RET_VAL(VSIFCloseL( fp ));
702 
703     return reinterpret_cast<GDALDataset *>(
704         GDALOpen( pszFilename, GA_Update ) );
705 }
706 
707 /************************************************************************/
708 /*                             FlushCache()                             */
709 /************************************************************************/
710 
FlushCache(void)711 void ROIPACDataset::FlushCache( void )
712 {
713     RawDataset::FlushCache();
714 
715     GDALRasterBand *band = (GetRasterCount() > 0) ? GetRasterBand(1) : nullptr;
716 
717     if ( eAccess == GA_ReadOnly || band == nullptr )
718         return;
719 
720     // If opening an existing file in Update mode (i.e. "r+") we need to make
721     // sure any existing content is cleared, otherwise the file may contain
722     // trailing content from the previous write.
723     CPL_IGNORE_RET_VAL(VSIFTruncateL( fpRsc, 0 ));
724 
725     CPL_IGNORE_RET_VAL(VSIFSeekL( fpRsc, 0, SEEK_SET ));
726 /* -------------------------------------------------------------------- */
727 /*      Rewrite out the header.                                         */
728 /* -------------------------------------------------------------------- */
729 /* -------------------------------------------------------------------- */
730 /*      Raster dimensions.                                              */
731 /* -------------------------------------------------------------------- */
732     CPL_IGNORE_RET_VAL(
733         VSIFPrintfL( fpRsc, "%-40s %d\n", "WIDTH", nRasterXSize ));
734     CPL_IGNORE_RET_VAL(
735         VSIFPrintfL( fpRsc, "%-40s %d\n", "FILE_LENGTH", nRasterYSize ));
736 
737 /* -------------------------------------------------------------------- */
738 /*      Georeferencing.                                                 */
739 /* -------------------------------------------------------------------- */
740     if ( pszProjection != nullptr )
741     {
742         OGRSpatialReference oSRS;
743         if( oSRS.importFromWkt( pszProjection ) == OGRERR_NONE )
744         {
745             int bNorth = FALSE;
746             int iUTMZone = oSRS.GetUTMZone( &bNorth );
747             if ( iUTMZone != 0 )
748             {
749                 CPL_IGNORE_RET_VAL(
750                     VSIFPrintfL( fpRsc, "%-40s %s%d\n", "PROJECTION", "UTM",
751                                  iUTMZone ) );
752             }
753             else if ( oSRS.IsGeographic() )
754             {
755                 CPL_IGNORE_RET_VAL(
756                     VSIFPrintfL( fpRsc, "%-40s %s\n", "PROJECTION", "LL" ) );
757             }
758             else
759             {
760                 CPLError( CE_Warning, CPLE_AppDefined,
761                           "ROI_PAC format only support Latitude/Longitude and "
762                           "UTM projections, discarding projection.");
763             }
764 
765             if ( oSRS.GetAttrValue( "DATUM" ) != nullptr )
766             {
767                 if ( strcmp( oSRS.GetAttrValue( "DATUM" ), "WGS_1984" ) == 0 )
768                 {
769                     CPL_IGNORE_RET_VAL(
770                         VSIFPrintfL( fpRsc, "%-40s %s\n", "DATUM", "WGS84" ) );
771                 }
772                 else
773                 {
774                     CPLError( CE_Warning, CPLE_AppDefined,
775                               "Datum \"%s\" probably not supported in the "
776                               "ROI_PAC format, saving it anyway",
777                                   oSRS.GetAttrValue( "DATUM" ) );
778                     CPL_IGNORE_RET_VAL(
779                         VSIFPrintfL( fpRsc, "%-40s %s\n", "DATUM",
780                                      oSRS.GetAttrValue( "DATUM" ) ) );
781                 }
782             }
783             if ( oSRS.GetAttrValue( "UNIT" ) != nullptr )
784             {
785                 CPL_IGNORE_RET_VAL(
786                     VSIFPrintfL( fpRsc, "%-40s %s\n", "X_UNIT",
787                                  oSRS.GetAttrValue( "UNIT" ) ));
788                 CPL_IGNORE_RET_VAL(
789                     VSIFPrintfL( fpRsc, "%-40s %s\n", "Y_UNIT",
790                                  oSRS.GetAttrValue( "UNIT" ) ));
791             }
792         }
793     }
794     if( bValidGeoTransform )
795     {
796         if ( adfGeoTransform[2] != 0 || adfGeoTransform[4] != 0 )
797         {
798             CPLError( CE_Warning, CPLE_AppDefined,
799                       "ROI_PAC format do not support geotransform with "
800                           "rotation, discarding info.");
801         }
802         else
803         {
804             CPL_IGNORE_RET_VAL(
805                 VSIFPrintfL( fpRsc, "%-40s %.16g\n",
806                              "X_FIRST", adfGeoTransform[0] ));
807             CPL_IGNORE_RET_VAL(
808                 VSIFPrintfL( fpRsc, "%-40s %.16g\n",
809                              "X_STEP", adfGeoTransform[1] ));
810             CPL_IGNORE_RET_VAL(VSIFPrintfL( fpRsc, "%-40s %.16g\n",
811                                             "Y_FIRST", adfGeoTransform[3] ));
812             CPL_IGNORE_RET_VAL(
813                 VSIFPrintfL( fpRsc, "%-40s %.16g\n",
814                              "Y_STEP", adfGeoTransform[5] ));
815             CPL_IGNORE_RET_VAL(
816                 VSIFPrintfL( fpRsc, "%-40s %.16g\n",
817                              "Z_OFFSET", band->GetOffset(nullptr) ));
818             CPL_IGNORE_RET_VAL(
819                 VSIFPrintfL( fpRsc, "%-40s %.16g\n",
820                              "Z_SCALE", band->GetScale(nullptr) ));
821         }
822     }
823 
824 /* -------------------------------------------------------------------- */
825 /*      Metadata stored in the ROI_PAC domain.                          */
826 /* -------------------------------------------------------------------- */
827     char** papszROIPACMetadata = GetMetadata( "ROI_PAC" );
828     for( int i = 0; i < CSLCount( papszROIPACMetadata ); i++ )
829     {
830         /* Get the tokens from the metadata item */
831         char **papszTokens = CSLTokenizeString2( papszROIPACMetadata[i],
832                                                  "=",
833                                                  CSLT_STRIPLEADSPACES
834                                                  | CSLT_STRIPENDSPACES);
835         if ( CSLCount( papszTokens ) != 2 )
836         {
837             CPLDebug( "ROI_PAC",
838                       "Line of header file could not be split at = "
839                       "into two elements: %s",
840                       papszROIPACMetadata[i] );
841             CSLDestroy( papszTokens );
842             continue;
843         }
844 
845         /* Don't write it out if it is one of the bits of metadata that is
846          * written out elsewhere in this routine */
847         if ( strcmp( papszTokens[0], "WIDTH" ) == 0
848               || strcmp( papszTokens[0], "FILE_LENGTH" ) == 0 )
849         {
850             CSLDestroy( papszTokens );
851             continue;
852         }
853         CPL_IGNORE_RET_VAL(
854             VSIFPrintfL( fpRsc, "%-40s %s\n",
855                          papszTokens[0], papszTokens[1] ));
856         CSLDestroy( papszTokens );
857     }
858 }
859 
860 /************************************************************************/
861 /*                         GetGeoTransform()                            */
862 /************************************************************************/
863 
GetGeoTransform(double * padfTransform)864 CPLErr ROIPACDataset::GetGeoTransform( double *padfTransform )
865 {
866     memcpy( padfTransform, adfGeoTransform, sizeof(adfGeoTransform) );
867     return bValidGeoTransform ? CE_None : CE_Failure;
868 }
869 
870 /************************************************************************/
871 /*                          SetGeoTransform()                           */
872 /************************************************************************/
873 
SetGeoTransform(double * padfTransform)874 CPLErr ROIPACDataset::SetGeoTransform( double *padfTransform )
875 {
876     memcpy( adfGeoTransform, padfTransform, sizeof(adfGeoTransform) );
877     bValidGeoTransform = true;
878     return CE_None;
879 }
880 
881 /************************************************************************/
882 /*                         GetProjectionRef()                           */
883 /************************************************************************/
884 
_GetProjectionRef(void)885 const char *ROIPACDataset::_GetProjectionRef( void )
886 {
887     return pszProjection != nullptr ? pszProjection : "";
888 }
889 
890 /************************************************************************/
891 /*                          SetProjection()                             */
892 /************************************************************************/
893 
_SetProjection(const char * pszNewProjection)894 CPLErr ROIPACDataset::_SetProjection( const char *pszNewProjection )
895 
896 {
897     CPLFree( pszProjection );
898     pszProjection = (pszNewProjection) ? CPLStrdup( pszNewProjection ) : nullptr;
899     return CE_None;
900 }
901 
902 /************************************************************************/
903 /*                            GetFileList()                             */
904 /************************************************************************/
905 
GetFileList()906 char **ROIPACDataset::GetFileList()
907 {
908     // Main data file, etc.
909     char **papszFileList = RawDataset::GetFileList();
910 
911     // RSC file.
912     papszFileList = CSLAddString( papszFileList, pszRscFilename );
913 
914     return papszFileList;
915 }
916 
917 /************************************************************************/
918 /*                         ROIPACRasterBand()                           */
919 /************************************************************************/
920 
ROIPACRasterBand(GDALDataset * poDSIn,int nBandIn,VSILFILE * fpRawIn,vsi_l_offset nImgOffsetIn,int nPixelOffsetIn,int nLineOffsetIn,GDALDataType eDataTypeIn,int bNativeOrderIn)921 ROIPACRasterBand::ROIPACRasterBand( GDALDataset *poDSIn, int nBandIn, VSILFILE *fpRawIn,
922                                     vsi_l_offset nImgOffsetIn, int nPixelOffsetIn,
923                                     int nLineOffsetIn,
924                                     GDALDataType eDataTypeIn, int bNativeOrderIn ) :
925     RawRasterBand(poDSIn, nBandIn, fpRawIn, nImgOffsetIn, nPixelOffsetIn,
926                   nLineOffsetIn, eDataTypeIn, bNativeOrderIn,
927                   RawRasterBand::OwnFP::NO)
928 {}
929 
930 /************************************************************************/
931 /*                        GDALRegister_ROIPAC()                         */
932 /************************************************************************/
933 
GDALRegister_ROIPAC()934 void GDALRegister_ROIPAC()
935 {
936     if( GDALGetDriverByName( "ROI_PAC" ) != nullptr )
937         return;
938 
939     GDALDriver *poDriver = new GDALDriver();
940 
941     poDriver->SetDescription( "ROI_PAC" );
942     poDriver->SetMetadataItem( GDAL_DMD_LONGNAME, "ROI_PAC raster" );
943     poDriver->SetMetadataItem( GDAL_DMD_HELPTOPIC,
944                                "drivers/raster/roi_pac.html" );
945     poDriver->SetMetadataItem( GDAL_DCAP_RASTER, "YES" );
946     poDriver->SetMetadataItem( GDAL_DCAP_VIRTUALIO, "YES" );
947 
948     poDriver->pfnOpen = ROIPACDataset::Open;
949     poDriver->pfnIdentify = ROIPACDataset::Identify;
950     poDriver->pfnCreate = ROIPACDataset::Create;
951 
952     GetGDALDriverManager()->RegisterDriver( poDriver );
953 }
954