Coverage Report

Created: 2025-12-31 08:30

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/frmts/idrisi/IdrisiDataset.cpp
Line
Count
Source
1
/*****************************************************************************
2
 *
3
 * Project:  Idrisi Raster Image File Driver
4
 * Purpose:  Read/write Idrisi Raster Image Format RST
5
 * Author:   Ivan Lucena, [lucena_ivan at hotmail.com]
6
 *
7
 * Revised by Hongmei Zhu, February, 2013
8
 * honzhu@clarku.edu
9
 * Clark Labs/Clark University
10
 *
11
 ******************************************************************************
12
 * Copyright( c ) 2006, Ivan Lucena
13
 * Copyright (c) 2007-2012, Even Rouault <even dot rouault at spatialys.com>
14
 *
15
 * SPDX-License-Identifier: MIT
16
 ****************************************************************************/
17
18
#include "cpl_conv.h"
19
#include "cpl_string.h"
20
#include "cpl_csv.h"
21
#include "gdal_frmts.h"
22
#include "gdal_pam.h"
23
#include "gdal_alg.h"
24
#include "gdal_priv.h"
25
#include "gdal_rat.h"
26
#include "ogr_spatialref.h"
27
#include "idrisi.h"
28
29
#include "proj_experimental.h"
30
#include "ogr_proj_p.h"
31
32
#include <cmath>
33
34
#ifdef _WIN32
35
#define PATHDELIM '\\'
36
#else
37
0
#define PATHDELIM '/'
38
#endif
39
40
//----- Safe numeric conversion, NULL as zero
41
0
#define atoi_nz(s) (s == nullptr ? (int)0 : atoi(s))
42
0
#define CPLAtof_nz(s) (s == nullptr ? (double)0.0 : CPLAtof(s))
43
44
//----- file extensions:
45
static const char *const extRST = "rst";
46
static const char *const extRDC = "rdc";
47
static const char *const extSMP = "smp";
48
static const char *const extREF = "ref";
49
// static const char * const extRSTu        = "RST";
50
static const char *const extRDCu = "RDC";
51
static const char *const extSMPu = "SMP";
52
static const char *const extREFu = "REF";
53
54
//----- field names on rdc file:
55
static const char *const rdcFILE_FORMAT = "file format ";
56
static const char *const rdcFILE_TITLE = "file title  ";
57
static const char *const rdcDATA_TYPE = "data type   ";
58
static const char *const rdcFILE_TYPE = "file type   ";
59
static const char *const rdcCOLUMNS = "columns     ";
60
static const char *const rdcROWS = "rows        ";
61
static const char *const rdcREF_SYSTEM = "ref. system ";
62
static const char *const rdcREF_UNITS = "ref. units  ";
63
static const char *const rdcUNIT_DIST = "unit dist.  ";
64
static const char *const rdcMIN_X = "min. X      ";
65
static const char *const rdcMAX_X = "max. X      ";
66
static const char *const rdcMIN_Y = "min. Y      ";
67
static const char *const rdcMAX_Y = "max. Y      ";
68
static const char *const rdcPOSN_ERROR = "pos'n error ";
69
static const char *const rdcRESOLUTION = "resolution  ";
70
static const char *const rdcMIN_VALUE = "min. value  ";
71
static const char *const rdcMAX_VALUE = "max. value  ";
72
static const char *const rdcDISPLAY_MIN = "display min ";
73
static const char *const rdcDISPLAY_MAX = "display max ";
74
static const char *const rdcVALUE_UNITS = "value units ";
75
static const char *const rdcVALUE_ERROR = "value error ";
76
static const char *const rdcFLAG_VALUE = "flag value  ";
77
static const char *const rdcFLAG_DEFN = "flag def'n  ";
78
static const char *const rdcFLAG_DEFN2 = "flag def`n  ";
79
static const char *const rdcLEGEND_CATS = "legend cats ";
80
static const char *const rdcLINEAGES = "lineage     ";
81
static const char *const rdcCOMMENTS = "comment     ";
82
static const char *const rdcCODE_N = "code %6d ";
83
84
//----- ".ref" file field names:
85
static const char *const refREF_SYSTEM = "ref. system ";
86
static const char *const refREF_SYSTEM2 = "ref.system  ";
87
static const char *const refPROJECTION = "projection  ";
88
static const char *const refDATUM = "datum       ";
89
static const char *const refDELTA_WGS84 = "delta WGS84 ";
90
static const char *const refELLIPSOID = "ellipsoid   ";
91
static const char *const refMAJOR_SAX = "major s-ax  ";
92
static const char *const refMINOR_SAX = "minor s-ax  ";
93
static const char *const refORIGIN_LONG = "origin long ";
94
static const char *const refORIGIN_LAT = "origin lat  ";
95
static const char *const refORIGIN_X = "origin X    ";
96
static const char *const refORIGIN_Y = "origin Y    ";
97
static const char *const refSCALE_FAC = "scale fac   ";
98
static const char *const refUNITS = "units       ";
99
static const char *const refPARAMETERS = "parameters  ";
100
static const char *const refSTANDL_1 = "stand ln 1  ";
101
static const char *const refSTANDL_2 = "stand ln 2  ";
102
103
//----- standard values:
104
static const char *const rstVERSION = "Idrisi Raster A.1";
105
static const char *const rstBYTE = "byte";
106
static const char *const rstINTEGER = "integer";
107
static const char *const rstREAL = "real";
108
static const char *const rstRGB24 = "rgb24";
109
static const char *const rstDEGREE = "deg";
110
static const char *const rstMETER = "m";
111
static const char *const rstLATLONG = "latlong";
112
static const char *const rstLATLONG2 = "lat/long";
113
static const char *const rstPLANE = "plane";
114
static const char *const rstUTM = "utm-%d%c";
115
static const char *const rstSPC = "spc%2d%2s%d";
116
117
//----- palette file( .smp ) header size:
118
constexpr int smpHEADERSIZE = 18;
119
120
//----- check if file exists:
121
bool FileExists(const char *pszPath);
122
123
//----- Reference Table
124
struct ReferenceTab
125
{
126
    int nCode;
127
    const char *pszName;
128
};
129
130
//----- USA State's reference table to USGS PCS Code
131
constexpr ReferenceTab aoUSStateTable[] = {
132
    {101, "al"},  {201, "az"},  {301, "ar"},  {401, "ca"},  {501, "co"},
133
    {600, "ct"},  {700, "de"},  {901, "fl"},  {1001, "ga"}, {1101, "id"},
134
    {1201, "il"}, {1301, "in"}, {1401, "ia"}, {1501, "ks"}, {1601, "ky"},
135
    {1701, "la"}, {1801, "me"}, {1900, "md"}, {2001, "ma"}, {2111, "mi"},
136
    {2201, "mn"}, {2301, "ms"}, {2401, "mo"}, {2500, "mt"}, {2600, "ne"},
137
    {2701, "nv"}, {2800, "nh"}, {2900, "nj"}, {3001, "nm"}, {3101, "ny"},
138
    {3200, "nc"}, {3301, "nd"}, {3401, "oh"}, {3501, "ok"}, {3601, "or"},
139
    {3701, "pa"}, {3800, "ri"}, {3900, "sc"}, {4001, "sd"}, {4100, "tn"},
140
    {4201, "tx"}, {4301, "ut"}, {4400, "vt"}, {4501, "va"}, {4601, "wa"},
141
    {4701, "wv"}, {4801, "wv"}, {4901, "wy"}, {5001, "ak"}, {5101, "hi"},
142
    {5200, "pr"}};
143
144
//---- The origin table for USA State Plane Systems
145
struct OriginTab83
146
{
147
    double longitude;
148
    double latitude;
149
    const char *spcs;
150
};
151
152
constexpr int ORIGIN_COUNT = 148;
153
154
//---- USA State plane coordinate system in IDRISI
155
static const OriginTab83 SPCS83Origin[] = {
156
    {85.83, 30.50, "SPC83AL1"},  {87.50, 30.00, "SPC83AL2"},
157
    {176.00, 51.00, "SPC83AK0"}, {142.00, 54.00, "SPC83AK2"},
158
    {146.00, 54.00, "SPC83AK3"}, {150.00, 54.00, "SPC83AK4"},
159
    {154.00, 54.00, "SPC83AK5"}, {158.00, 54.00, "SPC83AK6"},
160
    {162.00, 54.00, "SPC83AK7"}, {166.00, 54.00, "SPC83AK8"},
161
    {170.00, 54.00, "SPC83AK9"}, {110.17, 31.00, "SPC83AZ1"},
162
    {111.92, 31.00, "SPC83AZ2"}, {113.75, 31.00, "SPC83AZ3"},
163
    {92.00, 34.33, "SPC83AR1"},  {92.00, 32.67, "SPC83AR2"},
164
    {122.00, 39.33, "SPC83CA1"}, {122.00, 37.67, "SPC83CA2"},
165
    {120.50, 36.50, "SPC83CA3"}, {119.00, 35.33, "SPC83CA4"},
166
    {118.00, 33.50, "SPC83CA5"}, {116.25, 32.17, "SPC83CA6"},
167
    {105.50, 39.33, "SPC83CO1"}, {105.50, 37.83, "SPC83CO2"},
168
    {105.50, 36.67, "SPC83CO3"}, {72.75, 40.83, "SPC83CT1"},
169
    {75.42, 38.00, "SPC83DE1"},  {81.00, 24.33, "SPC83FL1"},
170
    {82.00, 24.33, "SPC83FL2"},  {84.50, 29.00, "SPC83FL3"},
171
    {82.17, 30.00, "SPC83GA1"},  {84.17, 30.00, "SPC83GA2"},
172
    {155.50, 18.83, "SPC83HI1"}, {156.67, 20.33, "SPC83HI2"},
173
    {158.00, 21.17, "SPC83HI3"}, {159.50, 21.83, "SPC83HI4"},
174
    {160.17, 21.67, "SPC83HI5"}, {112.17, 41.67, "SPC83ID1"},
175
    {114.00, 41.67, "SPC83ID2"}, {115.75, 41.67, "SPC83ID3"},
176
    {88.33, 36.67, "SPC83IL1"},  {90.17, 36.67, "SPC83IL1"},
177
    {85.67, 37.50, "SPC83IN1"},  {87.08, 37.50, "SPC83IN2"},
178
    {93.50, 41.50, "SPC83IA1"},  {93.50, 40.00, "SPC83IA1"},
179
    {98.00, 38.33, "SPC83KS1"},  {98.50, 36.67, "SPC83KS2"},
180
    {84.25, 37.50, "SPC83KY1"},  {85.75, 36.33, "SPC83KY2"},
181
    {92.50, 30.50, "SPC83LA1"},  {91.33, 28.50, "SPC83LA2"},
182
    {91.33, 25.50, "SPC83LA3"},  {92.50, 30.67, "SPC27LA1"},  // NAD27 system
183
    {91.33, 28.67, "SPC27LA2"},  {91.33, 25.67, "SPC27LA3"},  //
184
    {68.50, 43.67, "SPC83ME1"},  {68.50, 43.83, "SPC27ME1"},  // NAD27
185
    {70.17, 42.83, "SPC83ME2"},  {77.00, 37.67, "SPC83MD1"},  //
186
    {77.00, 37.83, "SPC27MD1"},                               // NAD27
187
    {71.50, 41.00, "SPC83MA1"},  {70.50, 41.00, "SPC83MA2"},
188
    {87.00, 44.78, "SPC83MI1"},  {84.37, 43.32, "SPC83MI2"},
189
    {84.37, 41.50, "SPC83MI3"},  {84.33, 43.32, "SPC27MI2"},  // NAD27 L
190
    {84.33, 41.50, "SPC27MI3"},                               // NAD27 L
191
    {83.67, 41.50, "SPC27MI4"},                               // NAD27 TM
192
    {85.75, 41.50, "SPC27MI5"},                               // NAD27 TM
193
    {88.75, 41.50, "SPC27MI6"},                               // NAD27 TM
194
    {93.10, 46.50, "SPC83MN1"},  {94.25, 45.00, "SPC83MN2"},
195
    {94.00, 43.00, "SPC83MN3"},  {88.83, 29.50, "SPC83MS1"},
196
    {90.33, 29.50, "SPC83MS2"},  {88.83, 29.67, "SPC83MS1"},  // NAD27
197
    {90.33, 30.50, "SPC83MS2"},                               //
198
    {90.50, 35.83, "SPC83MO1"},  {92.50, 35.83, "SPC83MO2"},
199
    {94.50, 36.17, "SPC83MO3"},  {109.50, 44.25, "SPC83MT1"},
200
    {109.50, 47.00, "SPC27MT1"},                               // NAD27
201
    {109.50, 45.83, "SPC27MT2"}, {109.50, 44.00, "SPC27MT3"},  //
202
    {100.00, 39.83, "SPC83NE1"}, {115.58, 34.75, "SPC83NV1"},
203
    {116.67, 34.75, "SPC83NV2"}, {118.58, 34.75, "SPC83NV3"},
204
    {71.67, 42.50, "SPC83NH1"},  {74.50, 38.83, "SPC83NJ1"},
205
    {74.67, 38.83, "SPC27NJ1"},  // NAD27
206
    {104.33, 31.00, "SPC83NM1"}, {106.25, 31.00, "SPC83NM2"},
207
    {107.83, 31.00, "SPC83NM3"}, {74.50, 38.83, "SPC83NY1"},
208
    {76.58, 40.00, "SPC83NY2"},  {78.58, 40.00, "SPC83NY3"},
209
    {74.00, 40.17, "SPC83NY4"},  {74.33, 40.00, "SPC27NY1"},  // NAD27
210
    {74.00, 40.50, "SPC27NY4"},                               //
211
    {79.00, 33.75, "SPC83NC1"},  {100.50, 47.00, "SPC83ND1"},
212
    {100.50, 45.67, "SPC83ND2"}, {82.50, 39.67, "SPC83OH1"},
213
    {82.50, 38.00, "SPC83OH2"},  {98.00, 35.00, "SPC83OK1"},
214
    {98.00, 33.33, "SPC83OK2"},  {120.50, 43.67, "SPC83OR1"},
215
    {120.50, 41.67, "SPC83OR2"}, {77.75, 40.17, "SPC83PA1"},
216
    {77.75, 39.33, "SPC83PA2"},  {71.50, 41.08, "SPC83RI1"},
217
    {81.00, 31.83, "SPC83SC1"},  {81.00, 33.00, "SPC27SC1"},  // NAD27
218
    {81.00, 31.83, "SPC27SC2"},                               // NAD27
219
    {100.00, 43.83, "SPC83SD1"}, {100.33, 42.33, "SPC83SD2"},
220
    {86.00, 34.33, "SPC83TN1"},  {86.00, 34.67, "SPC27TN1"},  // NAD27
221
    {101.50, 34.00, "SPC83TX1"},                              //
222
    {98.50, 31.67, "SPC83TX2"},  {100.33, 29.67, "SPC83TX3"},
223
    {99.00, 27.83, "SPC83TX4"},  {98.50, 25.67, "SPC83TX5"},
224
    {97.50, 31.67, "SPC27TX2"},  // NAD27
225
    {111.50, 40.33, "SPC83UT1"}, {111.50, 38.33, "SPC83UT2"},
226
    {111.50, 36.67, "SPC83UT3"}, {72.50, 42.50, "SPC83VT1"},
227
    {78.50, 37.67, "SPC83VA1"},  {78.50, 36.33, "SPC83VA2"},
228
    {120.83, 47.00, "SPC83WA1"}, {120.50, 45.33, "SPC83WA2"},
229
    {79.50, 38.50, "SPC83WV1"},  {81.00, 37.00, "SPC83WV2"},
230
    {90.00, 45.17, "SPC83WI1"},  {90.00, 43.83, "SPC83WI2"},
231
    {90.00, 42.00, "SPC83WI3"},  {105.17, 40.50, "SPC83WY1"},
232
    {107.33, 40.50, "SPC83WY2"}, {108.75, 40.50, "SPC83WY3"},
233
    {110.08, 40.50, "SPC83WY4"}, {105.17, 40.67, "SPC27WY1"},  // NAD27
234
    {105.33, 40.67, "SPC27WY2"}, {108.75, 40.67, "SPC27WY3"},
235
    {110.08, 40.67, "SPC27WY4"},  //
236
    {66.43, 17.83, "SPC83PR1"}};
237
238
// Get IDRISI State Plane name by origin
239
char *GetSpcs(double dfLon, double dfLat);
240
241
// change NAD from 83 to 27
242
void NAD83to27(char *pszOutRef, char *pszInRef);
243
244
0
#define US_STATE_COUNT (sizeof(aoUSStateTable) / sizeof(ReferenceTab))
245
246
//----- Get the Code of a US State
247
int GetStateCode(const char *pszState);
248
249
//----- Get the state name of a Code
250
const char *GetStateName(int nCode);
251
252
//----- Conversion Table definition
253
struct ConversionTab
254
{
255
    const char *pszName;
256
    int nDefaultI;
257
    int nDefaultG;
258
    double dfConv;
259
};
260
261
//----- Linear Unit Conversion Table
262
static const ConversionTab aoLinearUnitsConv[] = {
263
    {"m", /*  0 */ 0, 1, 1.0},
264
    {SRS_UL_METER, /*  1 */ 0, 1, 1.0},
265
    {"meters", /*  2 */ 0, 1, 1.0},
266
    {"metre", /*  3 */ 0, 1, 1.0},
267
268
    {"ft", /*  4 */ 4, 5, 0.3048},
269
    {SRS_UL_FOOT, /*  5 */ 4, 5, 0.3048},
270
    {"feet", /*  6 */ 4, 5, 0.3048},
271
    {"foot_us", /*  7 */ 4, 5, 0.3048006},
272
    {"u.s. foot", /*  8 */ 4, 5, 0.3048006},
273
274
    {"mi", /*  9 */ 9, 10, 1612.9},
275
    {"mile", /* 10 */ 9, 10, 1612.9},
276
    {"miles", /* 11 */ 9, 10, 1612.9},
277
278
    {"km", /* 12 */ 12, 13, 1000.0},
279
    {"kilometers", /* 13 */ 12, 13, 1000.0},
280
    {"kilometer", /* 14 */ 12, 13, 1000.0},
281
    {"kilometre", /* 15 */ 12, 13, 1000.0},
282
283
    {"deg", /* 16 */ 16, 17, 0.0},
284
    {SRS_UA_DEGREE, /* 17 */ 16, 17, 0.0},
285
    {"degrees", /* 18 */ 16, 17, 0.0},
286
287
    {"rad", /* 19 */ 19, 20, 0.0},
288
    {SRS_UA_RADIAN, /* 20 */ 19, 20, 0.0},
289
    {"radians", /* 21 */ 19, 20, 0.0}};
290
0
#define LINEAR_UNITS_COUNT (sizeof(aoLinearUnitsConv) / sizeof(ConversionTab))
291
292
//----- Get the index of a given linear unit
293
static int GetUnitIndex(const char *pszUnitName);
294
295
//----- Get the default name
296
static char *GetUnitDefault(const char *pszUnitName,
297
                            const char *pszToMeter = nullptr);
298
299
//----- Get the "to meter"
300
static int GetToMeterIndex(const char *pszToMeter);
301
302
//----- CSLSaveCRLF
303
static int SaveAsCRLF(char **papszStrList, const char *pszFname);
304
305
/************************************************************************/
306
/*                     myCSLFetchNameValue()                            */
307
/************************************************************************/
308
309
static const char *myCSLFetchNameValue(char **papszStrList, const char *pszName)
310
0
{
311
0
    if (papszStrList == nullptr || pszName == nullptr)
312
0
        return nullptr;
313
314
0
    size_t nLen = strlen(pszName);
315
0
    while (nLen > 0 && pszName[nLen - 1] == ' ')
316
0
        nLen--;
317
0
    while (*papszStrList != nullptr)
318
0
    {
319
0
        if (EQUALN(*papszStrList, pszName, nLen))
320
0
        {
321
0
            size_t i;
322
0
            for (i = nLen; (*papszStrList)[i] == ' '; ++i)
323
0
            {
324
0
            }
325
0
            if ((*papszStrList)[i] == '=' || (*papszStrList)[i] == ':')
326
0
            {
327
0
                return (*papszStrList) + i + 1;
328
0
            }
329
0
        }
330
0
        ++papszStrList;
331
0
    }
332
0
    return nullptr;
333
0
}
334
335
/************************************************************************/
336
/*                   myCSLSetNameValueSeparator()                       */
337
/************************************************************************/
338
339
static void myCSLSetNameValueSeparator(char **papszList,
340
                                       const char *pszSeparator)
341
0
{
342
0
    const int nLines = CSLCount(papszList);
343
344
0
    for (int iLine = 0; iLine < nLines; ++iLine)
345
0
    {
346
0
        char *pszSep = strchr(papszList[iLine], '=');
347
0
        if (pszSep == nullptr)
348
0
            pszSep = strchr(papszList[iLine], ':');
349
0
        if (pszSep == nullptr)
350
0
            continue;
351
0
        *pszSep = '\0';
352
0
        const char *pszKey = papszList[iLine];
353
0
        const char *pszValue = pszSep + 1;
354
0
        while (*pszValue == ' ')
355
0
            pszValue++;
356
357
0
        char *pszNewLine = static_cast<char *>(CPLMalloc(
358
0
            strlen(pszValue) + strlen(pszKey) + strlen(pszSeparator) + 1));
359
0
        strcpy(pszNewLine, pszKey);
360
0
        strcat(pszNewLine, pszSeparator);
361
0
        strcat(pszNewLine, pszValue);
362
0
        CPLFree(papszList[iLine]);
363
0
        papszList[iLine] = pszNewLine;
364
0
    }
365
0
}
366
367
//----- Classes pre-definition:
368
class IdrisiRasterBand;
369
370
//  ----------------------------------------------------------------------------
371
//        Idrisi GDALDataset
372
//  ----------------------------------------------------------------------------
373
374
class IdrisiDataset final : public GDALPamDataset
375
{
376
    friend class IdrisiRasterBand;
377
378
  private:
379
    VSILFILE *fp;
380
381
    char *pszFilename;
382
    char *pszDocFilename;
383
    char **papszRDC;
384
    GDALGeoTransform m_gt{};
385
386
    mutable OGRSpatialReference m_oSRS{};
387
    char **papszCategories;
388
    char *pszUnitType;
389
    // Move GeoReference2Wkt() into header file.
390
    CPLErr Wkt2GeoReference(const OGRSpatialReference &oSRS,
391
                            char **pszRefSystem, char **pszRefUnit);
392
393
  protected:
394
    GDALColorTable *poColorTable;
395
396
  public:
397
    IdrisiDataset();
398
    ~IdrisiDataset() override;
399
400
    static GDALDataset *Open(GDALOpenInfo *poOpenInfo);
401
    static GDALDataset *Create(const char *pszFilename, int nXSize, int nYSize,
402
                               int nBands, GDALDataType eType,
403
                               char **papszOptions);
404
    static GDALDataset *CreateCopy(const char *pszFilename,
405
                                   GDALDataset *poSrcDS, int bStrict,
406
                                   char **papszOptions,
407
                                   GDALProgressFunc pfnProgress,
408
                                   void *pProgressData);
409
    char **GetFileList(void) override;
410
    CPLErr GetGeoTransform(GDALGeoTransform &gt) const override;
411
    CPLErr SetGeoTransform(const GDALGeoTransform &gt) override;
412
413
    const OGRSpatialReference *GetSpatialRef() const override;
414
    CPLErr SetSpatialRef(const OGRSpatialReference *poSRS) override;
415
};
416
417
//  ----------------------------------------------------------------------------
418
//        Idrisi GDALPamRasterBand
419
//  ----------------------------------------------------------------------------
420
421
class IdrisiRasterBand final : public GDALPamRasterBand
422
{
423
    friend class IdrisiDataset;
424
425
    GDALRasterAttributeTable *poDefaultRAT;
426
427
  private:
428
    int nRecordSize;
429
    GByte *pabyScanLine;
430
431
  public:
432
    IdrisiRasterBand(IdrisiDataset *poDS, int nBand, GDALDataType eDataType);
433
    ~IdrisiRasterBand() override;
434
435
    double GetNoDataValue(int *pbSuccess = nullptr) override;
436
    double GetMinimum(int *pbSuccess = nullptr) override;
437
    double GetMaximum(int *pbSuccess = nullptr) override;
438
    CPLErr IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage) override;
439
    virtual CPLErr IWriteBlock(int nBlockXOff, int nBlockYOff,
440
                               void *pImage) override;
441
    GDALColorTable *GetColorTable() override;
442
    GDALColorInterp GetColorInterpretation() override;
443
    char **GetCategoryNames() override;
444
    const char *GetUnitType() override;
445
446
    CPLErr SetCategoryNames(char **papszCategoryNames) override;
447
    CPLErr SetNoDataValue(double dfNoDataValue) override;
448
    CPLErr SetColorTable(GDALColorTable *poColorTable) override;
449
    CPLErr SetUnitType(const char *pszUnitType) override;
450
    virtual CPLErr SetStatistics(double dfMin, double dfMax, double dfMean,
451
                                 double dfStdDev) override;
452
    CPLErr SetMinMax(double dfMin, double dfMax);
453
    GDALRasterAttributeTable *GetDefaultRAT() override;
454
    CPLErr SetDefaultRAT(const GDALRasterAttributeTable *) override;
455
456
    float fMaximum;
457
    float fMinimum;
458
    bool bFirstVal;
459
};
460
461
//  ------------------------------------------------------------------------  //
462
//                        Implementation of IdrisiDataset                     //
463
//  ------------------------------------------------------------------------  //
464
465
/************************************************************************/
466
/*                           IdrisiDataset()                            */
467
/************************************************************************/
468
469
IdrisiDataset::IdrisiDataset()
470
0
    : fp(nullptr), pszFilename(nullptr), pszDocFilename(nullptr),
471
0
      papszRDC(nullptr), papszCategories(nullptr), pszUnitType(nullptr),
472
0
      poColorTable(new GDALColorTable())
473
0
{
474
0
    m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
475
0
}
476
477
/************************************************************************/
478
/*                           ~IdrisiDataset()                           */
479
/************************************************************************/
480
481
IdrisiDataset::~IdrisiDataset()
482
0
{
483
0
    FlushCache(true);
484
485
0
    if (papszRDC != nullptr && eAccess == GA_Update)
486
0
    {
487
        // int bSuccessMin = FALSE;
488
        // int bSuccessMax = FALSE;
489
490
0
        double dfMin = 0.0;
491
0
        double dfMax = 0.0;
492
0
        double dfMean = 0.0;
493
0
        double dfStdDev = 0.0;
494
495
0
        for (int i = 0; i < nBands; i++)
496
0
        {
497
0
            IdrisiRasterBand *poBand =
498
0
                cpl::down_cast<IdrisiRasterBand *>(GetRasterBand(i + 1));
499
0
            poBand->ComputeStatistics(false, &dfMin, &dfMax, &dfMean, &dfStdDev,
500
0
                                      nullptr, nullptr);
501
            /*
502
              dfMin = poBand->GetMinimum( &bSuccessMin );
503
              dfMax = poBand->GetMaximum( &bSuccessMax );
504
505
              if( ! ( bSuccessMin && bSuccessMax ) )
506
              {
507
              poBand->GetStatistics( false, true, &dfMin, &dfMax, NULL, NULL );
508
              }
509
            */
510
0
            poBand->SetMinMax(dfMin, dfMax);
511
0
        }
512
513
0
        myCSLSetNameValueSeparator(papszRDC, ": ");
514
0
        SaveAsCRLF(papszRDC, pszDocFilename);
515
0
    }
516
0
    CSLDestroy(papszRDC);
517
518
0
    if (poColorTable)
519
0
    {
520
0
        delete poColorTable;
521
0
    }
522
0
    CPLFree(pszFilename);
523
0
    CPLFree(pszDocFilename);
524
0
    CSLDestroy(papszCategories);
525
0
    CPLFree(pszUnitType);
526
527
0
    if (fp != nullptr)
528
0
        VSIFCloseL(fp);
529
0
}
530
531
/************************************************************************/
532
/*                                Open()                                */
533
/************************************************************************/
534
535
GDALDataset *IdrisiDataset::Open(GDALOpenInfo *poOpenInfo)
536
536k
{
537
536k
    if ((poOpenInfo->fpL == nullptr) ||
538
82.3k
        (poOpenInfo->IsExtensionEqualToCI(extRST) == FALSE))  // modified
539
536k
        return nullptr;
540
541
    // --------------------------------------------------------------------
542
    //      Check the documentation file .rdc
543
    // --------------------------------------------------------------------
544
545
0
    std::string osLDocFilename =
546
0
        CPLResetExtensionSafe(poOpenInfo->pszFilename, extRDC);
547
548
0
    if (!FileExists(osLDocFilename.c_str()))
549
0
    {
550
0
        osLDocFilename =
551
0
            CPLResetExtensionSafe(poOpenInfo->pszFilename, extRDCu);
552
553
0
        if (!FileExists(osLDocFilename.c_str()))
554
0
        {
555
0
            return nullptr;
556
0
        }
557
0
    }
558
559
0
    char **papszLRDC = CSLLoad(osLDocFilename.c_str());
560
561
0
    myCSLSetNameValueSeparator(papszLRDC, ":");
562
563
0
    const char *pszVersion = myCSLFetchNameValue(papszLRDC, rdcFILE_FORMAT);
564
565
0
    if (pszVersion == nullptr || !EQUAL(pszVersion, rstVERSION))
566
0
    {
567
0
        CSLDestroy(papszLRDC);
568
0
        return nullptr;
569
0
    }
570
571
    // --------------------------------------------------------------------
572
    //      Create a corresponding GDALDataset
573
    // --------------------------------------------------------------------
574
575
0
    IdrisiDataset *poDS = new IdrisiDataset();
576
0
    poDS->eAccess = poOpenInfo->eAccess;
577
0
    poDS->pszFilename = CPLStrdup(poOpenInfo->pszFilename);
578
579
0
    if (poOpenInfo->eAccess == GA_ReadOnly)
580
0
    {
581
0
        poDS->fp = VSIFOpenL(poDS->pszFilename, "rb");
582
0
    }
583
0
    else
584
0
    {
585
0
        poDS->fp = VSIFOpenL(poDS->pszFilename, "r+b");
586
0
    }
587
588
0
    if (poDS->fp == nullptr)
589
0
    {
590
0
        CSLDestroy(papszLRDC);
591
0
        delete poDS;
592
0
        return nullptr;
593
0
    }
594
595
0
    poDS->pszDocFilename = CPLStrdup(osLDocFilename.c_str());
596
0
    poDS->papszRDC = CSLDuplicate(papszLRDC);
597
0
    CSLDestroy(papszLRDC);
598
599
    // --------------------------------------------------------------------
600
    //      Load information from rdc
601
    // --------------------------------------------------------------------
602
603
0
    poDS->nRasterXSize =
604
0
        atoi_nz(myCSLFetchNameValue(poDS->papszRDC, rdcCOLUMNS));
605
0
    poDS->nRasterYSize = atoi_nz(myCSLFetchNameValue(poDS->papszRDC, rdcROWS));
606
0
    if (!GDALCheckDatasetDimensions(poDS->nRasterXSize, poDS->nRasterYSize))
607
0
    {
608
0
        delete poDS;
609
0
        return nullptr;
610
0
    }
611
612
    // --------------------------------------------------------------------
613
    //      Create band information
614
    // --------------------------------------------------------------------
615
616
0
    const char *pszDataType = myCSLFetchNameValue(poDS->papszRDC, rdcDATA_TYPE);
617
0
    if (pszDataType == nullptr)
618
0
    {
619
0
        delete poDS;
620
0
        return nullptr;
621
0
    }
622
623
0
    if (EQUAL(pszDataType, rstBYTE))
624
0
    {
625
0
        poDS->nBands = 1;
626
0
        poDS->SetBand(1, new IdrisiRasterBand(poDS, 1, GDT_UInt8));
627
0
    }
628
0
    else if (EQUAL(pszDataType, rstINTEGER))
629
0
    {
630
0
        poDS->nBands = 1;
631
0
        poDS->SetBand(1, new IdrisiRasterBand(poDS, 1, GDT_Int16));
632
0
    }
633
0
    else if (EQUAL(pszDataType, rstREAL))
634
0
    {
635
0
        poDS->nBands = 1;
636
0
        poDS->SetBand(1, new IdrisiRasterBand(poDS, 1, GDT_Float32));
637
0
    }
638
0
    else if (EQUAL(pszDataType, rstRGB24))
639
0
    {
640
0
        poDS->nBands = 3;
641
0
        poDS->SetBand(1, new IdrisiRasterBand(poDS, 1, GDT_UInt8));
642
0
        poDS->SetBand(2, new IdrisiRasterBand(poDS, 2, GDT_UInt8));
643
0
        poDS->SetBand(3, new IdrisiRasterBand(poDS, 3, GDT_UInt8));
644
0
    }
645
0
    else
646
0
    {
647
0
        CPLError(CE_Failure, CPLE_AppDefined, "Unknown data type : %s",
648
0
                 pszDataType);
649
0
        delete poDS;
650
0
        return nullptr;
651
0
    }
652
653
0
    for (int i = 0; i < poDS->nBands; i++)
654
0
    {
655
0
        IdrisiRasterBand *band =
656
0
            cpl::down_cast<IdrisiRasterBand *>(poDS->GetRasterBand(i + 1));
657
0
        if (band->pabyScanLine == nullptr)
658
0
        {
659
0
            delete poDS;
660
0
            return nullptr;
661
0
        }
662
0
    }
663
664
    // --------------------------------------------------------------------
665
    //      Load the transformation matrix
666
    // --------------------------------------------------------------------
667
668
0
    const char *pszMinX = myCSLFetchNameValue(poDS->papszRDC, rdcMIN_X);
669
0
    const char *pszMaxX = myCSLFetchNameValue(poDS->papszRDC, rdcMAX_X);
670
0
    const char *pszMinY = myCSLFetchNameValue(poDS->papszRDC, rdcMIN_Y);
671
0
    const char *pszMaxY = myCSLFetchNameValue(poDS->papszRDC, rdcMAX_Y);
672
0
    const char *pszUnit = myCSLFetchNameValue(poDS->papszRDC, rdcUNIT_DIST);
673
674
0
    if (pszMinX != nullptr && strlen(pszMinX) > 0 && pszMaxX != nullptr &&
675
0
        strlen(pszMaxX) > 0 && pszMinY != nullptr && strlen(pszMinY) > 0 &&
676
0
        pszMaxY != nullptr && strlen(pszMaxY) > 0 && pszUnit != nullptr &&
677
0
        strlen(pszUnit) > 0)
678
0
    {
679
0
        double dfMinX, dfMaxX, dfMinY, dfMaxY, dfUnit, dfXPixSz, dfYPixSz;
680
681
0
        dfMinX = CPLAtof_nz(pszMinX);
682
0
        dfMaxX = CPLAtof_nz(pszMaxX);
683
0
        dfMinY = CPLAtof_nz(pszMinY);
684
0
        dfMaxY = CPLAtof_nz(pszMaxY);
685
0
        dfUnit = CPLAtof_nz(pszUnit);
686
687
0
        dfMinX = dfMinX * dfUnit;
688
0
        dfMaxX = dfMaxX * dfUnit;
689
0
        dfMinY = dfMinY * dfUnit;
690
0
        dfMaxY = dfMaxY * dfUnit;
691
692
0
        dfYPixSz = (dfMinY - dfMaxY) / poDS->nRasterYSize;
693
0
        dfXPixSz = (dfMaxX - dfMinX) / poDS->nRasterXSize;
694
695
0
        poDS->m_gt[0] = dfMinX;
696
0
        poDS->m_gt[1] = dfXPixSz;
697
0
        poDS->m_gt[2] = 0.0;
698
0
        poDS->m_gt[3] = dfMaxY;
699
0
        poDS->m_gt[4] = 0.0;
700
0
        poDS->m_gt[5] = dfYPixSz;
701
0
    }
702
703
    // --------------------------------------------------------------------
704
    //      Set Color Table in the presence of a smp file
705
    // --------------------------------------------------------------------
706
707
0
    if (poDS->nBands != 3)
708
0
    {
709
0
        const std::string osSMPFilename =
710
0
            CPLResetExtensionSafe(poDS->pszFilename, extSMP);
711
0
        VSILFILE *fpSMP = VSIFOpenL(osSMPFilename.c_str(), "rb");
712
0
        if (fpSMP != nullptr)
713
0
        {
714
0
            int dfMaxValue =
715
0
                atoi_nz(myCSLFetchNameValue(poDS->papszRDC, rdcMAX_VALUE));
716
0
            int nCatCount =
717
0
                atoi_nz(myCSLFetchNameValue(poDS->papszRDC, rdcLEGEND_CATS));
718
0
            if (nCatCount == 0)
719
0
                dfMaxValue = 255;
720
0
            VSIFSeekL(fpSMP, smpHEADERSIZE, SEEK_SET);
721
0
            GDALColorEntry oEntry;
722
0
            unsigned char aucRGB[3];
723
0
            int i = 0;
724
0
            while ((VSIFReadL(&aucRGB, sizeof(aucRGB), 1, fpSMP)) &&
725
0
                   (i <= dfMaxValue))
726
0
            {
727
0
                oEntry.c1 = (short)aucRGB[0];
728
0
                oEntry.c2 = (short)aucRGB[1];
729
0
                oEntry.c3 = (short)aucRGB[2];
730
0
                oEntry.c4 = (short)255;
731
0
                poDS->poColorTable->SetColorEntry(i, &oEntry);
732
0
                i++;
733
0
            }
734
0
            VSIFCloseL(fpSMP);
735
0
        }
736
0
    }
737
738
    // --------------------------------------------------------------------
739
    //      Check for Unit Type
740
    // --------------------------------------------------------------------
741
742
0
    const char *pszValueUnit =
743
0
        myCSLFetchNameValue(poDS->papszRDC, rdcVALUE_UNITS);
744
745
0
    if (pszValueUnit == nullptr)
746
0
        poDS->pszUnitType = CPLStrdup("unspecified");
747
0
    else
748
0
    {
749
0
        if (STARTS_WITH_CI(pszValueUnit, "meter"))
750
0
        {
751
0
            poDS->pszUnitType = CPLStrdup("m");
752
0
        }
753
0
        else if (STARTS_WITH_CI(pszValueUnit, "feet"))
754
0
        {
755
0
            poDS->pszUnitType = CPLStrdup("ft");
756
0
        }
757
0
        else
758
0
            poDS->pszUnitType = CPLStrdup(pszValueUnit);
759
0
    }
760
761
    // --------------------------------------------------------------------
762
    //      Check for category names.
763
    // --------------------------------------------------------------------
764
765
0
    int nCatCount =
766
0
        atoi_nz(myCSLFetchNameValue(poDS->papszRDC, rdcLEGEND_CATS));
767
768
0
    if (nCatCount > 0)
769
0
    {
770
        // ----------------------------------------------------------------
771
        //      Sequentialize categories names, from 0 to the last "code n"
772
        // ----------------------------------------------------------------
773
774
0
        int nLine = -1;
775
0
        for (int i = 0; (i < CSLCount(poDS->papszRDC)) && (nLine == -1); i++)
776
0
            if (EQUALN(poDS->papszRDC[i], rdcLEGEND_CATS, 11))
777
0
                nLine = i;  // get the line where legend cats is
778
779
0
        if (nLine > 0)
780
0
        {
781
0
            int nCode = 0;
782
0
            int nCount = 0;
783
0
            sscanf(poDS->papszRDC[++nLine], rdcCODE_N,
784
0
                   &nCode);  // assign legend cats to nCode
785
0
            for (int i = 0; (i < 255) && (nCount < nCatCount); i++)
786
0
            {
787
0
                if (i == nCode)
788
0
                {
789
0
                    poDS->papszCategories = CSLAddString(
790
0
                        poDS->papszCategories,
791
0
                        CPLParseNameValue(poDS->papszRDC[nLine], nullptr));
792
0
                    nCount++;
793
0
                    if (nCount < nCatCount)
794
0
                        sscanf(poDS->papszRDC[++nLine], rdcCODE_N, &nCode);
795
0
                }
796
0
                else
797
0
                    poDS->papszCategories =
798
0
                        CSLAddString(poDS->papszCategories, "");
799
0
            }
800
0
        }
801
0
    }
802
803
    /* -------------------------------------------------------------------- */
804
    /*      Automatic Generated Color Table                                 */
805
    /* -------------------------------------------------------------------- */
806
807
0
    if (poDS->papszCategories != nullptr &&
808
0
        (poDS->poColorTable->GetColorEntryCount() == 0))
809
0
    {
810
0
        int nEntryCount = CSLCount(poDS->papszCategories);
811
812
0
        GDALColorEntry sFromColor;
813
0
        sFromColor.c1 = (short)(255);
814
0
        sFromColor.c2 = (short)(0);
815
0
        sFromColor.c3 = (short)(0);
816
0
        sFromColor.c4 = (short)(255);
817
818
0
        GDALColorEntry sToColor;
819
0
        sToColor.c1 = (short)(0);
820
0
        sToColor.c2 = (short)(0);
821
0
        sToColor.c3 = (short)(255);
822
0
        sToColor.c4 = (short)(255);
823
824
0
        poDS->poColorTable->CreateColorRamp(0, &sFromColor, (nEntryCount - 1),
825
0
                                            &sToColor);
826
0
    }
827
828
    /* -------------------------------------------------------------------- */
829
    /*      Initialize any PAM information.                                 */
830
    /* -------------------------------------------------------------------- */
831
832
0
    poDS->SetDescription(poOpenInfo->pszFilename);
833
0
    poDS->TryLoadXML();
834
835
    /* -------------------------------------------------------------------- */
836
    /*      Check for external overviews.                                   */
837
    /* -------------------------------------------------------------------- */
838
839
0
    poDS->oOvManager.Initialize(poDS, poOpenInfo->pszFilename);
840
841
0
    return poDS;
842
0
}
843
844
/************************************************************************/
845
/*                               Create()                               */
846
/************************************************************************/
847
848
GDALDataset *IdrisiDataset::Create(const char *pszFilename, int nXSize,
849
                                   int nYSize, int nBandsIn, GDALDataType eType,
850
                                   char ** /* papszOptions */)
851
0
{
852
    // --------------------------------------------------------------------
853
    //      Check input options
854
    // --------------------------------------------------------------------
855
856
0
    if (nBandsIn != 1 && nBandsIn != 3)
857
0
    {
858
0
        CPLError(CE_Failure, CPLE_AppDefined,
859
0
                 "Attempt to create IDRISI dataset with an illegal number of "
860
0
                 "bands(%d)."
861
0
                 " Try again by selecting a specific band if possible. \n",
862
0
                 nBandsIn);
863
0
        return nullptr;
864
0
    }
865
866
0
    if (nBandsIn == 3 && eType != GDT_UInt8)
867
0
    {
868
0
        CPLError(
869
0
            CE_Failure, CPLE_AppDefined,
870
0
            "Attempt to create IDRISI dataset with an unsupported combination "
871
0
            "of the number of bands(%d) and data type(%s). \n",
872
0
            nBandsIn, GDALGetDataTypeName(eType));
873
0
        return nullptr;
874
0
    }
875
876
    // ----------------------------------------------------------------
877
    //  Create the header file with minimum information
878
    // ----------------------------------------------------------------
879
880
0
    const char *pszLDataType = nullptr;
881
882
0
    switch (eType)
883
0
    {
884
0
        case GDT_UInt8:
885
0
            if (nBandsIn == 1)
886
0
                pszLDataType = rstBYTE;
887
0
            else
888
0
                pszLDataType = rstRGB24;
889
0
            break;
890
0
        case GDT_Int16:
891
0
            pszLDataType = rstINTEGER;
892
0
            break;
893
0
        case GDT_Float32:
894
0
            pszLDataType = rstREAL;
895
0
            break;
896
        //--- process compatible data types
897
0
        case (GDT_UInt16):
898
0
            pszLDataType = rstINTEGER;
899
0
            CPLError(CE_Warning, CPLE_AppDefined,
900
0
                     "This process requires a conversion from %s to signed "
901
0
                     "16-bit %s, "
902
0
                     "which may cause data loss.\n",
903
0
                     GDALGetDataTypeName(eType), rstINTEGER);
904
0
            break;
905
0
        case GDT_UInt32:
906
0
            pszLDataType = rstINTEGER;
907
0
            CPLError(CE_Warning, CPLE_AppDefined,
908
0
                     "This process requires a conversion from %s to signed "
909
0
                     "16-bit %s, "
910
0
                     "which may cause data loss.\n",
911
0
                     GDALGetDataTypeName(eType), rstINTEGER);
912
0
            break;
913
0
        case GDT_Int32:
914
0
            pszLDataType = rstINTEGER;
915
0
            CPLError(CE_Warning, CPLE_AppDefined,
916
0
                     "This process requires a conversion from %s to signed "
917
0
                     "16-bit %s, "
918
0
                     "which may cause data loss.\n",
919
0
                     GDALGetDataTypeName(eType), rstINTEGER);
920
0
            break;
921
0
        case GDT_Float64:
922
0
            pszLDataType = rstREAL;
923
0
            CPLError(CE_Warning, CPLE_AppDefined,
924
0
                     "This process requires a conversion from %s to float "
925
0
                     "32-bit %s, "
926
0
                     "which may cause data loss.\n",
927
0
                     GDALGetDataTypeName(eType), rstREAL);
928
0
            break;
929
0
        default:
930
0
            CPLError(CE_Failure, CPLE_AppDefined,
931
0
                     "Attempt to create IDRISI dataset with an illegal "
932
0
                     "data type(%s).\n",
933
0
                     GDALGetDataTypeName(eType));
934
0
            return nullptr;
935
0
    };
936
937
0
    char **papszLRDC = nullptr;
938
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcFILE_FORMAT, rstVERSION);
939
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcFILE_TITLE, "");
940
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcDATA_TYPE, pszLDataType);
941
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcFILE_TYPE, "binary");
942
0
    papszLRDC =
943
0
        CSLAddNameValue(papszLRDC, rdcCOLUMNS, CPLSPrintf("%d", nXSize));
944
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcROWS, CPLSPrintf("%d", nYSize));
945
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcREF_SYSTEM, "plane");
946
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcREF_UNITS, "m");
947
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcUNIT_DIST, "1");
948
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMIN_X, "0");
949
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMAX_X, CPLSPrintf("%d", nXSize));
950
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMIN_Y, "0");
951
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMAX_Y, CPLSPrintf("%d", nYSize));
952
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcPOSN_ERROR, "unspecified");
953
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcRESOLUTION, "1.0");
954
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMIN_VALUE, "0");
955
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcMAX_VALUE, "0");
956
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcDISPLAY_MIN, "0");
957
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcDISPLAY_MAX, "0");
958
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcVALUE_UNITS, "unspecified");
959
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcVALUE_ERROR, "unspecified");
960
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcFLAG_VALUE, "none");
961
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcFLAG_DEFN, "none");
962
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcLEGEND_CATS, "0");
963
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcLINEAGES, "");
964
0
    papszLRDC = CSLAddNameValue(papszLRDC, rdcCOMMENTS, "");
965
966
0
    const std::string osLDocFilename =
967
0
        CPLResetExtensionSafe(pszFilename, extRDC);
968
969
0
    myCSLSetNameValueSeparator(papszLRDC, ": ");
970
0
    SaveAsCRLF(papszLRDC, osLDocFilename.c_str());
971
0
    CSLDestroy(papszLRDC);
972
973
    // ----------------------------------------------------------------
974
    //  Create an empty data file
975
    // ----------------------------------------------------------------
976
977
0
    VSILFILE *fp = VSIFOpenL(pszFilename, "wb+");
978
979
0
    if (fp == nullptr)
980
0
    {
981
0
        CPLError(CE_Failure, CPLE_OpenFailed,
982
0
                 "Attempt to create file %s' failed.\n", pszFilename);
983
0
        return nullptr;
984
0
    }
985
986
0
    const int nTargetDTSize = EQUAL(pszLDataType, rstBYTE)      ? 1
987
0
                              : EQUAL(pszLDataType, rstINTEGER) ? 2
988
0
                              : EQUAL(pszLDataType, rstRGB24)   ? 3
989
0
                                                                : 4;
990
0
    VSIFTruncateL(fp,
991
0
                  static_cast<vsi_l_offset>(nXSize) * nYSize * nTargetDTSize);
992
0
    VSIFCloseL(fp);
993
994
0
    return (IdrisiDataset *)GDALOpen(pszFilename, GA_Update);
995
0
}
996
997
/************************************************************************/
998
/*                             CreateCopy()                             */
999
/************************************************************************/
1000
1001
GDALDataset *IdrisiDataset::CreateCopy(const char *pszFilename,
1002
                                       GDALDataset *poSrcDS, int bStrict,
1003
                                       char **papszOptions,
1004
                                       GDALProgressFunc pfnProgress,
1005
                                       void *pProgressData)
1006
0
{
1007
0
    if (!pfnProgress(0.0, nullptr, pProgressData))
1008
0
        return nullptr;
1009
1010
    // ------------------------------------------------------------------------
1011
    //      Check number of bands
1012
    // ------------------------------------------------------------------------
1013
0
    if (!(poSrcDS->GetRasterCount() == 1) && !(poSrcDS->GetRasterCount() == 3))
1014
0
    {
1015
0
        CPLError(CE_Failure, CPLE_AppDefined,
1016
0
                 "Attempt to create IDRISI dataset with an illegal number of "
1017
0
                 "bands(%d)."
1018
0
                 " Try again by selecting a specific band if possible.\n",
1019
0
                 poSrcDS->GetRasterCount());
1020
0
        return nullptr;
1021
0
    }
1022
0
    if ((poSrcDS->GetRasterCount() == 3) &&
1023
0
        ((poSrcDS->GetRasterBand(1)->GetRasterDataType() != GDT_UInt8) ||
1024
0
         (poSrcDS->GetRasterBand(2)->GetRasterDataType() != GDT_UInt8) ||
1025
0
         (poSrcDS->GetRasterBand(3)->GetRasterDataType() != GDT_UInt8)))
1026
0
    {
1027
0
        CPLError(
1028
0
            CE_Failure, CPLE_AppDefined,
1029
0
            "Attempt to create IDRISI dataset with an unsupported "
1030
0
            "data type when there are three bands. Only BYTE allowed.\n"
1031
0
            "Try again by selecting a specific band to convert if possible.\n");
1032
0
        return nullptr;
1033
0
    }
1034
1035
    // ------------------------------------------------------------------------
1036
    //      Check Data types
1037
    // ------------------------------------------------------------------------
1038
1039
0
    for (int i = 1; i <= poSrcDS->GetRasterCount(); i++)
1040
0
    {
1041
0
        GDALDataType eType = poSrcDS->GetRasterBand(i)->GetRasterDataType();
1042
1043
0
        if (bStrict)
1044
0
        {
1045
0
            if (eType != GDT_UInt8 && eType != GDT_Int16 &&
1046
0
                eType != GDT_Float32)
1047
0
            {
1048
0
                CPLError(CE_Failure, CPLE_AppDefined,
1049
0
                         "Attempt to create IDRISI dataset in strict mode "
1050
0
                         "with an illegal data type(%s).\n",
1051
0
                         GDALGetDataTypeName(eType));
1052
0
                return nullptr;
1053
0
            }
1054
0
        }
1055
0
        else
1056
0
        {
1057
0
            if (eType != GDT_UInt8 && eType != GDT_Int16 &&
1058
0
                eType != GDT_UInt16 && eType != GDT_UInt32 &&
1059
0
                eType != GDT_Int32 && eType != GDT_Float32 &&
1060
0
                eType != GDT_Float64)
1061
0
            {
1062
0
                CPLError(CE_Failure, CPLE_AppDefined,
1063
0
                         "Attempt to create IDRISI dataset with an illegal "
1064
0
                         "data type(%s).\n",
1065
0
                         GDALGetDataTypeName(eType));
1066
0
                return nullptr;
1067
0
            }
1068
0
        }
1069
0
    }
1070
1071
    // --------------------------------------------------------------------
1072
    //      Define data type
1073
    // --------------------------------------------------------------------
1074
1075
0
    GDALRasterBand *poBand = poSrcDS->GetRasterBand(1);
1076
0
    GDALDataType eType = poBand->GetRasterDataType();
1077
1078
0
    int bSuccessMin = FALSE;
1079
0
    int bSuccessMax = FALSE;
1080
1081
0
    double dfMin = poBand->GetMinimum(&bSuccessMin);
1082
0
    double dfMax = poBand->GetMaximum(&bSuccessMax);
1083
1084
0
    if (!(bSuccessMin && bSuccessMax))
1085
0
    {
1086
0
        poBand->GetStatistics(false, true, &dfMin, &dfMax, nullptr, nullptr);
1087
0
    }
1088
1089
0
    if (!((eType == GDT_UInt8) || (eType == GDT_Int16) ||
1090
0
          (eType == GDT_Float32)))
1091
0
    {
1092
0
        if (eType == GDT_Float64)
1093
0
        {
1094
0
            eType = GDT_Float32;
1095
0
        }
1096
0
        else
1097
0
        {
1098
0
            if ((dfMin < (double)SHRT_MIN) || (dfMax > (double)SHRT_MAX))
1099
0
            {
1100
0
                eType = GDT_Float32;
1101
0
            }
1102
0
            else
1103
0
            {
1104
0
                eType = GDT_Int16;
1105
0
            }
1106
0
        }
1107
0
    }
1108
1109
    // --------------------------------------------------------------------
1110
    //      Create the dataset
1111
    // --------------------------------------------------------------------
1112
1113
0
    IdrisiDataset *poDS = cpl::down_cast<IdrisiDataset *>(IdrisiDataset::Create(
1114
0
        pszFilename, poSrcDS->GetRasterXSize(), poSrcDS->GetRasterYSize(),
1115
0
        poSrcDS->GetRasterCount(), eType, papszOptions));
1116
1117
0
    if (poDS == nullptr)
1118
0
        return nullptr;
1119
1120
    // --------------------------------------------------------------------
1121
    //      Copy information to the dataset
1122
    // --------------------------------------------------------------------
1123
1124
0
    GDALGeoTransform gt;
1125
0
    if (poSrcDS->GetGeoTransform(gt) == CE_None)
1126
0
    {
1127
0
        poDS->SetGeoTransform(gt);
1128
0
    }
1129
1130
0
    if (!EQUAL(poSrcDS->GetProjectionRef(), ""))
1131
0
    {
1132
0
        poDS->SetProjection(poSrcDS->GetProjectionRef());
1133
0
    }
1134
1135
    // --------------------------------------------------------------------
1136
    //      Copy information to the raster band(s)
1137
    // --------------------------------------------------------------------
1138
1139
0
    for (int i = 1; i <= poDS->nBands; i++)
1140
0
    {
1141
0
        GDALRasterBand *poSrcBand = poSrcDS->GetRasterBand(i);
1142
0
        IdrisiRasterBand *poDstBand =
1143
0
            cpl::down_cast<IdrisiRasterBand *>(poDS->GetRasterBand(i));
1144
1145
0
        if (poDS->nBands == 1)
1146
0
        {
1147
0
            poDstBand->SetUnitType(poSrcBand->GetUnitType());
1148
0
            poDstBand->SetColorTable(poSrcBand->GetColorTable());
1149
0
            poDstBand->SetCategoryNames(poSrcBand->GetCategoryNames());
1150
1151
0
            const GDALRasterAttributeTable *poRAT = poSrcBand->GetDefaultRAT();
1152
1153
0
            if (poRAT != nullptr)
1154
0
            {
1155
0
                poDstBand->SetDefaultRAT(poRAT);
1156
0
            }
1157
0
        }
1158
1159
0
        dfMin = poSrcBand->GetMinimum(nullptr);
1160
0
        dfMax = poSrcBand->GetMaximum(nullptr);
1161
0
        poDstBand->SetMinMax(dfMin, dfMax);
1162
0
        int bHasNoDataValue;
1163
0
        double dfNoDataValue = poSrcBand->GetNoDataValue(&bHasNoDataValue);
1164
0
        if (bHasNoDataValue)
1165
0
            poDstBand->SetNoDataValue(dfNoDataValue);
1166
0
    }
1167
1168
    // --------------------------------------------------------------------
1169
    //      Copy image data
1170
    // --------------------------------------------------------------------
1171
1172
0
    if (GDALDatasetCopyWholeRaster((GDALDatasetH)poSrcDS, (GDALDatasetH)poDS,
1173
0
                                   nullptr, pfnProgress,
1174
0
                                   pProgressData) != CE_None)
1175
0
    {
1176
0
        delete poDS;
1177
0
        return nullptr;
1178
0
    }
1179
1180
    // --------------------------------------------------------------------
1181
    //      Finalize
1182
    // --------------------------------------------------------------------
1183
1184
0
    poDS->FlushCache(false);
1185
1186
0
    return poDS;
1187
0
}
1188
1189
/************************************************************************/
1190
/*                            GetFileList()                             */
1191
/************************************************************************/
1192
1193
char **IdrisiDataset::GetFileList()
1194
0
{
1195
0
    char **papszFileList = GDALPamDataset::GetFileList();
1196
1197
    // --------------------------------------------------------------------
1198
    //      Symbol table file
1199
    // --------------------------------------------------------------------
1200
1201
0
    std::string osAssociated = CPLResetExtensionSafe(pszFilename, extSMP);
1202
1203
0
    if (FileExists(osAssociated.c_str()))
1204
0
    {
1205
0
        papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1206
0
    }
1207
0
    else
1208
0
    {
1209
0
        osAssociated = CPLResetExtensionSafe(pszFilename, extSMPu);
1210
1211
0
        if (FileExists(osAssociated.c_str()))
1212
0
        {
1213
0
            papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1214
0
        }
1215
0
    }
1216
1217
    // --------------------------------------------------------------------
1218
    //      Documentation file
1219
    // --------------------------------------------------------------------
1220
1221
0
    osAssociated = CPLResetExtensionSafe(pszFilename, extRDC);
1222
1223
0
    if (FileExists(osAssociated.c_str()))
1224
0
    {
1225
0
        papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1226
0
    }
1227
0
    else
1228
0
    {
1229
0
        osAssociated = CPLResetExtensionSafe(pszFilename, extRDCu);
1230
1231
0
        if (FileExists(osAssociated.c_str()))
1232
0
        {
1233
0
            papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1234
0
        }
1235
0
    }
1236
1237
    // --------------------------------------------------------------------
1238
    //      Reference file
1239
    // --------------------------------------------------------------------
1240
1241
0
    osAssociated = CPLResetExtensionSafe(pszFilename, extREF);
1242
1243
0
    if (FileExists(osAssociated.c_str()))
1244
0
    {
1245
0
        papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1246
0
    }
1247
0
    else
1248
0
    {
1249
0
        osAssociated = CPLResetExtensionSafe(pszFilename, extREFu);
1250
1251
0
        if (FileExists(osAssociated.c_str()))
1252
0
        {
1253
0
            papszFileList = CSLAddString(papszFileList, osAssociated.c_str());
1254
0
        }
1255
0
    }
1256
1257
0
    return papszFileList;
1258
0
}
1259
1260
/************************************************************************/
1261
/*                          GetGeoTransform()                           */
1262
/************************************************************************/
1263
1264
CPLErr IdrisiDataset::GetGeoTransform(GDALGeoTransform &gt) const
1265
0
{
1266
0
    if (GDALPamDataset::GetGeoTransform(gt) != CE_None)
1267
0
    {
1268
0
        gt = m_gt;
1269
0
    }
1270
1271
0
    return CE_None;
1272
0
}
1273
1274
/************************************************************************/
1275
/*                          SetGeoTransform()                           */
1276
/************************************************************************/
1277
1278
CPLErr IdrisiDataset::SetGeoTransform(const GDALGeoTransform &gt)
1279
0
{
1280
0
    if (gt[2] != 0.0 || gt[4] != 0.0)
1281
0
    {
1282
0
        CPLError(CE_Failure, CPLE_AppDefined,
1283
0
                 "Attempt to set rotated geotransform on Idrisi Raster file.\n"
1284
0
                 "Idrisi Raster does not support rotation.\n");
1285
0
        return CE_Failure;
1286
0
    }
1287
1288
    // --------------------------------------------------------------------
1289
    // Update the .rdc file
1290
    // --------------------------------------------------------------------
1291
1292
0
    double dfXPixSz = gt[1];
1293
0
    double dfYPixSz = gt[5];
1294
0
    double dfMinX = gt[0];
1295
0
    double dfMaxX = (dfXPixSz * nRasterXSize) + dfMinX;
1296
1297
0
    double dfMinY, dfMaxY;
1298
0
    if (dfYPixSz < 0)
1299
0
    {
1300
0
        dfMaxY = gt[3];
1301
0
        dfMinY = (dfYPixSz * nRasterYSize) + gt[3];
1302
0
    }
1303
0
    else
1304
0
    {
1305
0
        dfMaxY = (dfYPixSz * nRasterYSize) + gt[3];
1306
0
        dfMinY = gt[3];
1307
0
    }
1308
1309
0
    papszRDC = CSLSetNameValue(papszRDC, rdcMIN_X, CPLSPrintf("%.7f", dfMinX));
1310
0
    papszRDC = CSLSetNameValue(papszRDC, rdcMAX_X, CPLSPrintf("%.7f", dfMaxX));
1311
0
    papszRDC = CSLSetNameValue(papszRDC, rdcMIN_Y, CPLSPrintf("%.7f", dfMinY));
1312
0
    papszRDC = CSLSetNameValue(papszRDC, rdcMAX_Y, CPLSPrintf("%.7f", dfMaxY));
1313
0
    papszRDC = CSLSetNameValue(papszRDC, rdcRESOLUTION,
1314
0
                               CPLSPrintf("%.7f", fabs(dfYPixSz)));
1315
1316
    // --------------------------------------------------------------------
1317
    // Update the Dataset attribute
1318
    // --------------------------------------------------------------------
1319
1320
0
    m_gt = gt;
1321
1322
0
    return CE_None;
1323
0
}
1324
1325
/************************************************************************/
1326
/*                          GetSpatialRef()                             */
1327
/************************************************************************/
1328
1329
const OGRSpatialReference *IdrisiDataset::GetSpatialRef() const
1330
0
{
1331
0
    const auto poSRS = GDALPamDataset::GetSpatialRef();
1332
0
    if (poSRS)
1333
0
        return poSRS;
1334
1335
0
    if (m_oSRS.IsEmpty())
1336
0
    {
1337
0
        const char *pszRefSystem = myCSLFetchNameValue(papszRDC, rdcREF_SYSTEM);
1338
0
        const char *pszRefUnit = myCSLFetchNameValue(papszRDC, rdcREF_UNITS);
1339
0
        if (pszRefSystem != nullptr && pszRefUnit != nullptr)
1340
0
            IdrisiGeoReference2Wkt(pszFilename, pszRefSystem, pszRefUnit,
1341
0
                                   m_oSRS);
1342
0
    }
1343
0
    return m_oSRS.IsEmpty() ? nullptr : &m_oSRS;
1344
0
}
1345
1346
/************************************************************************/
1347
/*                           SetProjection()                            */
1348
/************************************************************************/
1349
1350
CPLErr IdrisiDataset::SetSpatialRef(const OGRSpatialReference *poSRS)
1351
0
{
1352
0
    m_oSRS.Clear();
1353
0
    if (poSRS)
1354
0
        m_oSRS = *poSRS;
1355
1356
0
    char *pszRefSystem = nullptr;
1357
0
    char *pszRefUnit = nullptr;
1358
1359
0
    CPLErr eResult = Wkt2GeoReference(m_oSRS, &pszRefSystem, &pszRefUnit);
1360
1361
0
    papszRDC = CSLSetNameValue(papszRDC, rdcREF_SYSTEM, pszRefSystem);
1362
0
    papszRDC = CSLSetNameValue(papszRDC, rdcREF_UNITS, pszRefUnit);
1363
1364
0
    CPLFree(pszRefSystem);
1365
0
    CPLFree(pszRefUnit);
1366
1367
0
    return eResult;
1368
0
}
1369
1370
/************************************************************************/
1371
/*                          IdrisiRasterBand()                          */
1372
/************************************************************************/
1373
1374
IdrisiRasterBand::IdrisiRasterBand(IdrisiDataset *poDSIn, int nBandIn,
1375
                                   GDALDataType eDataTypeIn)
1376
0
    : poDefaultRAT(nullptr),
1377
0
      nRecordSize(poDSIn->GetRasterXSize() * poDSIn->nBands *
1378
0
                  GDALGetDataTypeSizeBytes(eDataTypeIn)),
1379
0
      pabyScanLine(static_cast<GByte *>(VSI_MALLOC3_VERBOSE(
1380
0
          poDSIn->GetRasterXSize(), GDALGetDataTypeSizeBytes(eDataTypeIn),
1381
0
          poDSIn->nBands))),
1382
0
      fMaximum(0.0), fMinimum(0.0), bFirstVal(true)
1383
0
{
1384
0
    poDS = poDSIn;
1385
0
    nBand = nBandIn;
1386
0
    eDataType = eDataTypeIn;
1387
0
    nBlockYSize = 1;
1388
0
    nBlockXSize = poDS->GetRasterXSize();
1389
0
}
1390
1391
/************************************************************************/
1392
/*                         ~IdrisiRasterBand()                          */
1393
/************************************************************************/
1394
1395
IdrisiRasterBand::~IdrisiRasterBand()
1396
0
{
1397
0
    CPLFree(pabyScanLine);
1398
1399
0
    if (poDefaultRAT)
1400
0
    {
1401
0
        delete poDefaultRAT;
1402
0
    }
1403
0
}
1404
1405
/************************************************************************/
1406
/*                             IReadBlock()                             */
1407
/************************************************************************/
1408
1409
CPLErr IdrisiRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
1410
                                    void *pImage)
1411
0
{
1412
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1413
1414
0
    if (VSIFSeekL(poGDS->fp, vsi_l_offset(nRecordSize) * nBlockYOff, SEEK_SET) <
1415
0
        0)
1416
0
    {
1417
0
        CPLError(CE_Failure, CPLE_FileIO,
1418
0
                 "Can't seek(%s) block with X offset %d and Y offset %d.\n%s",
1419
0
                 poGDS->pszFilename, nBlockXOff, nBlockYOff,
1420
0
                 VSIStrerror(errno));
1421
0
        return CE_Failure;
1422
0
    }
1423
1424
0
    if ((int)VSIFReadL(pabyScanLine, 1, nRecordSize, poGDS->fp) < nRecordSize)
1425
0
    {
1426
0
        CPLError(CE_Failure, CPLE_FileIO,
1427
0
                 "Can't read(%s) block with X offset %d and Y offset %d.\n%s",
1428
0
                 poGDS->pszFilename, nBlockXOff, nBlockYOff,
1429
0
                 VSIStrerror(errno));
1430
0
        return CE_Failure;
1431
0
    }
1432
1433
0
    if (poGDS->nBands == 3)
1434
0
    {
1435
0
        for (int i = 0, j = (3 - nBand); i < nBlockXSize; i++, j += 3)
1436
0
        {
1437
0
            ((GByte *)pImage)[i] = pabyScanLine[j];
1438
0
        }
1439
0
    }
1440
0
    else
1441
0
    {
1442
0
        memcpy(pImage, pabyScanLine, nRecordSize);
1443
0
    }
1444
1445
#ifdef CPL_MSB
1446
    if (eDataType == GDT_Float32)
1447
        GDALSwapWords(pImage, 4, nBlockXSize * nBlockYSize, 4);
1448
#endif
1449
1450
0
    return CE_None;
1451
0
}
1452
1453
/************************************************************************/
1454
/*                            IWriteBlock()                             */
1455
/************************************************************************/
1456
1457
CPLErr IdrisiRasterBand::IWriteBlock(int nBlockXOff, int nBlockYOff,
1458
                                     void *pImage)
1459
0
{
1460
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1461
1462
#ifdef CPL_MSB
1463
    // Swap in input buffer if needed.
1464
    if (eDataType == GDT_Float32)
1465
        GDALSwapWords(pImage, 4, nBlockXSize * nBlockYSize, 4);
1466
#endif
1467
1468
0
    if (poGDS->nBands == 1)
1469
0
    {
1470
0
        memcpy(pabyScanLine, pImage, nRecordSize);
1471
0
    }
1472
0
    else
1473
0
    {
1474
0
        if (nBand > 1)
1475
0
        {
1476
0
            VSIFSeekL(poGDS->fp, vsi_l_offset(nRecordSize) * nBlockYOff,
1477
0
                      SEEK_SET);
1478
0
            VSIFReadL(pabyScanLine, 1, nRecordSize, poGDS->fp);
1479
0
        }
1480
0
        int i, j;
1481
0
        for (i = 0, j = (3 - nBand); i < nBlockXSize; i++, j += 3)
1482
0
        {
1483
0
            pabyScanLine[j] = ((GByte *)pImage)[i];
1484
0
        }
1485
0
    }
1486
1487
#ifdef CPL_MSB
1488
    // Swap input buffer back to original form.
1489
    if (eDataType == GDT_Float32)
1490
        GDALSwapWords(pImage, 4, nBlockXSize * nBlockYSize, 4);
1491
#endif
1492
1493
0
    VSIFSeekL(poGDS->fp, vsi_l_offset(nRecordSize) * nBlockYOff, SEEK_SET);
1494
1495
0
    if ((int)VSIFWriteL(pabyScanLine, 1, nRecordSize, poGDS->fp) < nRecordSize)
1496
0
    {
1497
0
        CPLError(CE_Failure, CPLE_FileIO,
1498
0
                 "Can't write(%s) block with X offset %d and Y offset %d.\n%s",
1499
0
                 poGDS->pszFilename, nBlockXOff, nBlockYOff,
1500
0
                 VSIStrerror(errno));
1501
0
        return CE_Failure;
1502
0
    }
1503
1504
0
    int bHasNoDataValue = FALSE;
1505
0
    float fNoDataValue = (float)GetNoDataValue(&bHasNoDataValue);
1506
1507
    // --------------------------------------------------------------------
1508
    //      Search for the minimum and maximum values
1509
    // --------------------------------------------------------------------
1510
1511
0
    if (eDataType == GDT_Float32)
1512
0
    {
1513
0
        for (int i = 0; i < nBlockXSize; i++)
1514
0
        {
1515
0
            float fVal =
1516
0
                reinterpret_cast<float *>(pabyScanLine)[i];  // this is fine
1517
0
            if (!bHasNoDataValue || fVal != fNoDataValue)
1518
0
            {
1519
0
                if (bFirstVal)
1520
0
                {
1521
0
                    fMinimum = fVal;
1522
0
                    fMaximum = fVal;
1523
0
                    bFirstVal = false;
1524
0
                }
1525
0
                else
1526
0
                {
1527
0
                    if (fVal < fMinimum)
1528
0
                        fMinimum = fVal;
1529
0
                    if (fVal > fMaximum)
1530
0
                        fMaximum = fVal;
1531
0
                }
1532
0
            }
1533
0
        }
1534
0
    }
1535
0
    else if (eDataType == GDT_Int16)
1536
0
    {
1537
0
        for (int i = 0; i < nBlockXSize; i++)
1538
0
        {
1539
0
            float fVal = (float)(reinterpret_cast<GInt16 *>(pabyScanLine))[i];
1540
0
            if (!bHasNoDataValue || fVal != fNoDataValue)
1541
0
            {
1542
0
                if (bFirstVal)
1543
0
                {
1544
0
                    fMinimum = fVal;
1545
0
                    fMaximum = fVal;
1546
0
                    bFirstVal = false;
1547
0
                }
1548
0
                else
1549
0
                {
1550
0
                    if (fVal < fMinimum)
1551
0
                        fMinimum = fVal;
1552
0
                    if (fVal > fMaximum)
1553
0
                        fMaximum = fVal;
1554
0
                }
1555
0
            }
1556
0
        }
1557
0
    }
1558
0
    else if (poGDS->nBands == 1)
1559
0
    {
1560
0
        for (int i = 0; i < nBlockXSize; i++)
1561
0
        {
1562
0
            float fVal = (float)((GByte *)pabyScanLine)[i];
1563
0
            if (!bHasNoDataValue || fVal != fNoDataValue)
1564
0
            {
1565
0
                if (bFirstVal)
1566
0
                {
1567
0
                    fMinimum = fVal;
1568
0
                    fMaximum = fVal;
1569
0
                    bFirstVal = false;
1570
0
                }
1571
0
                else
1572
0
                {
1573
0
                    if (fVal < fMinimum)
1574
0
                        fMinimum =
1575
0
                            fVal;  // I don't change this part, keep it as it is
1576
0
                    if (fVal > fMaximum)
1577
0
                        fMaximum = fVal;
1578
0
                }
1579
0
            }
1580
0
        }
1581
0
    }
1582
0
    else
1583
0
    {
1584
0
        for (int i = 0, j = (3 - nBand); i < nBlockXSize; i++, j += 3)
1585
0
        {
1586
0
            float fVal = (float)((GByte *)pabyScanLine)[j];
1587
0
            if (!bHasNoDataValue || fVal != fNoDataValue)
1588
0
            {
1589
0
                if (bFirstVal)
1590
0
                {
1591
0
                    fMinimum = fVal;
1592
0
                    fMaximum = fVal;
1593
0
                    bFirstVal = false;
1594
0
                }
1595
0
                else
1596
0
                {
1597
0
                    if (fVal < fMinimum)
1598
0
                        fMinimum = fVal;
1599
0
                    if (fVal > fMaximum)
1600
0
                        fMaximum = fVal;
1601
0
                }
1602
0
            }
1603
0
        }
1604
0
    }
1605
1606
0
    return CE_None;
1607
0
}
1608
1609
/************************************************************************/
1610
/*                             GetMinimum()                             */
1611
/************************************************************************/
1612
1613
double IdrisiRasterBand::GetMinimum(int *pbSuccess)
1614
0
{
1615
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1616
1617
0
    if (myCSLFetchNameValue(poGDS->papszRDC, rdcMIN_VALUE) == nullptr)
1618
0
        return GDALPamRasterBand::GetMinimum(pbSuccess);
1619
1620
0
    double adfMinValue[3];
1621
0
    CPLsscanf(myCSLFetchNameValue(poGDS->papszRDC, rdcMIN_VALUE), "%lf %lf %lf",
1622
0
              &adfMinValue[0], &adfMinValue[1], &adfMinValue[2]);
1623
1624
0
    if (pbSuccess)
1625
0
    {
1626
0
        *pbSuccess = true;
1627
0
    }
1628
1629
0
    return adfMinValue[this->nBand - 1];
1630
0
}
1631
1632
/************************************************************************/
1633
/*                             GetMaximum()                             */
1634
/************************************************************************/
1635
1636
double IdrisiRasterBand::GetMaximum(int *pbSuccess)
1637
0
{
1638
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1639
1640
0
    if (myCSLFetchNameValue(poGDS->papszRDC, rdcMAX_VALUE) == nullptr)
1641
0
        return GDALPamRasterBand::GetMaximum(pbSuccess);
1642
1643
0
    double adfMaxValue[3];
1644
0
    CPLsscanf(myCSLFetchNameValue(poGDS->papszRDC, rdcMAX_VALUE), "%lf %lf %lf",
1645
0
              &adfMaxValue[0], &adfMaxValue[1], &adfMaxValue[2]);
1646
1647
0
    if (pbSuccess)
1648
0
    {
1649
0
        *pbSuccess = true;
1650
0
    }
1651
1652
0
    return adfMaxValue[this->nBand - 1];
1653
0
}
1654
1655
/************************************************************************/
1656
/*                           GetNoDataValue()                           */
1657
/************************************************************************/
1658
1659
double IdrisiRasterBand::GetNoDataValue(int *pbSuccess)
1660
0
{
1661
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1662
1663
0
    const char *pszFlagDefn = nullptr;
1664
1665
0
    if (myCSLFetchNameValue(poGDS->papszRDC, rdcFLAG_DEFN) != nullptr)
1666
0
        pszFlagDefn = myCSLFetchNameValue(poGDS->papszRDC, rdcFLAG_DEFN);
1667
0
    else if (myCSLFetchNameValue(poGDS->papszRDC, rdcFLAG_DEFN2) != nullptr)
1668
0
        pszFlagDefn = myCSLFetchNameValue(poGDS->papszRDC, rdcFLAG_DEFN2);
1669
1670
    // ------------------------------------------------------------------------
1671
    // If Flag_Def is not "none", Flag_Value means "background"
1672
    // or "missing data"
1673
    // ------------------------------------------------------------------------
1674
1675
0
    double dfNoData;
1676
0
    if (pszFlagDefn != nullptr && !EQUAL(pszFlagDefn, "none"))
1677
0
    {
1678
0
        dfNoData =
1679
0
            CPLAtof_nz(myCSLFetchNameValue(poGDS->papszRDC, rdcFLAG_VALUE));
1680
0
        if (pbSuccess)
1681
0
            *pbSuccess = TRUE;
1682
0
    }
1683
0
    else
1684
0
    {
1685
0
        dfNoData = -9999.0; /* this value should be ignored */
1686
0
        if (pbSuccess)
1687
0
            *pbSuccess = FALSE;
1688
0
    }
1689
1690
0
    return dfNoData;
1691
0
}
1692
1693
/************************************************************************/
1694
/*                           SetNoDataValue()                           */
1695
/************************************************************************/
1696
1697
CPLErr IdrisiRasterBand::SetNoDataValue(double dfNoDataValue)
1698
0
{
1699
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1700
1701
0
    poGDS->papszRDC = CSLSetNameValue(poGDS->papszRDC, rdcFLAG_VALUE,
1702
0
                                      CPLSPrintf("%.7g", dfNoDataValue));
1703
0
    poGDS->papszRDC =
1704
0
        CSLSetNameValue(poGDS->papszRDC, rdcFLAG_DEFN, "missing data");
1705
1706
0
    return CE_None;
1707
0
}
1708
1709
/************************************************************************/
1710
/*                       GetColorInterpretation()                       */
1711
/************************************************************************/
1712
1713
GDALColorInterp IdrisiRasterBand::GetColorInterpretation()
1714
0
{
1715
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1716
1717
0
    if (poGDS->nBands == 3)
1718
0
    {
1719
0
        switch (nBand)
1720
0
        {
1721
0
            case 1:
1722
0
                return GCI_BlueBand;
1723
0
            case 2:
1724
0
                return GCI_GreenBand;
1725
0
            case 3:
1726
0
                return GCI_RedBand;
1727
0
        }
1728
0
    }
1729
0
    else if (poGDS->poColorTable->GetColorEntryCount() > 0)
1730
0
    {
1731
0
        return GCI_PaletteIndex;
1732
0
    }
1733
0
    return GCI_GrayIndex;
1734
0
}
1735
1736
/************************************************************************/
1737
/*                          GetCategoryNames()                          */
1738
/************************************************************************/
1739
1740
char **IdrisiRasterBand::GetCategoryNames()
1741
0
{
1742
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1743
1744
0
    return poGDS->papszCategories;
1745
0
}
1746
1747
/************************************************************************/
1748
/*                          SetCategoryNames()                          */
1749
/************************************************************************/
1750
1751
CPLErr IdrisiRasterBand::SetCategoryNames(char **papszCategoryNames)
1752
0
{
1753
0
    const int nCatCount = CSLCount(papszCategoryNames);
1754
1755
0
    if (nCatCount == 0)
1756
0
        return CE_None;
1757
1758
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1759
1760
0
    CSLDestroy(poGDS->papszCategories);
1761
0
    poGDS->papszCategories = CSLDuplicate(papszCategoryNames);
1762
1763
    // ------------------------------------------------------
1764
    //        Search for the "Legend cats  : N" line
1765
    // ------------------------------------------------------
1766
1767
0
    int nLine = -1;
1768
0
    for (int i = 0; (i < CSLCount(poGDS->papszRDC)) && (nLine == -1); i++)
1769
0
        if (EQUALN(poGDS->papszRDC[i], rdcLEGEND_CATS, 12))
1770
0
            nLine = i;
1771
1772
0
    if (nLine < 0)
1773
0
        return CE_None;
1774
1775
0
    int nCount = atoi_nz(myCSLFetchNameValue(poGDS->papszRDC, rdcLEGEND_CATS));
1776
1777
    // ------------------------------------------------------
1778
    //        Delete old instance of the category names
1779
    // ------------------------------------------------------
1780
1781
0
    if (nCount > 0)
1782
0
        poGDS->papszRDC =
1783
0
            CSLRemoveStrings(poGDS->papszRDC, nLine + 1, nCount, nullptr);
1784
1785
0
    nCount = 0;
1786
1787
0
    for (int i = 0; i < nCatCount; i++)
1788
0
    {
1789
0
        if ((strlen(papszCategoryNames[i]) > 0))
1790
0
        {
1791
0
            poGDS->papszRDC =
1792
0
                CSLInsertString(poGDS->papszRDC, (nLine + nCount + 1),
1793
0
                                CPLSPrintf("%s:%s", CPLSPrintf(rdcCODE_N, i),
1794
0
                                           papszCategoryNames[i]));
1795
0
            nCount++;
1796
0
        }
1797
0
    }
1798
1799
0
    poGDS->papszRDC =
1800
0
        CSLSetNameValue(poGDS->papszRDC, rdcLEGEND_CATS,
1801
0
                        CPLSPrintf("%d", nCount));  // this is fine
1802
1803
0
    return CE_None;
1804
0
}
1805
1806
/************************************************************************/
1807
/*                           GetColorTable()                            */
1808
/************************************************************************/
1809
1810
GDALColorTable *IdrisiRasterBand::GetColorTable()
1811
0
{
1812
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1813
1814
0
    if (poGDS->poColorTable->GetColorEntryCount() == 0)
1815
0
    {
1816
0
        return nullptr;
1817
0
    }
1818
1819
0
    return poGDS->poColorTable;
1820
0
}
1821
1822
/************************************************************************/
1823
/*                           SetColorTable()                            */
1824
/************************************************************************/
1825
1826
CPLErr IdrisiRasterBand::SetColorTable(GDALColorTable *poColorTable)
1827
0
{
1828
0
    if (poColorTable == nullptr)
1829
0
    {
1830
0
        return CE_None;
1831
0
    }
1832
1833
0
    if (poColorTable->GetColorEntryCount() == 0)
1834
0
    {
1835
0
        return CE_None;
1836
0
    }
1837
1838
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1839
1840
0
    delete poGDS->poColorTable;
1841
1842
0
    poGDS->poColorTable = poColorTable->Clone();
1843
1844
0
    const std::string osSMPFilename =
1845
0
        CPLResetExtensionSafe(poGDS->pszFilename, extSMP);
1846
0
    VSILFILE *fpSMP = VSIFOpenL(osSMPFilename.c_str(), "w");
1847
1848
0
    if (fpSMP != nullptr)
1849
0
    {
1850
0
        VSIFWriteL("[Idrisi]", 8, 1, fpSMP);
1851
0
        GByte nPlatform = 1;
1852
0
        VSIFWriteL(&nPlatform, 1, 1, fpSMP);
1853
0
        GByte nVersion = 11;
1854
0
        VSIFWriteL(&nVersion, 1, 1, fpSMP);
1855
0
        GByte nDepth = 8;
1856
0
        VSIFWriteL(&nDepth, 1, 1, fpSMP);
1857
0
        GByte nHeadSz = 18;
1858
0
        VSIFWriteL(&nHeadSz, 1, 1, fpSMP);
1859
0
        GUInt16 nCount = 255;
1860
0
        VSIFWriteL(&nCount, 2, 1, fpSMP);
1861
0
        GUInt16 nMix = 0;
1862
0
        VSIFWriteL(&nMix, 2, 1, fpSMP);
1863
0
        GUInt16 nMax = 255;
1864
0
        VSIFWriteL(&nMax, 2, 1, fpSMP);
1865
1866
0
        GDALColorEntry oEntry;
1867
0
        GByte aucRGB[3];
1868
1869
0
        for (int i = 0; i < poColorTable->GetColorEntryCount(); i++)
1870
0
        {
1871
0
            poColorTable->GetColorEntryAsRGB(i, &oEntry);
1872
0
            aucRGB[0] = (GByte)oEntry.c1;
1873
0
            aucRGB[1] = (GByte)oEntry.c2;
1874
0
            aucRGB[2] = (GByte)oEntry.c3;
1875
0
            VSIFWriteL(&aucRGB, 3, 1, fpSMP);
1876
0
        }
1877
        /* smp files always have 256 occurrences. */
1878
0
        for (int i = poColorTable->GetColorEntryCount(); i <= 255; i++)
1879
0
        {
1880
0
            poColorTable->GetColorEntryAsRGB(i, &oEntry);
1881
0
            aucRGB[0] = (GByte)0;
1882
0
            aucRGB[1] = (GByte)0;
1883
0
            aucRGB[2] = (GByte)0;
1884
0
            VSIFWriteL(&aucRGB, 3, 1, fpSMP);
1885
0
        }
1886
0
        VSIFCloseL(fpSMP);
1887
0
    }
1888
1889
0
    return CE_None;
1890
0
}
1891
1892
/************************************************************************/
1893
/*                           GetUnitType()                              */
1894
/************************************************************************/
1895
1896
const char *IdrisiRasterBand::GetUnitType()
1897
0
{
1898
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1899
1900
0
    return poGDS->pszUnitType;
1901
0
}
1902
1903
/************************************************************************/
1904
/*                            SetUnitType()                             */
1905
/************************************************************************/
1906
1907
CPLErr IdrisiRasterBand::SetUnitType(const char *pszUnitType)
1908
0
{
1909
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1910
1911
0
    if (strlen(pszUnitType) == 0)
1912
0
    {
1913
0
        poGDS->papszRDC =
1914
0
            CSLSetNameValue(poGDS->papszRDC, rdcVALUE_UNITS, "unspecified");
1915
0
    }
1916
0
    else
1917
0
    {
1918
0
        poGDS->papszRDC =
1919
0
            CSLSetNameValue(poGDS->papszRDC, rdcVALUE_UNITS, pszUnitType);
1920
0
    }
1921
1922
0
    return CE_None;
1923
0
}
1924
1925
/************************************************************************/
1926
/*                             SetMinMax()                              */
1927
/************************************************************************/
1928
1929
CPLErr IdrisiRasterBand::SetMinMax(double dfMin, double dfMax)
1930
0
{
1931
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
1932
1933
0
    fMinimum = (float)dfMin;
1934
0
    fMaximum = (float)dfMax;
1935
1936
0
    double adfMin[3] = {0.0, 0.0, 0.0};
1937
0
    double adfMax[3] = {0.0, 0.0, 0.0};
1938
1939
0
    if (myCSLFetchNameValue(poGDS->papszRDC, rdcMIN_VALUE) != nullptr)
1940
0
        CPLsscanf(myCSLFetchNameValue(poGDS->papszRDC, rdcMIN_VALUE),
1941
0
                  "%lf %lf %lf", &adfMin[0], &adfMin[1], &adfMin[2]);
1942
0
    if (myCSLFetchNameValue(poGDS->papszRDC, rdcMAX_VALUE) != nullptr)
1943
0
        CPLsscanf(myCSLFetchNameValue(poGDS->papszRDC, rdcMAX_VALUE),
1944
0
                  "%lf %lf %lf", &adfMax[0], &adfMax[1], &adfMax[2]);
1945
1946
0
    adfMin[nBand - 1] = dfMin;
1947
0
    adfMax[nBand - 1] = dfMax;
1948
1949
0
    if (poGDS->nBands == 3)
1950
0
    {
1951
0
        poGDS->papszRDC = CSLSetNameValue(
1952
0
            poGDS->papszRDC, rdcMIN_VALUE,
1953
0
            CPLSPrintf("%.8g %.8g %.8g", adfMin[0], adfMin[1], adfMin[2]));
1954
0
        poGDS->papszRDC = CSLSetNameValue(
1955
0
            poGDS->papszRDC, rdcMAX_VALUE,
1956
0
            CPLSPrintf("%.8g %.8g %.8g", adfMax[0], adfMax[1], adfMax[2]));
1957
0
        poGDS->papszRDC = CSLSetNameValue(
1958
0
            poGDS->papszRDC, rdcDISPLAY_MIN,
1959
0
            CPLSPrintf("%.8g %.8g %.8g", adfMin[0], adfMin[1], adfMin[2]));
1960
0
        poGDS->papszRDC = CSLSetNameValue(
1961
0
            poGDS->papszRDC, rdcDISPLAY_MAX,
1962
0
            CPLSPrintf("%.8g %.8g %.8g", adfMax[0], adfMax[1], adfMax[2]));
1963
0
    }
1964
0
    else
1965
0
    {
1966
0
        poGDS->papszRDC = CSLSetNameValue(poGDS->papszRDC, rdcMIN_VALUE,
1967
0
                                          CPLSPrintf("%.8g", adfMin[0]));
1968
0
        poGDS->papszRDC = CSLSetNameValue(poGDS->papszRDC, rdcMAX_VALUE,
1969
0
                                          CPLSPrintf("%.8g", adfMax[0]));
1970
0
        poGDS->papszRDC = CSLSetNameValue(poGDS->papszRDC, rdcDISPLAY_MIN,
1971
0
                                          CPLSPrintf("%.8g", adfMin[0]));
1972
0
        poGDS->papszRDC = CSLSetNameValue(poGDS->papszRDC, rdcDISPLAY_MAX,
1973
0
                                          CPLSPrintf("%.8g", adfMax[0]));
1974
0
    }
1975
1976
0
    return CE_None;
1977
0
}
1978
1979
/************************************************************************/
1980
/*                           SetStatistics()                            */
1981
/************************************************************************/
1982
1983
CPLErr IdrisiRasterBand::SetStatistics(double dfMin, double dfMax,
1984
                                       double dfMean, double dfStdDev)
1985
0
{
1986
0
    SetMinMax(dfMin, dfMax);
1987
1988
0
    return GDALPamRasterBand::SetStatistics(dfMin, dfMax, dfMean, dfStdDev);
1989
0
}
1990
1991
/************************************************************************/
1992
/*                           SetDefaultRAT()                            */
1993
/************************************************************************/
1994
1995
CPLErr IdrisiRasterBand::SetDefaultRAT(const GDALRasterAttributeTable *poRAT)
1996
0
{
1997
0
    if (!poRAT)
1998
0
    {
1999
0
        return CE_Failure;
2000
0
    }
2001
2002
    // ----------------------------------------------------------
2003
    // Get field indices
2004
    // ----------------------------------------------------------
2005
2006
0
    int iValue = -1;
2007
0
    int iRed = poRAT->GetColOfUsage(GFU_Red);
2008
0
    int iGreen = poRAT->GetColOfUsage(GFU_Green);
2009
0
    int iBlue = poRAT->GetColOfUsage(GFU_Blue);
2010
2011
0
    GDALColorTable *poCT = nullptr;
2012
0
    char **papszNames = nullptr;
2013
2014
0
    int nFact = 1;
2015
2016
    // ----------------------------------------------------------
2017
    // Seek for "Value" field index (AGIS standards field name)
2018
    // ----------------------------------------------------------
2019
2020
0
    if (GetColorTable() == nullptr ||
2021
0
        GetColorTable()->GetColorEntryCount() == 0)
2022
0
    {
2023
0
        for (int i = 0; i < poRAT->GetColumnCount(); i++)
2024
0
        {
2025
0
            if (STARTS_WITH_CI(poRAT->GetNameOfCol(i), "Value"))
2026
0
            {
2027
0
                iValue = i;
2028
0
                break;
2029
0
            }
2030
0
        }
2031
2032
0
        if (iRed != -1 && iGreen != -1 && iBlue != -1)
2033
0
        {
2034
0
            poCT = new GDALColorTable();
2035
0
            nFact = poRAT->GetTypeOfCol(iRed) == GFT_Real ? 255 : 1;
2036
0
        }
2037
0
    }
2038
2039
    // ----------------------------------------------------------
2040
    // Seek for Name field index
2041
    // ----------------------------------------------------------
2042
2043
0
    int iName = -1;
2044
0
    if (CSLCount(GetCategoryNames()) == 0)
2045
0
    {
2046
0
        iName = poRAT->GetColOfUsage(GFU_Name);
2047
0
        if (iName == -1)
2048
0
        {
2049
0
            for (int i = 0; i < poRAT->GetColumnCount(); i++)
2050
0
            {
2051
0
                if (STARTS_WITH_CI(poRAT->GetNameOfCol(i), "Class_Name"))
2052
0
                {
2053
0
                    iName = i;
2054
0
                    break;
2055
0
                }
2056
0
                else if (STARTS_WITH_CI(poRAT->GetNameOfCol(i), "Categor"))
2057
0
                {
2058
0
                    iName = i;
2059
0
                    break;
2060
0
                }
2061
0
                else if (STARTS_WITH_CI(poRAT->GetNameOfCol(i), "Name"))
2062
0
                {
2063
0
                    iName = i;
2064
0
                    break;
2065
0
                }
2066
0
            }
2067
0
        }
2068
2069
        /* if still can't find it use the first String column */
2070
2071
0
        if (iName == -1)
2072
0
        {
2073
0
            for (int i = 0; i < poRAT->GetColumnCount(); i++)
2074
0
            {
2075
0
                if (poRAT->GetTypeOfCol(i) == GFT_String)
2076
0
                {
2077
0
                    iName = i;
2078
0
                    break;
2079
0
                }
2080
0
            }
2081
0
        }
2082
2083
        // ----------------------------------------------------------
2084
        // Incomplete Attribute Table;
2085
        // ----------------------------------------------------------
2086
2087
0
        if (iName == -1)
2088
0
        {
2089
0
            iName = iValue;
2090
0
        }
2091
0
    }
2092
2093
    // ----------------------------------------------------------
2094
    // Load values
2095
    // ----------------------------------------------------------
2096
2097
0
    GDALColorEntry sColor;
2098
0
    int iEntry = 0;
2099
0
    int iOut = 0;
2100
0
    int nEntryCount = poRAT->GetRowCount();
2101
0
    int nValue = 0;
2102
2103
0
    if (iValue != -1)
2104
0
    {
2105
0
        nValue = poRAT->GetValueAsInt(iEntry, iValue);
2106
0
    }
2107
2108
0
    for (iOut = 0; iOut < 65535 && (iEntry < nEntryCount); iOut++)
2109
0
    {
2110
0
        if (iOut == nValue)
2111
0
        {
2112
0
            if (poCT)
2113
0
            {
2114
0
                const double dRed = poRAT->GetValueAsDouble(iEntry, iRed);
2115
0
                const double dGreen = poRAT->GetValueAsDouble(iEntry, iGreen);
2116
0
                const double dBlue = poRAT->GetValueAsDouble(iEntry, iBlue);
2117
0
                sColor.c1 = (short)(dRed * nFact);
2118
0
                sColor.c2 = (short)(dGreen * nFact);
2119
0
                sColor.c3 = (short)(dBlue * nFact);
2120
0
                sColor.c4 = (short)(255 / nFact);
2121
0
                poCT->SetColorEntry(iEntry, &sColor);
2122
0
            }
2123
2124
0
            if (iName != -1)
2125
0
            {
2126
0
                papszNames = CSLAddString(
2127
0
                    papszNames, poRAT->GetValueAsString(iEntry, iName));
2128
0
            }
2129
2130
            /* Advance on the table */
2131
2132
0
            if ((++iEntry) < nEntryCount)
2133
0
            {
2134
0
                if (iValue != -1)
2135
0
                    nValue = poRAT->GetValueAsInt(iEntry, iValue);
2136
0
                else
2137
0
                    nValue = iEntry;
2138
0
            }
2139
0
        }
2140
0
        else if (iOut < nValue)
2141
0
        {
2142
0
            if (poCT)
2143
0
            {
2144
0
                sColor.c1 = (short)0;
2145
0
                sColor.c2 = (short)0;
2146
0
                sColor.c3 = (short)0;
2147
0
                sColor.c4 = (short)255;
2148
0
                poCT->SetColorEntry(iEntry, &sColor);
2149
0
            }
2150
2151
0
            if (iName != -1)
2152
0
                papszNames = CSLAddString(papszNames, "");
2153
0
        }
2154
0
    }
2155
2156
    // ----------------------------------------------------------
2157
    // Set Color Table
2158
    // ----------------------------------------------------------
2159
2160
0
    if (poCT)
2161
0
    {
2162
0
        SetColorTable(poCT);
2163
0
        delete poCT;
2164
0
    }
2165
2166
    // ----------------------------------------------------------
2167
    // Update Category Names
2168
    // ----------------------------------------------------------
2169
2170
0
    if (papszNames)
2171
0
    {
2172
0
        SetCategoryNames(papszNames);
2173
0
        CSLDestroy(papszNames);
2174
0
    }
2175
2176
    // ----------------------------------------------------------
2177
    // Update Attribute Table
2178
    // ----------------------------------------------------------
2179
2180
0
    if (poDefaultRAT)
2181
0
    {
2182
0
        delete poDefaultRAT;
2183
0
    }
2184
2185
0
    poDefaultRAT = poRAT->Clone();
2186
2187
0
    return CE_None;
2188
0
}
2189
2190
/************************************************************************/
2191
/*                           GetDefaultRAT()                            */
2192
/************************************************************************/
2193
2194
GDALRasterAttributeTable *IdrisiRasterBand::GetDefaultRAT()
2195
0
{
2196
0
    IdrisiDataset *poGDS = cpl::down_cast<IdrisiDataset *>(poDS);
2197
2198
0
    if (poGDS->papszCategories == nullptr)
2199
0
    {
2200
0
        return nullptr;
2201
0
    }
2202
2203
0
    bool bHasColorTable = poGDS->poColorTable->GetColorEntryCount() > 0;
2204
2205
    // ----------------------------------------------------------
2206
    // Create the bands Attribute Table
2207
    // ----------------------------------------------------------
2208
2209
0
    if (poDefaultRAT)
2210
0
    {
2211
0
        delete poDefaultRAT;
2212
0
    }
2213
2214
0
    poDefaultRAT = new GDALDefaultRasterAttributeTable();
2215
2216
    // ----------------------------------------------------------
2217
    // Create (Value, Red, Green, Blue, Alpha, Class_Name) fields
2218
    // ----------------------------------------------------------
2219
2220
0
    poDefaultRAT->CreateColumn("Value", GFT_Integer, GFU_Generic);
2221
0
    poDefaultRAT->CreateColumn("Value_1", GFT_Integer, GFU_MinMax);
2222
2223
0
    if (bHasColorTable)
2224
0
    {
2225
0
        poDefaultRAT->CreateColumn("Red", GFT_Integer, GFU_Red);
2226
0
        poDefaultRAT->CreateColumn("Green", GFT_Integer, GFU_Green);
2227
0
        poDefaultRAT->CreateColumn("Blue", GFT_Integer, GFU_Blue);
2228
0
        poDefaultRAT->CreateColumn("Alpha", GFT_Integer, GFU_Alpha);
2229
0
    }
2230
0
    poDefaultRAT->CreateColumn("Class_name", GFT_String, GFU_Name);
2231
2232
    // ----------------------------------------------------------
2233
    // Loop through the Category Names.
2234
    // ----------------------------------------------------------
2235
2236
0
    GDALColorEntry sEntry;
2237
0
    int iName = poDefaultRAT->GetColOfUsage(GFU_Name);
2238
0
    int nEntryCount = CSLCount(poGDS->papszCategories);
2239
0
    int iRows = 0;
2240
2241
0
    for (int iEntry = 0; iEntry < nEntryCount; iEntry++)
2242
0
    {
2243
0
        if (EQUAL(poGDS->papszCategories[iEntry], ""))
2244
0
        {
2245
0
            continue;  // Eliminate the empty ones
2246
0
        }
2247
0
        poDefaultRAT->SetRowCount(poDefaultRAT->GetRowCount() + 1);
2248
0
        poDefaultRAT->SetValue(iRows, 0, iEntry);
2249
0
        poDefaultRAT->SetValue(iRows, 1, iEntry);
2250
0
        if (bHasColorTable)
2251
0
        {
2252
0
            poGDS->poColorTable->GetColorEntryAsRGB(iEntry, &sEntry);
2253
0
            poDefaultRAT->SetValue(iRows, 2, sEntry.c1);
2254
0
            poDefaultRAT->SetValue(iRows, 3, sEntry.c2);
2255
0
            poDefaultRAT->SetValue(iRows, 4, sEntry.c3);
2256
0
            poDefaultRAT->SetValue(iRows, 5, sEntry.c4);
2257
0
        }
2258
0
        poDefaultRAT->SetValue(iRows, iName, poGDS->papszCategories[iEntry]);
2259
0
        iRows++;
2260
0
    }
2261
2262
0
    return poDefaultRAT;
2263
0
}
2264
2265
/************************************************************************/
2266
/*                       IdrisiGeoReference2Wkt()                       */
2267
/************************************************************************/
2268
2269
/***
2270
 * Converts Idrisi geographic reference information to OpenGIS WKT.
2271
 *
2272
 * The Idrisi metadata file contain two fields that describe the
2273
 * geographic reference, RefSystem and RefUnit.
2274
 *
2275
 * RefSystem can contains the world "plane" or the name of a georeference
2276
 * file <refsystem>.ref that details the geographic reference
2277
 * system( coordinate system and projection parameters ). RefUnits
2278
 * indicates the unit of the image bounds.
2279
 *
2280
 * The georeference files are generally located in the product installation
2281
 * folder $IDRISIDIR\\Georef, but they are first looked for in the same
2282
 * folder as the data file.
2283
 *
2284
 * If a Reference system names can be recognized by a name convention
2285
 * it will be interpreted without the need to read the georeference file.
2286
 * That includes "latlong" and all the UTM and State Plane zones.
2287
 *
2288
 * RefSystem "latlong" means that the data is not project and the coordinate
2289
 * system is WGS84. RefSystem "plane" means that the there is no coordinate
2290
 * system but the it is possible to calculate areas and distance by looking
2291
 * at the RefUnits.
2292
 *
2293
 * If the environment variable IDRISIDIR is not set and the georeference file
2294
 * need to be read then the projection string will result as unknown.
2295
 ***/
2296
2297
CPLErr IdrisiGeoReference2Wkt(const char *pszFilename, const char *pszRefSystem,
2298
                              const char *pszRefUnits,
2299
                              OGRSpatialReference &oSRS)
2300
0
{
2301
    // ---------------------------------------------------------
2302
    //  Plane
2303
    // ---------------------------------------------------------
2304
2305
0
    if (EQUAL(pszRefSystem, rstPLANE))
2306
0
    {
2307
0
        oSRS.SetLocalCS("Plane");
2308
0
        int nUnit = GetUnitIndex(pszRefUnits);
2309
0
        if (nUnit > -1)
2310
0
        {
2311
0
            int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2312
0
            oSRS.SetLinearUnits(aoLinearUnitsConv[nDeft].pszName,
2313
0
                                aoLinearUnitsConv[nDeft].dfConv);
2314
0
        }
2315
0
        return CE_None;
2316
0
    }
2317
2318
    // ---------------------------------------------------------
2319
    //  Latlong
2320
    // ---------------------------------------------------------
2321
2322
0
    if (EQUAL(pszRefSystem, rstLATLONG) || EQUAL(pszRefSystem, rstLATLONG2))
2323
0
    {
2324
0
        oSRS.SetWellKnownGeogCS("WGS84");
2325
0
        return CE_None;
2326
0
    }
2327
2328
    // ---------------------------------------------------------
2329
    //  Prepare for scanning in lower case
2330
    // ---------------------------------------------------------
2331
2332
0
    char *pszRefSystemLower = CPLStrdup(pszRefSystem);
2333
0
    CPLStrlwr(pszRefSystemLower);
2334
2335
    // ---------------------------------------------------------
2336
    //  UTM naming convention( ex.: utm-30n )
2337
    // ---------------------------------------------------------
2338
2339
0
    if (EQUALN(pszRefSystem, rstUTM, 3))
2340
0
    {
2341
0
        int nZone;
2342
0
        char cNorth;
2343
0
        sscanf(pszRefSystemLower, rstUTM, &nZone, &cNorth);
2344
0
        oSRS.SetWellKnownGeogCS("WGS84");
2345
0
        oSRS.SetUTM(nZone, (cNorth == 'n'));
2346
0
        CPLFree(pszRefSystemLower);
2347
0
        return CE_None;
2348
0
    }
2349
2350
    // ---------------------------------------------------------
2351
    //  State Plane naming convention( ex.: spc83ma1 )
2352
    // ---------------------------------------------------------
2353
2354
0
    if (EQUALN(pszRefSystem, rstSPC, 3))
2355
0
    {
2356
0
        int nNAD;
2357
0
        int nZone;
2358
0
        char szState[3];
2359
0
        sscanf(pszRefSystemLower, rstSPC, &nNAD, szState, &nZone);
2360
0
        int nSPCode = GetStateCode(szState);
2361
0
        if (nSPCode != -1)
2362
0
        {
2363
0
            nZone = (nZone == 1 ? nSPCode : nSPCode + nZone - 1);
2364
2365
0
            if (oSRS.SetStatePlane(nZone, (nNAD == 83)) != OGRERR_FAILURE)
2366
0
            {
2367
0
                CPLFree(pszRefSystemLower);
2368
0
                return CE_None;
2369
0
            }
2370
2371
            // ----------------------------------------------------------
2372
            //  If SetStatePlane fails, set GeoCS as NAD Datum and let it
2373
            //  try to read the projection info from georeference file( * )
2374
            // ----------------------------------------------------------
2375
2376
0
            oSRS.SetWellKnownGeogCS(CPLSPrintf("NAD%d", nNAD));
2377
0
        }
2378
0
    }
2379
2380
0
    CPLFree(pszRefSystemLower);
2381
0
    pszRefSystemLower = nullptr;
2382
2383
    // ------------------------------------------------------------------
2384
    //  Search for georeference file <RefSystem>.ref
2385
    // ------------------------------------------------------------------
2386
2387
0
    const char *pszFName =
2388
0
        CPLSPrintf("%s%c%s.ref", CPLGetDirnameSafe(pszFilename).c_str(),
2389
0
                   PATHDELIM, pszRefSystem);
2390
2391
0
    if (!FileExists(pszFName))
2392
0
    {
2393
        // ------------------------------------------------------------------
2394
        //  Look at $IDRISIDIR\Georef\<RefSystem>.ref
2395
        // ------------------------------------------------------------------
2396
2397
0
        const char *pszIdrisiDir = CPLGetConfigOption("IDRISIDIR", nullptr);
2398
2399
0
        if ((pszIdrisiDir) != nullptr)
2400
0
        {
2401
0
            pszFName = CPLSPrintf("%s%cgeoref%c%s.ref", pszIdrisiDir, PATHDELIM,
2402
0
                                  PATHDELIM, pszRefSystem);
2403
0
        }
2404
0
    }
2405
2406
    // ------------------------------------------------------------------
2407
    //  Cannot find georeference file
2408
    // ------------------------------------------------------------------
2409
2410
0
    if (!FileExists(pszFName))
2411
0
    {
2412
0
        CPLDebug("RST", "Cannot find Idrisi georeference file %s",
2413
0
                 pszRefSystem);
2414
2415
0
        if (oSRS.IsGeographic() == FALSE) /* see State Plane remarks( * ) */
2416
0
        {
2417
0
            oSRS.SetLocalCS("Unknown");
2418
0
            int nUnit = GetUnitIndex(pszRefUnits);
2419
0
            if (nUnit > -1)
2420
0
            {
2421
0
                int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2422
0
                oSRS.SetLinearUnits(aoLinearUnitsConv[nDeft].pszName,
2423
0
                                    aoLinearUnitsConv[nDeft].dfConv);
2424
0
            }
2425
0
        }
2426
0
        return CE_Failure;
2427
0
    }
2428
2429
    // ------------------------------------------------------------------
2430
    //  Read values from georeference file
2431
    // ------------------------------------------------------------------
2432
2433
0
    char **papszRef = CSLLoad(pszFName);
2434
0
    myCSLSetNameValueSeparator(papszRef, ":");
2435
2436
0
    char *pszGeorefName = nullptr;
2437
2438
0
    const char *pszREF_SYSTEM = myCSLFetchNameValue(papszRef, refREF_SYSTEM);
2439
0
    if (pszREF_SYSTEM != nullptr && EQUAL(pszREF_SYSTEM, "") == FALSE)
2440
0
    {
2441
0
        pszGeorefName = CPLStrdup(pszREF_SYSTEM);
2442
0
    }
2443
0
    else
2444
0
    {
2445
0
        pszGeorefName =
2446
0
            CPLStrdup(myCSLFetchNameValue(papszRef, refREF_SYSTEM2));
2447
0
    }
2448
0
    char *pszProjName = CPLStrdup(myCSLFetchNameValue(papszRef, refPROJECTION));
2449
0
    char *pszDatum = CPLStrdup(myCSLFetchNameValue(papszRef, refDATUM));
2450
0
    char *pszEllipsoid = CPLStrdup(myCSLFetchNameValue(papszRef, refELLIPSOID));
2451
0
    const double dfCenterLat =
2452
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refORIGIN_LAT));
2453
0
    const double dfCenterLong =
2454
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refORIGIN_LONG));
2455
0
    const double dfSemiMajor =
2456
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refMAJOR_SAX));
2457
0
    const double dfSemiMinor =
2458
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refMINOR_SAX));
2459
0
    const double dfFalseEasting =
2460
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refORIGIN_X));
2461
0
    const double dfFalseNorthing =
2462
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refORIGIN_Y));
2463
0
    const double dfStdP1 =
2464
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refSTANDL_1));
2465
0
    const double dfStdP2 =
2466
0
        CPLAtof_nz(myCSLFetchNameValue(papszRef, refSTANDL_2));
2467
0
    double dfScale;
2468
0
    double adfToWGS84[3] = {0.0, 0.0, 0.0};
2469
2470
0
    const char *pszToWGS84 = myCSLFetchNameValue(papszRef, refDELTA_WGS84);
2471
0
    if (pszToWGS84)
2472
0
        CPLsscanf(pszToWGS84, "%lf %lf %lf", &adfToWGS84[0], &adfToWGS84[1],
2473
0
                  &adfToWGS84[2]);
2474
2475
0
    const char *pszSCALE_FAC = myCSLFetchNameValue(papszRef, refSCALE_FAC);
2476
0
    if (pszSCALE_FAC == nullptr || EQUAL(pszSCALE_FAC, "na"))
2477
0
        dfScale = 1.0;
2478
0
    else
2479
0
        dfScale = CPLAtof_nz(pszSCALE_FAC);
2480
2481
0
    CSLDestroy(papszRef);
2482
2483
    // ----------------------------------------------------------------------
2484
    //  Set the Geographic Coordinate System
2485
    // ----------------------------------------------------------------------
2486
2487
0
    if (oSRS.IsGeographic() == FALSE) /* see State Plane remarks(*) */
2488
0
    {
2489
0
        int nEPSG = 0;
2490
2491
        // ----------------------------------------------------------------------
2492
        //  Is it a WGS84 equivalent?
2493
        // ----------------------------------------------------------------------
2494
2495
0
        if ((STARTS_WITH_CI(pszEllipsoid, "WGS")) &&
2496
0
            (strstr(pszEllipsoid, "84")) && (STARTS_WITH_CI(pszDatum, "WGS")) &&
2497
0
            (strstr(pszDatum, "84")) && (adfToWGS84[0] == 0.0) &&
2498
0
            (adfToWGS84[1] == 0.0) && (adfToWGS84[2] == 0.0))
2499
0
        {
2500
0
            nEPSG = 4326;
2501
0
        }
2502
2503
        // ----------------------------------------------------------------------
2504
        //  Match GCS's DATUM_NAME by using 'ApproxString' over Datum
2505
        // ----------------------------------------------------------------------
2506
2507
0
        if (nEPSG == 0)
2508
0
        {
2509
0
            const PJ_TYPE nObjType = PJ_TYPE_GEODETIC_REFERENCE_FRAME;
2510
0
            auto datumList =
2511
0
                proj_create_from_name(OSRGetProjTLSContext(), "EPSG", pszDatum,
2512
0
                                      &nObjType, 1, true, 1, nullptr);
2513
0
            if (datumList && proj_list_get_count(datumList) == 1)
2514
0
            {
2515
0
                auto datum =
2516
0
                    proj_list_get(OSRGetProjTLSContext(), datumList, 0);
2517
0
                if (datum)
2518
0
                {
2519
0
                    const char *datumCode = proj_get_id_code(datum, 0);
2520
0
                    if (datumCode)
2521
0
                    {
2522
0
                        auto crsList = proj_query_geodetic_crs_from_datum(
2523
0
                            OSRGetProjTLSContext(), "EPSG", "EPSG", datumCode,
2524
0
                            "geographic 2D");
2525
0
                        if (crsList && proj_list_get_count(crsList) != 0)
2526
0
                        {
2527
0
                            auto crs = proj_list_get(OSRGetProjTLSContext(),
2528
0
                                                     crsList, 0);
2529
0
                            if (crs)
2530
0
                            {
2531
0
                                const char *crsCode = proj_get_id_code(crs, 0);
2532
0
                                if (crsCode)
2533
0
                                {
2534
0
                                    nEPSG = atoi(crsCode);
2535
0
                                }
2536
0
                                proj_destroy(crs);
2537
0
                            }
2538
0
                        }
2539
0
                        proj_list_destroy(crsList);
2540
0
                    }
2541
0
                    proj_destroy(datum);
2542
0
                }
2543
0
            }
2544
0
            proj_list_destroy(datumList);
2545
0
        }
2546
2547
        // ----------------------------------------------------------------------
2548
        //  Match GCS's COORD_REF_SYS_NAME by using 'ApproxString' over Datum
2549
        // ----------------------------------------------------------------------
2550
2551
0
        if (nEPSG == 0)
2552
0
        {
2553
0
            const PJ_TYPE nObjType = PJ_TYPE_GEOGRAPHIC_2D_CRS;
2554
0
            auto crsList =
2555
0
                proj_create_from_name(OSRGetProjTLSContext(), "EPSG", pszDatum,
2556
0
                                      &nObjType, 1, true, 1, nullptr);
2557
0
            if (crsList && proj_list_get_count(crsList) != 0)
2558
0
            {
2559
0
                auto crs = proj_list_get(OSRGetProjTLSContext(), crsList, 0);
2560
0
                if (crs)
2561
0
                {
2562
0
                    const char *crsCode = proj_get_id_code(crs, 0);
2563
0
                    if (crsCode)
2564
0
                    {
2565
0
                        nEPSG = atoi(crsCode);
2566
0
                    }
2567
0
                    proj_destroy(crs);
2568
0
                }
2569
0
            }
2570
0
            proj_list_destroy(crsList);
2571
0
        }
2572
2573
0
        if (nEPSG != 0)
2574
0
        {
2575
0
            oSRS.importFromEPSG(nEPSG);
2576
0
        }
2577
0
        else
2578
0
        {
2579
            // --------------------------------------------------
2580
            //  Create GeogCS based on the georeference file info
2581
            // --------------------------------------------------
2582
2583
0
            oSRS.SetGeogCS(pszRefSystem, pszDatum, pszEllipsoid, dfSemiMajor,
2584
0
                           (dfSemiMinor == dfSemiMajor)
2585
0
                               ? 0.0
2586
0
                               : (-1.0 / (dfSemiMinor / dfSemiMajor - 1.0)));
2587
0
        }
2588
2589
        // ----------------------------------------------------------------------
2590
        //  Note: That will override EPSG info:
2591
        // ----------------------------------------------------------------------
2592
2593
0
        oSRS.SetTOWGS84(adfToWGS84[0], adfToWGS84[1], adfToWGS84[2]);
2594
0
    }
2595
2596
    // ----------------------------------------------------------------------
2597
    //  If the georeference file tells that it is a non project system:
2598
    // ----------------------------------------------------------------------
2599
2600
0
    if (EQUAL(pszProjName, "none"))
2601
0
    {
2602
0
        CPLFree(pszGeorefName);
2603
0
        CPLFree(pszProjName);
2604
0
        CPLFree(pszDatum);
2605
0
        CPLFree(pszEllipsoid);
2606
2607
0
        return CE_None;
2608
0
    }
2609
2610
    // ----------------------------------------------------------------------
2611
    //  Create Projection information based on georeference file info
2612
    // ----------------------------------------------------------------------
2613
2614
    //  Idrisi user's Manual,   Supported Projection:
2615
    //
2616
    //      Mercator
2617
    //      Transverse Mercator
2618
    //      Gauss-Kruger
2619
    //      Lambert Conformal Conic
2620
    //      Plate Carree
2621
    //      Hammer Aitoff
2622
    //      Lambert North Polar Azimuthal Equal Area
2623
    //      Lambert South Polar Azimuthal Equal Area
2624
    //      Lambert Transverse Azimuthal Equal Area
2625
    //      Lambert Oblique Polar Azimuthal Equal Area
2626
    //      North Polar Stereographic
2627
    //      South Polar Stereographic
2628
    //      Transverse Stereographic
2629
    //      Oblique Stereographic
2630
    //      Albers Equal Area Conic
2631
    //      Sinusoidal
2632
    //      Cylindrical Equal Area
2633
2634
0
    if (EQUAL(pszProjName, "Mercator"))
2635
0
    {
2636
0
        oSRS.SetMercator(dfCenterLat, dfCenterLong, dfScale, dfFalseEasting,
2637
0
                         dfFalseNorthing);
2638
0
    }
2639
0
    else if (EQUAL(pszProjName, "Transverse Mercator"))
2640
0
    {
2641
0
        oSRS.SetTM(dfCenterLat, dfCenterLong, dfScale, dfFalseEasting,
2642
0
                   dfFalseNorthing);
2643
0
    }
2644
0
    else if (EQUAL(pszProjName, "Gauss-Kruger"))
2645
0
    {
2646
0
        oSRS.SetTM(dfCenterLat, dfCenterLong, dfScale, dfFalseEasting,
2647
0
                   dfFalseNorthing);
2648
0
    }
2649
0
    else if (EQUAL(pszProjName, "Lambert Conformal Conic"))
2650
0
    {
2651
0
        oSRS.SetLCC(dfStdP1, dfStdP2, dfCenterLat, dfCenterLong, dfFalseEasting,
2652
0
                    dfFalseNorthing);
2653
0
    }
2654
0
    else if (EQUAL(pszProjName, "Plate Carr"
2655
0
                                "\xE9"
2656
0
                                "e")) /* 'eacute' in ISO-8859-1 */
2657
0
    {
2658
0
        oSRS.SetEquirectangular(dfCenterLat, dfCenterLong, dfFalseEasting,
2659
0
                                dfFalseNorthing);
2660
0
    }
2661
0
    else if (EQUAL(pszProjName, "Hammer Aitoff"))
2662
0
    {
2663
0
        oSRS.SetProjection(pszProjName);
2664
0
        oSRS.SetProjParm(SRS_PP_LATITUDE_OF_ORIGIN, dfCenterLat);
2665
0
        oSRS.SetProjParm(SRS_PP_CENTRAL_MERIDIAN, dfCenterLong);
2666
0
        oSRS.SetProjParm(SRS_PP_FALSE_EASTING, dfFalseEasting);
2667
0
        oSRS.SetProjParm(SRS_PP_FALSE_NORTHING, dfFalseNorthing);
2668
0
    }
2669
0
    else if (EQUAL(pszProjName, "Lambert North Polar Azimuthal Equal Area") ||
2670
0
             EQUAL(pszProjName, "Lambert South Polar Azimuthal Equal Area") ||
2671
0
             EQUAL(pszProjName, "Lambert Transverse Azimuthal Equal Area") ||
2672
0
             EQUAL(pszProjName, "Lambert Oblique Polar Azimuthal Equal Area"))
2673
0
    {
2674
0
        oSRS.SetLAEA(dfCenterLat, dfCenterLong, dfFalseEasting,
2675
0
                     dfFalseNorthing);
2676
0
    }
2677
0
    else if (EQUAL(pszProjName, "North Polar Stereographic") ||
2678
0
             EQUAL(pszProjName, "South Polar Stereographic"))
2679
0
    {
2680
0
        oSRS.SetPS(dfCenterLat, dfCenterLong, dfScale, dfFalseEasting,
2681
0
                   dfFalseNorthing);
2682
0
    }
2683
0
    else if (EQUAL(pszProjName, "Transverse Stereographic"))
2684
0
    {
2685
0
        oSRS.SetStereographic(dfCenterLat, dfCenterLong, dfScale,
2686
0
                              dfFalseEasting, dfFalseNorthing);
2687
0
    }
2688
0
    else if (EQUAL(pszProjName, "Oblique Stereographic"))
2689
0
    {
2690
0
        oSRS.SetOS(dfCenterLat, dfCenterLong, dfScale, dfFalseEasting,
2691
0
                   dfFalseNorthing);
2692
0
    }
2693
0
    else if (EQUAL(pszProjName, "Alber's Equal Area Conic") ||
2694
0
             EQUAL(pszProjName, "Albers Equal Area Conic"))
2695
0
    {
2696
0
        oSRS.SetACEA(dfStdP1, dfStdP2, dfCenterLat, dfCenterLong,
2697
0
                     dfFalseEasting, dfFalseNorthing);
2698
0
    }
2699
0
    else if (EQUAL(pszProjName, "Sinusoidal"))
2700
0
    {
2701
0
        oSRS.SetSinusoidal(dfCenterLong, dfFalseEasting, dfFalseNorthing);
2702
0
    }
2703
0
    else if (EQUAL(pszProjName, "CylindricalEA") ||
2704
0
             EQUAL(pszProjName, "Cylindrical Equal Area"))
2705
0
    {
2706
0
        oSRS.SetCEA(dfStdP1, dfCenterLong, dfFalseEasting, dfFalseNorthing);
2707
0
    }
2708
0
    else
2709
0
    {
2710
0
        CPLError(
2711
0
            CE_Warning, CPLE_NotSupported,
2712
0
            "Projection not listed on Idrisi User's Manual( v.15.0/2005 ).\n\t"
2713
0
            "[\"%s\" in georeference file \"%s\"]",
2714
0
            pszProjName, pszFName);
2715
0
        oSRS.Clear();
2716
2717
0
        CPLFree(pszGeorefName);
2718
0
        CPLFree(pszProjName);
2719
0
        CPLFree(pszDatum);
2720
0
        CPLFree(pszEllipsoid);
2721
2722
0
        return CE_Warning;
2723
0
    }
2724
2725
    // ----------------------------------------------------------------------
2726
    //  Set the Linear Units
2727
    // ----------------------------------------------------------------------
2728
2729
0
    int nUnit = GetUnitIndex(pszRefUnits);
2730
0
    if (nUnit > -1)
2731
0
    {
2732
0
        int nDeft = aoLinearUnitsConv[nUnit].nDefaultG;
2733
0
        oSRS.SetLinearUnits(aoLinearUnitsConv[nDeft].pszName,
2734
0
                            aoLinearUnitsConv[nDeft].dfConv);
2735
0
    }
2736
0
    else
2737
0
    {
2738
0
        oSRS.SetLinearUnits("unknown", 1.0);
2739
0
    }
2740
2741
    // ----------------------------------------------------------------------
2742
    //  Name ProjCS with the name on the georeference file
2743
    // ----------------------------------------------------------------------
2744
2745
0
    oSRS.SetProjCS(pszGeorefName);
2746
2747
0
    CPLFree(pszGeorefName);
2748
0
    CPLFree(pszProjName);
2749
0
    CPLFree(pszDatum);
2750
0
    CPLFree(pszEllipsoid);
2751
2752
0
    return CE_None;
2753
0
}
2754
2755
/************************************************************************/
2756
/*                        Wkt2GeoReference()                            */
2757
/************************************************************************/
2758
2759
/***
2760
 * Converts OpenGIS WKT to Idrisi geographic reference information.
2761
 *
2762
 * That function will fill up the two parameters RefSystem and RefUnit
2763
 * that goes into the Idrisi metadata. But it could also create
2764
 * a accompanying georeference file to the output if necessary.
2765
 *
2766
 * First it will try to identify the ProjString as Local, WGS84 or
2767
 * one of the Idrisi name convention reference systems
2768
 * otherwise, if the projection system is supported by Idrisi,
2769
 * it will create a accompanying georeference files.
2770
 ***/
2771
2772
CPLErr IdrisiDataset::Wkt2GeoReference(const OGRSpatialReference &oSRS,
2773
                                       char **pszRefSystem, char **pszRefUnit)
2774
0
{
2775
    // -----------------------------------------------------
2776
    //  Plane with default "Meters"
2777
    // -----------------------------------------------------
2778
2779
0
    if (oSRS.IsEmpty())
2780
0
    {
2781
0
        *pszRefSystem = CPLStrdup(rstPLANE);
2782
0
        *pszRefUnit = CPLStrdup(rstMETER);
2783
0
        return CE_None;
2784
0
    }
2785
2786
    // -----------------------------------------------------
2787
    //  Local => Plane + Linear Unit
2788
    // -----------------------------------------------------
2789
2790
0
    if (oSRS.IsLocal())
2791
0
    {
2792
0
        *pszRefSystem = CPLStrdup(rstPLANE);
2793
0
        *pszRefUnit = GetUnitDefault(oSRS.GetAttrValue("UNIT"),
2794
0
                                     CPLSPrintf("%f", oSRS.GetLinearUnits()));
2795
0
        return CE_None;
2796
0
    }
2797
2798
    // -----------------------------------------------------
2799
    //  Test to identify WGS84 => Latlong + Angular Unit
2800
    // -----------------------------------------------------
2801
2802
0
    if (oSRS.IsGeographic())
2803
0
    {
2804
0
        char *pszSpheroid = CPLStrdup(oSRS.GetAttrValue("SPHEROID"));
2805
0
        char *pszAuthName = CPLStrdup(oSRS.GetAuthorityName("GEOGCS"));
2806
0
        char *pszDatum = CPLStrdup(oSRS.GetAttrValue("DATUM"));
2807
0
        int nGCSCode = -1;
2808
0
        if (EQUAL(pszAuthName, "EPSG"))
2809
0
        {
2810
0
            nGCSCode = atoi(oSRS.GetAuthorityCode("GEOGCS"));
2811
0
        }
2812
0
        if ((nGCSCode == 4326) ||
2813
0
            ((STARTS_WITH_CI(pszSpheroid, "WGS")) &&
2814
0
             (strstr(pszSpheroid, "84")) && (STARTS_WITH_CI(pszDatum, "WGS")) &&
2815
0
             (strstr(pszDatum, "84"))))
2816
0
        {
2817
0
            *pszRefSystem = CPLStrdup(rstLATLONG);
2818
0
            *pszRefUnit = CPLStrdup(rstDEGREE);
2819
2820
0
            CPLFree(pszSpheroid);
2821
0
            CPLFree(pszAuthName);
2822
0
            CPLFree(pszDatum);
2823
2824
0
            return CE_None;
2825
0
        }
2826
2827
0
        CPLFree(pszSpheroid);
2828
0
        CPLFree(pszAuthName);
2829
0
        CPLFree(pszDatum);
2830
0
    }
2831
2832
    // -----------------------------------------------------
2833
    //  Prepare to match some projections
2834
    // -----------------------------------------------------
2835
2836
0
    const char *pszProjName = oSRS.GetAttrValue("PROJECTION");
2837
2838
0
    if (pszProjName == nullptr)
2839
0
    {
2840
0
        pszProjName = "";
2841
0
    }
2842
2843
    // -----------------------------------------------------
2844
    //  Check for UTM zones
2845
    // -----------------------------------------------------
2846
2847
0
    if (EQUAL(pszProjName, SRS_PT_TRANSVERSE_MERCATOR))
2848
0
    {
2849
0
        int nZone = oSRS.GetUTMZone();
2850
2851
0
        if ((nZone != 0) && (EQUAL(oSRS.GetAttrValue("DATUM"), SRS_DN_WGS84)))
2852
0
        {
2853
0
            double dfNorth = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING);
2854
0
            *pszRefSystem = CPLStrdup(
2855
0
                CPLSPrintf(rstUTM, nZone, (dfNorth == 0.0 ? 'n' : 's')));
2856
0
            *pszRefUnit = CPLStrdup(rstMETER);
2857
0
            return CE_None;
2858
0
        }
2859
0
    }
2860
2861
    // -----------------------------------------------------
2862
    //  Check for State Plane
2863
    // -----------------------------------------------------
2864
2865
0
    if (EQUAL(pszProjName, SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP) ||
2866
0
        EQUAL(pszProjName, SRS_PT_TRANSVERSE_MERCATOR))
2867
0
    {
2868
0
        CPLString osPCSCode;
2869
0
        const char *pszID = oSRS.GetAuthorityCode("PROJCS");
2870
0
        if (pszID != nullptr && strlen(pszID) > 0)
2871
0
        {
2872
0
            const char *pszPCSCode =
2873
0
                CSVGetField(CSVFilename("stateplane.csv"), "EPSG_PCS_CODE",
2874
0
                            pszID, CC_Integer, "ID");
2875
0
            osPCSCode = (pszPCSCode) ? pszPCSCode : "";
2876
0
            if (!osPCSCode.empty())
2877
0
            {
2878
0
                int nZone = osPCSCode.back() - '0';
2879
0
                int nSPCode = atoi(osPCSCode);
2880
2881
0
                if (nZone == 0)
2882
0
                    nZone = 1;
2883
0
                else
2884
0
                    nSPCode = nSPCode - nZone + 1;
2885
2886
0
                int nNADYear = 83;
2887
0
                if (nSPCode > 10000)
2888
0
                {
2889
0
                    nNADYear = 27;
2890
0
                    nSPCode -= 10000;
2891
0
                }
2892
0
                char *pszState = CPLStrdup(GetStateName(nSPCode));
2893
0
                if (!EQUAL(pszState, ""))
2894
0
                {
2895
0
                    *pszRefSystem = CPLStrdup(
2896
0
                        CPLSPrintf(rstSPC, nNADYear, pszState, nZone));
2897
0
                    *pszRefUnit =
2898
0
                        GetUnitDefault(oSRS.GetAttrValue("UNIT"),
2899
0
                                       CPLSPrintf("%f", oSRS.GetLinearUnits()));
2900
0
                    CPLFree(pszState);
2901
0
                    return CE_None;
2902
0
                }
2903
0
                CPLFree(pszState);
2904
0
            }
2905
0
        }  //
2906
2907
        // if EPSG code is missing, go to following steps to work with origin
2908
2909
0
        const char *pszNAD83 = "83";
2910
0
        const char *pszNAD27 = "27";
2911
0
        bool bIsOldNAD = false;
2912
2913
0
        const char *pszDatumValue = oSRS.GetAttrValue("DATUM", 0);
2914
0
        if ((strstr(pszDatumValue, pszNAD83) == nullptr) &&
2915
0
            (strstr(pszDatumValue, pszNAD27) != nullptr))
2916
            // strcpy(pszNAD, "27");
2917
0
            bIsOldNAD = true;
2918
2919
0
        if ((oSRS.FindProjParm("central_meridian", nullptr) != -1) &&
2920
0
            (oSRS.FindProjParm("latitude_of_origin", nullptr) != -1))
2921
0
        {
2922
0
            double dfLon = oSRS.GetProjParm("central_meridian");
2923
0
            double dfLat = oSRS.GetProjParm("latitude_of_origin");
2924
0
            dfLon = (int)(fabs(dfLon) * 100 + 0.5) / 100.0;
2925
0
            dfLat = (int)(fabs(dfLat) * 100 + 0.5) / 100.0;
2926
0
            *pszRefSystem = CPLStrdup(GetSpcs(dfLon, dfLat));
2927
0
        }
2928
2929
0
        if (*pszRefSystem != nullptr)
2930
0
        {
2931
            // Convert 83 TO 27
2932
0
            if (bIsOldNAD)
2933
0
            {
2934
0
                char pszOutRefSystem[9];
2935
0
                NAD83to27(pszOutRefSystem, *pszRefSystem);
2936
0
                *pszRefSystem = CPLStrdup(pszOutRefSystem);
2937
0
            }
2938
0
            *pszRefUnit =
2939
0
                GetUnitDefault(oSRS.GetAttrValue("UNIT"),
2940
0
                               CPLSPrintf("%f", oSRS.GetLinearUnits()));
2941
0
            return CE_None;
2942
0
        }
2943
0
    }
2944
2945
0
    const char *pszProjectionOut = nullptr;
2946
2947
0
    if (oSRS.IsProjected())
2948
0
    {
2949
        // ---------------------------------------------------------
2950
        //  Check for supported projections
2951
        // ---------------------------------------------------------
2952
2953
0
        if (EQUAL(pszProjName, SRS_PT_MERCATOR_1SP))
2954
0
        {
2955
0
            pszProjectionOut = "Mercator";
2956
0
        }
2957
0
        else if (EQUAL(pszProjName, SRS_PT_TRANSVERSE_MERCATOR))
2958
0
        {
2959
0
            pszProjectionOut = "Transverse Mercator";
2960
0
        }
2961
0
        else if (EQUAL(pszProjName, SRS_PT_LAMBERT_CONFORMAL_CONIC_2SP))
2962
0
        {
2963
0
            pszProjectionOut = "Lambert Conformal Conic";
2964
0
        }
2965
0
        else if (EQUAL(pszProjName, SRS_PT_EQUIRECTANGULAR))
2966
0
        {
2967
0
            pszProjectionOut = "Plate Carr"
2968
0
                               "\xE9"
2969
0
                               "e"; /* 'eacute' in ISO-8859-1 */
2970
0
        }
2971
0
        else if (EQUAL(pszProjName, SRS_PT_LAMBERT_AZIMUTHAL_EQUAL_AREA))
2972
0
        {
2973
0
            double dfCenterLat =
2974
0
                oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN, 0.0, nullptr);
2975
0
            if (dfCenterLat == 0.0)
2976
0
                pszProjectionOut = "Lambert Transverse Azimuthal Equal Area";
2977
0
            else if (fabs(dfCenterLat) == 90.0)
2978
0
                pszProjectionOut = "Lambert Oblique Polar Azimuthal Equal Area";
2979
0
            else if (dfCenterLat > 0.0)
2980
0
                pszProjectionOut = "Lambert North Oblique Azimuthal Equal Area";
2981
0
            else
2982
0
                pszProjectionOut = "Lambert South Oblique Azimuthal Equal Area";
2983
0
        }
2984
0
        else if (EQUAL(pszProjName, SRS_PT_POLAR_STEREOGRAPHIC))
2985
0
        {
2986
0
            if (oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN, 0.0, nullptr) > 0)
2987
0
                pszProjectionOut = "North Polar Stereographic";
2988
0
            else
2989
0
                pszProjectionOut = "South Polar Stereographic";
2990
0
        }
2991
0
        else if (EQUAL(pszProjName, SRS_PT_STEREOGRAPHIC))
2992
0
        {
2993
0
            pszProjectionOut = "Transverse Stereographic";
2994
0
        }
2995
0
        else if (EQUAL(pszProjName, SRS_PT_OBLIQUE_STEREOGRAPHIC))
2996
0
        {
2997
0
            pszProjectionOut = "Oblique Stereographic";
2998
0
        }
2999
0
        else if (EQUAL(pszProjName, SRS_PT_SINUSOIDAL))
3000
0
        {
3001
0
            pszProjectionOut = "Sinusoidal";
3002
0
        }
3003
0
        else if (EQUAL(pszProjName, SRS_PT_ALBERS_CONIC_EQUAL_AREA))
3004
0
        {
3005
0
            pszProjectionOut = "Alber's Equal Area Conic";
3006
0
        }
3007
0
        else if (EQUAL(pszProjName, SRS_PT_CYLINDRICAL_EQUAL_AREA))
3008
0
        {
3009
0
            pszProjectionOut = "Cylindrical Equal Area";
3010
0
        }
3011
3012
        // ---------------------------------------------------------
3013
        //  Failure, Projection system not supported
3014
        // ---------------------------------------------------------
3015
3016
0
        if (pszProjectionOut == nullptr)
3017
0
        {
3018
0
            CPLDebug("RST", "Not supported by RST driver: PROJECTION[\"%s\"]",
3019
0
                     pszProjName);
3020
3021
0
            *pszRefSystem = CPLStrdup(rstPLANE);
3022
0
            *pszRefUnit = CPLStrdup(rstMETER);
3023
0
            return CE_Failure;
3024
0
        }
3025
0
    }
3026
0
    else
3027
0
    {
3028
0
        pszProjectionOut = "none";
3029
0
    }
3030
3031
    // ---------------------------------------------------------
3032
    //  Prepare to write ref file
3033
    // ---------------------------------------------------------
3034
3035
0
    char *pszGeorefName = CPLStrdup("Unknown");
3036
0
    char *pszDatum = CPLStrdup(oSRS.GetAttrValue("DATUM"));
3037
0
    char *pszEllipsoid = CPLStrdup(oSRS.GetAttrValue("SPHEROID"));
3038
0
    double dfSemiMajor = oSRS.GetSemiMajor();
3039
0
    double dfSemiMinor = oSRS.GetSemiMinor();
3040
0
    double adfToWGS84[3];
3041
0
    oSRS.GetTOWGS84(adfToWGS84, 3);
3042
3043
0
    double dfCenterLat = 0.0;
3044
0
    double dfCenterLong = 0.0;
3045
0
    double dfFalseNorthing = 0.0;
3046
0
    double dfFalseEasting = 0.0;
3047
0
    double dfScale = 1.0;
3048
0
    int nParameters = 0;
3049
0
    double dfStdP1 = 0.0;
3050
0
    double dfStdP2 = 0.0;
3051
0
    char *pszAngularUnit = CPLStrdup(oSRS.GetAttrValue("GEOGCS|UNIT"));
3052
0
    char *pszLinearUnit = nullptr;
3053
3054
0
    if (oSRS.IsProjected())
3055
0
    {
3056
0
        CPLFree(pszGeorefName);
3057
0
        pszGeorefName = CPLStrdup(oSRS.GetAttrValue("PROJCS"));
3058
0
        dfCenterLat = oSRS.GetProjParm(SRS_PP_LATITUDE_OF_ORIGIN, 0.0, nullptr);
3059
0
        dfCenterLong = oSRS.GetProjParm(SRS_PP_CENTRAL_MERIDIAN, 0.0, nullptr);
3060
0
        dfFalseNorthing = oSRS.GetProjParm(SRS_PP_FALSE_NORTHING, 0.0, nullptr);
3061
0
        dfFalseEasting = oSRS.GetProjParm(SRS_PP_FALSE_EASTING, 0.0, nullptr);
3062
0
        dfScale = oSRS.GetProjParm(SRS_PP_SCALE_FACTOR, 0.0, nullptr);
3063
0
        dfStdP1 = oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_1, -0.1, nullptr);
3064
0
        if (EQUAL(pszProjectionOut, "Cylindrical Equal Area"))
3065
0
        {
3066
0
            dfStdP2 = -dfStdP1;
3067
0
            dfScale = 1.0;
3068
0
        }
3069
0
        else
3070
0
        {
3071
0
            dfStdP2 =
3072
0
                oSRS.GetProjParm(SRS_PP_STANDARD_PARALLEL_2, -0.1, nullptr);
3073
0
        }
3074
0
        if (dfStdP1 != -0.1)
3075
0
        {
3076
0
            nParameters = 1;
3077
0
            if (dfStdP2 != -0.1)
3078
0
                nParameters = 2;
3079
0
        }
3080
0
        pszLinearUnit = GetUnitDefault(oSRS.GetAttrValue("PROJCS|UNIT"),
3081
0
                                       CPLSPrintf("%f", oSRS.GetLinearUnits()));
3082
0
    }
3083
0
    else
3084
0
    {
3085
0
        pszLinearUnit = GetUnitDefault(pszAngularUnit);
3086
0
    }
3087
3088
    // ---------------------------------------------------------
3089
    //  Create a companion georeference file for this dataset
3090
    // ---------------------------------------------------------
3091
3092
0
    char **papszRef = nullptr;
3093
0
    papszRef = CSLAddNameValue(papszRef, refREF_SYSTEM, pszGeorefName);
3094
0
    papszRef = CSLAddNameValue(papszRef, refPROJECTION, pszProjectionOut);
3095
0
    papszRef = CSLAddNameValue(papszRef, refDATUM, pszDatum);
3096
0
    papszRef = CSLAddNameValue(papszRef, refDELTA_WGS84,
3097
0
                               CPLSPrintf("%.3g %.3g %.3g", adfToWGS84[0],
3098
0
                                          adfToWGS84[1], adfToWGS84[2]));
3099
0
    papszRef = CSLAddNameValue(papszRef, refELLIPSOID, pszEllipsoid);
3100
0
    papszRef = CSLAddNameValue(papszRef, refMAJOR_SAX,
3101
0
                               CPLSPrintf("%.3f", dfSemiMajor));
3102
0
    papszRef = CSLAddNameValue(papszRef, refMINOR_SAX,
3103
0
                               CPLSPrintf("%.3f", dfSemiMinor));
3104
0
    papszRef = CSLAddNameValue(papszRef, refORIGIN_LONG,
3105
0
                               CPLSPrintf("%.9g", dfCenterLong));
3106
0
    papszRef = CSLAddNameValue(papszRef, refORIGIN_LAT,
3107
0
                               CPLSPrintf("%.9g", dfCenterLat));
3108
0
    papszRef = CSLAddNameValue(papszRef, refORIGIN_X,
3109
0
                               CPLSPrintf("%.9g", dfFalseEasting));
3110
0
    papszRef = CSLAddNameValue(papszRef, refORIGIN_Y,
3111
0
                               CPLSPrintf("%.9g", dfFalseNorthing));
3112
0
    papszRef =
3113
0
        CSLAddNameValue(papszRef, refSCALE_FAC, CPLSPrintf("%.9g", dfScale));
3114
0
    papszRef = CSLAddNameValue(papszRef, refUNITS, pszLinearUnit);
3115
0
    papszRef = CSLAddNameValue(papszRef, refPARAMETERS,
3116
0
                               CPLSPrintf("%1d", nParameters));
3117
0
    if (nParameters > 0)
3118
0
        papszRef =
3119
0
            CSLAddNameValue(papszRef, refSTANDL_1, CPLSPrintf("%.9g", dfStdP1));
3120
0
    if (nParameters > 1)
3121
0
        papszRef =
3122
0
            CSLAddNameValue(papszRef, refSTANDL_2, CPLSPrintf("%.9g", dfStdP2));
3123
0
    myCSLSetNameValueSeparator(papszRef, ": ");
3124
0
    SaveAsCRLF(papszRef, CPLResetExtensionSafe(pszFilename, extREF).c_str());
3125
0
    CSLDestroy(papszRef);
3126
3127
0
    *pszRefSystem = CPLStrdup(CPLGetBasenameSafe(pszFilename).c_str());
3128
0
    *pszRefUnit = CPLStrdup(pszLinearUnit);
3129
3130
0
    CPLFree(pszGeorefName);
3131
0
    CPLFree(pszDatum);
3132
0
    CPLFree(pszEllipsoid);
3133
0
    CPLFree(pszLinearUnit);
3134
0
    CPLFree(pszAngularUnit);
3135
3136
0
    return CE_None;
3137
0
}
3138
3139
/************************************************************************/
3140
/*                             FileExists()                             */
3141
/************************************************************************/
3142
3143
bool FileExists(const char *pszPath)
3144
0
{
3145
0
    VSIStatBufL sStat;
3146
3147
0
    return VSIStatL(pszPath, &sStat) == 0;
3148
0
}
3149
3150
/************************************************************************/
3151
/*                            GetStateCode()                            */
3152
/************************************************************************/
3153
3154
int GetStateCode(const char *pszState)
3155
0
{
3156
0
    for (unsigned int i = 0; i < US_STATE_COUNT; i++)
3157
0
    {
3158
0
        if (EQUAL(pszState, aoUSStateTable[i].pszName))
3159
0
        {
3160
0
            return aoUSStateTable[i].nCode;
3161
0
        }
3162
0
    }
3163
0
    return -1;
3164
0
}
3165
3166
/************************************************************************/
3167
/*                            GetStateName()                            */
3168
/************************************************************************/
3169
3170
const char *GetStateName(int nCode)
3171
0
{
3172
0
    for (unsigned int i = 0; i < US_STATE_COUNT; i++)
3173
0
    {
3174
0
        if (nCode == aoUSStateTable[i].nCode)
3175
0
        {
3176
0
            return aoUSStateTable[i].pszName;
3177
0
        }
3178
0
    }
3179
0
    return nullptr;
3180
0
}
3181
3182
/************************************************************************/
3183
/*                            GetSpcs()                                 */
3184
/************************************************************************/
3185
3186
char *GetSpcs(double dfLon, double dfLat)
3187
0
{
3188
0
    for (int i = 0; i < ORIGIN_COUNT; i++)
3189
0
    {
3190
0
        if ((dfLon == SPCS83Origin[i].longitude) &&
3191
0
            (dfLat == SPCS83Origin[i].latitude))
3192
0
        {
3193
0
            return (char *)SPCS83Origin[i].spcs;
3194
0
        }
3195
0
    }
3196
0
    return nullptr;
3197
0
}
3198
3199
/************************************************************************/
3200
/*                            NAD83to27()                               */
3201
/************************************************************************/
3202
void NAD83to27(char *pszOutRef, char *pszInRef)
3203
0
{
3204
0
    char *pOutput = pszOutRef;
3205
0
    char *pInput = pszInRef;
3206
0
    strncpy(pOutput, pInput, 3);
3207
3208
0
    pOutput = pOutput + 3;
3209
0
    pInput = pInput + 3;
3210
3211
0
    memcpy(pOutput, "27", 2);
3212
0
    pOutput = pOutput + 2;
3213
0
    pInput = pInput + 2;
3214
0
    strcpy(pOutput, pInput);
3215
0
}
3216
3217
/************************************************************************/
3218
/*                            GetUnitIndex()                            */
3219
/************************************************************************/
3220
3221
int GetUnitIndex(const char *pszUnitName)
3222
0
{
3223
0
    for (int i = 0; i < (int)LINEAR_UNITS_COUNT; i++)
3224
0
    {
3225
0
        if (EQUAL(pszUnitName, aoLinearUnitsConv[i].pszName))
3226
0
        {
3227
0
            return i;
3228
0
        }
3229
0
    }
3230
0
    return -1;
3231
0
}
3232
3233
/************************************************************************/
3234
/*                            GetToMeterIndex()                         */
3235
/************************************************************************/
3236
3237
int GetToMeterIndex(const char *pszToMeter)
3238
0
{
3239
0
    const double dfToMeter = CPLAtof_nz(pszToMeter);
3240
3241
0
    if (dfToMeter != 0.0)
3242
0
    {
3243
0
        for (int i = 0; i < (int)LINEAR_UNITS_COUNT; i++)
3244
0
        {
3245
0
            if (std::abs(aoLinearUnitsConv[i].dfConv - dfToMeter) < 0.00001)
3246
0
            {
3247
0
                return i;
3248
0
            }
3249
0
        }
3250
0
    }
3251
3252
0
    return -1;
3253
0
}
3254
3255
/************************************************************************/
3256
/*                            GetUnitDefault()                          */
3257
/************************************************************************/
3258
3259
char *GetUnitDefault(const char *pszUnitName, const char *pszToMeter)
3260
0
{
3261
0
    int nIndex = GetUnitIndex(pszUnitName);
3262
3263
0
    if (nIndex == -1 && pszToMeter != nullptr)
3264
0
    {
3265
0
        nIndex = GetToMeterIndex(pszToMeter);
3266
0
    }
3267
3268
0
    if (nIndex == -1)
3269
0
    {
3270
0
        return CPLStrdup("Unknown");
3271
0
    }
3272
3273
0
    return CPLStrdup(
3274
0
        aoLinearUnitsConv[aoLinearUnitsConv[nIndex].nDefaultI].pszName);
3275
0
}
3276
3277
/************************************************************************/
3278
/*                               CSLSaveCRLF()                          */
3279
/************************************************************************/
3280
3281
/***
3282
 * Write a stringlist to a CR + LF terminated text file.
3283
 *
3284
 * Returns the number of lines written, or 0 if the file could not
3285
 * be written.
3286
 */
3287
3288
int SaveAsCRLF(char **papszStrList, const char *pszFname)
3289
0
{
3290
0
    int nLines = 0;
3291
3292
0
    if (papszStrList)
3293
0
    {
3294
0
        VSILFILE *fp = VSIFOpenL(pszFname, "wt");
3295
0
        if (fp != nullptr)
3296
0
        {
3297
0
            while (*papszStrList != nullptr)
3298
0
            {
3299
0
                if (VSIFPrintfL(fp, "%s\r\n", *papszStrList) < 1)
3300
0
                {
3301
0
                    CPLError(CE_Failure, CPLE_FileIO,
3302
0
                             "CSLSaveCRLF(\"%s\") failed: unable to write to "
3303
0
                             "output file.",
3304
0
                             pszFname);
3305
0
                    break;
3306
0
                }
3307
3308
0
                nLines++;
3309
0
                papszStrList++;
3310
0
            }
3311
3312
0
            VSIFCloseL(fp);
3313
0
        }
3314
0
        else
3315
0
        {
3316
0
            CPLError(CE_Failure, CPLE_OpenFailed,
3317
0
                     "CSLSaveCRLF(\"%s\") failed: unable to open output file.",
3318
0
                     pszFname);
3319
0
        }
3320
0
    }
3321
3322
0
    return nLines;
3323
0
}
3324
3325
/************************************************************************/
3326
/*                        GDALRegister_IDRISI()                         */
3327
/************************************************************************/
3328
3329
void GDALRegister_IDRISI()
3330
22
{
3331
22
    if (GDALGetDriverByName("RST") != nullptr)
3332
0
        return;
3333
3334
22
    GDALDriver *poDriver = new GDALDriver();
3335
3336
22
    poDriver->SetDescription("RST");
3337
22
    poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
3338
22
    poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, rstVERSION);
3339
22
    poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/Idrisi.html");
3340
22
    poDriver->SetMetadataItem(GDAL_DMD_EXTENSION, extRST);
3341
22
    poDriver->SetMetadataItem(GDAL_DMD_CREATIONDATATYPES, "Byte Int16 Float32");
3342
3343
22
    poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
3344
3345
22
    poDriver->pfnOpen = IdrisiDataset::Open;
3346
22
    poDriver->pfnCreate = IdrisiDataset::Create;
3347
22
    poDriver->pfnCreateCopy = IdrisiDataset::CreateCopy;
3348
3349
22
    GetGDALDriverManager()->RegisterDriver(poDriver);
3350
22
}