Coverage Report

Created: 2026-09-26 08:22

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/frmts/pds/vicardataset.cpp
Line
Count
Source
1
/******************************************************************************
2
 *
3
 * Project:  VICAR Driver; JPL/MIPL VICAR Format
4
 * Purpose:  Implementation of VICARDataset
5
 * Author:   Sebastian Walter <sebastian dot walter at fu-berlin dot de>
6
 *
7
 * NOTE: This driver code is loosely based on the ISIS and PDS drivers.
8
 * It is not intended to diminish the contribution of the original authors
9
 ******************************************************************************
10
 * Copyright (c) 2014, Sebastian Walter <sebastian dot walter at fu-berlin dot
11
 *de>
12
 *
13
 * SPDX-License-Identifier: MIT
14
 ****************************************************************************/
15
16
constexpr int VICAR_NULL1 = 0;
17
constexpr int VICAR_NULL2 = -32768;
18
constexpr double VICAR_NULL3 = -32768.0;
19
20
#include "cpl_port.h"
21
22
#include "cpl_safemaths.hpp"
23
#include "cpl_vax.h"
24
#include "cpl_vsi_error.h"
25
#include "vicardataset.h"
26
#include "nasakeywordhandler.h"
27
#include "vicarkeywordhandler.h"
28
#include "pdsdrivercore.h"
29
#include "json_utils.h"
30
31
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
32
#include "gtiff.h"
33
#include "geotiff.h"
34
#include "tifvsi.h"
35
#include "xtiffio.h"
36
#include "gt_wkt_srs_priv.h"
37
#endif
38
39
#include <exception>
40
#include <limits>
41
#include <string>
42
43
#ifdef EMBED_RESOURCE_FILES
44
#include "embedded_resources.h"
45
#endif
46
47
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
48
/* GeoTIFF 1.0 geokeys */
49
50
static const geokey_t GTiffAsciiKeys[] = {GTCitationGeoKey, GeogCitationGeoKey,
51
                                          PCSCitationGeoKey,
52
                                          VerticalCitationGeoKey};
53
54
static const geokey_t GTiffDoubleKeys[] = {
55
    GeogInvFlatteningGeoKey,      GeogSemiMajorAxisGeoKey,
56
    GeogSemiMinorAxisGeoKey,      ProjAzimuthAngleGeoKey,
57
    ProjCenterLatGeoKey,          ProjCenterLongGeoKey,
58
    ProjFalseEastingGeoKey,       ProjFalseNorthingGeoKey,
59
    ProjFalseOriginEastingGeoKey, ProjFalseOriginLatGeoKey,
60
    ProjFalseOriginLongGeoKey,    ProjFalseOriginNorthingGeoKey,
61
    ProjLinearUnitSizeGeoKey,     ProjNatOriginLatGeoKey,
62
    ProjNatOriginLongGeoKey,      ProjOriginLatGeoKey,
63
    ProjOriginLongGeoKey,         ProjRectifiedGridAngleGeoKey,
64
    ProjScaleAtNatOriginGeoKey,   ProjScaleAtOriginGeoKey,
65
    ProjStdParallel1GeoKey,       ProjStdParallel2GeoKey,
66
    ProjStdParallelGeoKey,        ProjStraightVertPoleLongGeoKey,
67
    GeogLinearUnitSizeGeoKey,     GeogAngularUnitSizeGeoKey,
68
    GeogPrimeMeridianLongGeoKey,  ProjCenterEastingGeoKey,
69
    ProjCenterNorthingGeoKey,     ProjScaleAtCenterGeoKey};
70
71
static const geokey_t GTiffShortKeys[] = {
72
    GTModelTypeGeoKey,      GTRasterTypeGeoKey,      GeogAngularUnitsGeoKey,
73
    GeogEllipsoidGeoKey,    GeogGeodeticDatumGeoKey, GeographicTypeGeoKey,
74
    ProjCoordTransGeoKey,   ProjLinearUnitsGeoKey,   ProjectedCSTypeGeoKey,
75
    ProjectionGeoKey,       GeogPrimeMeridianGeoKey, GeogLinearUnitsGeoKey,
76
    GeogAzimuthUnitsGeoKey, VerticalCSTypeGeoKey,    VerticalDatumGeoKey,
77
    VerticalUnitsGeoKey};
78
#endif
79
80
/************************************************************************/
81
/*                     OGRVICARBinaryPrefixesLayer                      */
82
/************************************************************************/
83
84
class OGRVICARBinaryPrefixesLayer final : public OGRLayer
85
{
86
    VSILFILE *m_fp = nullptr;
87
    OGRFeatureDefn *m_poFeatureDefn = nullptr;
88
    int m_iRecord = 0;
89
    int m_nRecords = 0;
90
    vsi_l_offset m_nFileOffset = 0;
91
    vsi_l_offset m_nStride = 0;
92
    bool m_bError = false;
93
    bool m_bByteSwapIntegers = false;
94
    RawRasterBand::ByteOrder m_eBREALByteOrder{};
95
96
    enum Type
97
    {
98
        FIELD_UNKNOWN,
99
        FIELD_UNSIGNED_CHAR,
100
        FIELD_UNSIGNED_SHORT,
101
        FIELD_UNSIGNED_INT,
102
        FIELD_SHORT,
103
        FIELD_INT,
104
        FIELD_FLOAT,
105
        FIELD_DOUBLE,
106
    };
107
108
    static Type GetTypeFromString(const char *pszStr);
109
110
    struct Field
111
    {
112
        int nOffset;
113
        Type eType;
114
    };
115
116
    std::vector<Field> m_aoFields{};
117
    std::vector<GByte> m_abyRecord{};
118
119
    CPL_DISALLOW_COPY_ASSIGN(OGRVICARBinaryPrefixesLayer)
120
121
    OGRFeature *GetNextRawFeature();
122
123
  public:
124
    OGRVICARBinaryPrefixesLayer(VSILFILE *fp, int nRecords,
125
                                const CPLJSONObject &oDef,
126
                                vsi_l_offset nFileOffset, vsi_l_offset nStride,
127
                                RawRasterBand::ByteOrder eBINTByteOrder,
128
                                RawRasterBand::ByteOrder eBREALByteOrder);
129
    ~OGRVICARBinaryPrefixesLayer() override;
130
131
    bool HasError() const
132
0
    {
133
0
        return m_bError;
134
0
    }
135
136
    void ResetReading() override
137
0
    {
138
0
        m_iRecord = 0;
139
0
    }
140
141
    const OGRFeatureDefn *GetLayerDefn() const override
142
0
    {
143
0
        return m_poFeatureDefn;
144
0
    }
145
146
    OGRFeature *GetNextFeature() override;
147
148
    bool TestCapability(const char *) const override
149
0
    {
150
0
        return false;
151
0
    }
152
};
153
154
/************************************************************************/
155
/*                         GetTypeFromString()                          */
156
/************************************************************************/
157
158
OGRVICARBinaryPrefixesLayer::Type
159
OGRVICARBinaryPrefixesLayer::GetTypeFromString(const char *pszStr)
160
0
{
161
0
    if (EQUAL(pszStr, "unsigned char") || EQUAL(pszStr, "unsigned byte"))
162
0
        return FIELD_UNSIGNED_CHAR;
163
0
    if (EQUAL(pszStr, "unsigned short"))
164
0
        return FIELD_UNSIGNED_SHORT;
165
0
    if (EQUAL(pszStr, "unsigned int"))
166
0
        return FIELD_UNSIGNED_INT;
167
0
    if (EQUAL(pszStr, "short"))
168
0
        return FIELD_SHORT;
169
0
    if (EQUAL(pszStr, "int"))
170
0
        return FIELD_INT;
171
0
    if (EQUAL(pszStr, "float"))
172
0
        return FIELD_FLOAT;
173
0
    if (EQUAL(pszStr, "double"))
174
0
        return FIELD_DOUBLE;
175
0
    return FIELD_UNKNOWN;
176
0
}
177
178
/************************************************************************/
179
/*                    OGRVICARBinaryPrefixesLayer()                     */
180
/************************************************************************/
181
182
OGRVICARBinaryPrefixesLayer::OGRVICARBinaryPrefixesLayer(
183
    VSILFILE *fp, int nRecords, const CPLJSONObject &oDef,
184
    vsi_l_offset nFileOffset, vsi_l_offset nStride,
185
    RawRasterBand::ByteOrder eBINTByteOrder,
186
    RawRasterBand::ByteOrder eBREALByteOrder)
187
0
    : m_fp(fp), m_nRecords(nRecords), m_nFileOffset(nFileOffset),
188
0
      m_nStride(nStride),
189
#ifdef CPL_LSB
190
0
      m_bByteSwapIntegers(eBINTByteOrder !=
191
0
                          RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN),
192
#else
193
      m_bByteSwapIntegers(eBINTByteOrder !=
194
                          RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN),
195
#endif
196
0
      m_eBREALByteOrder(eBREALByteOrder)
197
0
{
198
0
    m_poFeatureDefn = new OGRFeatureDefn("binary_prefixes");
199
0
    SetDescription(m_poFeatureDefn->GetName());
200
0
    m_poFeatureDefn->Reference();
201
0
    m_poFeatureDefn->SetGeomType(wkbNone);
202
0
    int nRecordSize = oDef.GetInteger("size");
203
0
    const auto oFields = oDef.GetObj("fields");
204
0
    if (oFields.IsValid() && oFields.GetType() == CPLJSONObject::Type::Array)
205
0
    {
206
0
        auto oFieldsArray = oFields.ToArray();
207
0
        int nOffset = 0;
208
0
        for (int i = 0; i < oFieldsArray.Size(); i++)
209
0
        {
210
0
            auto oField = oFieldsArray[i];
211
0
            if (oField.GetType() == CPLJSONObject::Type::Object)
212
0
            {
213
0
                auto osName = oField.GetString("name");
214
0
                auto osType = oField.GetString("type");
215
0
                auto bHidden = oField.GetBool("hidden");
216
0
                auto eType = GetTypeFromString(osType.c_str());
217
0
                if (eType == FIELD_UNKNOWN)
218
0
                {
219
0
                    CPLError(CE_Failure, CPLE_NotSupported,
220
0
                             "Field %s of type %s not supported",
221
0
                             osName.c_str(), osType.c_str());
222
0
                    m_bError = true;
223
0
                    return;
224
0
                }
225
0
                else if (!osName.empty())
226
0
                {
227
0
                    OGRFieldType eFieldType(OFTMaxType);
228
0
                    Field f;
229
0
                    f.nOffset = nOffset;
230
0
                    f.eType = eType;
231
0
                    switch (eType)
232
0
                    {
233
0
                        case FIELD_UNSIGNED_CHAR:
234
0
                            nOffset += 1;
235
0
                            eFieldType = OFTInteger;
236
0
                            break;
237
0
                        case FIELD_UNSIGNED_SHORT:
238
0
                            nOffset += 2;
239
0
                            eFieldType = OFTInteger;
240
0
                            break;
241
0
                        case FIELD_UNSIGNED_INT:
242
0
                            nOffset += 4;
243
0
                            eFieldType = OFTInteger64;
244
0
                            break;
245
0
                        case FIELD_SHORT:
246
0
                            nOffset += 2;
247
0
                            eFieldType = OFTInteger;
248
0
                            break;
249
0
                        case FIELD_INT:
250
0
                            nOffset += 4;
251
0
                            eFieldType = OFTInteger;
252
0
                            break;
253
0
                        case FIELD_FLOAT:
254
0
                            nOffset += 4;
255
0
                            eFieldType = OFTReal;
256
0
                            break;
257
0
                        case FIELD_DOUBLE:
258
0
                            nOffset += 8;
259
0
                            eFieldType = OFTReal;
260
0
                            break;
261
0
                        default:
262
0
                            CPLAssert(false);
263
0
                            break;
264
0
                    }
265
0
                    if (nOffset > nRecordSize)
266
0
                    {
267
0
                        CPLError(CE_Failure, CPLE_AppDefined,
268
0
                                 "Field definitions not consistent with "
269
0
                                 "declared record size");
270
0
                        m_bError = true;
271
0
                        return;
272
0
                    }
273
0
                    if (!bHidden)
274
0
                    {
275
0
                        m_aoFields.push_back(f);
276
0
                        OGRFieldDefn oFieldDefn(osName.c_str(), eFieldType);
277
0
                        m_poFeatureDefn->AddFieldDefn(&oFieldDefn);
278
0
                    }
279
0
                }
280
0
                else
281
0
                {
282
0
                    m_bError = true;
283
0
                }
284
0
            }
285
0
            else
286
0
            {
287
0
                m_bError = true;
288
0
            }
289
0
            if (m_bError)
290
0
            {
291
0
                CPLError(CE_Failure, CPLE_AppDefined,
292
0
                         "Error while reading binary prefix definition");
293
0
                return;
294
0
            }
295
0
        }
296
0
    }
297
0
    m_abyRecord.resize(nRecordSize);
298
0
}
299
300
/************************************************************************/
301
/*                    ~OGRVICARBinaryPrefixesLayer()                    */
302
/************************************************************************/
303
304
OGRVICARBinaryPrefixesLayer::~OGRVICARBinaryPrefixesLayer()
305
0
{
306
0
    m_poFeatureDefn->Release();
307
0
}
308
309
/************************************************************************/
310
/*                         GetNextRawFeature()                          */
311
/************************************************************************/
312
313
OGRFeature *OGRVICARBinaryPrefixesLayer::GetNextRawFeature()
314
0
{
315
0
    if (m_iRecord >= m_nRecords)
316
0
        return nullptr;
317
318
0
    if (VSIFSeekL(m_fp, m_nFileOffset + m_iRecord * m_nStride, SEEK_SET) != 0 ||
319
0
        VSIFReadL(&m_abyRecord[0], m_abyRecord.size(), 1, m_fp) != 1)
320
0
    {
321
0
        return nullptr;
322
0
    }
323
324
0
    OGRFeature *poFeature = new OGRFeature(m_poFeatureDefn);
325
0
    for (int i = 0; i < poFeature->GetFieldCount(); i++)
326
0
    {
327
0
        int nOffset = m_aoFields[i].nOffset;
328
0
        switch (m_aoFields[i].eType)
329
0
        {
330
0
            case FIELD_UNSIGNED_CHAR:
331
0
                poFeature->SetField(i, m_abyRecord[nOffset]);
332
0
                break;
333
0
            case FIELD_UNSIGNED_SHORT:
334
0
            {
335
0
                unsigned short v;
336
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
337
0
                if (m_bByteSwapIntegers)
338
0
                {
339
0
                    CPL_SWAP16PTR(&v);
340
0
                }
341
0
                poFeature->SetField(i, v);
342
0
                break;
343
0
            }
344
0
            case FIELD_UNSIGNED_INT:
345
0
            {
346
0
                unsigned int v;
347
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
348
0
                if (m_bByteSwapIntegers)
349
0
                {
350
0
                    CPL_SWAP32PTR(&v);
351
0
                }
352
0
                poFeature->SetField(i, static_cast<GIntBig>(v));
353
0
                break;
354
0
            }
355
0
            case FIELD_SHORT:
356
0
            {
357
0
                short v;
358
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
359
0
                if (m_bByteSwapIntegers)
360
0
                {
361
0
                    CPL_SWAP16PTR(&v);
362
0
                }
363
0
                poFeature->SetField(i, v);
364
0
                break;
365
0
            }
366
0
            case FIELD_INT:
367
0
            {
368
0
                int v;
369
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
370
0
                if (m_bByteSwapIntegers)
371
0
                {
372
0
                    CPL_SWAP32PTR(&v);
373
0
                }
374
0
                poFeature->SetField(i, v);
375
0
                break;
376
0
            }
377
0
            case FIELD_FLOAT:
378
0
            {
379
0
                float v;
380
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
381
0
                if (m_eBREALByteOrder == RawRasterBand::ByteOrder::ORDER_VAX)
382
0
                {
383
0
                    CPLVaxToIEEEFloat(&v);
384
0
                }
385
0
                else if (m_eBREALByteOrder !=
386
0
#ifdef CPL_LSB
387
0
                         RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN
388
#else
389
                         RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN
390
#endif
391
0
                )
392
0
                {
393
0
                    CPL_SWAP32PTR(&v);
394
0
                }
395
0
                poFeature->SetField(i, v);
396
0
                break;
397
0
            }
398
0
            case FIELD_DOUBLE:
399
0
            {
400
0
                double v;
401
0
                memcpy(&v, &m_abyRecord[nOffset], sizeof(v));
402
0
                if (m_eBREALByteOrder == RawRasterBand::ByteOrder::ORDER_VAX)
403
0
                {
404
0
                    CPLVaxToIEEEDouble(&v);
405
0
                }
406
0
                else if (m_eBREALByteOrder !=
407
0
#ifdef CPL_LSB
408
0
                         RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN
409
#else
410
                         RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN
411
#endif
412
0
                )
413
0
                {
414
0
                    CPL_SWAP64PTR(&v);
415
0
                }
416
0
                poFeature->SetField(i, v);
417
0
                break;
418
0
            }
419
0
            default:
420
0
                CPLAssert(false);
421
0
        }
422
0
    }
423
0
    poFeature->SetFID(m_iRecord);
424
0
    m_iRecord++;
425
0
    return poFeature;
426
0
}
427
428
/************************************************************************/
429
/*                           GetNextFeature()                           */
430
/************************************************************************/
431
432
OGRFeature *OGRVICARBinaryPrefixesLayer::GetNextFeature()
433
0
{
434
0
    while (true)
435
0
    {
436
0
        auto poFeature = GetNextRawFeature();
437
0
        if (poFeature == nullptr)
438
0
            return nullptr;
439
440
0
        if ((m_poFilterGeom == nullptr ||
441
0
             FilterGeometry(poFeature->GetGeometryRef())) &&
442
0
            (m_poAttrQuery == nullptr || m_poAttrQuery->Evaluate(poFeature)))
443
0
        {
444
0
            return poFeature;
445
0
        }
446
0
        else
447
0
            delete poFeature;
448
0
    }
449
0
}
450
451
/************************************************************************/
452
/*                          VICARRawRasterBand                          */
453
/************************************************************************/
454
455
class VICARRawRasterBand final : public RawRasterBand
456
{
457
  protected:
458
    friend class VICARDataset;
459
460
  public:
461
    VICARRawRasterBand(VICARDataset *poDSIn, int nBandIn, VSILFILE *fpRawIn,
462
                       vsi_l_offset nImgOffsetIn, int nPixelOffsetIn,
463
                       int nLineOffsetIn, GDALDataType eDataTypeIn,
464
                       ByteOrder eByteOrderIn);
465
466
    CPLErr IReadBlock(int, int, void *) override;
467
    CPLErr IWriteBlock(int, int, void *) override;
468
469
    CPLErr IRasterIO(GDALRWFlag, int, int, int, int, void *, int, int,
470
                     GDALDataType, GSpacing nPixelSpace, GSpacing nLineSpace,
471
                     GDALRasterIOExtraArg *psExtraArg) override;
472
};
473
474
/************************************************************************/
475
/*                         VICARRawRasterBand()                         */
476
/************************************************************************/
477
478
VICARRawRasterBand::VICARRawRasterBand(VICARDataset *poDSIn, int nBandIn,
479
                                       VSILFILE *fpRawIn,
480
                                       vsi_l_offset nImgOffsetIn,
481
                                       int nPixelOffsetIn, int nLineOffsetIn,
482
                                       GDALDataType eDataTypeIn,
483
                                       ByteOrder eByteOrderIn)
484
19
    : RawRasterBand(poDSIn, nBandIn, fpRawIn, nImgOffsetIn, nPixelOffsetIn,
485
19
                    nLineOffsetIn, eDataTypeIn, eByteOrderIn,
486
19
                    RawRasterBand::OwnFP::NO)
487
19
{
488
19
}
489
490
/************************************************************************/
491
/*                             IReadBlock()                             */
492
/************************************************************************/
493
494
CPLErr VICARRawRasterBand::IReadBlock(int nXBlock, int nYBlock, void *pImage)
495
496
116
{
497
116
    VICARDataset *poGDS = cpl::down_cast<VICARDataset *>(poDS);
498
116
    if (!poGDS->m_bIsLabelWritten)
499
0
        poGDS->WriteLabel();
500
116
    return RawRasterBand::IReadBlock(nXBlock, nYBlock, pImage);
501
116
}
502
503
/************************************************************************/
504
/*                            IWriteBlock()                             */
505
/************************************************************************/
506
507
CPLErr VICARRawRasterBand::IWriteBlock(int nXBlock, int nYBlock, void *pImage)
508
509
0
{
510
0
    VICARDataset *poGDS = cpl::down_cast<VICARDataset *>(poDS);
511
0
    if (!poGDS->m_bIsLabelWritten)
512
0
        poGDS->WriteLabel();
513
0
    return RawRasterBand::IWriteBlock(nXBlock, nYBlock, pImage);
514
0
}
515
516
/************************************************************************/
517
/*                             IRasterIO()                              */
518
/************************************************************************/
519
520
CPLErr VICARRawRasterBand::IRasterIO(GDALRWFlag eRWFlag, int nXOff, int nYOff,
521
                                     int nXSize, int nYSize, void *pData,
522
                                     int nBufXSize, int nBufYSize,
523
                                     GDALDataType eBufType,
524
                                     GSpacing nPixelSpace, GSpacing nLineSpace,
525
                                     GDALRasterIOExtraArg *psExtraArg)
526
527
116
{
528
116
    VICARDataset *poGDS = cpl::down_cast<VICARDataset *>(poDS);
529
116
    if (!poGDS->m_bIsLabelWritten)
530
0
        poGDS->WriteLabel();
531
116
    return RawRasterBand::IRasterIO(eRWFlag, nXOff, nYOff, nXSize, nYSize,
532
116
                                    pData, nBufXSize, nBufYSize, eBufType,
533
116
                                    nPixelSpace, nLineSpace, psExtraArg);
534
116
}
535
536
/************************************************************************/
537
/*                         VICARBASICRasterBand                         */
538
/************************************************************************/
539
540
class VICARBASICRasterBand final : public GDALPamRasterBand
541
{
542
  public:
543
    VICARBASICRasterBand(VICARDataset *poDSIn, int nBandIn, GDALDataType eType);
544
545
    CPLErr IReadBlock(int, int, void *) override;
546
    CPLErr IWriteBlock(int, int, void *) override;
547
};
548
549
/************************************************************************/
550
/*                        VICARBASICRasterBand()                        */
551
/************************************************************************/
552
553
VICARBASICRasterBand::VICARBASICRasterBand(VICARDataset *poDSIn, int nBandIn,
554
                                           GDALDataType eType)
555
22
{
556
22
    poDS = poDSIn;
557
22
    nBand = nBandIn;
558
22
    nBlockXSize = poDSIn->GetRasterXSize();
559
22
    nBlockYSize = 1;
560
22
    eDataType = eType;
561
22
}
562
563
namespace
564
{
565
class DecodeEncodeException final : public std::exception
566
{
567
  public:
568
4
    DecodeEncodeException() = default;
569
};
570
}  // namespace
571
572
//////////////////////////////////////////////////////////////////////////
573
/// Below functions are adapted from Public Domain VICAR project
574
/// from
575
/// https://github.com/nasa/VICAR/blob/master/vos/rtl/source/basic_compression.c
576
//////////////////////////////////////////////////////////////////////////
577
578
/* masking array used in the algorithm to take out bits in memory */
579
const unsigned int cod1mask[25] = {
580
    0x0,      0x1,      0x3,      0x7,     0xf,     0x1f,    0x3f,
581
    0x7f,     0xff,     0x1ff,    0x3ff,   0x7ff,   0xfff,   0x1fff,
582
    0x3fff,   0x7fff,   0xffff,   0x1ffff, 0x3ffff, 0x7ffff, 0xfffff,
583
    0x1fffff, 0x3fffff, 0x7fffff, 0xffffff};
584
585
/*****************************************************/
586
/* This function is a helper function for the BASIC  */
587
/* compression algorithm to get a specified number   */
588
/* of bits from buffer, convert it to a number and   */
589
/* return it.                                        */
590
/*****************************************************/
591
static unsigned char grab1(int nbit, const unsigned char *buffer,
592
                           size_t buffer_size, size_t &buffer_pos,
593
                           int &bit1ptr) /* bit position in the current byte of
594
                                            encrypted or decrypted buffer*/
595
295k
{
596
295k
    unsigned char val;
597
295k
    int shift = 8 - nbit - (bit1ptr);
598
599
295k
    if (buffer_pos >= buffer_size)
600
1
    {
601
1
        CPLError(CE_Failure, CPLE_AppDefined, "Out of decoding buffer");
602
1
        throw DecodeEncodeException();
603
1
    }
604
605
295k
    if (shift > 0)
606
182k
    {
607
182k
        val = (buffer[buffer_pos] >> shift) & cod1mask[nbit];
608
182k
        bit1ptr += nbit;
609
610
182k
        return val;
611
182k
    }
612
113k
    if (shift < 0)
613
74.7k
    {
614
74.7k
        unsigned v1 = buffer[buffer_pos] & cod1mask[nbit + shift];
615
74.7k
        buffer_pos++;
616
617
74.7k
        if (buffer_pos >= buffer_size)
618
3
        {
619
3
            CPLError(CE_Failure, CPLE_AppDefined, "Out of decoding buffer");
620
3
            throw DecodeEncodeException();
621
3
        }
622
623
74.7k
        unsigned v2 = (buffer[buffer_pos] >> (8 + shift)) & cod1mask[-shift];
624
625
74.7k
        val = static_cast<unsigned char>((v1 << (-shift)) + v2);
626
627
74.7k
        bit1ptr = -shift;
628
629
74.7k
        return val;
630
74.7k
    }
631
38.3k
    val = buffer[buffer_pos] & cod1mask[nbit];
632
38.3k
    buffer_pos++;
633
38.3k
    bit1ptr = 0;
634
635
38.3k
    return val;
636
113k
}
637
638
/*****************************************************/
639
/* This function is the decoding algorithm for BASIC */
640
/* compression.  The encoded buffer is passed into   */
641
/* code and the decoded buffer is passed out in buf. */
642
/*****************************************************/
643
static void basic_decode(const unsigned char *code, size_t code_size,
644
                         unsigned char *buf, int ns, int wid)
645
282
{
646
282
    int runInt = -3;
647
282
    unsigned char runChar;
648
282
    unsigned int nval = 999999;
649
282
    static const int cmprtrns1[7] = {-3, -2, -1, 0, 1, 2, 3};
650
282
    size_t buffer_pos = 0;
651
282
    int bit1ptr = 0;
652
282
    unsigned int old = 0;
653
282
    const int ptop = ns * wid;
654
655
564
    for (int iw = 0; iw < wid; iw++)
656
282
    {
657
556k
        for (int ip = iw; ip < ptop; ip += wid)
658
556k
        {
659
556k
            if (runInt > (-3))
660
281k
            {
661
281k
                buf[ip] = static_cast<unsigned char>(nval);
662
281k
                runInt--;
663
281k
                continue;
664
281k
            }
665
274k
            unsigned char val = grab1(3, code, code_size, buffer_pos, bit1ptr);
666
667
274k
            if (val < 7)
668
266k
            {
669
266k
                nval = CPLUnsanitizedAdd<unsigned>(old, cmprtrns1[val]);
670
266k
                buf[ip] = static_cast<unsigned char>(nval);
671
266k
                old = nval;
672
266k
                continue;
673
266k
            }
674
8.22k
            val = grab1(1, code, code_size, buffer_pos, bit1ptr);
675
676
8.22k
            if (val)
677
2.92k
            {
678
2.92k
                runChar = grab1(4, code, code_size, buffer_pos, bit1ptr);
679
2.92k
                if (runChar == 15)
680
224
                {
681
224
                    runChar = grab1(8, code, code_size, buffer_pos, bit1ptr);
682
683
224
                    if (runChar == 255)
684
161
                    {
685
161
                        unsigned char part0 =
686
161
                            grab1(8, code, code_size, buffer_pos, bit1ptr);
687
161
                        unsigned char part1 =
688
161
                            grab1(8, code, code_size, buffer_pos, bit1ptr);
689
161
                        unsigned char part2 =
690
161
                            grab1(8, code, code_size, buffer_pos, bit1ptr);
691
161
                        runInt = part0 | (part1 << 8) | (part2 << 16);
692
161
                    }
693
63
                    else
694
63
                        runInt = runChar + 15;
695
224
                }
696
2.70k
                else
697
2.70k
                    runInt = runChar;
698
699
2.92k
                val = grab1(3, code, code_size, buffer_pos, bit1ptr);
700
2.92k
                if (val < 7)
701
2.32k
                    nval = CPLUnsanitizedAdd<unsigned>(old, cmprtrns1[val]);
702
602
                else
703
602
                    nval = grab1(8, code, code_size, buffer_pos, bit1ptr);
704
2.92k
                buf[ip] = static_cast<unsigned char>(nval);
705
2.92k
                old = nval;
706
2.92k
            }
707
5.30k
            else
708
5.30k
            {
709
5.30k
                val = grab1(8, code, code_size, buffer_pos, bit1ptr);
710
5.30k
                buf[ip] = val;
711
5.30k
                old = val;
712
5.30k
            }
713
8.22k
        }
714
282
    }
715
282
}
716
717
/*****************************************************/
718
/* This function is a helper function for the BASIC  */
719
/* encoding operation.  It puts the value in val into*/
720
/* memory location pointed by pcode1+reg1 into nbit  */
721
/* number of bits.                                   */
722
/*****************************************************/
723
static void emit1(unsigned char val, int nbit, unsigned char *reg1,
724
                  int &bit1ptr, unsigned char *coded_buffer,
725
                  size_t &coded_buffer_pos, size_t coded_buffer_size)
726
0
{
727
0
    int shift;
728
729
0
    shift = 8 - nbit - bit1ptr;
730
0
    if (shift > 0)
731
0
    {
732
0
        *reg1 = static_cast<unsigned char>(*reg1 | (val << shift));
733
0
        bit1ptr += nbit;
734
0
        return;
735
0
    }
736
0
    if (shift < 0)
737
0
    {
738
0
        if (coded_buffer_pos >= coded_buffer_size)
739
0
        {
740
0
            CPLError(CE_Failure, CPLE_AppDefined, "Out of encoding buffer");
741
0
            throw DecodeEncodeException();
742
0
        }
743
0
        coded_buffer[coded_buffer_pos] =
744
0
            static_cast<unsigned char>(*reg1 | (val >> (-shift)));
745
0
        coded_buffer_pos++;
746
0
        *reg1 = static_cast<unsigned char>(val << (8 + shift));
747
0
        bit1ptr = -shift;
748
0
        return;
749
0
    }
750
0
    if (coded_buffer_pos >= coded_buffer_size)
751
0
    {
752
0
        CPLError(CE_Failure, CPLE_AppDefined, "Out of encoding buffer");
753
0
        throw DecodeEncodeException();
754
0
    }
755
0
    coded_buffer[coded_buffer_pos] = static_cast<unsigned char>(*reg1 | val);
756
0
    coded_buffer_pos++;
757
758
0
    *reg1 = 0;
759
0
    bit1ptr = 0;
760
0
}
761
762
/*****************************************************/
763
/* This function is meat of the BASIC encoding       */
764
/* algorithm.  This function is called repeatedly by */
765
/* the basic_encode function to compress the data    */
766
/* according to its run length (run), last 2         */
767
/* different values (vold and old), and the current  */
768
/* value (val), into the memory location pointed by  */
769
/* pcode1+reg1.                                      */
770
/*****************************************************/
771
static void basic_encrypt(int *run, int *old, int *vold, int val,
772
                          unsigned char *reg1, int &bit1ptr,
773
                          unsigned char *coded_buffer, size_t &coded_buffer_pos,
774
                          size_t coded_buffer_size)
775
0
{
776
0
    if (*run < 4)
777
0
    {
778
0
        if (abs(*old - *vold) < 4)
779
0
            emit1(static_cast<unsigned char>(*old - *vold + 3), 3, reg1,
780
0
                  bit1ptr, coded_buffer, coded_buffer_pos, coded_buffer_size);
781
0
        else
782
0
        {
783
0
            emit1(static_cast<unsigned char>(14), 4, reg1, bit1ptr,
784
0
                  coded_buffer, coded_buffer_pos, coded_buffer_size);
785
0
            emit1(static_cast<unsigned char>(*old), 8, reg1, bit1ptr,
786
0
                  coded_buffer, coded_buffer_pos, coded_buffer_size);
787
0
        }
788
789
0
        while (*run > 1)
790
0
        {
791
0
            emit1(static_cast<unsigned char>(3), 3, reg1, bit1ptr, coded_buffer,
792
0
                  coded_buffer_pos, coded_buffer_size);
793
0
            (*run)--;
794
0
        }
795
796
0
        *vold = *old;
797
0
        *old = val;
798
0
    }
799
0
    else
800
0
    {
801
0
        emit1(static_cast<unsigned char>(15), 4, reg1, bit1ptr, coded_buffer,
802
0
              coded_buffer_pos, coded_buffer_size);
803
0
        if (*run < 19)
804
0
        {
805
0
            emit1(static_cast<unsigned char>(*run - 4), 4, reg1, bit1ptr,
806
0
                  coded_buffer, coded_buffer_pos, coded_buffer_size);
807
0
        }
808
0
        else
809
0
        {
810
0
            emit1(static_cast<unsigned char>(15), 4, reg1, bit1ptr,
811
0
                  coded_buffer, coded_buffer_pos, coded_buffer_size);
812
0
            if (*run < 274)
813
0
            {
814
0
                emit1(static_cast<char>(*run - 19), 8, reg1, bit1ptr,
815
0
                      coded_buffer, coded_buffer_pos, coded_buffer_size);
816
0
            }
817
0
            else
818
0
            {
819
0
                emit1(static_cast<unsigned char>(255), 8, reg1, bit1ptr,
820
0
                      coded_buffer, coded_buffer_pos, coded_buffer_size);
821
822
0
                unsigned char part0 =
823
0
                    static_cast<unsigned char>((*run - 4) & 0xff);
824
0
                unsigned char part1 =
825
0
                    static_cast<unsigned char>(((*run - 4) >> 8) & 0xff);
826
0
                unsigned char part2 =
827
0
                    static_cast<unsigned char>(((*run - 4) >> 16) & 0xff);
828
0
                emit1(part0, 8, reg1, bit1ptr, coded_buffer, coded_buffer_pos,
829
0
                      coded_buffer_size);
830
0
                emit1(part1, 8, reg1, bit1ptr, coded_buffer, coded_buffer_pos,
831
0
                      coded_buffer_size);
832
0
                emit1(part2, 8, reg1, bit1ptr, coded_buffer, coded_buffer_pos,
833
0
                      coded_buffer_size);
834
0
            }
835
0
        }
836
0
        if (abs(*old - *vold) < 4)
837
0
        {
838
0
            emit1(static_cast<unsigned char>(*old - *vold + 3), 3, reg1,
839
0
                  bit1ptr, coded_buffer, coded_buffer_pos, coded_buffer_size);
840
0
        }
841
0
        else
842
0
        {
843
0
            emit1(static_cast<unsigned char>(7), 3, reg1, bit1ptr, coded_buffer,
844
0
                  coded_buffer_pos, coded_buffer_size);
845
0
            emit1(static_cast<unsigned char>(*old), 8, reg1, bit1ptr,
846
0
                  coded_buffer, coded_buffer_pos, coded_buffer_size);
847
0
        }
848
0
        *vold = *old;
849
0
        *old = val;
850
0
        *run = 1;
851
0
    }
852
0
}
853
854
/*****************************************************/
855
/* This function loops through the data given by     */
856
/* unencodedBuf, keeping track of run length.  When  */
857
/* the value of the data changes, it passes the run  */
858
/* length, last 2 differing values, the current      */
859
/* value, and the pointer in pcode1 buffer to encode */
860
/* the data into, to basic_encrypt function.         */
861
/*****************************************************/
862
static void basic_encode(const unsigned char *unencodedBuf,
863
                         unsigned char *coded_buffer, size_t coded_buffer_size,
864
                         int ns, int wid, size_t *totBytes)
865
0
{
866
0
    int val = 0;
867
0
    int bit1ptr = 0;
868
0
    const int ptop = ns * wid;
869
0
    unsigned char reg1 = 0;
870
0
    int run = 0;
871
0
    int old = unencodedBuf[0];
872
0
    int vold = 999999;
873
874
0
    size_t coded_buffer_pos = 0;
875
876
0
    for (int iw = 0; iw < wid; iw++)
877
0
    {
878
0
        for (int ip = iw; ip < ptop; ip += wid)
879
0
        {
880
0
            val = unencodedBuf[ip];
881
882
0
            if (val == old)
883
0
                run++;
884
0
            else
885
0
                basic_encrypt(&run, &old, &vold, val, &reg1, bit1ptr,
886
0
                              coded_buffer, coded_buffer_pos,
887
0
                              coded_buffer_size);
888
0
        }
889
0
    }
890
891
    /* purge of last code */
892
0
    basic_encrypt(&run, &old, &vold, val, &reg1, bit1ptr, coded_buffer,
893
0
                  coded_buffer_pos, coded_buffer_size);
894
895
0
    if (coded_buffer_pos >= coded_buffer_size)
896
0
    {
897
0
        CPLError(CE_Failure, CPLE_AppDefined, "Out of encoding buffer");
898
0
        throw DecodeEncodeException();
899
0
    }
900
0
    coded_buffer[coded_buffer_pos] = reg1;
901
902
0
    *totBytes = coded_buffer_pos;
903
0
    if (bit1ptr > 0)
904
0
        (*totBytes)++;
905
0
}
906
907
//////////////////////////////////////////////////////////////////////////
908
/// End of VICAR code
909
//////////////////////////////////////////////////////////////////////////
910
911
/************************************************************************/
912
/*                             IReadBlock()                             */
913
/************************************************************************/
914
915
CPLErr VICARBASICRasterBand::IReadBlock(int /*nXBlock*/, int nYBlock,
916
                                        void *pImage)
917
918
299
{
919
299
    VICARDataset *poGDS = cpl::down_cast<VICARDataset *>(poDS);
920
921
299
    const int nRecord = (nBand - 1) * nRasterYSize + nYBlock;
922
299
    const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
923
924
299
    if (poGDS->eAccess == GA_Update &&
925
0
        poGDS->m_anRecordOffsets[nRecord + 1] == 0)
926
0
    {
927
0
        memset(pImage, 0, static_cast<size_t>(nDTSize) * nRasterXSize);
928
0
        return CE_None;
929
0
    }
930
931
    // Find at which offset the compressed record is.
932
    // For BASIC compression, each compressed run is preceded by a uint32 value
933
    // given its size, including the size of this uint32 value
934
    // For BASIC2 compression, the uint32 sizes of all records are put
935
    // immediately after the label.
936
598
    for (; poGDS->m_nLastRecordOffset <= nRecord; poGDS->m_nLastRecordOffset++)
937
299
    {
938
299
        CPLAssert(poGDS->m_anRecordOffsets[poGDS->m_nLastRecordOffset + 1] ==
939
299
                  0);
940
941
299
        int nRet;
942
299
        if (poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC)
943
0
        {
944
0
            nRet =
945
0
                VSIFSeekL(poGDS->fpImage,
946
0
                          poGDS->m_anRecordOffsets[poGDS->m_nLastRecordOffset] -
947
0
                              sizeof(GUInt32),
948
0
                          SEEK_SET);
949
0
        }
950
299
        else
951
299
        {
952
299
            nRet = VSIFSeekL(poGDS->fpImage,
953
299
                             poGDS->m_nImageOffsetWithoutNBB +
954
299
                                 static_cast<vsi_l_offset>(sizeof(GUInt32)) *
955
299
                                     poGDS->m_nLastRecordOffset,
956
299
                             SEEK_SET);
957
299
        }
958
299
        GUInt32 nSize;
959
299
        if (nRet != 0 ||
960
299
            VSIFReadL(&nSize, sizeof(nSize), 1, poGDS->fpImage) != 1)
961
0
        {
962
0
            CPLError(CE_Failure, CPLE_AppDefined, "Cannot read record %d size",
963
0
                     poGDS->m_nLastRecordOffset);
964
0
            return CE_Failure;
965
0
        }
966
299
        CPL_LSBPTR32(&nSize);
967
299
        if ((poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC &&
968
0
             nSize <= sizeof(GUInt32)) ||
969
299
            (poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC2 &&
970
299
             nSize == 0) ||
971
299
            poGDS->m_anRecordOffsets[poGDS->m_nLastRecordOffset] >
972
299
                std::numeric_limits<uint64_t>::max() - nSize)
973
0
        {
974
0
            CPLError(CE_Failure, CPLE_AppDefined, "Wrong size at record %d",
975
0
                     poGDS->m_nLastRecordOffset);
976
0
            return CE_Failure;
977
0
        }
978
979
299
        poGDS->m_anRecordOffsets[poGDS->m_nLastRecordOffset + 1] =
980
299
            poGDS->m_anRecordOffsets[poGDS->m_nLastRecordOffset] + nSize;
981
299
    }
982
983
299
    unsigned int nSize;
984
299
    if (poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC)
985
0
    {
986
0
        nSize = static_cast<unsigned>(poGDS->m_anRecordOffsets[nRecord + 1] -
987
0
                                      poGDS->m_anRecordOffsets[nRecord] -
988
0
                                      sizeof(GUInt32));
989
0
    }
990
299
    else
991
299
    {
992
299
        nSize = static_cast<unsigned>(poGDS->m_anRecordOffsets[nRecord + 1] -
993
299
                                      poGDS->m_anRecordOffsets[nRecord]);
994
299
    }
995
299
    if (nSize > 100 * 1000 * 1000 ||
996
298
        (nSize > 1000 &&
997
199
         (nSize - 11) / 4 > static_cast<unsigned>(nRasterXSize) * nDTSize))
998
1
    {
999
1
        return CE_Failure;
1000
1
    }
1001
298
    if (poGDS->m_abyCodedBuffer.size() < nSize)
1002
22
    {
1003
22
        try
1004
22
        {
1005
22
            poGDS->m_abyCodedBuffer.resize(nSize);
1006
22
        }
1007
22
        catch (const std::exception &e)
1008
22
        {
1009
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
1010
0
            return CE_Failure;
1011
0
        }
1012
22
    }
1013
298
    if (VSIFSeekL(poGDS->fpImage, poGDS->m_anRecordOffsets[nRecord],
1014
298
                  SEEK_SET) != 0 ||
1015
298
        VSIFReadL(&poGDS->m_abyCodedBuffer[0], nSize, 1, poGDS->fpImage) != 1)
1016
16
    {
1017
16
        CPLError(CE_Failure, CPLE_FileIO, "Cannot read record %d", nRecord);
1018
16
        return CE_Failure;
1019
16
    }
1020
1021
282
    try
1022
282
    {
1023
282
        basic_decode(poGDS->m_abyCodedBuffer.data(), nSize,
1024
282
                     static_cast<unsigned char *>(pImage), nRasterXSize,
1025
282
                     nDTSize);
1026
282
    }
1027
282
    catch (const DecodeEncodeException &)
1028
282
    {
1029
4
        return CE_Failure;
1030
4
    }
1031
#ifdef CPL_MSB
1032
    if (nDTSize > 1)
1033
    {
1034
        GDALSwapWords(pImage, nDTSize, nRasterXSize, nDTSize);
1035
    }
1036
#endif
1037
278
    return CE_None;
1038
282
}
1039
1040
/************************************************************************/
1041
/*                            IWriteBlock()                             */
1042
/************************************************************************/
1043
1044
CPLErr VICARBASICRasterBand::IWriteBlock(int /*nXBlock*/, int nYBlock,
1045
                                         void *pImage)
1046
1047
0
{
1048
0
    VICARDataset *poGDS = cpl::down_cast<VICARDataset *>(poDS);
1049
0
    if (poGDS->eAccess == GA_ReadOnly)
1050
0
        return CE_Failure;
1051
0
    if (!poGDS->m_bIsLabelWritten)
1052
0
    {
1053
0
        poGDS->WriteLabel();
1054
0
        poGDS->m_nLabelSize = VSIFTellL(poGDS->fpImage);
1055
0
        poGDS->m_anRecordOffsets[0] = poGDS->m_nLabelSize;
1056
0
        if (poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC)
1057
0
        {
1058
0
            poGDS->m_anRecordOffsets[0] += sizeof(GUInt32);
1059
0
        }
1060
0
        else
1061
0
        {
1062
0
            poGDS->m_anRecordOffsets[0] +=
1063
0
                static_cast<vsi_l_offset>(sizeof(GUInt32)) * nRasterYSize;
1064
0
        }
1065
0
    }
1066
0
    if (nYBlock != poGDS->m_nLastRecordOffset)
1067
0
    {
1068
0
        CPLError(CE_Failure, CPLE_NotSupported,
1069
0
                 "Lines must be written in sequential order");
1070
0
        return CE_Failure;
1071
0
    }
1072
1073
0
    const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
1074
0
    const size_t nMaxEncodedSize =
1075
0
        static_cast<size_t>(nRasterXSize) * nDTSize +
1076
0
        static_cast<size_t>(nRasterXSize) * nDTSize / 2 + 11;
1077
0
    if (poGDS->m_abyCodedBuffer.size() < nMaxEncodedSize)
1078
0
    {
1079
0
        try
1080
0
        {
1081
0
            poGDS->m_abyCodedBuffer.resize(nMaxEncodedSize);
1082
0
        }
1083
0
        catch (const std::exception &e)
1084
0
        {
1085
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
1086
0
            return CE_Failure;
1087
0
        }
1088
0
    }
1089
1090
#ifdef CPL_MSB
1091
    if (nDTSize > 1)
1092
    {
1093
        GDALSwapWords(pImage, nDTSize, nRasterXSize, nDTSize);
1094
    }
1095
#endif
1096
1097
0
    size_t nCodedSize = 0;
1098
0
    try
1099
0
    {
1100
0
        basic_encode(
1101
0
            static_cast<unsigned char *>(pImage), &poGDS->m_abyCodedBuffer[0],
1102
0
            poGDS->m_abyCodedBuffer.size(), nRasterXSize, nDTSize, &nCodedSize);
1103
0
    }
1104
0
    catch (const DecodeEncodeException &)
1105
0
    {
1106
0
        return CE_Failure;
1107
0
    }
1108
1109
#ifdef CPL_MSB
1110
    if (nDTSize > 1)
1111
    {
1112
        GDALSwapWords(pImage, nDTSize, nRasterXSize, nDTSize);
1113
    }
1114
#endif
1115
1116
0
    if (poGDS->m_eCompress == VICARDataset::COMPRESS_BASIC)
1117
0
    {
1118
0
        VSIFSeekL(poGDS->fpImage,
1119
0
                  poGDS->m_anRecordOffsets[nYBlock] - sizeof(GUInt32),
1120
0
                  SEEK_SET);
1121
0
        GUInt32 nSizeToWrite =
1122
0
            static_cast<GUInt32>(nCodedSize + sizeof(GUInt32));
1123
0
        CPL_LSBPTR32(&nSizeToWrite);
1124
0
        VSIFWriteL(&nSizeToWrite, sizeof(GUInt32), 1, poGDS->fpImage);
1125
0
        VSIFWriteL(poGDS->m_abyCodedBuffer.data(), nCodedSize, 1,
1126
0
                   poGDS->fpImage);
1127
0
        poGDS->m_anRecordOffsets[nYBlock + 1] =
1128
0
            poGDS->m_anRecordOffsets[nYBlock] + nCodedSize + sizeof(GUInt32);
1129
0
    }
1130
0
    else
1131
0
    {
1132
0
        VSIFSeekL(poGDS->fpImage,
1133
0
                  poGDS->m_nLabelSize +
1134
0
                      static_cast<vsi_l_offset>(nYBlock) * sizeof(GUInt32),
1135
0
                  SEEK_SET);
1136
0
        GUInt32 nSizeToWrite = static_cast<GUInt32>(nCodedSize);
1137
0
        CPL_LSBPTR32(&nSizeToWrite);
1138
0
        VSIFWriteL(&nSizeToWrite, sizeof(GUInt32), 1, poGDS->fpImage);
1139
0
        VSIFSeekL(poGDS->fpImage, poGDS->m_anRecordOffsets[nYBlock], SEEK_SET);
1140
0
        VSIFWriteL(poGDS->m_abyCodedBuffer.data(), nCodedSize, 1,
1141
0
                   poGDS->fpImage);
1142
0
        poGDS->m_anRecordOffsets[nYBlock + 1] =
1143
0
            poGDS->m_anRecordOffsets[nYBlock] + nCodedSize;
1144
0
    }
1145
1146
0
    poGDS->m_nLastRecordOffset++;
1147
1148
0
    return CE_None;
1149
0
}
1150
1151
/************************************************************************/
1152
/*                            VICARDataset()                            */
1153
/************************************************************************/
1154
1155
VICARDataset::VICARDataset()
1156
1157
271
{
1158
271
    m_oJSonLabel.Deinit();
1159
271
    m_oSrcJSonLabel.Deinit();
1160
271
}
1161
1162
/************************************************************************/
1163
/*                           ~VICARDataset()                            */
1164
/************************************************************************/
1165
1166
VICARDataset::~VICARDataset()
1167
1168
271
{
1169
271
    VICARDataset::Close();
1170
271
}
1171
1172
/************************************************************************/
1173
/*                               Close()                                */
1174
/************************************************************************/
1175
1176
CPLErr VICARDataset::Close(GDALProgressFunc, void *)
1177
312
{
1178
312
    CPLErr eErr = CE_None;
1179
312
    if (nOpenFlags != OPEN_FLAGS_CLOSED)
1180
271
    {
1181
271
        if (!m_bIsLabelWritten)
1182
0
            WriteLabel();
1183
1184
271
        if (VICARDataset::FlushCache(true) != CE_None)
1185
0
            eErr = CE_Failure;
1186
1187
271
        PatchLabel();
1188
271
        if (fpImage)
1189
271
            VSIFCloseL(fpImage);
1190
1191
271
        if (GDALPamDataset::Close() != CE_None)
1192
0
            eErr = CE_Failure;
1193
271
    }
1194
312
    return eErr;
1195
312
}
1196
1197
/************************************************************************/
1198
/*                           GetSpatialRef()                            */
1199
/************************************************************************/
1200
1201
const OGRSpatialReference *VICARDataset::GetSpatialRef() const
1202
1203
38
{
1204
38
    if (!m_oSRS.IsEmpty())
1205
6
        return &m_oSRS;
1206
1207
32
    return GDALPamDataset::GetSpatialRef();
1208
38
}
1209
1210
/************************************************************************/
1211
/*                           SetSpatialRef()                            */
1212
/************************************************************************/
1213
1214
CPLErr VICARDataset::SetSpatialRef(const OGRSpatialReference *poSRS)
1215
0
{
1216
0
    if (eAccess == GA_ReadOnly)
1217
0
        return GDALPamDataset::SetSpatialRef(poSRS);
1218
0
    if (poSRS)
1219
0
        m_oSRS = *poSRS;
1220
0
    else
1221
0
        m_oSRS.Clear();
1222
0
    InvalidateLabel();
1223
0
    return CE_None;
1224
0
}
1225
1226
/************************************************************************/
1227
/*                          GetGeoTransform()                           */
1228
/************************************************************************/
1229
1230
CPLErr VICARDataset::GetGeoTransform(GDALGeoTransform &gt) const
1231
38
{
1232
38
    if (m_bGotTransform)
1233
30
    {
1234
30
        gt = m_gt;
1235
30
        return CE_None;
1236
30
    }
1237
1238
8
    return GDALPamDataset::GetGeoTransform(gt);
1239
38
}
1240
1241
/************************************************************************/
1242
/*                          SetGeoTransform()                           */
1243
/************************************************************************/
1244
1245
CPLErr VICARDataset::SetGeoTransform(const GDALGeoTransform &gt)
1246
1247
0
{
1248
0
    if (eAccess == GA_ReadOnly)
1249
0
        return GDALPamDataset::SetGeoTransform(gt);
1250
0
    if (gt.xscale <= 0.0 || gt.xscale != -gt.yscale || gt.xrot != 0.0 ||
1251
0
        gt.yrot != 0.0)
1252
0
    {
1253
0
        CPLError(CE_Failure, CPLE_NotSupported,
1254
0
                 "Only north-up geotransform with square pixels supported");
1255
0
        return CE_Failure;
1256
0
    }
1257
0
    m_bGotTransform = true;
1258
0
    m_gt = gt;
1259
0
    InvalidateLabel();
1260
0
    return CE_None;
1261
0
}
1262
1263
/************************************************************************/
1264
/*                         GetRawBinaryLayout()                         */
1265
/************************************************************************/
1266
1267
bool VICARDataset::GetRawBinaryLayout(GDALDataset::RawBinaryLayout &sLayout)
1268
0
{
1269
0
    if (!RawDataset::GetRawBinaryLayout(sLayout))
1270
0
        return false;
1271
0
    sLayout.osRawFilename = GetDescription();
1272
0
    return true;
1273
0
}
1274
1275
/************************************************************************/
1276
/*                       GetMetadataDomainList()                        */
1277
/************************************************************************/
1278
1279
char **VICARDataset::GetMetadataDomainList()
1280
0
{
1281
0
    return BuildMetadataDomainList(nullptr, FALSE, "", "json:VICAR", nullptr);
1282
0
}
1283
1284
/************************************************************************/
1285
/*                            GetMetadata()                             */
1286
/************************************************************************/
1287
1288
CSLConstList VICARDataset::GetMetadata(const char *pszDomain)
1289
38
{
1290
38
    if (pszDomain != nullptr && EQUAL(pszDomain, "json:VICAR"))
1291
0
    {
1292
0
        if (m_aosVICARMD.empty())
1293
0
        {
1294
0
            if (eAccess == GA_Update && !m_oJSonLabel.IsValid())
1295
0
            {
1296
0
                BuildLabel();
1297
0
            }
1298
0
            CPLAssert(m_oJSonLabel.IsValid());
1299
0
            const CPLString osJson =
1300
0
                m_oJSonLabel.Format(CPLJSONObject::PrettyFormat::Pretty);
1301
0
            m_aosVICARMD.InsertString(0, osJson.c_str());
1302
0
        }
1303
0
        return m_aosVICARMD.List();
1304
0
    }
1305
38
    return GDALPamDataset::GetMetadata(pszDomain);
1306
38
}
1307
1308
/************************************************************************/
1309
/*                          InvalidateLabel()                           */
1310
/************************************************************************/
1311
1312
void VICARDataset::InvalidateLabel()
1313
0
{
1314
0
    m_oJSonLabel.Deinit();
1315
0
    m_aosVICARMD.Clear();
1316
0
}
1317
1318
/************************************************************************/
1319
/*                            SetMetadata()                             */
1320
/************************************************************************/
1321
1322
CPLErr VICARDataset::SetMetadata(CSLConstList papszMD, const char *pszDomain)
1323
0
{
1324
0
    if (m_bUseSrcLabel && eAccess == GA_Update && pszDomain != nullptr &&
1325
0
        EQUAL(pszDomain, "json:VICAR"))
1326
0
    {
1327
0
        m_oSrcJSonLabel.Deinit();
1328
0
        InvalidateLabel();
1329
0
        if (papszMD != nullptr && papszMD[0] != nullptr)
1330
0
        {
1331
0
            CPLJSONDocument oJSONDocument;
1332
0
            const GByte *pabyData = reinterpret_cast<const GByte *>(papszMD[0]);
1333
0
            if (!oJSONDocument.LoadMemory(pabyData))
1334
0
            {
1335
0
                return CE_Failure;
1336
0
            }
1337
1338
0
            m_oSrcJSonLabel = oJSONDocument.GetRoot();
1339
0
            if (!m_oSrcJSonLabel.IsValid())
1340
0
            {
1341
0
                return CE_Failure;
1342
0
            }
1343
0
        }
1344
0
        return CE_None;
1345
0
    }
1346
0
    return GDALPamDataset::SetMetadata(papszMD, pszDomain);
1347
0
}
1348
1349
/************************************************************************/
1350
/*                          SerializeString()                           */
1351
/************************************************************************/
1352
1353
static std::string SerializeString(const std::string &s)
1354
0
{
1355
0
    return '\'' + CPLString(s).replaceAll('\'', "''").replaceAll('\n', "\\n") +
1356
0
           '\'';
1357
0
}
1358
1359
/************************************************************************/
1360
/*                        WriteLabelItemValue()                         */
1361
/************************************************************************/
1362
1363
static void WriteLabelItemValue(std::string &osLabel, const CPLJSONObject &obj)
1364
0
{
1365
0
    const auto eType(obj.GetType());
1366
0
    if (eType == CPLJSONObject::Type::Boolean)
1367
0
    {
1368
0
        osLabel += CPLSPrintf("%d", obj.ToBool() ? 1 : 0);
1369
0
    }
1370
0
    else if (eType == CPLJSONObject::Type::Integer)
1371
0
    {
1372
0
        osLabel += CPLSPrintf("%d", obj.ToInteger());
1373
0
    }
1374
0
    else if (eType == CPLJSONObject::Type::Long)
1375
0
    {
1376
0
        std::string osVal(
1377
0
            CPLSPrintf("%.17g", static_cast<double>(obj.ToLong())));
1378
0
        if (osVal.find('.') == std::string::npos)
1379
0
            osVal += ".0";
1380
0
        osLabel += osVal;
1381
0
    }
1382
0
    else if (eType == CPLJSONObject::Type::Double)
1383
0
    {
1384
0
        double dfVal = obj.ToDouble();
1385
0
        if (dfVal >= static_cast<double>(std::numeric_limits<GIntBig>::min()) &&
1386
0
            dfVal <= static_cast<double>(std::numeric_limits<GIntBig>::max()) &&
1387
0
            static_cast<double>(static_cast<GIntBig>(dfVal)) == dfVal)
1388
0
        {
1389
0
            std::string osVal(CPLSPrintf("%.17g", dfVal));
1390
0
            if (osVal.find('.') == std::string::npos)
1391
0
                osVal += ".0";
1392
0
            osLabel += osVal;
1393
0
        }
1394
0
        else
1395
0
        {
1396
0
            osLabel += CPLSPrintf("%.15g", dfVal);
1397
0
        }
1398
0
    }
1399
0
    else if (eType == CPLJSONObject::Type::String)
1400
0
    {
1401
0
        osLabel += SerializeString(obj.ToString());
1402
0
    }
1403
0
    else if (eType == CPLJSONObject::Type::Array)
1404
0
    {
1405
0
        const auto oArray = obj.ToArray();
1406
0
        osLabel += '(';
1407
0
        for (int i = 0; i < oArray.Size(); i++)
1408
0
        {
1409
0
            if (i > 0)
1410
0
                osLabel += ',';
1411
0
            WriteLabelItemValue(osLabel, oArray[i]);
1412
0
        }
1413
0
        osLabel += ')';
1414
0
    }
1415
0
    else if (eType == CPLJSONObject::Type::Null)
1416
0
    {
1417
0
        osLabel += "'NULL'";
1418
0
    }
1419
0
    else
1420
0
    {
1421
0
        osLabel +=
1422
0
            SerializeString(obj.Format(CPLJSONObject::PrettyFormat::Plain));
1423
0
    }
1424
0
}
1425
1426
/************************************************************************/
1427
/*                          SanitizeItemName()                          */
1428
/************************************************************************/
1429
1430
static std::string SanitizeItemName(const std::string &osItemName)
1431
0
{
1432
0
    std::string osRet(osItemName);
1433
0
    if (osRet.size() > 32)
1434
0
        osRet.resize(32);
1435
0
    if (osRet.empty())
1436
0
        return "UNNAMED";
1437
0
    if (osRet[0] < 'A' || osRet[0] > 'Z')
1438
0
        osRet[0] = 'X';  // item name must start with a letter
1439
0
    for (size_t i = 1; i < osRet.size(); i++)
1440
0
    {
1441
0
        char ch = osRet[i];
1442
0
        if (ch >= 'a' && ch <= 'z')
1443
0
            osRet[i] = ch - 'a' + 'A';
1444
0
        else if (!((ch >= 'A' && ch <= 'Z') || (ch >= '0' && ch <= '9') ||
1445
0
                   ch == '_'))
1446
0
        {
1447
0
            osRet[i] = '_';
1448
0
        }
1449
0
    }
1450
0
    if (osRet != osItemName)
1451
0
    {
1452
0
        CPLError(CE_Warning, CPLE_AppDefined,
1453
0
                 "Label item name %s has been sanitized to %s",
1454
0
                 osItemName.c_str(), osRet.c_str());
1455
0
    }
1456
0
    return osRet;
1457
0
}
1458
1459
/************************************************************************/
1460
/*                           WriteLabelItem()                           */
1461
/************************************************************************/
1462
1463
static void WriteLabelItem(std::string &osLabel, const CPLJSONObject &obj,
1464
                           const std::string &osItemName = std::string())
1465
0
{
1466
0
    osLabel += ' ';
1467
0
    osLabel +=
1468
0
        SanitizeItemName(osItemName.empty() ? obj.GetName() : osItemName);
1469
0
    osLabel += '=';
1470
0
    WriteLabelItemValue(osLabel, obj);
1471
0
}
1472
1473
/************************************************************************/
1474
/*                             WriteLabel()                             */
1475
/************************************************************************/
1476
1477
void VICARDataset::WriteLabel()
1478
0
{
1479
0
    m_bIsLabelWritten = true;
1480
1481
0
    if (!m_oJSonLabel.IsValid())
1482
0
        BuildLabel();
1483
1484
0
    std::string osLabel;
1485
0
    auto children = m_oJSonLabel.GetChildren();
1486
0
    for (const auto &child : children)
1487
0
    {
1488
0
        const auto osName(child.GetName());
1489
0
        if (osName == "LBLSIZE" || osName == "PROPERTY" || osName == "TASK")
1490
0
            continue;
1491
0
        std::string osNameSubst;
1492
0
        if (osName == "DAT_TIM" || osName == "USER")
1493
0
        {
1494
0
            osNameSubst = osName + '_';
1495
0
        }
1496
0
        WriteLabelItem(osLabel, child, osNameSubst);
1497
0
    }
1498
1499
0
    auto property = m_oJSonLabel.GetObj("PROPERTY");
1500
0
    if (property.IsValid() && property.GetType() == CPLJSONObject::Type::Object)
1501
0
    {
1502
0
        children = property.GetChildren();
1503
0
        for (const auto &child : children)
1504
0
        {
1505
0
            if (child.GetType() == CPLJSONObject::Type::Object)
1506
0
            {
1507
0
                osLabel += " PROPERTY=" + SerializeString(child.GetName());
1508
0
                auto childrenProperty = child.GetChildren();
1509
0
                for (const auto &childProperty : childrenProperty)
1510
0
                {
1511
0
                    const auto osName(child.GetName());
1512
0
                    std::string osNameSubst;
1513
0
                    if (osName == "LBLSIZE" || osName == "PROPERTY" ||
1514
0
                        osName == "TASK" || osName == "DAT_TIM" ||
1515
0
                        osName == "USER")
1516
0
                    {
1517
0
                        osNameSubst = osName + '_';
1518
0
                    }
1519
0
                    WriteLabelItem(osLabel, childProperty, osNameSubst);
1520
0
                }
1521
0
            }
1522
0
        }
1523
0
    }
1524
1525
0
    auto task = m_oJSonLabel.GetObj("TASK");
1526
0
    if (task.IsValid() && task.GetType() == CPLJSONObject::Type::Object)
1527
0
    {
1528
0
        children = task.GetChildren();
1529
0
        for (const auto &child : children)
1530
0
        {
1531
0
            if (child.GetType() == CPLJSONObject::Type::Object)
1532
0
            {
1533
0
                osLabel += " TASK=" + SerializeString(child.GetName());
1534
0
                auto oUser = child.GetObj("USER");
1535
0
                if (oUser.IsValid())
1536
0
                    WriteLabelItem(osLabel, oUser);
1537
0
                auto oDatTim = child.GetObj("DAT_TIM");
1538
0
                if (oDatTim.IsValid())
1539
0
                    WriteLabelItem(osLabel, oDatTim);
1540
0
                auto childrenProperty = child.GetChildren();
1541
0
                for (const auto &childProperty : childrenProperty)
1542
0
                {
1543
0
                    const auto osName(child.GetName());
1544
0
                    if (osName == "USER" || osName == "DAT_TIM")
1545
0
                        continue;
1546
0
                    std::string osNameSubst;
1547
0
                    if (osName == "LBLSIZE" || osName == "PROPERTY" ||
1548
0
                        osName == "TASK")
1549
0
                    {
1550
0
                        osNameSubst = osName + '_';
1551
0
                    }
1552
0
                    WriteLabelItem(osLabel, childProperty, osNameSubst);
1553
0
                }
1554
0
            }
1555
0
        }
1556
0
    }
1557
1558
    // Figure out label size, round it to the next multiple of RECSIZE
1559
0
    constexpr size_t MAX_LOG10_LBLSIZE = 10;
1560
0
    size_t nLabelSize = strlen("LBLSIZE=") + MAX_LOG10_LBLSIZE + osLabel.size();
1561
0
    nLabelSize = DIV_ROUND_UP(nLabelSize, m_nRecordSize) * m_nRecordSize;
1562
0
    std::string osLabelSize(
1563
0
        CPLSPrintf("LBLSIZE=%d", static_cast<int>(nLabelSize)));
1564
0
    while (osLabelSize.size() < strlen("LBLSIZE=") + MAX_LOG10_LBLSIZE)
1565
0
        osLabelSize += ' ';
1566
0
    osLabel = osLabelSize + osLabel;
1567
0
    CPLAssert(osLabel.size() <= nLabelSize);
1568
1569
    // Write label
1570
0
    VSIFSeekL(fpImage, 0, SEEK_SET);
1571
0
    VSIFWriteL(osLabel.data(), 1, osLabel.size(), fpImage);
1572
0
    const size_t nZeroPadding = nLabelSize - osLabel.size();
1573
0
    if (nZeroPadding)
1574
0
    {
1575
0
        VSIFWriteL(std::string(nZeroPadding, '\0').data(), 1, nZeroPadding,
1576
0
                   fpImage);
1577
0
    }
1578
1579
0
    if (m_bInitToNodata && m_eCompress == COMPRESS_NONE)
1580
0
    {
1581
0
        const int nDTSize =
1582
0
            GDALGetDataTypeSizeBytes(GetRasterBand(1)->GetRasterDataType());
1583
0
        VSIFTruncateL(fpImage, VSIFTellL(fpImage) +
1584
0
                                   static_cast<vsi_l_offset>(nRasterXSize) *
1585
0
                                       nRasterYSize * nBands * nDTSize);
1586
0
    }
1587
1588
    // Patch band offsets to take into account label
1589
0
    for (int i = 0; i < nBands; i++)
1590
0
    {
1591
0
        auto poBand = dynamic_cast<VICARRawRasterBand *>(GetRasterBand(i + 1));
1592
0
        if (poBand)
1593
0
            poBand->nImgOffset += nLabelSize;
1594
0
    }
1595
0
}
1596
1597
/************************************************************************/
1598
/*                             PatchLabel()                             */
1599
/************************************************************************/
1600
1601
void VICARDataset::PatchLabel()
1602
271
{
1603
271
    if (eAccess == GA_ReadOnly || m_eCompress == COMPRESS_NONE)
1604
271
        return;
1605
1606
0
    VSIFSeekL(fpImage, 0, SEEK_END);
1607
0
    const vsi_l_offset nFileSize = VSIFTellL(fpImage);
1608
0
    VSIFSeekL(fpImage, 0, SEEK_SET);
1609
0
    std::string osBuffer;
1610
0
    osBuffer.resize(1024);
1611
0
    size_t nRead = VSIFReadL(&osBuffer[0], 1, 1024, fpImage);
1612
1613
0
    {
1614
0
        CPLString osEOCI1;
1615
0
        osEOCI1.Printf("%u", static_cast<unsigned>(nFileSize));
1616
0
        while (osEOCI1.size() < 10)
1617
0
            osEOCI1 += ' ';
1618
0
        size_t nPos = osBuffer.find("EOCI1=");
1619
0
        CPLAssert(nPos <= nRead - (strlen("EOCI1=") + 10));
1620
0
        memcpy(&osBuffer[nPos + strlen("EOCI1=")], osEOCI1.data(), 10);
1621
0
    }
1622
1623
0
    {
1624
0
        CPLString osEOCI2;
1625
0
        osEOCI2.Printf("%u", static_cast<unsigned>(nFileSize >> 32));
1626
0
        while (osEOCI2.size() < 10)
1627
0
            osEOCI2 += ' ';
1628
0
        size_t nPos = osBuffer.find("EOCI2=");
1629
0
        CPLAssert(nPos <= nRead - (strlen("EOCI2=") + 10));
1630
0
        memcpy(&osBuffer[nPos + strlen("EOCI2=")], osEOCI2.data(), 10);
1631
0
    }
1632
0
    VSIFSeekL(fpImage, 0, SEEK_SET);
1633
0
    VSIFWriteL(&osBuffer[0], 1, nRead, fpImage);
1634
0
}
1635
1636
/************************************************************************/
1637
/*                             BuildLabel()                             */
1638
/************************************************************************/
1639
1640
void VICARDataset::BuildLabel()
1641
0
{
1642
0
    CPLJSONObject oLabel = m_oSrcJSonLabel;
1643
0
    if (!oLabel.IsValid())
1644
0
    {
1645
0
        oLabel = CPLJSONObject();
1646
0
    }
1647
1648
0
    oLabel.Set("LBLSIZE", 0);  // to be overridden later
1649
1650
0
    if (!oLabel.GetObj("TYPE").IsValid())
1651
0
        oLabel.Set("TYPE", "IMAGE");
1652
1653
0
    const auto eType = GetRasterBand(1)->GetRasterDataType();
1654
0
    const char *pszFormat = "";
1655
0
    CPL_IGNORE_RET_VAL(pszFormat);  // Make CSA happy
1656
0
    switch (eType)
1657
0
    {
1658
0
        case GDT_UInt8:
1659
0
            pszFormat = "BYTE";
1660
0
            break;
1661
0
        case GDT_Int16:
1662
0
            pszFormat = "HALF";
1663
0
            break;
1664
0
        case GDT_Int32:
1665
0
            pszFormat = "FULL";
1666
0
            break;
1667
0
        case GDT_Float32:
1668
0
            pszFormat = "REAL";
1669
0
            break;
1670
0
        case GDT_Float64:
1671
0
            pszFormat = "DOUB";
1672
0
            break;
1673
0
        case GDT_CFloat32:
1674
0
            pszFormat = "COMP";
1675
0
            break;
1676
0
        default:
1677
0
            CPLAssert(false);
1678
0
            break;
1679
0
    }
1680
0
    oLabel.Set("FORMAT", pszFormat);
1681
1682
0
    oLabel.Set("BUFSIZ", m_nRecordSize);  // arbitrary value
1683
0
    oLabel.Set("DIM", 3);
1684
0
    oLabel.Set("EOL", 0);
1685
0
    oLabel.Set("RECSIZE", m_nRecordSize);
1686
0
    oLabel.Set("ORG", "BSQ");
1687
0
    oLabel.Set("NL", nRasterYSize);
1688
0
    oLabel.Set("NS", nRasterXSize);
1689
0
    oLabel.Set("NB", nBands);
1690
0
    oLabel.Set("N1", nRasterXSize);
1691
0
    oLabel.Set("N2", nRasterYSize);
1692
0
    oLabel.Set("N3", nBands);
1693
0
    oLabel.Set("N4", 0);
1694
0
    oLabel.Set("NBB", 0);
1695
0
    oLabel.Set("NLB", 0);
1696
0
    oLabel.Set("HOST", "X86-64-LINX");
1697
0
    oLabel.Set("INTFMT", "LOW");
1698
0
    oLabel.Set("REALFMT", "RIEEE");
1699
0
    oLabel.Set("BHOST", "X86-64-LINX");
1700
0
    oLabel.Set("BINTFMT", "LOW");
1701
0
    if (!oLabel.GetObj("BLTYPE").IsValid())
1702
0
        oLabel.Set("BLTYPE", "");
1703
0
    oLabel.Set("COMPRESS", m_eCompress == COMPRESS_BASIC    ? "BASIC"
1704
0
                           : m_eCompress == COMPRESS_BASIC2 ? "BASIC2"
1705
0
                                                            : "NONE");
1706
0
    if (m_eCompress == COMPRESS_NONE)
1707
0
    {
1708
0
        oLabel.Set("EOCI1", 0);
1709
0
        oLabel.Set("EOCI2", 0);
1710
0
    }
1711
0
    else
1712
0
    {
1713
        // To be later patched. Those fake values must take 10 bytes
1714
        // (8 + 2 single quotes) so that they can be later replaced by a
1715
        // integer of maximum value 4294967295 (10 digits)
1716
0
        oLabel.Set("EOCI1", "XXXXXXXX");
1717
0
        oLabel.Set("EOCI2", "XXXXXXXX");
1718
0
    }
1719
1720
0
    if (m_bUseSrcMap)
1721
0
    {
1722
0
        auto oMap = oLabel.GetObj("PROPERTY/MAP");
1723
0
        if (oMap.IsValid() && oMap.GetType() == CPLJSONObject::Type::Object)
1724
0
        {
1725
0
            if (!m_osTargetName.empty())
1726
0
                oMap.Set("TARGET_NAME", m_osTargetName);
1727
0
            if (!m_osLatitudeType.empty())
1728
0
                oMap.Set("COORDINATE_SYSTEM_NAME", m_osLatitudeType);
1729
0
            if (!m_osLongitudeDirection.empty())
1730
0
                oMap.Set("POSITIVE_LONGITUDE_DIRECTION",
1731
0
                         m_osLongitudeDirection);
1732
0
        }
1733
0
    }
1734
0
    else if (m_bGeoRefFormatIsMIPL)
1735
0
    {
1736
0
        auto oProperty = oLabel.GetObj("PROPERTY");
1737
0
        if (oProperty.IsValid())
1738
0
        {
1739
0
            oProperty.Delete("MAP");
1740
0
            oProperty.Delete("GEOTIFF");
1741
0
        }
1742
0
        if (!m_oSRS.IsEmpty())
1743
0
        {
1744
0
            BuildLabelPropertyMap(oLabel);
1745
0
        }
1746
0
    }
1747
0
    else
1748
0
    {
1749
0
        auto oProperty = oLabel.GetObj("PROPERTY");
1750
0
        if (oProperty.IsValid())
1751
0
        {
1752
0
            oProperty.Delete("MAP");
1753
0
            oProperty.Delete("GEOTIFF");
1754
0
        }
1755
0
        if (!m_oSRS.IsEmpty())
1756
0
        {
1757
0
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
1758
0
            BuildLabelPropertyGeoTIFF(oLabel);
1759
0
#endif
1760
0
        }
1761
0
    }
1762
1763
0
    m_oJSonLabel = std::move(oLabel);
1764
0
}
1765
1766
/************************************************************************/
1767
/*                       BuildLabelPropertyMap()                        */
1768
/************************************************************************/
1769
1770
void VICARDataset::BuildLabelPropertyMap(CPLJSONObject &oLabel)
1771
0
{
1772
0
    if (m_oSRS.IsProjected() || m_oSRS.IsGeographic())
1773
0
    {
1774
0
        auto oProperty = GetOrCreateJSONObject(oLabel, "PROPERTY");
1775
0
        auto oMap = GetOrCreateJSONObject(oProperty, "MAP");
1776
1777
0
        const char *pszDatum = m_oSRS.GetAttrValue("DATUM");
1778
0
        CPLString osTargetName(m_osTargetName);
1779
0
        if (osTargetName.empty())
1780
0
        {
1781
0
            if (pszDatum && STARTS_WITH(pszDatum, "D_"))
1782
0
            {
1783
0
                osTargetName = pszDatum + 2;
1784
0
            }
1785
0
            else if (pszDatum)
1786
0
            {
1787
0
                osTargetName = pszDatum;
1788
0
            }
1789
0
        }
1790
0
        if (!osTargetName.empty())
1791
0
            oMap.Add("TARGET_NAME", osTargetName);
1792
1793
0
        oMap.Add("A_AXIS_RADIUS", m_oSRS.GetSemiMajor() / 1000.0);
1794
0
        oMap.Add("B_AXIS_RADIUS", m_oSRS.GetSemiMajor() / 1000.0);
1795
0
        oMap.Add("C_AXIS_RADIUS", m_oSRS.GetSemiMinor() / 1000.0);
1796
1797
0
        if (!m_osLatitudeType.empty())
1798
0
            oMap.Add("COORDINATE_SYSTEM_NAME", m_osLatitudeType);
1799
0
        else
1800
0
            oMap.Add("COORDINATE_SYSTEM_NAME", "PLANETOCENTRIC");
1801
1802
0
        if (!m_osLongitudeDirection.empty())
1803
0
            oMap.Add("POSITIVE_LONGITUDE_DIRECTION", m_osLongitudeDirection);
1804
0
        else
1805
0
            oMap.Add("POSITIVE_LONGITUDE_DIRECTION", "EAST");
1806
1807
0
        const char *pszProjection = m_oSRS.GetAttrValue("PROJECTION");
1808
0
        if (pszProjection == nullptr)
1809
0
        {
1810
0
            oMap.Add("MAP_PROJECTION_TYPE", "SIMPLE_CYLINDRICAL");
1811
0
            oMap.Add("CENTER_LONGITUDE", 0.0);
1812
0
            oMap.Add("CENTER_LATITUDE", 0.0);
1813
0
        }
1814
0
        else if (EQUAL(pszProjection, SRS_PT_EQUIRECTANGULAR))
1815
0
        {
1816
0
            oMap.Add("MAP_PROJECTION_TYPE", "EQUIRECTANGULAR");
1817
0
            if (m_oSRS.GetNormProjParm(SRS_PP_LATITUDE_OF_ORIGIN, 0.0) != 0.0)
1818
0
            {
1819
0
                CPLError(CE_Warning, CPLE_NotSupported,
1820
0
                         "Ignoring %s. Only 0 value supported",
1821
0
                         SRS_PP_LATITUDE_OF_ORIGIN);
1822
0
            }
1823
0
            oMap.Add("CENTER_LONGITUDE",
1824
0
                     m_oSRS.GetNormProjParm(SRS_PP_CENTRAL_MERIDIAN, 0.0));
1825
0
            const double dfCenterLat =
1826
0
                m_oSRS.GetNormProjParm(SRS_PP_STANDARD_PARALLEL_1, 0.0);
1827
0
            oMap.Add("CENTER_LATITUDE", dfCenterLat);
1828
0
        }
1829
0
        else if (EQUAL(pszProjection, SRS_PT_SINUSOIDAL))
1830
0
        {
1831
0
            oMap.Add("MAP_PROJECTION_TYPE", "SINUSOIDAL");
1832
0
            oMap.Add("CENTER_LONGITUDE",
1833
0
                     m_oSRS.GetNormProjParm(SRS_PP_LONGITUDE_OF_CENTER, 0.0));
1834
0
            oMap.Add("CENTER_LATITUDE", 0.0);
1835
0
        }
1836
0
        else
1837
0
        {
1838
0
            CPLError(CE_Warning, CPLE_NotSupported,
1839
0
                     "Projection %s not supported", pszProjection);
1840
0
        }
1841
1842
0
        if (oMap["MAP_PROJECTION_TYPE"].IsValid())
1843
0
        {
1844
0
            if (m_oSRS.GetNormProjParm(SRS_PP_FALSE_EASTING, 0.0) != 0.0)
1845
0
            {
1846
0
                CPLError(CE_Warning, CPLE_NotSupported,
1847
0
                         "Ignoring %s. Only 0 value supported",
1848
0
                         SRS_PP_FALSE_EASTING);
1849
0
            }
1850
0
            if (m_oSRS.GetNormProjParm(SRS_PP_FALSE_NORTHING, 0.0) != 0.0)
1851
0
            {
1852
0
                CPLError(CE_Warning, CPLE_AppDefined,
1853
0
                         "Ignoring %s. Only 0 value supported",
1854
0
                         SRS_PP_FALSE_NORTHING);
1855
0
            }
1856
1857
0
            if (m_bGotTransform)
1858
0
            {
1859
0
                const double dfDegToMeter =
1860
0
                    m_oSRS.GetSemiMajor() * M_PI / 180.0;
1861
0
                if (m_oSRS.IsProjected())
1862
0
                {
1863
0
                    const double dfLinearUnits = m_oSRS.GetLinearUnits();
1864
0
                    const double dfScale = m_gt.xscale * dfLinearUnits;
1865
0
                    oMap.Add("SAMPLE_PROJECTION_OFFSET",
1866
0
                             -m_gt.xorig * dfLinearUnits / dfScale - 0.5);
1867
0
                    oMap.Add("LINE_PROJECTION_OFFSET",
1868
0
                             m_gt.yorig * dfLinearUnits / dfScale - 0.5);
1869
0
                    oMap.Add("MAP_SCALE", dfScale / 1000.0);
1870
0
                }
1871
0
                else if (m_oSRS.IsGeographic())
1872
0
                {
1873
0
                    const double dfScale = m_gt.xscale * dfDegToMeter;
1874
0
                    oMap.Add("SAMPLE_PROJECTION_OFFSET",
1875
0
                             -m_gt.xorig * dfDegToMeter / dfScale - 0.5);
1876
0
                    oMap.Add("LINE_PROJECTION_OFFSET",
1877
0
                             m_gt.yorig * dfDegToMeter / dfScale - 0.5);
1878
0
                    oMap.Add("MAP_SCALE", dfScale / 1000.0);
1879
0
                }
1880
0
            }
1881
0
        }
1882
0
    }
1883
0
    else
1884
0
    {
1885
0
        CPLError(CE_Warning, CPLE_NotSupported, "SRS not supported");
1886
0
    }
1887
0
}
1888
1889
/************************************************************************/
1890
/*                     BuildLabelPropertyGeoTIFF()                      */
1891
/************************************************************************/
1892
1893
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
1894
void VICARDataset::BuildLabelPropertyGeoTIFF(CPLJSONObject &oLabel)
1895
0
{
1896
0
    auto oProperty = GetOrCreateJSONObject(oLabel, "PROPERTY");
1897
0
    auto oGeoTIFF = GetOrCreateJSONObject(oProperty, "GEOTIFF");
1898
1899
    // Ported from Vicar Open Source: Afids expects to be able to read
1900
    // NITF_NROWS and NITF_NCOLS
1901
1902
0
    oGeoTIFF.Add("NITF_NROWS", nRasterYSize);
1903
0
    oGeoTIFF.Add("NITF_NCOLS", nRasterXSize);
1904
1905
    // Create a in-memory GeoTIFF file
1906
1907
0
    const std::string osTmpFilename(
1908
0
        VSIMemGenerateHiddenFilename("vicar_tmp.tif"));
1909
0
    GDALDriver *poGTiffDriver =
1910
0
        GDALDriver::FromHandle(GDALGetDriverByName("GTiff"));
1911
0
    if (poGTiffDriver == nullptr)
1912
0
    {
1913
0
        CPLError(CE_Failure, CPLE_AppDefined, "GTiff driver not available");
1914
0
        return;
1915
0
    }
1916
0
    const char *const apszOptions[] = {"GEOTIFF_VERSION=1.0", nullptr};
1917
0
    auto poDS = std::unique_ptr<GDALDataset>(poGTiffDriver->Create(
1918
0
        osTmpFilename.c_str(), 1, 1, 1, GDT_UInt8, apszOptions));
1919
0
    if (!poDS)
1920
0
        return;
1921
0
    poDS->SetSpatialRef(&m_oSRS);
1922
0
    if (m_bGotTransform)
1923
0
        poDS->SetGeoTransform(m_gt);
1924
0
    poDS->SetMetadataItem(GDALMD_AREA_OR_POINT,
1925
0
                          GetMetadataItem(GDALMD_AREA_OR_POINT));
1926
0
    poDS.reset();
1927
1928
    // Open it with libtiff/libgeotiff
1929
0
    VSILFILE *fpL = VSIFOpenL(osTmpFilename.c_str(), "r");
1930
0
    if (fpL == nullptr)
1931
0
    {
1932
0
        VSIUnlink(osTmpFilename.c_str());
1933
0
        return;
1934
0
    }
1935
1936
0
    TIFF *hTIFF = VSI_TIFFOpen(osTmpFilename.c_str(), "r", fpL);
1937
0
    CPLAssert(hTIFF);
1938
1939
0
    GTIF *hGTIF = GTIFNew(hTIFF);
1940
0
    CPLAssert(hGTIF);
1941
1942
    // Get geotiff keys and write them as VICAR metadata
1943
0
    for (const auto &gkey : GTiffShortKeys)
1944
0
    {
1945
0
        unsigned short val = 0;
1946
0
        if (GDALGTIFKeyGetSHORT(hGTIF, gkey, &val, 0, 1))
1947
0
        {
1948
0
            oGeoTIFF.Add(
1949
0
                CPLString(GTIFKeyName(gkey)).toupper(),
1950
0
                CPLSPrintf("%d(%s)", val, GTIFValueNameEx(hGTIF, gkey, val)));
1951
0
        }
1952
0
    }
1953
1954
0
    for (const auto &gkey : GTiffDoubleKeys)
1955
0
    {
1956
0
        double val = 0;
1957
0
        if (GDALGTIFKeyGetDOUBLE(hGTIF, gkey, &val, 0, 1))
1958
0
        {
1959
0
            oGeoTIFF.Add(CPLString(GTIFKeyName(gkey)).toupper(),
1960
0
                         CPLSPrintf("%.17g", val));
1961
0
        }
1962
0
    }
1963
1964
0
    for (const auto &gkey : GTiffAsciiKeys)
1965
0
    {
1966
0
        char szAscii[1024];
1967
0
        if (GDALGTIFKeyGetASCII(hGTIF, gkey, szAscii,
1968
0
                                static_cast<int>(sizeof(szAscii))))
1969
0
        {
1970
0
            oGeoTIFF.Add(CPLString(GTIFKeyName(gkey)).toupper(), szAscii);
1971
0
        }
1972
0
    }
1973
1974
0
    GTIFFree(hGTIF);
1975
1976
    // Get geotiff tags and write them as VICAR metadata
1977
0
    const std::map<int, const char *> oMapTagCodeToName = {
1978
0
        {TIFFTAG_GEOPIXELSCALE, "MODELPIXELSCALETAG"},
1979
0
        {TIFFTAG_GEOTIEPOINTS, "MODELTIEPOINTTAG"},
1980
0
        {TIFFTAG_GEOTRANSMATRIX, "MODELTRANSFORMATIONTAG"}};
1981
1982
0
    for (const auto &kv : oMapTagCodeToName)
1983
0
    {
1984
0
        uint16_t nCount = 0;
1985
0
        double *padfValues = nullptr;
1986
0
        if (TIFFGetField(hTIFF, kv.first, &nCount, &padfValues))
1987
0
        {
1988
0
            std::string osVal("(");
1989
0
            for (uint16_t i = 0; i < nCount; ++i)
1990
0
            {
1991
0
                if (i > 0)
1992
0
                    osVal += ',';
1993
0
                osVal += CPLSPrintf("%.17g", padfValues[i]);
1994
0
            }
1995
0
            osVal += ')';
1996
0
            oGeoTIFF.Add(kv.second, osVal);
1997
0
        }
1998
0
    }
1999
2000
0
    XTIFFClose(hTIFF);
2001
0
    CPL_IGNORE_RET_VAL(VSIFCloseL(fpL));
2002
0
    VSIUnlink(osTmpFilename.c_str());
2003
0
}
2004
#endif
2005
2006
/************************************************************************/
2007
/*                     ReadProjectionFromMapGroup()                     */
2008
/************************************************************************/
2009
2010
void VICARDataset::ReadProjectionFromMapGroup()
2011
9
{
2012
9
    double dfXDim = 1.0;
2013
9
    double dfYDim = 1.0;
2014
2015
9
    const char *value = GetKeyword("MAP.MAP_SCALE");
2016
9
    if (strlen(value) > 0)
2017
7
    {
2018
7
        dfXDim = CPLAtof(value) * 1000.0;
2019
7
        dfYDim = CPLAtof(value) * -1 * 1000.0;
2020
7
    }
2021
2022
9
    const double dfSampleOffset_Shift =
2023
9
        CPLAtof(CPLGetConfigOption("PDS_SampleProjOffset_Shift", "0.5"));
2024
2025
9
    const double dfLineOffset_Shift =
2026
9
        CPLAtof(CPLGetConfigOption("PDS_LineProjOffset_Shift", "0.5"));
2027
2028
9
    const double dfSampleOffset_Mult =
2029
9
        CPLAtof(CPLGetConfigOption("PDS_SampleProjOffset_Mult", "-1.0"));
2030
2031
9
    const double dfLineOffset_Mult =
2032
9
        CPLAtof(CPLGetConfigOption("PDS_LineProjOffset_Mult", "1.0"));
2033
2034
    /***********   Grab LINE_PROJECTION_OFFSET ************/
2035
9
    double dfULYMap = 0.5;
2036
2037
9
    value = GetKeyword("MAP.LINE_PROJECTION_OFFSET");
2038
9
    if (strlen(value) > 0)
2039
7
    {
2040
7
        const double yulcenter = CPLAtof(value);
2041
7
        dfULYMap =
2042
7
            ((yulcenter + dfLineOffset_Shift) * -dfYDim * dfLineOffset_Mult);
2043
7
    }
2044
    /***********   Grab SAMPLE_PROJECTION_OFFSET ************/
2045
9
    double dfULXMap = 0.5;
2046
2047
9
    value = GetKeyword("MAP.SAMPLE_PROJECTION_OFFSET");
2048
9
    if (strlen(value) > 0)
2049
6
    {
2050
6
        const double xulcenter = CPLAtof(value);
2051
6
        dfULXMap =
2052
6
            ((xulcenter + dfSampleOffset_Shift) * dfXDim * dfSampleOffset_Mult);
2053
6
    }
2054
2055
    /* ==================================================================== */
2056
    /*      Get the coordinate system.                                      */
2057
    /* ==================================================================== */
2058
9
    bool bProjectionSet = true;
2059
2060
    /***********  Grab TARGET_NAME  ************/
2061
    /**** This is the planets name i.e. MARS ***/
2062
9
    CPLString target_name = GetKeyword("MAP.TARGET_NAME");
2063
2064
    /**********   Grab MAP_PROJECTION_TYPE *****/
2065
9
    const CPLString map_proj_name = GetKeyword("MAP.MAP_PROJECTION_TYPE");
2066
2067
    /******  Grab semi_major & convert to KM ******/
2068
9
    const double semi_major = CPLAtof(GetKeyword("MAP.A_AXIS_RADIUS")) * 1000.0;
2069
2070
    /******  Grab semi-minor & convert to KM ******/
2071
9
    const double semi_minor = CPLAtof(GetKeyword("MAP.C_AXIS_RADIUS")) * 1000.0;
2072
2073
    /***********   Grab CENTER_LAT ************/
2074
9
    const double center_lat = CPLAtof(GetKeyword("MAP.CENTER_LATITUDE"));
2075
2076
    /***********   Grab CENTER_LON ************/
2077
9
    const double center_lon = CPLAtof(GetKeyword("MAP.CENTER_LONGITUDE"));
2078
2079
    /**********   Grab 1st std parallel *******/
2080
9
    const double first_std_parallel =
2081
9
        CPLAtof(GetKeyword("MAP.FIRST_STANDARD_PARALLEL"));
2082
2083
    /**********   Grab 2nd std parallel *******/
2084
9
    const double second_std_parallel =
2085
9
        CPLAtof(GetKeyword("MAP.SECOND_STANDARD_PARALLEL"));
2086
2087
    /*** grab  PROJECTION_LATITUDE_TYPE = "PLANETOCENTRIC" ****/
2088
    // Need to further study how ocentric/ographic will effect the gdal library.
2089
    // So far we will use this fact to define a sphere or ellipse for some
2090
    // projections Frank - may need to talk this over
2091
9
    bool bIsGeographic = true;
2092
9
    value = GetKeyword("MAP.COORDINATE_SYSTEM_NAME");
2093
9
    if (EQUAL(value, "PLANETOCENTRIC"))
2094
6
        bIsGeographic = false;
2095
2096
    /**   Set oSRS projection and parameters --- all PDS supported types added
2097
    if apparently supported in oSRS "AITOFF",  ** Not supported in GDAL??
2098
          "ALBERS",
2099
          "BONNE",
2100
          "BRIESEMEISTER",   ** Not supported in GDAL??
2101
          "CYLINDRICAL EQUAL AREA",
2102
          "EQUIDISTANT",
2103
          "EQUIRECTANGULAR",
2104
          "GNOMONIC",
2105
          "HAMMER",    ** Not supported in GDAL??
2106
          "HENDU",     ** Not supported in GDAL??
2107
          "LAMBERT AZIMUTHAL EQUAL AREA",
2108
          "LAMBERT CONFORMAL",
2109
          "MERCATOR",
2110
          "MOLLWEIDE",
2111
          "OBLIQUE CYLINDRICAL",
2112
          "ORTHOGRAPHIC",
2113
          "SIMPLE CYLINDRICAL",
2114
          "SINUSOIDAL",
2115
          "STEREOGRAPHIC",
2116
          "TRANSVERSE MERCATOR",
2117
          "VAN DER GRINTEN",     ** Not supported in GDAL??
2118
          "WERNER"     ** Not supported in GDAL??
2119
    **/
2120
9
    CPLDebug("PDS", "using projection %s\n\n", map_proj_name.c_str());
2121
2122
9
    OGRSpatialReference oSRS;
2123
2124
9
    if ((EQUAL(map_proj_name, "EQUIRECTANGULAR")) ||
2125
9
        (EQUAL(map_proj_name, "SIMPLE_CYLINDRICAL")) ||
2126
9
        (EQUAL(map_proj_name, "EQUIDISTANT")))
2127
0
    {
2128
0
        oSRS.SetEquirectangular2(0.0, center_lon, center_lat, 0, 0);
2129
0
    }
2130
9
    else if (EQUAL(map_proj_name, "ORTHOGRAPHIC"))
2131
0
    {
2132
0
        oSRS.SetOrthographic(center_lat, center_lon, 0, 0);
2133
0
    }
2134
9
    else if (EQUAL(map_proj_name, "SINUSOIDAL"))
2135
8
    {
2136
8
        oSRS.SetSinusoidal(center_lon, 0, 0);
2137
8
    }
2138
1
    else if (EQUAL(map_proj_name, "MERCATOR"))
2139
0
    {
2140
0
        oSRS.SetMercator(center_lat, center_lon, 1, 0, 0);
2141
0
    }
2142
1
    else if (EQUAL(map_proj_name, "STEREOGRAPHIC"))
2143
0
    {
2144
0
        if ((fabs(center_lat) - 90) < 0.0000001)
2145
0
        {
2146
0
            oSRS.SetPS(center_lat, center_lon, 1, 0, 0);
2147
0
        }
2148
0
        else
2149
0
        {
2150
0
            oSRS.SetStereographic(center_lat, center_lon, 1, 0, 0);
2151
0
        }
2152
0
    }
2153
1
    else if (EQUAL(map_proj_name, "POLAR_STEREOGRAPHIC"))
2154
0
    {
2155
0
        oSRS.SetPS(center_lat, center_lon, 1, 0, 0);
2156
0
    }
2157
1
    else if (EQUAL(map_proj_name, "TRANSVERSE_MERCATOR"))
2158
0
    {
2159
0
        oSRS.SetTM(center_lat, center_lon, 1, 0, 0);
2160
0
    }
2161
1
    else if (EQUAL(map_proj_name, "LAMBERT_CONFORMAL_CONIC"))
2162
0
    {
2163
0
        oSRS.SetLCC(first_std_parallel, second_std_parallel, center_lat,
2164
0
                    center_lon, 0, 0);
2165
0
    }
2166
1
    else if (EQUAL(map_proj_name, "LAMBERT_AZIMUTHAL_EQUAL_AREA"))
2167
0
    {
2168
0
        oSRS.SetLAEA(center_lat, center_lon, 0, 0);
2169
0
    }
2170
1
    else if (EQUAL(map_proj_name, "CYLINDRICAL_EQUAL_AREA"))
2171
0
    {
2172
0
        oSRS.SetCEA(first_std_parallel, center_lon, 0, 0);
2173
0
    }
2174
1
    else if (EQUAL(map_proj_name, "MOLLWEIDE"))
2175
0
    {
2176
0
        oSRS.SetMollweide(center_lon, 0, 0);
2177
0
    }
2178
1
    else if (EQUAL(map_proj_name, "ALBERS"))
2179
0
    {
2180
0
        oSRS.SetACEA(first_std_parallel, second_std_parallel, center_lat,
2181
0
                     center_lon, 0, 0);
2182
0
    }
2183
1
    else if (EQUAL(map_proj_name, "BONNE"))
2184
0
    {
2185
0
        oSRS.SetBonne(first_std_parallel, center_lon, 0, 0);
2186
0
    }
2187
1
    else if (EQUAL(map_proj_name, "GNOMONIC"))
2188
0
    {
2189
0
        oSRS.SetGnomonic(center_lat, center_lon, 0, 0);
2190
#ifdef FIXME
2191
    }
2192
    else if (EQUAL(map_proj_name, "OBLIQUE_CYLINDRICAL"))
2193
    {
2194
        // hope Swiss Oblique Cylindrical is the same
2195
        oSRS.SetSOC(center_lat, center_lon, 0, 0);
2196
#endif
2197
0
    }
2198
1
    else
2199
1
    {
2200
1
        CPLDebug("VICAR",
2201
1
                 "Dataset projection %s is not supported. Continuing...",
2202
1
                 map_proj_name.c_str());
2203
1
        bProjectionSet = false;
2204
1
    }
2205
2206
9
    if (bProjectionSet)
2207
8
    {
2208
        // Create projection name, i.e. MERCATOR MARS and set as ProjCS keyword
2209
8
        const CPLString proj_target_name = map_proj_name + " " + target_name;
2210
8
        oSRS.SetProjCS(proj_target_name);  // set ProjCS keyword
2211
2212
        // The geographic/geocentric name will be the same basic name as the
2213
        // body name 'GCS' = Geographic/Geocentric Coordinate System
2214
8
        const CPLString geog_name = "GCS_" + target_name;
2215
2216
        // The datum and sphere names will be the same basic name as the planet
2217
8
        const CPLString datum_name = "D_" + target_name;
2218
2219
8
        CPLString sphere_name = std::move(target_name);
2220
2221
        // calculate inverse flattening from major and minor axis: 1/f = a/(a-b)
2222
8
        double iflattening = 0.0;
2223
8
        if ((semi_major - semi_minor) < 0.0000001)
2224
8
            iflattening = 0;
2225
0
        else
2226
0
            iflattening = semi_major / (semi_major - semi_minor);
2227
2228
        // Set the body size but take into consideration which proj is being
2229
        // used to help w/ compatibility Notice that most PDS projections are
2230
        // spherical based on the fact that ISIS/PICS are spherical Set the body
2231
        // size but take into consideration which proj is being used to help w/
2232
        // proj4 compatibility The use of a Sphere, polar radius or ellipse here
2233
        // is based on how ISIS does it internally
2234
8
        if (((EQUAL(map_proj_name, "STEREOGRAPHIC") &&
2235
0
              (fabs(center_lat) == 90))) ||
2236
8
            (EQUAL(map_proj_name, "POLAR_STEREOGRAPHIC")))
2237
0
        {
2238
0
            if (bIsGeographic)
2239
0
            {
2240
                // Geograpraphic, so set an ellipse
2241
0
                oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major,
2242
0
                               iflattening, "Reference_Meridian", 0.0);
2243
0
            }
2244
0
            else
2245
0
            {
2246
                // Geocentric, so force a sphere using the semi-minor axis. I
2247
                // hope...
2248
0
                sphere_name += "_polarRadius";
2249
0
                oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_minor,
2250
0
                               0.0, "Reference_Meridian", 0.0);
2251
0
            }
2252
0
        }
2253
8
        else if ((EQUAL(map_proj_name, "SIMPLE_CYLINDRICAL")) ||
2254
8
                 (EQUAL(map_proj_name, "EQUIDISTANT")) ||
2255
8
                 (EQUAL(map_proj_name, "ORTHOGRAPHIC")) ||
2256
8
                 (EQUAL(map_proj_name, "STEREOGRAPHIC")) ||
2257
8
                 (EQUAL(map_proj_name, "SINUSOIDAL")))
2258
8
        {
2259
            // isis uses the spherical equation for these projections so force a
2260
            // sphere
2261
8
            oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, 0.0,
2262
8
                           "Reference_Meridian", 0.0);
2263
8
        }
2264
0
        else if (EQUAL(map_proj_name, "EQUIRECTANGULAR"))
2265
0
        {
2266
            // isis uses local radius as a sphere, which is pre-calculated in
2267
            // the PDS label as the semi-major
2268
0
            sphere_name += "_localRadius";
2269
0
            oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major, 0.0,
2270
0
                           "Reference_Meridian", 0.0);
2271
0
        }
2272
0
        else
2273
0
        {
2274
            // All other projections: Mercator, Transverse Mercator, Lambert
2275
            // Conformal, etc. Geographic, so set an ellipse
2276
0
            if (bIsGeographic)
2277
0
            {
2278
0
                oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major,
2279
0
                               iflattening, "Reference_Meridian", 0.0);
2280
0
            }
2281
0
            else
2282
0
            {
2283
                // Geocentric, so force a sphere. I hope...
2284
0
                oSRS.SetGeogCS(geog_name, datum_name, sphere_name, semi_major,
2285
0
                               0.0, "Reference_Meridian", 0.0);
2286
0
            }
2287
0
        }
2288
2289
8
        m_oSRS = std::move(oSRS);
2290
8
        m_oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
2291
8
    }
2292
9
    if (bProjectionSet)
2293
8
    {
2294
8
        m_bGotTransform = true;
2295
8
        m_gt.xorig = dfULXMap;
2296
8
        m_gt.xscale = dfXDim;
2297
8
        m_gt.xrot = 0.0;
2298
8
        m_gt.yorig = dfULYMap;
2299
8
        m_gt.yrot = 0.0;
2300
8
        m_gt.yscale = dfYDim;
2301
8
    }
2302
9
}
2303
2304
/************************************************************************/
2305
/*                   ReadProjectionFromGeoTIFFGroup()                   */
2306
/************************************************************************/
2307
2308
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
2309
void VICARDataset::ReadProjectionFromGeoTIFFGroup()
2310
22
{
2311
22
    m_bGeoRefFormatIsMIPL = true;
2312
2313
    // We will build a in-memory temporary GeoTIFF file from the VICAR GEOTIFF
2314
    // metadata items.
2315
2316
22
    const std::string osTmpFilename(
2317
22
        VSIMemGenerateHiddenFilename("vicar_tmp.tif"));
2318
2319
    /* -------------------------------------------------------------------- */
2320
    /*      Initialization of libtiff and libgeotiff.                       */
2321
    /* -------------------------------------------------------------------- */
2322
22
    GTiffOneTimeInit();
2323
22
    LibgeotiffOneTimeInit();
2324
2325
    /* -------------------------------------------------------------------- */
2326
    /*      Initialize access to the memory geotiff structure.              */
2327
    /* -------------------------------------------------------------------- */
2328
22
    VSILFILE *fpL = VSIFOpenL(osTmpFilename.c_str(), "w");
2329
22
    if (fpL == nullptr)
2330
0
        return;
2331
2332
22
    TIFF *hTIFF = VSI_TIFFOpen(osTmpFilename.c_str(), "w", fpL);
2333
2334
22
    if (hTIFF == nullptr)
2335
0
    {
2336
0
        CPLError(CE_Failure, CPLE_AppDefined,
2337
0
                 "TIFF/GeoTIFF structure is corrupt.");
2338
0
        CPL_IGNORE_RET_VAL(VSIFCloseL(fpL));
2339
0
        return;
2340
0
    }
2341
2342
    /* -------------------------------------------------------------------- */
2343
    /*      Write some minimal set of image parameters.                     */
2344
    /* -------------------------------------------------------------------- */
2345
22
    TIFFSetField(hTIFF, TIFFTAG_IMAGEWIDTH, 1);
2346
22
    TIFFSetField(hTIFF, TIFFTAG_IMAGELENGTH, 1);
2347
22
    TIFFSetField(hTIFF, TIFFTAG_BITSPERSAMPLE, 8);
2348
22
    TIFFSetField(hTIFF, TIFFTAG_SAMPLESPERPIXEL, 1);
2349
22
    TIFFSetField(hTIFF, TIFFTAG_ROWSPERSTRIP, 1);
2350
22
    TIFFSetField(hTIFF, TIFFTAG_PLANARCONFIG, PLANARCONFIG_CONTIG);
2351
22
    TIFFSetField(hTIFF, TIFFTAG_PHOTOMETRIC, PHOTOMETRIC_MINISBLACK);
2352
2353
    /* -------------------------------------------------------------------- */
2354
    /*      Write geotiff keys from VICAR metadata                          */
2355
    /* -------------------------------------------------------------------- */
2356
22
    GTIF *hGTIF = GTIFNew(hTIFF);
2357
22
    CPLAssert(hGTIF);
2358
2359
22
    for (const auto &gkey : GTiffAsciiKeys)
2360
88
    {
2361
88
        const char *pszValue = GetKeyword(
2362
88
            ("GEOTIFF." + CPLString(GTIFKeyName(gkey)).toupper()).c_str(),
2363
88
            nullptr);
2364
88
        if (pszValue)
2365
0
        {
2366
0
            GTIFKeySet(hGTIF, gkey, TYPE_ASCII,
2367
0
                       static_cast<int>(strlen(pszValue)), pszValue);
2368
0
        }
2369
88
    }
2370
2371
22
    for (const auto &gkey : GTiffDoubleKeys)
2372
660
    {
2373
660
        const char *pszValue = GetKeyword(
2374
660
            ("GEOTIFF." + CPLString(GTIFKeyName(gkey)).toupper()).c_str(),
2375
660
            nullptr);
2376
660
        if (pszValue)
2377
0
        {
2378
0
            GTIFKeySet(hGTIF, gkey, TYPE_DOUBLE, 1, CPLAtof(pszValue));
2379
0
        }
2380
660
    }
2381
2382
22
    for (const auto &gkey : GTiffShortKeys)
2383
352
    {
2384
352
        const char *pszValue = GetKeyword(
2385
352
            ("GEOTIFF." + CPLString(GTIFKeyName(gkey)).toupper()).c_str(),
2386
352
            nullptr);
2387
352
        if (pszValue)
2388
0
        {
2389
0
            GTIFKeySet(hGTIF, gkey, TYPE_SHORT, 1, atoi(pszValue));
2390
0
        }
2391
352
    }
2392
2393
22
    GTIFWriteKeys(hGTIF);
2394
22
    GTIFFree(hGTIF);
2395
2396
    /* -------------------------------------------------------------------- */
2397
    /*      Write geotiff tags from VICAR metadata                          */
2398
    /* -------------------------------------------------------------------- */
2399
2400
22
    const std::map<const char *, int> oMapTagNameToCode = {
2401
22
        {"MODELPIXELSCALETAG", TIFFTAG_GEOPIXELSCALE},
2402
22
        {"MODELTIEPOINTTAG", TIFFTAG_GEOTIEPOINTS},
2403
22
        {"MODELTRANSFORMATIONTAG", TIFFTAG_GEOTRANSMATRIX},
2404
22
    };
2405
2406
22
    for (const auto &kv : oMapTagNameToCode)
2407
66
    {
2408
66
        const char *pszValue =
2409
66
            GetKeyword((std::string("GEOTIFF.") + kv.first).c_str(), nullptr);
2410
66
        if (pszValue)
2411
44
        {
2412
            // Remove leading ( and trailing ), and replace comma by space
2413
            // to separate on it.
2414
44
            const CPLStringList aosTokens(
2415
44
                CSLTokenizeString2(CPLString(pszValue)
2416
44
                                       .replaceAll('(', "")
2417
44
                                       .replaceAll(')', "")
2418
44
                                       .replaceAll(',', ' ')
2419
44
                                       .c_str(),
2420
44
                                   " ", 0));
2421
44
            if (!aosTokens.empty())
2422
44
            {
2423
44
                std::vector<double> adfValues;
2424
242
                for (int i = 0; i < aosTokens.size(); ++i)
2425
198
                    adfValues.push_back(CPLAtof(aosTokens[i]));
2426
44
                TIFFSetField(hTIFF, kv.second, aosTokens.size(), &adfValues[0]);
2427
44
            }
2428
44
        }
2429
66
    }
2430
2431
    /* -------------------------------------------------------------------- */
2432
    /*      Finalize the geotiff file.                                      */
2433
    /* -------------------------------------------------------------------- */
2434
2435
22
    char bySmallImage = 0;
2436
2437
22
    TIFFWriteEncodedStrip(hTIFF, 0, &bySmallImage, 1);
2438
22
    TIFFWriteDirectory(hTIFF);
2439
2440
22
    XTIFFClose(hTIFF);
2441
22
    CPL_IGNORE_RET_VAL(VSIFCloseL(fpL));
2442
2443
    /* -------------------------------------------------------------------- */
2444
    /*      Get georeferencing from file.                                   */
2445
    /* -------------------------------------------------------------------- */
2446
22
    auto poGTiffDS =
2447
22
        std::unique_ptr<GDALDataset>(GDALDataset::Open(osTmpFilename.c_str()));
2448
22
    if (poGTiffDS)
2449
22
    {
2450
22
        auto poSRS = poGTiffDS->GetSpatialRef();
2451
22
        if (poSRS)
2452
0
            m_oSRS = *poSRS;
2453
2454
22
        if (poGTiffDS->GetGeoTransform(m_gt) == CE_None)
2455
22
        {
2456
22
            m_bGotTransform = true;
2457
22
        }
2458
2459
22
        const char *pszAreaOrPoint =
2460
22
            poGTiffDS->GetMetadataItem(GDALMD_AREA_OR_POINT);
2461
22
        if (pszAreaOrPoint)
2462
0
            GDALDataset::SetMetadataItem(GDALMD_AREA_OR_POINT, pszAreaOrPoint);
2463
22
    }
2464
2465
22
    VSIUnlink(osTmpFilename.c_str());
2466
22
}
2467
#endif
2468
2469
/************************************************************************/
2470
/*                                Open()                                */
2471
/************************************************************************/
2472
2473
GDALDataset *VICARDataset::Open(GDALOpenInfo *poOpenInfo)
2474
271
{
2475
    /* -------------------------------------------------------------------- */
2476
    /*      Does this look like a VICAR dataset?                            */
2477
    /* -------------------------------------------------------------------- */
2478
271
    const vsi_l_offset nLabelOffset = VICARGetLabelOffset(poOpenInfo);
2479
271
    if (nLabelOffset == static_cast<vsi_l_offset>(-1))
2480
0
        return nullptr;
2481
271
    if (nLabelOffset > 0)
2482
0
    {
2483
0
        CPLString osSubFilename;
2484
0
        osSubFilename.Printf("/vsisubfile/" CPL_FRMT_GUIB ",%s",
2485
0
                             static_cast<GUIntBig>(nLabelOffset),
2486
0
                             poOpenInfo->pszFilename);
2487
0
        GDALOpenInfo oOpenInfo(osSubFilename.c_str(), poOpenInfo->eAccess);
2488
0
        return Open(&oOpenInfo);
2489
0
    }
2490
2491
271
    auto poDS = std::make_unique<VICARDataset>();
2492
271
    poDS->fpImage = poOpenInfo->fpL;
2493
271
    poOpenInfo->fpL = nullptr;
2494
271
    if (!poDS->oKeywords.Ingest(poDS->fpImage, poOpenInfo->pabyHeader))
2495
65
    {
2496
65
        return nullptr;
2497
65
    }
2498
2499
    /************ CHECK INSTRUMENT/DATA *****************/
2500
2501
206
    bool bIsDTM = false;
2502
206
    const char *value = poDS->GetKeyword("DTM.DTM_OFFSET");
2503
206
    if (!EQUAL(value, ""))
2504
0
    {
2505
0
        bIsDTM = true;
2506
0
    }
2507
2508
206
    bool bInstKnown = false;
2509
    // Check for HRSC
2510
206
    if (EQUAL(poDS->GetKeyword("BLTYPE"), "M94_HRSC"))
2511
16
        bInstKnown = true;
2512
    // Check for Framing Camera on Dawn
2513
190
    else if (EQUAL(poDS->GetKeyword("INSTRUMENT_ID"), "FC2"))
2514
0
        bInstKnown = true;
2515
2516
    /************ Grab dimensions *****************/
2517
2518
206
    const int nCols = atoi(poDS->GetKeyword("NS"));
2519
206
    const int nRows = atoi(poDS->GetKeyword("NL"));
2520
206
    const int nBands = atoi(poDS->GetKeyword("NB"));
2521
2522
206
    if (!GDALCheckDatasetDimensions(nCols, nRows) ||
2523
106
        !GDALCheckBandCount(nBands, false))
2524
159
    {
2525
159
        CPLError(CE_Failure, CPLE_AppDefined,
2526
159
                 "File %s appears to be a VICAR file, but failed to find some "
2527
159
                 "required keywords.",
2528
159
                 poOpenInfo->pszFilename);
2529
159
        return nullptr;
2530
159
    }
2531
2532
47
    const GDALDataType eDataType =
2533
47
        GetDataTypeFromFormat(poDS->GetKeyword("FORMAT"));
2534
47
    if (eDataType == GDT_Unknown)
2535
6
    {
2536
6
        CPLError(CE_Failure, CPLE_AppDefined,
2537
6
                 "Could not find known VICAR label entries!");
2538
6
        return nullptr;
2539
6
    }
2540
41
    double dfNoData = 0.0;
2541
41
    if (eDataType == GDT_UInt8)
2542
41
    {
2543
41
        dfNoData = VICAR_NULL1;
2544
41
    }
2545
0
    else if (eDataType == GDT_Int16)
2546
0
    {
2547
0
        dfNoData = VICAR_NULL2;
2548
0
    }
2549
0
    else if (eDataType == GDT_Float32)
2550
0
    {
2551
0
        dfNoData = VICAR_NULL3;
2552
0
    }
2553
2554
    /***** CHECK ENDIANNESS **************/
2555
2556
41
    RawRasterBand::ByteOrder eByteOrder =
2557
41
        RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2558
41
    if (GDALDataTypeIsInteger(eDataType))
2559
41
    {
2560
41
        value = poDS->GetKeyword("INTFMT", "LOW");
2561
41
        if (EQUAL(value, "LOW"))
2562
41
        {
2563
41
            eByteOrder = RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2564
41
        }
2565
0
        else if (EQUAL(value, "HIGH"))
2566
0
        {
2567
0
            eByteOrder = RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN;
2568
0
        }
2569
0
        else
2570
0
        {
2571
0
            CPLError(CE_Failure, CPLE_NotSupported,
2572
0
                     "INTFMT=%s layout not supported.", value);
2573
0
            return nullptr;
2574
0
        }
2575
41
    }
2576
0
    else
2577
0
    {
2578
0
        value = poDS->GetKeyword("REALFMT", "VAX");
2579
0
        if (EQUAL(value, "RIEEE"))
2580
0
        {
2581
0
            eByteOrder = RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2582
0
        }
2583
0
        else if (EQUAL(value, "IEEE"))
2584
0
        {
2585
0
            eByteOrder = RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN;
2586
0
        }
2587
0
        else if (EQUAL(value, "VAX"))
2588
0
        {
2589
0
            eByteOrder = RawRasterBand::ByteOrder::ORDER_VAX;
2590
0
        }
2591
0
        else
2592
0
        {
2593
0
            CPLError(CE_Failure, CPLE_NotSupported,
2594
0
                     "REALFMT=%s layout not supported.", value);
2595
0
            return nullptr;
2596
0
        }
2597
0
    }
2598
2599
    /* -------------------------------------------------------------------- */
2600
    /*      Capture some information from the file that is of interest.     */
2601
    /* -------------------------------------------------------------------- */
2602
41
    poDS->nRasterXSize = nCols;
2603
41
    poDS->nRasterYSize = nRows;
2604
2605
41
    if (poDS->GetKeyword("MAP.MAP_PROJECTION_TYPE")[0] != '\0')
2606
9
    {
2607
9
        poDS->ReadProjectionFromMapGroup();
2608
9
    }
2609
32
#if defined(HAVE_TIFF) && defined(HAVE_GEOTIFF)
2610
32
    else if (poDS->GetKeyword("GEOTIFF.GTMODELTYPEGEOKEY")[0] != '\0' ||
2611
32
             poDS->GetKeyword("GEOTIFF.MODELTIEPOINTTAG")[0] != '\0')
2612
22
    {
2613
22
        poDS->ReadProjectionFromGeoTIFFGroup();
2614
22
    }
2615
41
#endif
2616
2617
41
    if (!poDS->m_bGotTransform)
2618
11
        poDS->m_bGotTransform = CPL_TO_BOOL(GDALReadWorldFile(
2619
11
            poOpenInfo->pszFilename, "wld", poDS->m_gt.data()));
2620
2621
41
    poDS->eAccess = poOpenInfo->eAccess;
2622
41
    poDS->m_oJSonLabel = poDS->oKeywords.GetJsonObject();
2623
2624
    /* -------------------------------------------------------------------- */
2625
    /*      Compute the line offsets.                                        */
2626
    /* -------------------------------------------------------------------- */
2627
2628
41
    uint64_t nPixelOffset;
2629
41
    uint64_t nLineOffset;
2630
41
    uint64_t nBandOffset;
2631
41
    uint64_t nImageOffsetWithoutNBB;
2632
41
    uint64_t nNBB;
2633
41
    uint64_t nImageSize;
2634
41
    if (!GetSpacings(poDS->oKeywords, nPixelOffset, nLineOffset, nBandOffset,
2635
41
                     nImageOffsetWithoutNBB, nNBB, nImageSize) ||
2636
41
        nImageOffsetWithoutNBB > std::numeric_limits<uint64_t>::max() -
2637
41
                                     (nNBB + nBandOffset * (nBands - 1)))
2638
0
    {
2639
0
        CPLDebug("VICAR", "Invalid spacings found");
2640
0
        return nullptr;
2641
0
    }
2642
2643
41
    poDS->m_nRecordSize = atoi(poDS->GetKeyword("RECSIZE", ""));
2644
2645
41
    if (nNBB != 0)
2646
3
    {
2647
3
        const char *pszBLType = poDS->GetKeyword("BLTYPE", nullptr);
2648
#ifdef USE_ONLY_EMBEDDED_RESOURCE_FILES
2649
        const char *pszVicarConf = nullptr;
2650
#else
2651
3
        const char *pszVicarConf = CPLFindFile("gdal", "vicar.json");
2652
3
#endif
2653
3
        CPLJSONDocument oDoc;
2654
3
        if (!pszVicarConf || EQUAL(pszVicarConf, "vicar.json"))
2655
3
        {
2656
3
#ifdef EMBED_RESOURCE_FILES
2657
3
            oDoc.LoadMemory(VICARGetEmbeddedConf());
2658
3
            pszVicarConf = "__embedded__";
2659
3
#endif
2660
3
        }
2661
2662
3
        if (pszBLType && pszVicarConf && poDS->m_nRecordSize > 0)
2663
1
        {
2664
2665
1
            RawRasterBand::ByteOrder eBINTByteOrder =
2666
1
                RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2667
1
            value = poDS->GetKeyword("BINTFMT", "LOW");
2668
1
            if (EQUAL(value, "LOW"))
2669
1
            {
2670
1
                eBINTByteOrder = RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2671
1
            }
2672
0
            else if (EQUAL(value, "HIGH"))
2673
0
            {
2674
0
                eBINTByteOrder = RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN;
2675
0
            }
2676
0
            else
2677
0
            {
2678
0
                CPLError(CE_Failure, CPLE_NotSupported,
2679
0
                         "BINTFMT=%s layout not supported.", value);
2680
0
            }
2681
2682
1
            RawRasterBand::ByteOrder eBREALByteOrder =
2683
1
                RawRasterBand::ByteOrder::ORDER_VAX;
2684
1
            value = poDS->GetKeyword("BREALFMT", "VAX");
2685
1
            if (EQUAL(value, "RIEEE"))
2686
1
            {
2687
1
                eBREALByteOrder = RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN;
2688
1
            }
2689
0
            else if (EQUAL(value, "IEEE"))
2690
0
            {
2691
0
                eBREALByteOrder = RawRasterBand::ByteOrder::ORDER_BIG_ENDIAN;
2692
0
            }
2693
0
            else if (EQUAL(value, "VAX"))
2694
0
            {
2695
0
                eBREALByteOrder = RawRasterBand::ByteOrder::ORDER_VAX;
2696
0
            }
2697
0
            else
2698
0
            {
2699
0
                CPLError(CE_Failure, CPLE_NotSupported,
2700
0
                         "BREALFMT=%s layout not supported.", value);
2701
0
            }
2702
2703
1
            if (EQUAL(pszVicarConf, "__embedded__") || oDoc.Load(pszVicarConf))
2704
1
            {
2705
1
                const auto oRoot = oDoc.GetRoot();
2706
1
                if (oRoot.GetType() == CPLJSONObject::Type::Object)
2707
1
                {
2708
1
                    auto oDef = oRoot.GetObj(pszBLType);
2709
1
                    if (oDef.IsValid() &&
2710
0
                        oDef.GetType() == CPLJSONObject::Type::Object &&
2711
0
                        static_cast<GUInt64>(oDef.GetInteger("size")) == nNBB)
2712
0
                    {
2713
0
                        auto poLayer =
2714
0
                            std::unique_ptr<OGRVICARBinaryPrefixesLayer>(
2715
0
                                new OGRVICARBinaryPrefixesLayer(
2716
0
                                    poDS->fpImage,
2717
0
                                    static_cast<int>(nImageSize /
2718
0
                                                     poDS->m_nRecordSize),
2719
0
                                    oDef, nImageOffsetWithoutNBB,
2720
0
                                    poDS->m_nRecordSize, eBINTByteOrder,
2721
0
                                    eBREALByteOrder));
2722
0
                        if (!poLayer->HasError())
2723
0
                        {
2724
0
                            poDS->m_poLayer = std::move(poLayer);
2725
0
                        }
2726
0
                    }
2727
1
                }
2728
1
            }
2729
1
        }
2730
3
    }
2731
2732
41
    poDS->m_nImageOffsetWithoutNBB =
2733
41
        static_cast<vsi_l_offset>(nImageOffsetWithoutNBB);
2734
2735
41
    CPLString osCompress = poDS->GetKeyword("COMPRESS", "NONE");
2736
41
    if (EQUAL(osCompress, "BASIC") || EQUAL(osCompress, "BASIC2"))
2737
22
    {
2738
22
        if (poOpenInfo->eAccess == GA_Update)
2739
0
        {
2740
0
            CPLError(CE_Failure, CPLE_NotSupported,
2741
0
                     "Update of compressed VICAR file not supported");
2742
0
            return nullptr;
2743
0
        }
2744
22
        poDS->SetMetadataItem("COMPRESS", osCompress, GDAL_MDD_IMAGE_STRUCTURE);
2745
22
        poDS->m_eCompress =
2746
22
            EQUAL(osCompress, "BASIC") ? COMPRESS_BASIC : COMPRESS_BASIC2;
2747
22
        if (poDS->nRasterYSize > 100 * 1000 * 1000 / nBands)
2748
0
        {
2749
0
            CPLError(CE_Failure, CPLE_NotSupported,
2750
0
                     "Too many records for compressed dataset");
2751
0
            return nullptr;
2752
0
        }
2753
22
        if (!GDALDataTypeIsInteger(eDataType))
2754
0
        {
2755
0
            CPLError(CE_Failure, CPLE_NotSupported,
2756
0
                     "Data type incompatible with compression");
2757
0
            return nullptr;
2758
0
        }
2759
        // To avoid potential issues in basic_decode()
2760
22
        const int nDTSize = GDALGetDataTypeSizeBytes(eDataType);
2761
22
        if (nDTSize == 0 || poDS->nRasterXSize > INT_MAX / nDTSize)
2762
0
        {
2763
0
            CPLError(CE_Failure, CPLE_NotSupported, "Too large scanline");
2764
0
            return nullptr;
2765
0
        }
2766
22
        const int nRecords = poDS->nRasterYSize * nBands;
2767
22
        try
2768
22
        {
2769
            // + 1 to store implicitly the size of the last record
2770
22
            poDS->m_anRecordOffsets.resize(nRecords + 1);
2771
22
        }
2772
22
        catch (const std::exception &e)
2773
22
        {
2774
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
2775
0
            return nullptr;
2776
0
        }
2777
22
        if (poDS->m_eCompress == COMPRESS_BASIC)
2778
0
        {
2779
0
            poDS->m_anRecordOffsets[0] =
2780
0
                poDS->m_nImageOffsetWithoutNBB + sizeof(GUInt32);
2781
0
        }
2782
22
        else
2783
22
        {
2784
22
            poDS->m_anRecordOffsets[0] =
2785
22
                poDS->m_nImageOffsetWithoutNBB + sizeof(GUInt32) * nRecords;
2786
22
        }
2787
22
    }
2788
19
    else if (!EQUAL(osCompress, "NONE"))
2789
0
    {
2790
0
        CPLError(CE_Failure, CPLE_NotSupported, "COMPRESS=%s not supported",
2791
0
                 osCompress.c_str());
2792
0
        return nullptr;
2793
0
    }
2794
2795
    /* -------------------------------------------------------------------- */
2796
    /*      Create band information objects.                                */
2797
    /* -------------------------------------------------------------------- */
2798
82
    for (int i = 0; i < nBands; i++)
2799
41
    {
2800
41
        std::unique_ptr<GDALRasterBand> poBand;
2801
2802
41
        if (poDS->m_eCompress == COMPRESS_BASIC ||
2803
41
            poDS->m_eCompress == COMPRESS_BASIC2)
2804
22
        {
2805
22
            poBand = std::make_unique<VICARBASICRasterBand>(poDS.get(), i + 1,
2806
22
                                                            eDataType);
2807
22
        }
2808
19
        else
2809
19
        {
2810
19
            auto poRawBand = std::make_unique<VICARRawRasterBand>(
2811
19
                poDS.get(), i + 1, poDS->fpImage,
2812
19
                static_cast<vsi_l_offset>(nImageOffsetWithoutNBB + nNBB +
2813
19
                                          nBandOffset * i),
2814
19
                static_cast<int>(nPixelOffset), static_cast<int>(nLineOffset),
2815
19
                eDataType, eByteOrder);
2816
19
            if (!poRawBand->IsValid())
2817
0
            {
2818
0
                return nullptr;
2819
0
            }
2820
19
            poBand = std::move(poRawBand);
2821
19
        }
2822
2823
        // only set NoData if instrument is supported
2824
41
        if (bInstKnown)
2825
15
            poBand->SetNoDataValue(dfNoData);
2826
41
        if (bIsDTM)
2827
0
        {
2828
0
            poBand->SetScale(static_cast<double>(
2829
0
                CPLAtof(poDS->GetKeyword("DTM.DTM_SCALING_FACTOR"))));
2830
0
            poBand->SetOffset(static_cast<double>(
2831
0
                CPLAtof(poDS->GetKeyword("DTM.DTM_OFFSET"))));
2832
0
            const char *pszMin =
2833
0
                poDS->GetKeyword("DTM.DTM_MINIMUM_DN", nullptr);
2834
0
            const char *pszMax =
2835
0
                poDS->GetKeyword("DTM.DTM_MAXIMUM_DN", nullptr);
2836
0
            if (pszMin != nullptr && pszMax != nullptr)
2837
0
                poBand->SetStatistics(CPLAtofM(pszMin), CPLAtofM(pszMax), 0, 0);
2838
0
            const char *pszNoData =
2839
0
                poDS->GetKeyword("DTM.DTM_MISSING_DN", nullptr);
2840
0
            if (pszNoData != nullptr)
2841
0
                poBand->SetNoDataValue(CPLAtofM(pszNoData));
2842
0
        }
2843
41
        else if (EQUAL(poDS->GetKeyword("BLTYPE"), "M94_HRSC"))
2844
15
        {
2845
15
            double scale = CPLAtof(
2846
15
                poDS->GetKeyword("DLRTO8.REFLECTANCE_SCALING_FACTOR", "-1."));
2847
15
            if (scale < 0.)
2848
9
            {
2849
9
                scale = CPLAtof(
2850
9
                    poDS->GetKeyword("HRCAL.REFLECTANCE_SCALING_FACTOR", "1."));
2851
9
            }
2852
15
            poBand->SetScale(scale);
2853
15
            double offset =
2854
15
                CPLAtof(poDS->GetKeyword("DLRTO8.REFLECTANCE_OFFSET", "-1."));
2855
15
            if (offset < 0.)
2856
9
            {
2857
9
                offset =
2858
9
                    CPLAtof(poDS->GetKeyword("HRCAL.REFLECTANCE_OFFSET", "0."));
2859
9
            }
2860
15
            poBand->SetOffset(offset);
2861
15
        }
2862
41
        const char *pszMin = poDS->GetKeyword("STATISTICS.MINIMUM", nullptr);
2863
41
        const char *pszMax = poDS->GetKeyword("STATISTICS.MAXIMUM", nullptr);
2864
41
        const char *pszMean = poDS->GetKeyword("STATISTICS.MEAN", nullptr);
2865
41
        const char *pszStdDev =
2866
41
            poDS->GetKeyword("STATISTICS.STANDARD_DEVIATION", nullptr);
2867
41
        if (pszMin != nullptr && pszMax != nullptr && pszMean != nullptr &&
2868
0
            pszStdDev != nullptr)
2869
0
            poBand->SetStatistics(CPLAtofM(pszMin), CPLAtofM(pszMax),
2870
0
                                  CPLAtofM(pszMean), CPLAtofM(pszStdDev));
2871
2872
41
        poDS->SetBand(i + 1, std::move(poBand));
2873
41
    }
2874
2875
    /* -------------------------------------------------------------------- */
2876
    /*      Instrument-specific keywords as metadata.                       */
2877
    /* -------------------------------------------------------------------- */
2878
2879
    /******************   HRSC    ******************************/
2880
2881
41
    if (EQUAL(poDS->GetKeyword("BLTYPE"), "M94_HRSC"))
2882
15
    {
2883
15
        poDS->SetMetadataItem(
2884
15
            "SPACECRAFT_NAME",
2885
15
            poDS->GetKeyword("M94_INSTRUMENT.INSTRUMENT_HOST_NAME"));
2886
15
        poDS->SetMetadataItem("PRODUCT_TYPE", poDS->GetKeyword("TYPE"));
2887
2888
15
        if (EQUAL(poDS->GetKeyword("M94_INSTRUMENT.DETECTOR_ID"),
2889
15
                  "MEX_HRSC_SRC"))
2890
0
        {
2891
0
            static const char *const apszKeywords[] = {
2892
0
                "M94_ORBIT.IMAGE_TIME",
2893
0
                "FILE.EVENT_TYPE",
2894
0
                "FILE.PROCESSING_LEVEL_ID",
2895
0
                "M94_INSTRUMENT.DETECTOR_ID",
2896
0
                "M94_CAMERAS.EXPOSURE_DURATION",
2897
0
                "HRCONVER.INSTRUMENT_TEMPERATURE",
2898
0
                nullptr};
2899
0
            for (int i = 0; apszKeywords[i] != nullptr; i++)
2900
0
            {
2901
0
                const char *pszKeywordValue = poDS->GetKeyword(apszKeywords[i]);
2902
0
                if (pszKeywordValue != nullptr)
2903
0
                    poDS->SetMetadataItem(apszKeywords[i], pszKeywordValue);
2904
0
            }
2905
0
        }
2906
15
        else
2907
15
        {
2908
15
            static const char *const apszKeywords[] = {
2909
15
                "M94_ORBIT.START_TIME",
2910
15
                "M94_ORBIT.STOP_TIME",
2911
15
                "M94_INSTRUMENT.DETECTOR_ID",
2912
15
                "M94_CAMERAS.MACROPIXEL_SIZE",
2913
15
                "FILE.EVENT_TYPE",
2914
15
                "M94_INSTRUMENT.MISSION_PHASE_NAME",
2915
15
                "HRORTHO.SPICE_FILE_NAME",
2916
15
                "HRCONVER.MISSING_FRAMES",
2917
15
                "HRCONVER.OVERFLOW_FRAMES",
2918
15
                "HRCONVER.ERROR_FRAMES",
2919
15
                "HRFOOT.BEST_GROUND_SAMPLING_DISTANCE",
2920
15
                "DLRTO8.RADIANCE_SCALING_FACTOR",
2921
15
                "DLRTO8.RADIANCE_OFFSET",
2922
15
                "DLRTO8.REFLECTANCE_SCALING_FACTOR",
2923
15
                "DLRTO8.REFLECTANCE_OFFSET",
2924
15
                "HRCAL.RADIANCE_SCALING_FACTOR",
2925
15
                "HRCAL.RADIANCE_OFFSET",
2926
15
                "HRCAL.REFLECTANCE_SCALING_FACTOR",
2927
15
                "HRCAL.REFLECTANCE_OFFSET",
2928
15
                "HRORTHO.DTM_NAME",
2929
15
                "HRORTHO.EXTORI_FILE_NAME",
2930
15
                "HRORTHO.GEOMETRIC_CALIB_FILE_NAME",
2931
15
                nullptr};
2932
345
            for (int i = 0; apszKeywords[i] != nullptr; i++)
2933
330
            {
2934
330
                const char *pszKeywordValue =
2935
330
                    poDS->GetKeyword(apszKeywords[i], nullptr);
2936
330
                if (pszKeywordValue != nullptr)
2937
120
                    poDS->SetMetadataItem(apszKeywords[i], pszKeywordValue);
2938
330
            }
2939
15
        }
2940
15
    }
2941
41
    if (bIsDTM && EQUAL(poDS->GetKeyword("MAP.TARGET_NAME"), "MARS"))
2942
0
    {
2943
0
        poDS->SetMetadataItem("SPACECRAFT_NAME", "MARS_EXPRESS");
2944
0
        poDS->SetMetadataItem("PRODUCT_TYPE", "DTM");
2945
0
        static const char *const apszKeywords[] = {
2946
0
            "DTM.DTM_MISSING_DN",     "DTM.DTM_OFFSET",
2947
0
            "DTM.DTM_SCALING_FACTOR", "DTM.DTM_A_AXIS_RADIUS",
2948
0
            "DTM.DTM_B_AXIS_RADIUS",  "DTM.DTM_C_AXIS_RADIUS",
2949
0
            "DTM.DTM_DESC",           "DTM.DTM_MINIMUM_DN",
2950
0
            "DTM.DTM_MAXIMUM_DN",     nullptr};
2951
0
        for (int i = 0; apszKeywords[i] != nullptr; i++)
2952
0
        {
2953
0
            const char *pszKeywordValue = poDS->GetKeyword(apszKeywords[i]);
2954
0
            if (pszKeywordValue != nullptr)
2955
0
                poDS->SetMetadataItem(apszKeywords[i], pszKeywordValue);
2956
0
        }
2957
0
    }
2958
2959
    /******************   DAWN   ******************************/
2960
41
    else if (EQUAL(poDS->GetKeyword("INSTRUMENT_ID"), "FC2"))
2961
0
    {
2962
0
        poDS->SetMetadataItem("SPACECRAFT_NAME", "DAWN");
2963
0
        static const char *const apszKeywords[] = {
2964
0
            "ORBIT_NUMBER",
2965
0
            "FILTER_NUMBER",
2966
0
            "FRONT_DOOR_STATUS",
2967
0
            "FIRST_LINE",
2968
0
            "FIRST_LINE_SAMPLE",
2969
0
            "PRODUCER_INSTITUTION_NAME",
2970
0
            "SOURCE_FILE_NAME",
2971
0
            "PROCESSING_LEVEL_ID",
2972
0
            "TARGET_NAME",
2973
0
            "LIMB_IN_IMAGE",
2974
0
            "POLE_IN_IMAGE",
2975
0
            "REFLECTANCE_SCALING_FACTOR",
2976
0
            "SPICE_FILE_NAME",
2977
0
            "SPACECRAFT_CENTRIC_LATITUDE",
2978
0
            "SPACECRAFT_EASTERN_LONGITUDE",
2979
0
            "FOOTPRINT_POSITIVE_LONGITUDE",
2980
0
            nullptr};
2981
0
        for (int i = 0; apszKeywords[i] != nullptr; i++)
2982
0
        {
2983
0
            const char *pszKeywordValue = poDS->GetKeyword(apszKeywords[i]);
2984
0
            if (pszKeywordValue != nullptr)
2985
0
                poDS->SetMetadataItem(apszKeywords[i], pszKeywordValue);
2986
0
        }
2987
0
    }
2988
41
    else if (bIsDTM && (EQUAL(poDS->GetKeyword("TARGET_NAME"), "VESTA") ||
2989
0
                        EQUAL(poDS->GetKeyword("TARGET_NAME"), "CERES")))
2990
0
    {
2991
0
        poDS->SetMetadataItem("SPACECRAFT_NAME", "DAWN");
2992
0
        poDS->SetMetadataItem("PRODUCT_TYPE", "DTM");
2993
0
        static const char *const apszKeywords[] = {
2994
0
            "DTM_MISSING_DN",
2995
0
            "DTM_OFFSET",
2996
0
            "DTM_SCALING_FACTOR",
2997
0
            "DTM_A_AXIS_RADIUS",
2998
0
            "DTM_B_AXIS_RADIUS",
2999
0
            "DTM_C_AXIS_RADIUS",
3000
0
            "DTM_MINIMUM_DN",
3001
0
            "DTM_MAXIMUM_DN",
3002
0
            "MAP_PROJECTION_TYPE",
3003
0
            "COORDINATE_SYSTEM_NAME",
3004
0
            "POSITIVE_LONGITUDE_DIRECTION",
3005
0
            "MAP_SCALE",
3006
0
            "CENTER_LONGITUDE",
3007
0
            "LINE_PROJECTION_OFFSET",
3008
0
            "SAMPLE_PROJECTION_OFFSET",
3009
0
            nullptr};
3010
0
        for (int i = 0; apszKeywords[i] != nullptr; i++)
3011
0
        {
3012
0
            const char *pszKeywordValue = poDS->GetKeyword(apszKeywords[i]);
3013
0
            if (pszKeywordValue != nullptr)
3014
0
                poDS->SetMetadataItem(apszKeywords[i], pszKeywordValue);
3015
0
        }
3016
0
    }
3017
3018
    /* -------------------------------------------------------------------- */
3019
    /*      Initialize any PAM information.                                 */
3020
    /* -------------------------------------------------------------------- */
3021
41
    poDS->TryLoadXML();
3022
3023
    /* -------------------------------------------------------------------- */
3024
    /*      Check for overviews.                                            */
3025
    /* -------------------------------------------------------------------- */
3026
41
    poDS->oOvManager.Initialize(poDS.get(), poOpenInfo->pszFilename);
3027
3028
41
    return poDS.release();
3029
41
}
3030
3031
/************************************************************************/
3032
/*                             GetKeyword()                             */
3033
/************************************************************************/
3034
3035
const char *VICARDataset::GetKeyword(const char *pszPath,
3036
                                     const char *pszDefault)
3037
3038
3.48k
{
3039
3.48k
    return oKeywords.GetKeyword(pszPath, pszDefault);
3040
3.48k
}
3041
3042
/************************************************************************/
3043
/*                       GetDataTypeFromFormat()                        */
3044
/************************************************************************/
3045
3046
GDALDataType VICARDataset::GetDataTypeFromFormat(const char *pszFormat)
3047
88
{
3048
88
    if (EQUAL(pszFormat, "BYTE"))
3049
82
        return GDT_UInt8;
3050
3051
6
    if (EQUAL(pszFormat, "HALF") || EQUAL(pszFormat, "WORD"))
3052
0
        return GDT_Int16;
3053
3054
6
    if (EQUAL(pszFormat, "FULL") || EQUAL(pszFormat, "LONG"))
3055
0
        return GDT_Int32;
3056
3057
6
    if (EQUAL(pszFormat, "REAL"))
3058
0
        return GDT_Float32;
3059
3060
6
    if (EQUAL(pszFormat, "DOUB"))
3061
0
        return GDT_Float64;
3062
3063
6
    if (EQUAL(pszFormat, "COMP") || EQUAL(pszFormat, "COMPLEX"))
3064
0
        return GDT_CFloat32;
3065
3066
6
    return GDT_Unknown;
3067
6
}
3068
3069
/************************************************************************/
3070
/*                            GetSpacings()                             */
3071
/************************************************************************/
3072
3073
bool VICARDataset::GetSpacings(const VICARKeywordHandler &keywords,
3074
                               uint64_t &nPixelOffset, uint64_t &nLineOffset,
3075
                               uint64_t &nBandOffset,
3076
                               uint64_t &nImageOffsetWithoutNBB, uint64_t &nNBB,
3077
                               uint64_t &nImageSize)
3078
41
{
3079
41
    const GDALDataType eDataType =
3080
41
        GetDataTypeFromFormat(keywords.GetKeyword("FORMAT", ""));
3081
41
    if (eDataType == GDT_Unknown)
3082
0
        return false;
3083
41
    const uint64_t nItemSize = GDALGetDataTypeSizeBytes(eDataType);
3084
41
    const char *value = keywords.GetKeyword("ORG", "BSQ");
3085
    // number of bytes of binary prefix before each record
3086
41
    nNBB = atoi(keywords.GetKeyword("NBB", ""));
3087
41
    const uint64_t nCols64 = atoi(keywords.GetKeyword("NS", ""));
3088
41
    const uint64_t nRows64 = atoi(keywords.GetKeyword("NL", ""));
3089
41
    const uint64_t nBands64 = atoi(keywords.GetKeyword("NB", ""));
3090
41
    try
3091
41
    {
3092
41
        if (EQUAL(value, "BIP"))
3093
0
        {
3094
0
            nPixelOffset = (CPLSM(nItemSize) * CPLSM(nBands64)).v();
3095
0
            nBandOffset = nItemSize;
3096
0
            nLineOffset =
3097
0
                (CPLSM(nNBB) + CPLSM(nPixelOffset) * CPLSM(nCols64)).v();
3098
0
            nImageSize = (CPLSM(nLineOffset) * CPLSM(nRows64)).v();
3099
0
        }
3100
41
        else if (EQUAL(value, "BIL"))
3101
0
        {
3102
0
            nPixelOffset = nItemSize;
3103
0
            nBandOffset = (CPLSM(nItemSize) * CPLSM(nCols64)).v();
3104
0
            nLineOffset =
3105
0
                (CPLSM(nNBB) + CPLSM(nBandOffset) * CPLSM(nBands64)).v();
3106
0
            nImageSize = (CPLSM(nLineOffset) * CPLSM(nRows64)).v();
3107
0
        }
3108
41
        else if (EQUAL(value, "BSQ"))
3109
41
        {
3110
41
            nPixelOffset = nItemSize;
3111
41
            nLineOffset =
3112
41
                (CPLSM(nNBB) + CPLSM(nPixelOffset) * CPLSM(nCols64)).v();
3113
41
            nBandOffset = (CPLSM(nLineOffset) * CPLSM(nRows64)).v();
3114
41
            nImageSize = (CPLSM(nBandOffset) * CPLSM(nBands64)).v();
3115
41
        }
3116
0
        else
3117
0
        {
3118
0
            CPLError(CE_Failure, CPLE_NotSupported,
3119
0
                     "ORG=%s layout not supported.", value);
3120
0
            return false;
3121
0
        }
3122
41
    }
3123
41
    catch (const CPLSafeIntOverflow &)
3124
41
    {
3125
0
        return false;
3126
0
    }
3127
3128
41
    const uint64_t nLabelSize = atoi(keywords.GetKeyword("LBLSIZE", ""));
3129
41
    const uint64_t nRecordSize = atoi(keywords.GetKeyword("RECSIZE", ""));
3130
41
    const uint64_t nNLB = atoi(keywords.GetKeyword("NLB", ""));
3131
41
    try
3132
41
    {
3133
41
        nImageOffsetWithoutNBB =
3134
41
            (CPLSM(nLabelSize) + CPLSM(nRecordSize) * CPLSM(nNLB) + CPLSM(nNBB))
3135
41
                .v();
3136
41
        nImageOffsetWithoutNBB -= nNBB;
3137
41
    }
3138
41
    catch (const CPLSafeIntOverflow &)
3139
41
    {
3140
0
        return false;
3141
0
    }
3142
41
    return true;
3143
41
}
3144
3145
/************************************************************************/
3146
/*                               Create()                               */
3147
/************************************************************************/
3148
3149
GDALDataset *VICARDataset::Create(const char *pszFilename, int nXSize,
3150
                                  int nYSize, int nBandsIn, GDALDataType eType,
3151
                                  CSLConstList papszOptions)
3152
0
{
3153
0
    return CreateInternal(pszFilename, nXSize, nYSize, nBandsIn, eType,
3154
0
                          papszOptions);
3155
0
}
3156
3157
VICARDataset *VICARDataset::CreateInternal(const char *pszFilename, int nXSize,
3158
                                           int nYSize, int nBandsIn,
3159
                                           GDALDataType eType,
3160
                                           CSLConstList papszOptions)
3161
0
{
3162
0
    if (eType != GDT_UInt8 && eType != GDT_Int16 && eType != GDT_Int32 &&
3163
0
        eType != GDT_Float32 && eType != GDT_Float64 && eType != GDT_CFloat32)
3164
0
    {
3165
0
        CPLError(CE_Failure, CPLE_NotSupported, "Unsupported data type");
3166
0
        return nullptr;
3167
0
    }
3168
3169
0
    const int nPixelOffset = GDALGetDataTypeSizeBytes(eType);
3170
0
    if (nXSize == 0 || nYSize == 0 || nPixelOffset > INT_MAX / nXSize)
3171
0
    {
3172
0
        CPLError(CE_Failure, CPLE_NotSupported,
3173
0
                 "Unsupported raster dimensions");
3174
0
        return nullptr;
3175
0
    }
3176
0
    const int nLineOffset = nXSize * nPixelOffset;
3177
3178
0
    if (nBandsIn == 0 || nBandsIn > 32767)
3179
0
    {
3180
0
        CPLError(CE_Failure, CPLE_NotSupported, "Unsupported band count");
3181
0
        return nullptr;
3182
0
    }
3183
3184
0
    const char *pszCompress =
3185
0
        CSLFetchNameValueDef(papszOptions, "COMPRESS", "NONE");
3186
0
    CompressMethod eCompress = COMPRESS_NONE;
3187
0
    if (EQUAL(pszCompress, "NONE"))
3188
0
    {
3189
0
        eCompress = COMPRESS_NONE;
3190
0
    }
3191
0
    else if (EQUAL(pszCompress, "BASIC"))
3192
0
    {
3193
0
        eCompress = COMPRESS_BASIC;
3194
0
    }
3195
0
    else if (EQUAL(pszCompress, "BASIC2"))
3196
0
    {
3197
0
        eCompress = COMPRESS_BASIC2;
3198
0
    }
3199
0
    else
3200
0
    {
3201
0
        CPLError(CE_Failure, CPLE_NotSupported, "Unsupported COMPRESS value");
3202
0
        return nullptr;
3203
0
    }
3204
0
    if (eCompress != COMPRESS_NONE &&
3205
0
        (!GDALDataTypeIsInteger(eType) || nBandsIn != 1))
3206
0
    {
3207
0
        CPLError(
3208
0
            CE_Failure, CPLE_NotSupported,
3209
0
            "BASIC/BASIC2 compression only supports one-band integer datasets");
3210
0
        return nullptr;
3211
0
    }
3212
3213
0
    std::vector<vsi_l_offset> anRecordOffsets;
3214
0
    if (eCompress != COMPRESS_NONE)
3215
0
    {
3216
0
        const GUInt64 nMaxEncodedSize =
3217
0
            static_cast<GUInt64>(nXSize) * nPixelOffset +
3218
0
            static_cast<GUInt64>(nXSize) * nPixelOffset / 2 + 11;
3219
        // To avoid potential later int overflows
3220
0
        if (nMaxEncodedSize > static_cast<GUInt64>(INT_MAX))
3221
0
        {
3222
0
            CPLError(CE_Failure, CPLE_NotSupported, "Too large scanline");
3223
0
            return nullptr;
3224
0
        }
3225
0
        if (nYSize > 100 * 1000 * 1000)
3226
0
        {
3227
0
            CPLError(CE_Failure, CPLE_NotSupported,
3228
0
                     "Too many records for compressed dataset");
3229
0
            return nullptr;
3230
0
        }
3231
0
        try
3232
0
        {
3233
            // + 1 to store implicitly the size of the last record
3234
0
            anRecordOffsets.resize(nYSize + 1);
3235
0
        }
3236
0
        catch (const std::exception &e)
3237
0
        {
3238
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "%s", e.what());
3239
0
            return nullptr;
3240
0
        }
3241
0
    }
3242
3243
0
    CPLJSONObject oSrcJSonLabel;
3244
0
    oSrcJSonLabel.Deinit();
3245
3246
0
    const char *pszLabel = CSLFetchNameValue(papszOptions, "LABEL");
3247
0
    if (pszLabel)
3248
0
    {
3249
0
        CPLJSONDocument oJSONDocument;
3250
0
        if (pszLabel[0] == '{')
3251
0
        {
3252
0
            const GByte *pabyData = reinterpret_cast<const GByte *>(pszLabel);
3253
0
            if (!oJSONDocument.LoadMemory(pabyData))
3254
0
            {
3255
0
                return nullptr;
3256
0
            }
3257
0
        }
3258
0
        else
3259
0
        {
3260
0
            if (!oJSONDocument.Load(pszLabel))
3261
0
            {
3262
0
                return nullptr;
3263
0
            }
3264
0
        }
3265
3266
0
        oSrcJSonLabel = oJSONDocument.GetRoot();
3267
0
        if (!oSrcJSonLabel.IsValid())
3268
0
        {
3269
0
            return nullptr;
3270
0
        }
3271
0
    }
3272
3273
0
    VSILFILE *fp = VSIFOpenExL(pszFilename, "wb+", true);
3274
0
    if (fp == nullptr)
3275
0
    {
3276
0
        CPLError(CE_Failure, CPLE_FileIO, "Cannot create %s: %s", pszFilename,
3277
0
                 VSIGetLastErrorMsg());
3278
0
        return nullptr;
3279
0
    }
3280
3281
0
    VICARDataset *poDS = new VICARDataset();
3282
0
    poDS->fpImage = fp;
3283
0
    poDS->nRasterXSize = nXSize;
3284
0
    poDS->nRasterYSize = nYSize;
3285
0
    poDS->m_nRecordSize = nLineOffset;
3286
0
    poDS->m_bIsLabelWritten = false;
3287
0
    poDS->m_bGeoRefFormatIsMIPL = EQUAL(
3288
0
        CSLFetchNameValueDef(papszOptions, "GEOREF_FORMAT", "MIPL"), "MIPL");
3289
0
    poDS->m_bUseSrcLabel = CPLFetchBool(papszOptions, "USE_SRC_LABEL", true);
3290
0
    poDS->m_bUseSrcMap = CPLFetchBool(papszOptions, "USE_SRC_MAP", false);
3291
0
    poDS->m_osLatitudeType =
3292
0
        CSLFetchNameValueDef(papszOptions, "COORDINATE_SYSTEM_NAME", "");
3293
0
    poDS->m_osLongitudeDirection =
3294
0
        CSLFetchNameValueDef(papszOptions, "POSITIVE_LONGITUDE_DIRECTION", "");
3295
0
    poDS->m_osTargetName =
3296
0
        CSLFetchNameValueDef(papszOptions, "TARGET_NAME", "");
3297
0
    poDS->m_bInitToNodata = true;
3298
0
    poDS->m_oSrcJSonLabel = std::move(oSrcJSonLabel);
3299
0
    poDS->m_eCompress = eCompress;
3300
0
    poDS->m_anRecordOffsets = std::move(anRecordOffsets);
3301
0
    poDS->eAccess = GA_Update;
3302
3303
    /* -------------------------------------------------------------------- */
3304
    /*      Create band information objects.                                */
3305
    /* -------------------------------------------------------------------- */
3306
0
    const vsi_l_offset nBandOffset =
3307
0
        static_cast<vsi_l_offset>(nLineOffset) * nYSize;
3308
0
    for (int i = 0; i < nBandsIn; i++)
3309
0
    {
3310
0
        GDALRasterBand *poBand;
3311
0
        if (eCompress != COMPRESS_NONE)
3312
0
        {
3313
0
            poBand = new VICARBASICRasterBand(poDS, i + 1, eType);
3314
0
        }
3315
0
        else
3316
0
        {
3317
0
            poBand = new VICARRawRasterBand(
3318
0
                poDS, i + 1, poDS->fpImage,
3319
0
                i * nBandOffset,  // will be set later to final value since we
3320
                                  // need to include the label size
3321
0
                nPixelOffset, nLineOffset, eType,
3322
0
                RawRasterBand::ByteOrder::ORDER_LITTLE_ENDIAN);
3323
0
        }
3324
0
        poDS->SetBand(i + 1, poBand);
3325
0
    }
3326
3327
0
    return poDS;
3328
0
}
3329
3330
/************************************************************************/
3331
/*                             CreateCopy()                             */
3332
/************************************************************************/
3333
3334
GDALDataset *VICARDataset::CreateCopy(const char *pszFilename,
3335
                                      GDALDataset *poSrcDS, int /*bStrict*/,
3336
                                      CSLConstList papszOptions,
3337
                                      GDALProgressFunc pfnProgress,
3338
                                      void *pProgressData)
3339
0
{
3340
0
    if (poSrcDS->GetRasterCount() == 0)
3341
0
    {
3342
0
        CPLError(CE_Failure, CPLE_NotSupported, "Unsupported band count");
3343
0
        return nullptr;
3344
0
    }
3345
3346
0
    const int nXSize = poSrcDS->GetRasterXSize();
3347
0
    const int nYSize = poSrcDS->GetRasterYSize();
3348
0
    const int nBands = poSrcDS->GetRasterCount();
3349
0
    GDALDataType eType = poSrcDS->GetRasterBand(1)->GetRasterDataType();
3350
0
    auto poDS = std::unique_ptr<VICARDataset>(CreateInternal(
3351
0
        pszFilename, nXSize, nYSize, nBands, eType, papszOptions));
3352
0
    if (poDS == nullptr)
3353
0
        return nullptr;
3354
3355
0
    GDALGeoTransform gt;
3356
0
    if (poSrcDS->GetGeoTransform(gt) == CE_None && gt != GDALGeoTransform())
3357
0
    {
3358
0
        poDS->SetGeoTransform(gt);
3359
0
    }
3360
3361
0
    auto poSrcSRS = poSrcDS->GetSpatialRef();
3362
0
    if (poSrcSRS)
3363
0
    {
3364
0
        poDS->SetSpatialRef(poSrcSRS);
3365
0
    }
3366
3367
0
    if (poDS->m_bUseSrcLabel && !poDS->m_oSrcJSonLabel.IsValid())
3368
0
    {
3369
0
        CSLConstList papszMD_VICAR = poSrcDS->GetMetadata("json:VICAR");
3370
0
        if (papszMD_VICAR != nullptr)
3371
0
        {
3372
0
            poDS->SetMetadata(papszMD_VICAR, "json:VICAR");
3373
0
        }
3374
0
    }
3375
3376
0
    poDS->m_bInitToNodata = false;
3377
0
    CPLErr eErr = GDALDatasetCopyWholeRaster(poSrcDS, poDS.get(), nullptr,
3378
0
                                             pfnProgress, pProgressData);
3379
0
    poDS->FlushCache(false);
3380
0
    if (eErr != CE_None)
3381
0
    {
3382
0
        return nullptr;
3383
0
    }
3384
3385
0
    return poDS.release();
3386
0
}
3387
3388
/************************************************************************/
3389
/*                         GDALRegister_VICAR()                         */
3390
/************************************************************************/
3391
3392
void GDALRegister_VICAR()
3393
3394
24
{
3395
24
    if (GDALGetDriverByName(VICAR_DRIVER_NAME) != nullptr)
3396
0
        return;
3397
3398
24
    GDALDriver *poDriver = new GDALDriver();
3399
24
    VICARDriverSetCommonMetadata(poDriver);
3400
3401
24
    poDriver->pfnOpen = VICARDataset::Open;
3402
24
    poDriver->pfnCreate = VICARDataset::Create;
3403
24
    poDriver->pfnCreateCopy = VICARDataset::CreateCopy;
3404
3405
24
    GetGDALDriverManager()->RegisterDriver(poDriver);
3406
24
}