Coverage Report

Created: 2026-02-14 09:00

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/frmts/safe/safedataset.cpp
Line
Count
Source
1
/******************************************************************************
2
 *
3
 * Project:  Sentinel SAFE products
4
 * Purpose:  Sentinel Products (manifest.safe) driver
5
 * Author:   Delfim Rego, delfimrego@gmail.com
6
 *
7
 ******************************************************************************
8
 * Copyright (c) 2015, Delfim Rego <delfimrego@gmail.com>
9
 *
10
 * SPDX-License-Identifier: MIT
11
 ****************************************************************************/
12
13
#include "safedataset.h"
14
15
#include "gdal_driver.h"
16
#include "gdal_drivermanager.h"
17
#include "gdal_openinfo.h"
18
#include "gdal_cpp_functions.h"
19
20
#include "cpl_time.h"
21
22
#include <algorithm>
23
24
#ifdef USE_OMP
25
#include <omp.h>
26
#endif
27
28
#ifdef USE_OMP
29
static int GetNumThreadsToUse()
30
{
31
    unsigned int nCores = std::thread::hardware_concurrency();
32
    return (nCores / 2 > 1) ? nCores / 2 : 1;
33
}
34
#endif
35
36
/************************************************************************/
37
/*                            SAFERasterBand                            */
38
/************************************************************************/
39
40
SAFERasterBand::SAFERasterBand(SAFEDataset *poDSIn, GDALDataType eDataTypeIn,
41
                               const CPLString &osSwath,
42
                               const CPLString &osPolarization,
43
                               std::unique_ptr<GDALDataset> &&poBandFileIn)
44
0
    : poBandFile(std::move(poBandFileIn))
45
0
{
46
0
    poDS = poDSIn;
47
0
    GDALRasterBand *poSrcBand = poBandFile->GetRasterBand(1);
48
0
    poSrcBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
49
0
    eDataType = eDataTypeIn;
50
51
0
    if (!osSwath.empty())
52
0
        SetMetadataItem("SWATH", osSwath.c_str());
53
54
0
    if (!osPolarization.empty())
55
0
        SetMetadataItem("POLARIZATION", osPolarization.c_str());
56
0
}
57
58
/************************************************************************/
59
/*                             IReadBlock()                             */
60
/************************************************************************/
61
62
CPLErr SAFERasterBand::IReadBlock(int nBlockXOff, int nBlockYOff, void *pImage)
63
64
0
{
65
    /* -------------------------------------------------------------------- */
66
    /*      If the last strip is partial, we need to avoid                  */
67
    /*      over-requesting.  We also need to initialize the extra part     */
68
    /*      of the block to zero.                                           */
69
    /* -------------------------------------------------------------------- */
70
0
    const int nYOff = nBlockYOff * nBlockYSize;
71
0
    const int nRequestYSize = std::min(nBlockYSize, nRasterYSize - nYOff);
72
73
    /*-------------------------------------------------------------------- */
74
    /*      If the input imagery is tiled, also need to avoid over-        */
75
    /*      requesting in the X-direction.                                 */
76
    /* ------------------------------------------------------------------- */
77
0
    const int nXOff = nBlockXOff * nBlockXSize;
78
0
    const int nRequestXSize = std::min(nBlockXSize, nRasterXSize - nXOff);
79
80
0
    if (nRequestXSize < nBlockXSize || nRequestYSize < nBlockYSize)
81
0
    {
82
0
        memset(pImage, 0,
83
0
               static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
84
0
                   nBlockXSize * nBlockYSize);
85
0
    }
86
87
0
    if (eDataType == GDT_CInt16 && poBandFile->GetRasterCount() == 2)
88
0
        return poBandFile->RasterIO(GF_Read, nXOff, nYOff, nRequestXSize,
89
0
                                    nRequestYSize, pImage, nRequestXSize,
90
0
                                    nRequestYSize, GDT_Int16, 2, nullptr, 4,
91
0
                                    nBlockXSize * 4, 2, nullptr);
92
93
    /* -------------------------------------------------------------------- */
94
    /*      File has one sample marked as sample format void, a 32bits.     */
95
    /* -------------------------------------------------------------------- */
96
0
    else if (eDataType == GDT_CInt16 && poBandFile->GetRasterCount() == 1)
97
0
    {
98
0
        CPLErr eErr = poBandFile->RasterIO(
99
0
            GF_Read, nXOff, nYOff, nRequestXSize, nRequestYSize, pImage,
100
0
            nRequestXSize, nRequestYSize, GDT_CInt16, 1, nullptr, 4,
101
0
            nBlockXSize * 4, 0, nullptr);
102
103
0
        return eErr;
104
0
    }
105
106
    /* -------------------------------------------------------------------- */
107
    /*      The 16bit case is straight forward.  The underlying file        */
108
    /*      looks like a 16bit unsigned data too.                           */
109
    /* -------------------------------------------------------------------- */
110
0
    else if (eDataType == GDT_UInt16)
111
0
        return poBandFile->RasterIO(GF_Read, nXOff, nYOff, nRequestXSize,
112
0
                                    nRequestYSize, pImage, nRequestXSize,
113
0
                                    nRequestYSize, GDT_UInt16, 1, nullptr, 2,
114
0
                                    nBlockXSize * 2, 0, nullptr);
115
116
0
    else if (eDataType == GDT_UInt8)
117
0
        return poBandFile->RasterIO(GF_Read, nXOff, nYOff, nRequestXSize,
118
0
                                    nRequestYSize, pImage, nRequestXSize,
119
0
                                    nRequestYSize, GDT_UInt8, 1, nullptr, 1,
120
0
                                    nBlockXSize, 0, nullptr);
121
122
0
    CPLAssert(false);
123
0
    return CE_Failure;
124
0
}
125
126
/************************************************************************/
127
/*                          SAFESLCRasterBand                           */
128
/************************************************************************/
129
130
SAFESLCRasterBand::SAFESLCRasterBand(
131
    SAFEDataset *poDSIn, GDALDataType eDataTypeIn, const CPLString &osSwath,
132
    const CPLString &osPolarization,
133
    std::unique_ptr<GDALDataset> &&poBandFileIn, BandType eBandType)
134
0
    : poBandFile(std::move(poBandFileIn))
135
0
{
136
0
    poDS = poDSIn;
137
0
    eDataType = eDataTypeIn;
138
0
    m_eInputDataType = eDataTypeIn;
139
0
    GDALRasterBand *poSrcBand = poBandFile->GetRasterBand(1);
140
0
    poSrcBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
141
0
    m_eBandType = eBandType;
142
143
0
    if (!osSwath.empty())
144
0
        SetMetadataItem("SWATH", osSwath.c_str());
145
146
0
    if (!osPolarization.empty())
147
0
        SetMetadataItem("POLARIZATION", osPolarization.c_str());
148
149
    // For intensity band
150
0
    if (m_eBandType == INTENSITY)
151
0
        eDataType = GDT_Float32;
152
0
    else
153
        // For complex bands
154
0
        eDataType = GDT_CInt16;
155
0
}
156
157
/************************************************************************/
158
/*                             IReadBlock()                             */
159
/************************************************************************/
160
161
CPLErr SAFESLCRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
162
                                     void *pImage)
163
164
0
{
165
    /* -------------------------------------------------------------------- */
166
    /*      If the last strip is partial, we need to avoid                  */
167
    /*      over-requesting.  We also need to initialize the extra part     */
168
    /*      of the block to zero.                                           */
169
    /* -------------------------------------------------------------------- */
170
0
    const int nXOff = nBlockXOff * nBlockXSize;
171
0
    const int nRequestXSize = std::min(nBlockXSize, nRasterXSize - nXOff);
172
0
    const int nYOff = nBlockYOff * nBlockYSize;
173
0
    const int nRequestYSize = std::min(nBlockYSize, nRasterYSize - nYOff);
174
175
0
    if (nRequestYSize < nBlockYSize)
176
0
    {
177
0
        memset(pImage, 0,
178
0
               static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
179
0
                   nBlockXSize * nBlockYSize);
180
0
    }
181
182
    /*-------------------------------------------------------------------- */
183
    /*      If the input imagery is tiled, also need to avoid over-        */
184
    /*      requesting in the X-direction.                                 */
185
    /* ------------------------------------------------------------------- */
186
0
    if (nRequestXSize < nBlockXSize)
187
0
    {
188
0
        memset(pImage, 0,
189
0
               static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
190
0
                   nBlockXSize * nBlockYSize);
191
0
    }
192
193
0
    if (m_eInputDataType == GDT_CInt16 && poBandFile->GetRasterCount() == 2)
194
0
    {
195
0
        return poBandFile->RasterIO(GF_Read, nXOff, nYOff, nRequestXSize,
196
0
                                    nRequestYSize, pImage, nRequestXSize,
197
0
                                    nRequestYSize, GDT_Int16, 2, nullptr, 4,
198
0
                                    nBlockXSize * 4, 2, nullptr);
199
0
    }
200
    // File has one sample marked as sample format void, a 32bits.
201
0
    else if (m_eInputDataType == GDT_CInt16 &&
202
0
             poBandFile->GetRasterCount() == 1)
203
0
    {
204
0
        if (m_eBandType == COMPLEX)
205
0
        {
206
0
            CPLErr eErr = poBandFile->RasterIO(
207
0
                GF_Read, nXOff, nYOff, nRequestXSize, nRequestYSize, pImage,
208
0
                nRequestXSize, nRequestYSize, GDT_CInt16, 1, nullptr, 4,
209
0
                nBlockXSize * 4, 0, nullptr);
210
0
            if (eErr != CE_None)
211
0
            {
212
0
                return eErr;
213
0
            }
214
0
        }
215
0
        else if (m_eBandType == INTENSITY)
216
0
        {
217
0
            GInt16 *pnImageTmp = static_cast<GInt16 *>(VSI_MALLOC3_VERBOSE(
218
0
                2 * sizeof(int16_t), nBlockXSize, nBlockYSize));
219
0
            if (!pnImageTmp)
220
0
            {
221
0
                return CE_Failure;
222
0
            }
223
224
0
            CPLErr eErr = poBandFile->RasterIO(
225
0
                GF_Read, nXOff, nYOff, nRequestXSize, nRequestYSize, pnImageTmp,
226
0
                nRequestXSize, nRequestYSize, GDT_CInt16, 1, nullptr, 4,
227
0
                nBlockXSize * 4, 0, nullptr);
228
0
            if (eErr != CE_None)
229
0
            {
230
0
                CPLFree(pnImageTmp);
231
0
                return eErr;
232
0
            }
233
234
0
            float *pfBuffer = static_cast<float *>(pImage);
235
#ifdef USE_OMP
236
            omp_set_num_threads(GetNumThreadsToUse());
237
#pragma omp parallel
238
#endif
239
0
            for (int i = 0; i < nBlockYSize; i++)
240
0
            {
241
#ifdef USE_OMP
242
#pragma omp for nowait
243
#endif
244
0
                for (int j = 0; j < nBlockXSize; j++)
245
0
                {
246
0
                    int nPixOff = (2 * (i * nBlockXSize)) + (j * 2);
247
0
                    int nOutPixOff = (i * nBlockXSize) + j;
248
0
                    pfBuffer[nOutPixOff] = static_cast<float>(
249
0
                        static_cast<double>(pnImageTmp[nPixOff] *
250
0
                                            pnImageTmp[nPixOff]) +
251
0
                        static_cast<double>(pnImageTmp[nPixOff + 1] *
252
0
                                            pnImageTmp[nPixOff + 1]));
253
0
                }
254
0
            }
255
0
            CPLFree(pnImageTmp);
256
0
        }
257
0
        return CE_None;
258
0
    }
259
260
0
    CPLAssert(false);
261
0
    return CE_Failure;
262
0
}
263
264
/************************************************************************/
265
/*                         SAFECalibRasterBand                          */
266
/************************************************************************/
267
268
SAFECalibratedRasterBand::SAFECalibratedRasterBand(
269
    SAFEDataset *poDSIn, GDALDataType eDataTypeIn, const CPLString &osSwath,
270
    const CPLString &osPolarization,
271
    std::unique_ptr<GDALDataset> &&poBandDatasetIn,
272
    const char *pszCalibrationFilename, CalibrationType eCalibrationType)
273
0
    : poBandDataset(std::move(poBandDatasetIn))
274
0
{
275
0
    poDS = poDSIn;
276
0
    GDALRasterBand *poSrcBand = poBandDataset->GetRasterBand(1);
277
0
    poSrcBand->GetBlockSize(&nBlockXSize, &nBlockYSize);
278
0
    eDataType = eDataTypeIn;
279
280
0
    if (!osSwath.empty())
281
0
        SetMetadataItem("SWATH", osSwath.c_str());
282
283
0
    if (!osPolarization.empty())
284
0
        SetMetadataItem("POLARIZATION", osPolarization.c_str());
285
286
0
    m_osCalibrationFilename = pszCalibrationFilename;
287
0
    m_eInputDataType = eDataTypeIn;
288
0
    eDataType = GDT_Float32;
289
0
    m_eCalibrationType = eCalibrationType;
290
0
}
291
292
/************************************************************************/
293
/*                              ReadLUT()                               */
294
/************************************************************************/
295
/* Read the provided LUT in to m_ndTable                                */
296
/************************************************************************/
297
bool SAFECalibratedRasterBand::ReadLUT()
298
0
{
299
0
    const char *const papszCalibrationNodes[3] = {
300
0
        "=calibrationVector.sigmaNought", "=calibrationVector.betaNought",
301
0
        "=calibrationVector.gamma"};
302
0
    CPLString osCalibrationNodeName = papszCalibrationNodes[m_eCalibrationType];
303
0
    const char *pszEndSpace = " ";
304
0
    CPLString osStartTime, osEndTime;
305
0
    CPLXMLNode *poLUT = CPLParseXMLFile(m_osCalibrationFilename);
306
0
    if (!poLUT)
307
0
        return false;
308
309
0
    CPLString osParseLUT;
310
0
    CPLString osParsePixelLUT;
311
0
    CPLString osParseAzimuthLUT;
312
0
    CPLString osParseLineNoLUT;
313
0
    for (CPLXMLNode *psNode = poLUT; psNode != nullptr;)
314
0
    {
315
0
        if (psNode->psNext != nullptr)
316
0
            psNode = psNode->psNext;
317
318
0
        if (EQUAL(psNode->pszValue, "calibration"))
319
0
        {
320
0
            for (psNode = psNode->psChild; psNode != nullptr;)
321
0
            {
322
0
                if (EQUAL(psNode->pszValue, "adsHeader"))
323
0
                {
324
0
                    osStartTime =
325
0
                        CPLGetXMLValue(psNode, "=adsHeader.startTime", " ");
326
0
                    osEndTime =
327
0
                        CPLGetXMLValue(psNode, "=adsHeader.stopTime", " ");
328
0
                }
329
330
0
                if (psNode->psNext != nullptr)
331
0
                    psNode = psNode->psNext;
332
333
0
                if (EQUAL(psNode->pszValue, "calibrationVectorList"))
334
0
                {
335
0
                    for (psNode = psNode->psChild; psNode != nullptr;
336
0
                         psNode = psNode->psNext)
337
0
                    {
338
0
                        if (EQUAL(psNode->pszValue, "calibrationVector"))
339
0
                        {
340
0
                            osParseAzimuthLUT += CPLGetXMLValue(
341
0
                                psNode, "=calibrationVector.azimuthTime", " ");
342
0
                            osParseAzimuthLUT += pszEndSpace;
343
0
                            osParseLineNoLUT += CPLGetXMLValue(
344
0
                                psNode, "=calibrationVector.line", " ");
345
0
                            osParseLineNoLUT += pszEndSpace;
346
0
                            osParsePixelLUT =
347
0
                                CPLGetXMLValue(psNode, "pixel", " ");
348
0
                            m_nNumPixels = static_cast<int>(CPLAtof(
349
0
                                CPLGetXMLValue(psNode, "pixel.count", " ")));
350
0
                            osParseLUT += CPLGetXMLValue(
351
0
                                psNode, osCalibrationNodeName, " ");
352
0
                            osParseLUT += pszEndSpace;
353
0
                        }
354
0
                    }
355
0
                }
356
0
            }
357
0
        }
358
0
    }
359
0
    CPLDestroyXMLNode(poLUT);
360
361
0
    osParsePixelLUT += pszEndSpace;
362
363
0
    CPLStringList oStartTimeList(
364
0
        CSLTokenizeString2(osStartTime, " ", CSLT_HONOURSTRINGS));
365
0
    if (!oStartTimeList.size())
366
0
        return false;
367
0
    m_oStartTimePoint = getTimePoint(oStartTimeList[0]);
368
0
    CPLStringList oEndTimeList(
369
0
        CSLTokenizeString2(osEndTime, " ", CSLT_HONOURSTRINGS));
370
0
    if (!oEndTimeList.size())
371
0
        return false;
372
0
    m_oStopTimePoint = getTimePoint(oEndTimeList[0]);
373
0
    m_oAzimuthList.Assign(
374
0
        CSLTokenizeString2(osParseAzimuthLUT, " ", CSLT_HONOURSTRINGS));
375
0
    CPLStringList oLUTList(
376
0
        CSLTokenizeString2(osParseLUT, " ", CSLT_HONOURSTRINGS));
377
0
    CPLStringList oPixelList(
378
0
        CSLTokenizeString2(osParsePixelLUT, " ", CSLT_HONOURSTRINGS));
379
0
    CPLStringList oLineNoList(
380
0
        CSLTokenizeString2(osParseLineNoLUT, " ", CSLT_HONOURSTRINGS));
381
382
0
    m_anPixelLUT.resize(m_nNumPixels);
383
0
    for (int i = 0; i < m_nNumPixels; i++)
384
0
        m_anPixelLUT[i] = static_cast<int>(CPLAtof(oPixelList[i]));
385
386
0
    int nTableSize = oLUTList.size();
387
0
    m_afTable.resize(nTableSize);
388
0
    for (int i = 0; i < nTableSize; i++)
389
0
        m_afTable[i] = static_cast<float>(CPLAtof(oLUTList[i]));
390
391
0
    int nLineListSize = oLineNoList.size();
392
0
    m_anLineLUT.resize(nLineListSize);
393
0
    for (int i = 0; i < nLineListSize; i++)
394
0
        m_anLineLUT[i] = static_cast<int>(CPLAtof(oLineNoList[i]));
395
396
0
    return true;
397
0
}
398
399
/************************************************************************/
400
/*                             IReadBlock()                             */
401
/************************************************************************/
402
403
CPLErr SAFECalibratedRasterBand::IReadBlock(int nBlockXOff, int nBlockYOff,
404
                                            void *pImage)
405
406
0
{
407
    /* -------------------------------------------------------------------- */
408
    /*      If the last strip is partial, we need to avoid                  */
409
    /*      over-requesting.  We also need to initialize the extra part     */
410
    /*      of the block to zero.                                           */
411
    /* -------------------------------------------------------------------- */
412
0
    const int nXOff = nBlockXOff * nBlockXSize;
413
0
    const int nRequestXSize = std::min(nBlockXSize, nRasterXSize - nXOff);
414
0
    const int nYOff = nBlockYOff * nBlockYSize;
415
0
    const int nRequestYSize = std::min(nBlockYSize, nRasterYSize - nYOff);
416
417
0
    if (nRequestYSize < nBlockYSize)
418
0
    {
419
0
        memset(pImage, 0,
420
0
               static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
421
0
                   nBlockXSize * nBlockYSize);
422
0
    }
423
424
    /*-------------------------------------------------------------------- */
425
    /*      If the input imagery is tiled, also need to avoid over-        */
426
    /*      requesting in the X-direction.                                 */
427
    /* ------------------------------------------------------------------- */
428
0
    if (nRequestXSize < nBlockXSize)
429
0
    {
430
0
        memset(pImage, 0,
431
0
               static_cast<size_t>(GDALGetDataTypeSizeBytes(eDataType)) *
432
0
                   nBlockXSize * nBlockYSize);
433
0
    }
434
435
    // Check LUT values and fail before reading
436
0
    int nLineCalVecIdx = getCalibrationVectorIndex(nBlockYOff);
437
0
    const char *pszVec0Str = m_oAzimuthList[nLineCalVecIdx];
438
0
    const char *pszVec1Str = m_oAzimuthList[nLineCalVecIdx + 1];
439
0
    if (((m_eInputDataType == GDT_CInt16) || (m_eInputDataType == GDT_Int16)) &&
440
0
        (!pszVec0Str || !pszVec1Str))
441
0
        return CE_Failure;
442
443
0
    TimePoint azTime = getazTime(m_oStartTimePoint, m_oStopTimePoint,
444
0
                                 nRasterYSize, nBlockYOff);
445
0
    TimePoint oVec0Time = getTimePoint(pszVec0Str);
446
0
    TimePoint oVec1Time = getTimePoint(pszVec1Str);
447
0
    double dfMuY =
448
0
        getTimeDiff(oVec0Time, azTime) / getTimeDiff(oVec0Time, oVec1Time);
449
450
0
    if (m_eInputDataType == GDT_CInt16)
451
0
    {
452
0
        CPLErr eErr = CE_None;
453
0
        GInt16 *pnImageTmp = static_cast<GInt16 *>(
454
0
            VSI_MALLOC3_VERBOSE(2 * sizeof(int16_t), nBlockXSize, nBlockYSize));
455
0
        if (!pnImageTmp)
456
0
            return CE_Failure;
457
458
0
        if (poBandDataset->GetRasterCount() == 2)
459
0
        {
460
0
            eErr = poBandDataset->RasterIO(
461
0
                GF_Read, nBlockXOff * nBlockXSize, nBlockYOff * nBlockYSize,
462
0
                nRequestXSize, nRequestYSize, pnImageTmp, nRequestXSize,
463
0
                nRequestYSize, GDT_Int16, 2, nullptr, 4, nBlockXSize * 4, 2,
464
0
                nullptr);
465
0
        }
466
        /* --------------------------------------------------------------------
467
         */
468
        /*      File has one sample marked as sample format void, a 32bits. */
469
        /* --------------------------------------------------------------------
470
         */
471
0
        else if (poBandDataset->GetRasterCount() == 1)
472
0
        {
473
0
            eErr = poBandDataset->RasterIO(
474
0
                GF_Read, nBlockXOff * nBlockXSize, nBlockYOff * nBlockYSize,
475
0
                nRequestXSize, nRequestYSize, pnImageTmp, nRequestXSize,
476
0
                nRequestYSize, GDT_CInt16, 1, nullptr, 4, nBlockXSize * 4, 0,
477
0
                nullptr);
478
0
        }
479
480
        // Interpolation of LUT Value
481
#ifdef USE_OMP
482
        omp_set_num_threads(GetNumThreadsToUse());
483
#pragma omp parallel
484
#endif
485
0
        for (int i = 0; i < nBlockYSize; i++)
486
0
        {
487
#ifdef USE_OMP
488
#pragma omp for nowait
489
#endif
490
0
            for (int j = 0; j < nBlockXSize; j++)
491
0
            {
492
0
                int nPixOff = (2 * (i * nBlockXSize)) + (j * 2);
493
0
                int nOutPixOff = (i * nBlockXSize) + j;
494
0
                int nPixelCalvecIdx = getPixelIndex(j);
495
0
                double dfMuX = (double)(j - m_anPixelLUT[nPixelCalvecIdx]) /
496
0
                               (double)(m_anPixelLUT[nPixelCalvecIdx + 1] -
497
0
                                        m_anPixelLUT[nPixelCalvecIdx]);
498
0
                int lutIdx1 = (nLineCalVecIdx * m_nNumPixels) + nPixelCalvecIdx;
499
0
                int lutIdx2 =
500
0
                    (nLineCalVecIdx * m_nNumPixels) + (nPixelCalvecIdx + 1);
501
0
                int lutIdx3 =
502
0
                    ((nLineCalVecIdx + 1) * m_nNumPixels) + nPixelCalvecIdx;
503
0
                int lutIdx4 = ((nLineCalVecIdx + 1) * m_nNumPixels) +
504
0
                              (nPixelCalvecIdx + 1);
505
0
                double dfLutValue =
506
0
                    ((1 - dfMuY) * (((1 - dfMuX) * m_afTable[lutIdx1]) +
507
0
                                    (dfMuX * m_afTable[lutIdx2]))) +
508
0
                    (dfMuY * (((1 - dfMuX) * m_afTable[lutIdx3]) +
509
0
                              (dfMuX * m_afTable[lutIdx4])));
510
0
                double dfNum = static_cast<double>(
511
0
                    (pnImageTmp[nPixOff] * pnImageTmp[nPixOff]) +
512
0
                    (pnImageTmp[nPixOff + 1] * pnImageTmp[nPixOff + 1]));
513
0
                double dfCalibValue = dfNum / (dfLutValue * dfLutValue);
514
0
                ((float *)pImage)[nOutPixOff] = (float)dfCalibValue;
515
0
            }
516
0
        }
517
0
        CPLFree(pnImageTmp);
518
0
        return eErr;
519
0
    }
520
0
    else if (m_eInputDataType == GDT_UInt16)
521
0
    {
522
0
        CPLErr eErr = CE_None;
523
0
        GUInt16 *pnImageTmp = static_cast<GUInt16 *>(VSI_MALLOC3_VERBOSE(
524
0
            nBlockXSize, nBlockYSize, GDALGetDataTypeSizeBytes(GDT_UInt16)));
525
0
        if (!pnImageTmp)
526
0
            return CE_Failure;
527
0
        eErr = poBandDataset->RasterIO(GF_Read, nBlockXOff * nBlockXSize,
528
0
                                       nBlockYOff * nBlockYSize, nRequestXSize,
529
0
                                       nRequestYSize, pnImageTmp, nRequestXSize,
530
0
                                       nRequestYSize, GDT_UInt16, 1, nullptr, 2,
531
0
                                       nBlockXSize * 2, 0, nullptr);
532
533
#ifdef USE_OMP
534
        omp_set_num_threads(GetNumThreadsToUse());
535
#pragma omp parallel
536
#endif
537
0
        for (int i = 0; i < nBlockYSize; i++)
538
0
        {
539
#ifdef USE_OMP
540
#pragma omp for nowait
541
#endif
542
0
            for (int j = 0; j < nBlockXSize; j++)
543
0
            {
544
0
                int nPixOff = (i * nBlockXSize) + j;
545
0
                int nPixelCalvecIdx = getPixelIndex(j);
546
0
                double dfMuX = (double)(j - m_anPixelLUT[nPixelCalvecIdx]) /
547
0
                               (double)(m_anPixelLUT[nPixelCalvecIdx + 1] -
548
0
                                        m_anPixelLUT[nPixelCalvecIdx]);
549
0
                int lutIdx1 = (nLineCalVecIdx * m_nNumPixels) + nPixelCalvecIdx;
550
0
                int lutIdx2 =
551
0
                    (nLineCalVecIdx * m_nNumPixels) + (nPixelCalvecIdx + 1);
552
0
                int lutIdx3 =
553
0
                    ((nLineCalVecIdx + 1) * m_nNumPixels) + (nPixelCalvecIdx);
554
0
                int lutIdx4 = ((nLineCalVecIdx + 1) * m_nNumPixels) +
555
0
                              (nPixelCalvecIdx + 1);
556
0
                double dfLutValue =
557
0
                    ((1 - dfMuY) * (((1 - dfMuX) * m_afTable[lutIdx1]) +
558
0
                                    (dfMuX * m_afTable[lutIdx2]))) +
559
0
                    (dfMuY * (((1 - dfMuX) * m_afTable[lutIdx3]) +
560
0
                              (dfMuX * m_afTable[lutIdx4])));
561
0
                double dfCalibValue =
562
0
                    (double)(pnImageTmp[nPixOff] * pnImageTmp[nPixOff]) /
563
0
                    (dfLutValue * dfLutValue);
564
0
                ((float *)pImage)[nPixOff] = (float)dfCalibValue;
565
0
            }
566
0
        }
567
0
        CPLFree(pnImageTmp);
568
0
        return eErr;
569
0
    }
570
0
    else if (eDataType == GDT_UInt8)  // Check if this is required.
571
0
        return poBandDataset->RasterIO(
572
0
            GF_Read, nBlockXOff * nBlockXSize, nBlockYOff * nBlockYSize,
573
0
            nRequestXSize, nRequestYSize, pImage, nRequestXSize, nRequestYSize,
574
0
            GDT_UInt8, 1, nullptr, 1, nBlockXSize, 0, nullptr);
575
576
0
    CPLAssert(false);
577
0
    return CE_Failure;
578
0
}
579
580
/************************************************************************/
581
/* ==================================================================== */
582
/*                              SAFEDataset                              */
583
/* ==================================================================== */
584
/************************************************************************/
585
586
/************************************************************************/
587
/*                            SAFEDataset()                             */
588
/************************************************************************/
589
590
SAFEDataset::SAFEDataset()
591
0
{
592
0
    m_oGCPSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
593
0
}
594
595
/************************************************************************/
596
/*                            ~SAFEDataset()                            */
597
/************************************************************************/
598
599
SAFEDataset::~SAFEDataset()
600
601
0
{
602
0
    SAFEDataset::FlushCache(true);
603
604
0
    if (nGCPCount > 0)
605
0
    {
606
0
        GDALDeinitGCPs(nGCPCount, pasGCPList);
607
0
        CPLFree(pasGCPList);
608
0
    }
609
610
0
    SAFEDataset::CloseDependentDatasets();
611
612
0
    CSLDestroy(papszSubDatasets);
613
0
    CSLDestroy(papszExtraFiles);
614
0
}
615
616
/************************************************************************/
617
/*                       CloseDependentDatasets()                       */
618
/************************************************************************/
619
620
int SAFEDataset::CloseDependentDatasets()
621
0
{
622
0
    int bHasDroppedRef = GDALPamDataset::CloseDependentDatasets();
623
624
0
    if (nBands != 0)
625
0
        bHasDroppedRef = TRUE;
626
627
0
    for (int iBand = 0; iBand < nBands; iBand++)
628
0
    {
629
0
        delete papoBands[iBand];
630
0
    }
631
0
    nBands = 0;
632
633
0
    return bHasDroppedRef;
634
0
}
635
636
/************************************************************************/
637
/*                         GetMetaDataObject()                          */
638
/************************************************************************/
639
640
const CPLXMLNode *
641
SAFEDataset::GetMetaDataObject(const CPLXMLNode *psMetaDataObjects,
642
                               const char *metadataObjectId)
643
0
{
644
    /* -------------------------------------------------------------------- */
645
    /*      Look for DataObject Element by ID.                              */
646
    /* -------------------------------------------------------------------- */
647
0
    for (const CPLXMLNode *psMDO = psMetaDataObjects->psChild; psMDO != nullptr;
648
0
         psMDO = psMDO->psNext)
649
0
    {
650
0
        if (psMDO->eType != CXT_Element ||
651
0
            !(EQUAL(psMDO->pszValue, "metadataObject")))
652
0
        {
653
0
            continue;
654
0
        }
655
656
0
        const char *pszElementID = CPLGetXMLValue(psMDO, "ID", "");
657
658
0
        if (EQUAL(pszElementID, metadataObjectId))
659
0
        {
660
0
            return psMDO;
661
0
        }
662
0
    }
663
664
0
    CPLError(CE_Warning, CPLE_AppDefined, "MetadataObject not found with ID=%s",
665
0
             metadataObjectId);
666
667
0
    return nullptr;
668
0
}
669
670
/************************************************************************/
671
/*                           GetDataObject()                            */
672
/************************************************************************/
673
674
const CPLXMLNode *SAFEDataset::GetDataObject(const CPLXMLNode *psDataObjects,
675
                                             const char *dataObjectId)
676
0
{
677
    /* -------------------------------------------------------------------- */
678
    /*      Look for DataObject Element by ID.                              */
679
    /* -------------------------------------------------------------------- */
680
0
    for (const CPLXMLNode *psDO = psDataObjects->psChild; psDO != nullptr;
681
0
         psDO = psDO->psNext)
682
0
    {
683
0
        if (psDO->eType != CXT_Element ||
684
0
            !(EQUAL(psDO->pszValue, "dataObject")))
685
0
        {
686
0
            continue;
687
0
        }
688
689
0
        const char *pszElementID = CPLGetXMLValue(psDO, "ID", "");
690
691
0
        if (EQUAL(pszElementID, dataObjectId))
692
0
        {
693
0
            return psDO;
694
0
        }
695
0
    }
696
697
0
    CPLError(CE_Warning, CPLE_AppDefined, "DataObject not found with ID=%s",
698
0
             dataObjectId);
699
700
0
    return nullptr;
701
0
}
702
703
const CPLXMLNode *
704
SAFEDataset::GetDataObject(const CPLXMLNode *psMetaDataObjects,
705
                           const CPLXMLNode *psDataObjects,
706
                           const char *metadataObjectId)
707
0
{
708
    /* -------------------------------------------------------------------- */
709
    /*      Look for MetadataObject Element by ID.                          */
710
    /* -------------------------------------------------------------------- */
711
0
    const CPLXMLNode *psMDO =
712
0
        SAFEDataset::GetMetaDataObject(psMetaDataObjects, metadataObjectId);
713
714
0
    if (psMDO != nullptr)
715
0
    {
716
0
        const char *dataObjectId =
717
0
            CPLGetXMLValue(psMDO, "dataObjectPointer.dataObjectID", "");
718
0
        if (*dataObjectId != '\0')
719
0
        {
720
0
            return SAFEDataset::GetDataObject(psDataObjects, dataObjectId);
721
0
        }
722
0
    }
723
724
0
    CPLError(CE_Warning, CPLE_AppDefined, "DataObject not found with MetaID=%s",
725
0
             metadataObjectId);
726
727
0
    return nullptr;
728
0
}
729
730
/************************************************************************/
731
/*                            GetFileList()                             */
732
/************************************************************************/
733
734
char **SAFEDataset::GetFileList()
735
736
0
{
737
0
    char **papszFileList = GDALPamDataset::GetFileList();
738
739
0
    papszFileList = CSLInsertStrings(papszFileList, -1, papszExtraFiles);
740
741
0
    return papszFileList;
742
0
}
743
744
/************************************************************************/
745
/*                              Identify()                              */
746
/************************************************************************/
747
748
int SAFEDataset::Identify(GDALOpenInfo *poOpenInfo)
749
358k
{
750
    /* Check for the case where we're trying to read the calibrated data: */
751
358k
    if (STARTS_WITH_CI(poOpenInfo->pszFilename, "SENTINEL1_CALIB:"))
752
0
    {
753
0
        return TRUE;
754
0
    }
755
756
    /* Check for directory access when there is a manifest.safe file in the
757
       directory. */
758
358k
    if (poOpenInfo->bIsDirectory)
759
3.38k
    {
760
3.38k
        VSIStatBufL sStat;
761
3.38k
        CPLString osMDFilename = CPLFormCIFilenameSafe(
762
3.38k
            poOpenInfo->pszFilename, "manifest.safe", nullptr);
763
764
3.38k
        if (VSIStatL(osMDFilename, &sStat) == 0 && VSI_ISREG(sStat.st_mode))
765
0
        {
766
0
            GDALOpenInfo oOpenInfo(osMDFilename, GA_ReadOnly, nullptr);
767
0
            return Identify(&oOpenInfo);
768
0
        }
769
3.38k
        return FALSE;
770
3.38k
    }
771
772
    /* otherwise, do our normal stuff */
773
355k
    if (!EQUAL(CPLGetFilename(poOpenInfo->pszFilename), "manifest.safe"))
774
355k
        return FALSE;
775
776
0
    if (poOpenInfo->nHeaderBytes < 100)
777
0
        return FALSE;
778
779
0
    const char *pszHeader =
780
0
        reinterpret_cast<const char *>(poOpenInfo->pabyHeader);
781
0
    if (!strstr(pszHeader, "<xfdu:XFDU"))
782
0
        return FALSE;
783
784
    // This driver doesn't handle Sentinel-2 or RCM (RADARSAT Constellation Mission) data
785
0
    if (strstr(pszHeader, "sentinel-2") ||
786
0
        strstr(pszHeader, "rcm_prod_manifest.xsd"))
787
0
        return FALSE;
788
789
0
    return TRUE;
790
0
}
791
792
/************************************************************************/
793
/*                                Open()                                */
794
/************************************************************************/
795
796
GDALDataset *SAFEDataset::Open(GDALOpenInfo *poOpenInfo)
797
798
0
{
799
800
    // Is this a SENTINEL-1 manifest.safe definition?
801
0
    if (!SAFEDataset::Identify(poOpenInfo))
802
0
    {
803
0
        return nullptr;
804
0
    }
805
806
    /* -------------------------------------------------------------------- */
807
    /*        Get subdataset information, if relevant                       */
808
    /* -------------------------------------------------------------------- */
809
0
    CPLString osMDFilename;
810
0
    bool bIsSubDS = false;
811
0
    bool bIsSLC = false;
812
    // Subdataset 1st level selection (ex: for swath selection)
813
0
    CPLString osSelectedSubDS1;
814
    // Subdataset 2nd level selection (ex: for polarization selection)
815
0
    CPLString osSelectedSubDS2;
816
0
    CPLString osSelectedSubDS3;
817
    // 0 for SIGMA , 1 for BETA, 2 for GAMMA and # for UNCALIB dataset
818
0
    SAFECalibratedRasterBand::CalibrationType eCalibrationType =
819
0
        SAFECalibratedRasterBand::SIGMA_NOUGHT;
820
0
    bool bCalibrated = false;
821
822
    // 0 for amplitude, 1 for complex (2 band : I , Q) and 2 for INTENSITY
823
0
    typedef enum
824
0
    {
825
0
        UNKNOWN = -1,
826
0
        AMPLITUDE,
827
0
        COMPLEX,
828
0
        INTENSITY
829
0
    } RequestDataType;
830
831
0
    RequestDataType eRequestType = UNKNOWN;
832
    // Calibration Information selection
833
0
    CPLString osSubdatasetName;
834
835
0
    if (STARTS_WITH_CI(poOpenInfo->pszFilename, "SENTINEL1_CALIB:"))
836
0
    {
837
0
        bIsSubDS = true;
838
0
        osMDFilename = poOpenInfo->pszFilename + strlen("SENTINEL1_CALIB:");
839
0
        const char *pszSelectionCalib = strchr(osMDFilename.c_str(), ':');
840
0
        if (pszSelectionCalib == nullptr ||
841
0
            pszSelectionCalib == osMDFilename.c_str())
842
0
        {
843
0
            CPLError(CE_Failure, CPLE_AppDefined,
844
0
                     "Invalid syntax for SENTINEL1_CALIB:");
845
0
            return nullptr;
846
0
        }
847
848
0
        CPLString osCalibrationValue = osMDFilename;
849
0
        osCalibrationValue.resize(pszSelectionCalib - osMDFilename.c_str());
850
0
        osMDFilename = pszSelectionCalib + strlen(":");
851
0
        if (EQUAL(osCalibrationValue.c_str(), "UNCALIB"))
852
0
        {
853
0
            bCalibrated = false;
854
0
        }
855
0
        else if (EQUAL(osCalibrationValue.c_str(), "SIGMA0"))
856
0
        {
857
0
            bCalibrated = true;
858
0
            eCalibrationType = SAFECalibratedRasterBand::SIGMA_NOUGHT;
859
0
        }
860
0
        else if (EQUAL(osCalibrationValue.c_str(), "BETA0"))
861
0
        {
862
0
            bCalibrated = true;
863
0
            eCalibrationType = SAFECalibratedRasterBand::BETA_NOUGHT;
864
0
        }
865
0
        else if (EQUAL(osCalibrationValue.c_str(), "GAMMA"))
866
0
        {
867
0
            bCalibrated = true;
868
0
            eCalibrationType = SAFECalibratedRasterBand::GAMMA;
869
0
        }
870
0
        else
871
0
        {
872
0
            CPLError(CE_Failure, CPLE_AppDefined,
873
0
                     "Invalid syntax for SENTINEL1_CALIB:");
874
0
            return nullptr;
875
0
        }
876
877
0
        const auto nSelectionUnitPos = osMDFilename.rfind(':');
878
0
        if (nSelectionUnitPos == std::string::npos || nSelectionUnitPos == 0)
879
0
        {
880
0
            CPLError(CE_Failure, CPLE_AppDefined,
881
0
                     "Invalid syntax for SENTINEL1_CALIB:");
882
0
            return nullptr;
883
0
        }
884
885
0
        CPLString osUnitValue = osMDFilename.substr(nSelectionUnitPos + 1);
886
0
        osMDFilename.resize(nSelectionUnitPos);
887
0
        if (EQUAL(osUnitValue.c_str(), "AMPLITUDE"))
888
0
            eRequestType = AMPLITUDE;
889
0
        else if (EQUAL(osUnitValue.c_str(), "COMPLEX"))
890
0
            eRequestType = COMPLEX;
891
0
        else if (EQUAL(osUnitValue.c_str(), "INTENSITY"))
892
0
            eRequestType = INTENSITY;
893
0
        else
894
0
        {
895
0
            CPLError(CE_Failure, CPLE_AppDefined,
896
0
                     "Invalid syntax for SENTINEL1_CALIB:");
897
0
            return nullptr;
898
0
        }
899
900
0
        const auto nSelection1Pos = osMDFilename.rfind(':');
901
0
        if (nSelection1Pos == std::string::npos || nSelection1Pos == 0)
902
0
        {
903
0
            CPLError(CE_Failure, CPLE_AppDefined,
904
0
                     "Invalid syntax for SENTINEL1_CALIB:");
905
0
            return nullptr;
906
0
        }
907
0
        osSelectedSubDS1 = osMDFilename.substr(nSelection1Pos + 1);
908
0
        osMDFilename.resize(nSelection1Pos);
909
910
0
        const auto nSelection2Pos = osSelectedSubDS1.find('_');
911
0
        if (nSelection2Pos != std::string::npos && nSelection2Pos != 0)
912
0
        {
913
0
            osSelectedSubDS2 = osSelectedSubDS1.substr(nSelection2Pos + 1);
914
0
            osSelectedSubDS1.resize(nSelection2Pos);
915
0
            const auto nSelection3Pos = osSelectedSubDS2.find('_');
916
0
            if (nSelection3Pos != std::string::npos && nSelection3Pos != 0)
917
0
            {
918
0
                osSelectedSubDS3 = osSelectedSubDS2.substr(nSelection3Pos + 1);
919
0
                osSelectedSubDS2.resize(nSelection3Pos);
920
0
            }
921
0
        }
922
923
        // update directory check:
924
0
        VSIStatBufL sStat;
925
0
        if (VSIStatL(osMDFilename.c_str(), &sStat) == 0)
926
0
            poOpenInfo->bIsDirectory = VSI_ISDIR(sStat.st_mode);
927
0
        if (!bCalibrated)
928
0
            osSubdatasetName = "UNCALIB";
929
0
        else if (eCalibrationType == SAFECalibratedRasterBand::SIGMA_NOUGHT)
930
0
            osSubdatasetName = "SIGMA0";
931
0
        else if (eCalibrationType == SAFECalibratedRasterBand::BETA_NOUGHT)
932
0
            osSubdatasetName = "BETA0";
933
0
        else if (eCalibrationType == SAFECalibratedRasterBand::GAMMA)
934
0
            osSubdatasetName = "GAMMA";
935
0
        osSubdatasetName += ":";
936
0
        if (!osUnitValue.empty())
937
0
        {
938
0
            osSubdatasetName += osUnitValue;
939
0
            osSubdatasetName += ":";
940
0
        }
941
0
        if (!osSelectedSubDS1.empty())
942
0
        {
943
0
            osSubdatasetName += osSelectedSubDS1;
944
0
            osSubdatasetName += ":";
945
0
        }
946
0
        if (!osSelectedSubDS2.empty())
947
0
        {
948
0
            osSubdatasetName += osSelectedSubDS2;
949
0
            osSubdatasetName += ":";
950
0
        }
951
0
        if (!osSelectedSubDS3.empty())
952
0
        {
953
0
            osSubdatasetName += osSelectedSubDS3;
954
0
            osSubdatasetName += ":";
955
0
        }
956
0
        if (!osSubdatasetName.empty())
957
0
        {
958
0
            if (osSubdatasetName.back() == ':')
959
0
                osSubdatasetName.pop_back();
960
0
        }
961
0
    }
962
0
    else
963
0
    {
964
0
        osMDFilename = poOpenInfo->pszFilename;
965
0
    }
966
967
0
    if (poOpenInfo->bIsDirectory)
968
0
    {
969
0
        osMDFilename = CPLFormCIFilenameSafe(osMDFilename.c_str(),
970
0
                                             "manifest.safe", nullptr);
971
0
    }
972
973
    /* -------------------------------------------------------------------- */
974
    /*      Ingest the manifest.safe file.                                  */
975
    /* -------------------------------------------------------------------- */
976
977
0
    auto psManifest = CPLXMLTreeCloser(CPLParseXMLFile(osMDFilename));
978
0
    if (psManifest == nullptr)
979
0
        return nullptr;
980
981
0
    CPLString osPath(CPLGetPathSafe(osMDFilename));
982
983
    /* -------------------------------------------------------------------- */
984
    /*      Confirm the requested access is supported.                      */
985
    /* -------------------------------------------------------------------- */
986
0
    if (poOpenInfo->eAccess == GA_Update)
987
0
    {
988
0
        ReportUpdateNotSupportedByDriver("SAFE");
989
0
        return nullptr;
990
0
    }
991
992
    /* -------------------------------------------------------------------- */
993
    /*      Get contentUnit parent element.                                 */
994
    /* -------------------------------------------------------------------- */
995
0
    const CPLXMLNode *psContentUnits = CPLGetXMLNode(
996
0
        psManifest.get(), "=xfdu:XFDU.informationPackageMap.xfdu:contentUnit");
997
0
    if (psContentUnits == nullptr)
998
0
    {
999
0
        CPLError(CE_Failure, CPLE_OpenFailed,
1000
0
                 "Failed to find <xfdu:XFDU><informationPackageMap>"
1001
0
                 "<xfdu:contentUnit> in manifest file.");
1002
0
        return nullptr;
1003
0
    }
1004
1005
    /* -------------------------------------------------------------------- */
1006
    /*      Get Metadata Objects element.                                   */
1007
    /* -------------------------------------------------------------------- */
1008
0
    const CPLXMLNode *psMetaDataObjects =
1009
0
        CPLGetXMLNode(psManifest.get(), "=xfdu:XFDU.metadataSection");
1010
0
    if (psMetaDataObjects == nullptr)
1011
0
    {
1012
0
        CPLError(CE_Failure, CPLE_OpenFailed,
1013
0
                 "Failed to find <xfdu:XFDU><metadataSection>"
1014
0
                 "in manifest file.");
1015
0
        return nullptr;
1016
0
    }
1017
1018
    /* -------------------------------------------------------------------- */
1019
    /*      Get Data Objects element.                                       */
1020
    /* -------------------------------------------------------------------- */
1021
0
    const CPLXMLNode *psDataObjects =
1022
0
        CPLGetXMLNode(psManifest.get(), "=xfdu:XFDU.dataObjectSection");
1023
0
    if (psDataObjects == nullptr)
1024
0
    {
1025
0
        CPLError(CE_Failure, CPLE_OpenFailed,
1026
0
                 "Failed to find <xfdu:XFDU><dataObjectSection> in document.");
1027
0
        return nullptr;
1028
0
    }
1029
1030
    /* -------------------------------------------------------------------- */
1031
    /*      Create the dataset.                                             */
1032
    /* -------------------------------------------------------------------- */
1033
0
    auto poDS = std::make_unique<SAFEDataset>();
1034
1035
0
    poDS->psManifest = std::move(psManifest);
1036
1037
    /* -------------------------------------------------------------------- */
1038
    /*      Look for "Measurement Data Unit" contentUnit elements.          */
1039
    /* -------------------------------------------------------------------- */
1040
    // Map with all measures aggregated by swath
1041
0
    std::map<CPLString, std::set<CPLString>> oMapSwaths2Pols;
1042
0
    std::vector<CPLString> oImageNumberSwPol;
1043
0
    bool isWave = false;
1044
1045
0
    for (CPLXMLNode *psContentUnit = psContentUnits->psChild;
1046
0
         psContentUnit != nullptr; psContentUnit = psContentUnit->psNext)
1047
0
    {
1048
0
        if (psContentUnit->eType != CXT_Element ||
1049
0
            !(EQUAL(psContentUnit->pszValue, "xfdu:contentUnit")))
1050
0
        {
1051
0
            continue;
1052
0
        }
1053
1054
0
        const char *pszUnitType = CPLGetXMLValue(psContentUnit, "unitType", "");
1055
1056
0
        const char *pszAnnotation = nullptr;
1057
0
        const char *pszCalibration = nullptr;
1058
0
        const char *pszMeasurement = nullptr;
1059
1060
0
        if (EQUAL(pszUnitType, "Measurement Data Unit"))
1061
0
        {
1062
            /* Get dmdID and dataObjectID */
1063
0
            const char *pszDmdID = CPLGetXMLValue(psContentUnit, "dmdID", "");
1064
0
            const char *pszDataObjectID = CPLGetXMLValue(
1065
0
                psContentUnit, "dataObjectPointer.dataObjectID", "");
1066
0
            if (*pszDataObjectID == '\0' || *pszDmdID == '\0')
1067
0
                continue;
1068
1069
0
            const CPLXMLNode *psDataObject =
1070
0
                SAFEDataset::GetDataObject(psDataObjects, pszDataObjectID);
1071
1072
0
            const char *pszRepId = CPLGetXMLValue(psDataObject, "repID", "");
1073
0
            if (!EQUAL(pszRepId, "s1Level1MeasurementSchema"))
1074
0
                continue;
1075
1076
0
            pszMeasurement = CPLGetXMLValue(psDataObject,
1077
0
                                            "byteStream.fileLocation.href", "");
1078
0
            if (*pszMeasurement == '\0')
1079
0
                continue;
1080
1081
0
            if (CPLHasPathTraversal(pszMeasurement))
1082
0
            {
1083
0
                CPLError(CE_Failure, CPLE_AppDefined,
1084
0
                         "Path traversal detected in %s", pszMeasurement);
1085
0
                return nullptr;
1086
0
            }
1087
1088
0
            const CPLStringList aosTokens(CSLTokenizeString2(
1089
0
                pszDmdID, " ",
1090
0
                CSLT_ALLOWEMPTYTOKENS | CSLT_STRIPLEADSPACES |
1091
0
                    CSLT_STRIPENDSPACES));
1092
1093
0
            for (const char *pszId : aosTokens)
1094
0
            {
1095
0
                if (*pszId == '\0')
1096
0
                    continue;
1097
1098
                // Map the metadata ID to the object element
1099
0
                const CPLXMLNode *psDO = SAFEDataset::GetDataObject(
1100
0
                    psMetaDataObjects, psDataObjects, pszId);
1101
0
                if (psDO == nullptr)
1102
0
                    continue;
1103
1104
                // check object type
1105
0
                pszRepId = CPLGetXMLValue(psDO, "repID", "");
1106
0
                if (EQUAL(pszRepId, "s1Level1ProductSchema"))
1107
0
                {
1108
                    /* Get annotation filename */
1109
0
                    pszAnnotation = CPLGetXMLValue(
1110
0
                        psDO, "byteStream.fileLocation.href", "");
1111
0
                    if (*pszAnnotation == '\0')
1112
0
                        continue;
1113
0
                }
1114
0
                else if (EQUAL(pszRepId, "s1Level1CalibrationSchema"))
1115
0
                {
1116
0
                    pszCalibration = CPLGetXMLValue(
1117
0
                        psDO, "byteStream.fileLocation.href", "");
1118
0
                    if (*pszCalibration == '\0')
1119
0
                        continue;
1120
0
                }
1121
0
                else
1122
0
                {
1123
0
                    continue;
1124
0
                }
1125
0
            }
1126
1127
0
            if (pszAnnotation == nullptr || pszCalibration == nullptr)
1128
0
                continue;
1129
1130
0
            if (CPLHasPathTraversal(pszAnnotation))
1131
0
            {
1132
0
                CPLError(CE_Failure, CPLE_AppDefined,
1133
0
                         "Path traversal detected in %s", pszAnnotation);
1134
0
                return nullptr;
1135
0
            }
1136
1137
0
            if (CPLHasPathTraversal(pszCalibration))
1138
0
            {
1139
0
                CPLError(CE_Failure, CPLE_AppDefined,
1140
0
                         "Path traversal detected in %s", pszCalibration);
1141
0
                return nullptr;
1142
0
            }
1143
1144
            // open Annotation XML file
1145
0
            const CPLString osAnnotationFilePath =
1146
0
                CPLFormFilenameSafe(osPath, pszAnnotation, nullptr);
1147
0
            const CPLString osCalibrationFilePath =
1148
0
                CPLFormFilenameSafe(osPath, pszCalibration, nullptr);
1149
1150
0
            CPLXMLTreeCloser psAnnotation(
1151
0
                CPLParseXMLFile(osAnnotationFilePath));
1152
0
            if (psAnnotation.get() == nullptr)
1153
0
                continue;
1154
0
            CPLXMLTreeCloser psCalibration(
1155
0
                CPLParseXMLFile(osCalibrationFilePath));
1156
0
            if (psCalibration.get() == nullptr)
1157
0
                continue;
1158
1159
            /* --------------------------------------------------------------------
1160
             */
1161
            /*      Get overall image information. */
1162
            /* --------------------------------------------------------------------
1163
             */
1164
0
            const CPLString osProductType = CPLGetXMLValue(
1165
0
                psAnnotation.get(), "=product.adsHeader.productType", "UNK");
1166
0
            const CPLString osMissionId = CPLGetXMLValue(
1167
0
                psAnnotation.get(), "=product.adsHeader.missionId", "UNK");
1168
0
            const CPLString osPolarization = CPLGetXMLValue(
1169
0
                psAnnotation.get(), "=product.adsHeader.polarisation", "UNK");
1170
0
            const CPLString osMode = CPLGetXMLValue(
1171
0
                psAnnotation.get(), "=product.adsHeader.mode", "UNK");
1172
0
            const CPLString osSwath = CPLGetXMLValue(
1173
0
                psAnnotation.get(), "=product.adsHeader.swath", "UNK");
1174
0
            const CPLString osImageNumber = CPLGetXMLValue(
1175
0
                psAnnotation.get(), "=product.adsHeader.imageNumber", "UNK");
1176
1177
0
            oMapSwaths2Pols[osSwath].insert(osPolarization);
1178
0
            oImageNumberSwPol.push_back(osImageNumber + " " + osSwath + " " +
1179
0
                                        osPolarization);
1180
0
            if (EQUAL(osMode.c_str(), "WV"))
1181
0
                isWave = true;
1182
1183
0
            if (EQUAL(osProductType.c_str(), "SLC"))
1184
0
                bIsSLC = true;
1185
1186
            // if the dataunit is amplitude or complex and there is calibration
1187
            // applied it's not possible as calibrated datasets are intensity.
1188
0
            if (eRequestType != INTENSITY && bCalibrated)
1189
0
                continue;
1190
1191
0
            if (osSelectedSubDS1.empty())
1192
0
            {
1193
                // If not subdataset was selected,
1194
                // open the first one we can find.
1195
0
                osSelectedSubDS1 = osSwath;
1196
0
            }
1197
1198
0
            if (osSelectedSubDS3.empty() && isWave)
1199
0
            {
1200
                // If the selected mode is Wave mode (different file structure)
1201
                // open the first vignette in the dataset.
1202
0
                osSelectedSubDS3 = osImageNumber;
1203
0
            }
1204
0
            if (!EQUAL(osSelectedSubDS1.c_str(), osSwath.c_str()))
1205
0
            {
1206
                // do not mix swath, otherwise it does not work for SLC products
1207
0
                continue;
1208
0
            }
1209
1210
0
            if (!osSelectedSubDS2.empty() &&
1211
0
                (osSelectedSubDS2.find(osPolarization) == std::string::npos))
1212
0
            {
1213
                // Add only selected polarizations.
1214
0
                continue;
1215
0
            }
1216
1217
0
            if (!osSelectedSubDS3.empty() &&
1218
0
                !EQUAL(osSelectedSubDS3.c_str(), osImageNumber.c_str()))
1219
0
            {
1220
                // Add only selected image number (for Wave)
1221
0
                continue;
1222
0
            }
1223
1224
            // Changed the position of this code segment till nullptr
1225
0
            poDS->nRasterXSize = atoi(CPLGetXMLValue(
1226
0
                psAnnotation.get(),
1227
0
                "=product.imageAnnotation.imageInformation.numberOfSamples",
1228
0
                "-1"));
1229
0
            poDS->nRasterYSize = atoi(CPLGetXMLValue(
1230
0
                psAnnotation.get(),
1231
0
                "=product.imageAnnotation.imageInformation.numberOfLines",
1232
0
                "-1"));
1233
0
            if (poDS->nRasterXSize <= 1 || poDS->nRasterYSize <= 1)
1234
0
            {
1235
0
                CPLError(
1236
0
                    CE_Failure, CPLE_OpenFailed,
1237
0
                    "Non-sane raster dimensions provided in manifest.safe. "
1238
0
                    "If this is a valid SENTINEL-1 scene, please contact your "
1239
0
                    "data provider for a corrected dataset.");
1240
0
                return nullptr;
1241
0
            }
1242
1243
0
            poDS->SetMetadataItem("PRODUCT_TYPE", osProductType.c_str());
1244
0
            poDS->SetMetadataItem("MISSION_ID", osMissionId.c_str());
1245
0
            poDS->SetMetadataItem("MODE", osMode.c_str());
1246
0
            poDS->SetMetadataItem("SWATH", osSwath.c_str());
1247
1248
            /* --------------------------------------------------------------------
1249
             */
1250
            /*      Get dataType (so we can recognize complex data), and the */
1251
            /*      bitsPerSample. */
1252
            /* --------------------------------------------------------------------
1253
             */
1254
1255
0
            const char *pszDataType = CPLGetXMLValue(
1256
0
                psAnnotation.get(),
1257
0
                "=product.imageAnnotation.imageInformation.outputPixels", "");
1258
1259
0
            GDALDataType eDataType;
1260
0
            if (EQUAL(pszDataType, "16 bit Signed Integer"))
1261
0
                eDataType = GDT_CInt16;
1262
0
            else if (EQUAL(pszDataType, "16 bit Unsigned Integer"))
1263
0
                eDataType = GDT_UInt16;
1264
0
            else
1265
0
            {
1266
0
                CPLError(CE_Failure, CPLE_AppDefined,
1267
0
                         "dataType=%s: not a supported configuration.",
1268
0
                         pszDataType);
1269
0
                return nullptr;
1270
0
            }
1271
1272
            /* Extract pixel spacing information */
1273
0
            const char *pszPixelSpacing = CPLGetXMLValue(
1274
0
                psAnnotation.get(),
1275
0
                "=product.imageAnnotation.imageInformation.rangePixelSpacing",
1276
0
                "UNK");
1277
0
            poDS->SetMetadataItem("PIXEL_SPACING", pszPixelSpacing);
1278
1279
0
            const char *pszLineSpacing = CPLGetXMLValue(
1280
0
                psAnnotation.get(),
1281
0
                "=product.imageAnnotation.imageInformation.azimuthPixelSpacing",
1282
0
                "UNK");
1283
0
            poDS->SetMetadataItem("LINE_SPACING", pszLineSpacing);
1284
1285
            /* --------------------------------------------------------------------
1286
             */
1287
            /*      Form full filename (path of manifest.safe + measurement
1288
             * file).  */
1289
            /* --------------------------------------------------------------------
1290
             */
1291
0
            const std::string osFullFilename(
1292
0
                CPLFormFilenameSafe(osPath, pszMeasurement, nullptr));
1293
1294
            /* --------------------------------------------------------------------
1295
             */
1296
            /*      Try and open the file. */
1297
            /* --------------------------------------------------------------------
1298
             */
1299
0
            std::unique_ptr<GDALDataset> poBandFile;
1300
0
            {
1301
0
                CPLTurnFailureIntoWarningBackuper oErrorsToWarnings{};
1302
0
                poBandFile = std::unique_ptr<GDALDataset>(
1303
0
                    GDALDataset::Open(osFullFilename.c_str(),
1304
0
                                      GDAL_OF_RASTER | GDAL_OF_VERBOSE_ERROR));
1305
0
            }
1306
1307
0
            if (poBandFile == nullptr)
1308
0
            {
1309
                // NOP
1310
0
            }
1311
0
            else if (poBandFile->GetRasterCount() == 0)
1312
0
            {
1313
0
                poBandFile.reset();
1314
0
            }
1315
0
            else
1316
0
            {
1317
0
                poDS->papszExtraFiles =
1318
0
                    CSLAddString(poDS->papszExtraFiles, osAnnotationFilePath);
1319
0
                poDS->papszExtraFiles =
1320
0
                    CSLAddString(poDS->papszExtraFiles, osCalibrationFilePath);
1321
0
                poDS->papszExtraFiles =
1322
0
                    CSLAddString(poDS->papszExtraFiles, osFullFilename.c_str());
1323
                /* --------------------------------------------------------------------
1324
                 */
1325
                /*      Collect Annotation Processing Information */
1326
                /* --------------------------------------------------------------------
1327
                 */
1328
0
                CPLXMLNode *psProcessingInfo = CPLGetXMLNode(
1329
0
                    psAnnotation.get(),
1330
0
                    "=product.imageAnnotation.processingInformation");
1331
1332
0
                if (psProcessingInfo != nullptr)
1333
0
                {
1334
0
                    OGRSpatialReference oLL, oPrj;
1335
1336
0
                    const char *pszEllipsoidName =
1337
0
                        CPLGetXMLValue(psProcessingInfo, "ellipsoidName", "");
1338
0
                    const double minor_axis = CPLAtof(CPLGetXMLValue(
1339
0
                        psProcessingInfo, "ellipsoidSemiMinorAxis", "0.0"));
1340
0
                    const double major_axis = CPLAtof(CPLGetXMLValue(
1341
0
                        psProcessingInfo, "ellipsoidSemiMajorAxis", "0.0"));
1342
1343
0
                    if (EQUAL(pszEllipsoidName, "") || (minor_axis == 0.0) ||
1344
0
                        (major_axis == 0.0))
1345
0
                    {
1346
0
                        CPLError(CE_Warning, CPLE_AppDefined,
1347
0
                                 "Warning- incomplete"
1348
0
                                 " ellipsoid information.  Using wgs-84 "
1349
0
                                 "parameters.\n");
1350
0
                        oLL.SetWellKnownGeogCS("WGS84");
1351
0
                        oPrj.SetWellKnownGeogCS("WGS84");
1352
0
                    }
1353
0
                    else if (EQUAL(pszEllipsoidName, "WGS84"))
1354
0
                    {
1355
0
                        oLL.SetWellKnownGeogCS("WGS84");
1356
0
                        oPrj.SetWellKnownGeogCS("WGS84");
1357
0
                    }
1358
0
                    else
1359
0
                    {
1360
0
                        const double inv_flattening =
1361
0
                            major_axis / (major_axis - minor_axis);
1362
0
                        oLL.SetGeogCS("", "", pszEllipsoidName, major_axis,
1363
0
                                      inv_flattening);
1364
0
                        oPrj.SetGeogCS("", "", pszEllipsoidName, major_axis,
1365
0
                                       inv_flattening);
1366
0
                    }
1367
1368
0
                    oLL.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
1369
0
                    poDS->m_oGCPSRS = std::move(oLL);
1370
0
                }
1371
1372
                /* --------------------------------------------------------------------
1373
                 */
1374
                /*      Collect GCPs. */
1375
                /* --------------------------------------------------------------------
1376
                 */
1377
0
                CPLXMLNode *psGeoGrid = CPLGetXMLNode(
1378
0
                    psAnnotation.get(),
1379
0
                    "=product.geolocationGrid.geolocationGridPointList");
1380
1381
0
                if (psGeoGrid != nullptr)
1382
0
                {
1383
0
                    if (poDS->nGCPCount > 0)
1384
0
                    {
1385
0
                        GDALDeinitGCPs(poDS->nGCPCount, poDS->pasGCPList);
1386
0
                        CPLFree(poDS->pasGCPList);
1387
0
                    }
1388
1389
                    /* count GCPs */
1390
0
                    poDS->nGCPCount = 0;
1391
1392
0
                    for (CPLXMLNode *psNode = psGeoGrid->psChild;
1393
0
                         psNode != nullptr; psNode = psNode->psNext)
1394
0
                    {
1395
0
                        if (EQUAL(psNode->pszValue, "geolocationGridPoint"))
1396
0
                            poDS->nGCPCount++;
1397
0
                    }
1398
1399
0
                    if (poDS->nGCPCount > 0)
1400
0
                    {
1401
0
                        poDS->pasGCPList = static_cast<GDAL_GCP *>(
1402
0
                            CPLCalloc(sizeof(GDAL_GCP), poDS->nGCPCount));
1403
1404
0
                        poDS->nGCPCount = 0;
1405
1406
0
                        for (CPLXMLNode *psNode = psGeoGrid->psChild;
1407
0
                             psNode != nullptr; psNode = psNode->psNext)
1408
0
                        {
1409
0
                            GDAL_GCP *psGCP =
1410
0
                                poDS->pasGCPList + poDS->nGCPCount;
1411
1412
0
                            if (!EQUAL(psNode->pszValue,
1413
0
                                       "geolocationGridPoint"))
1414
0
                                continue;
1415
1416
0
                            poDS->nGCPCount++;
1417
1418
0
                            char szID[32];
1419
0
                            snprintf(szID, sizeof(szID), "%d", poDS->nGCPCount);
1420
0
                            psGCP->pszId = CPLStrdup(szID);
1421
0
                            psGCP->pszInfo = CPLStrdup("");
1422
0
                            psGCP->dfGCPPixel =
1423
0
                                CPLAtof(CPLGetXMLValue(psNode, "pixel", "0"));
1424
0
                            psGCP->dfGCPLine =
1425
0
                                CPLAtof(CPLGetXMLValue(psNode, "line", "0"));
1426
0
                            psGCP->dfGCPX = CPLAtof(
1427
0
                                CPLGetXMLValue(psNode, "longitude", ""));
1428
0
                            psGCP->dfGCPY =
1429
0
                                CPLAtof(CPLGetXMLValue(psNode, "latitude", ""));
1430
0
                            psGCP->dfGCPZ =
1431
0
                                CPLAtof(CPLGetXMLValue(psNode, "height", ""));
1432
0
                        }
1433
0
                    }
1434
0
                }
1435
1436
                // Create bands
1437
0
                if (EQUAL(osProductType.c_str(), "SLC"))
1438
0
                {
1439
                    // only add bands if no subDS or uncalibrated subdataset
1440
                    // with complex data. (Calibrated will always be intensity
1441
                    // only)
1442
0
                    if (!bCalibrated &&
1443
0
                        (eRequestType == UNKNOWN || eRequestType == COMPLEX))
1444
0
                    {
1445
0
                        poBandFile->MarkAsShared();
1446
0
                        SAFESLCRasterBand *poBand0 = new SAFESLCRasterBand(
1447
0
                            poDS.get(), eDataType, osSwath, osPolarization,
1448
0
                            std::move(poBandFile), SAFESLCRasterBand::COMPLEX);
1449
0
                        poDS->SetBand(poDS->GetRasterCount() + 1, poBand0);
1450
0
                    }
1451
0
                    else if (eRequestType == INTENSITY)  // Intensity
1452
0
                    {
1453
0
                        SAFESLCRasterBand *poBand1 = new SAFESLCRasterBand(
1454
0
                            poDS.get(), eDataType, osSwath, osPolarization,
1455
0
                            std::move(poBandFile),
1456
0
                            SAFESLCRasterBand::INTENSITY);
1457
0
                        poDS->SetBand(poDS->GetRasterCount() + 1, poBand1);
1458
0
                    }
1459
0
                }
1460
0
                else if (!bCalibrated &&
1461
0
                         (eRequestType == UNKNOWN || eRequestType == AMPLITUDE))
1462
0
                {
1463
0
                    SAFERasterBand *poBand = new SAFERasterBand(
1464
0
                        poDS.get(), eDataType, osSwath, osPolarization,
1465
0
                        std::move(poBandFile));
1466
0
                    poDS->SetBand(poDS->GetRasterCount() + 1, poBand);
1467
0
                }
1468
0
                else if (bCalibrated &&
1469
0
                         (eRequestType == UNKNOWN || eRequestType == COMPLEX))
1470
0
                {
1471
0
                    auto poBand = std::make_unique<SAFECalibratedRasterBand>(
1472
0
                        poDS.get(), eDataType, osSwath, osPolarization,
1473
0
                        std::move(poBandFile), osCalibrationFilePath,
1474
0
                        eCalibrationType);
1475
0
                    if (!poBand->ReadLUT())
1476
0
                    {
1477
0
                        CPLError(CE_Failure, CPLE_OpenFailed,
1478
0
                                 "Reading calibration LUT(s) failed: %s.",
1479
0
                                 osCalibrationFilePath.c_str());
1480
0
                        return nullptr;
1481
0
                    }
1482
0
                    poDS->SetBand(poDS->GetRasterCount() + 1,
1483
0
                                  std::move(poBand));
1484
0
                }
1485
0
            }
1486
0
        }
1487
0
    }
1488
1489
    // loop through all Swath/pols to add subdatasets
1490
0
    if (!bIsSubDS)
1491
0
    {
1492
0
        const CPLString aosCalibrationValues[4] = {"SIGMA0", "BETA0", "GAMMA",
1493
0
                                                   "UNCALIB"};
1494
0
        const CPLString aosDataUnitValues[3] = {"AMPLITUDE", "COMPLEX",
1495
0
                                                "INTENSITY"};
1496
0
        if (!isWave)
1497
0
        {
1498
0
            for (const auto &iterSwath : oMapSwaths2Pols)
1499
0
            {
1500
0
                CPLString osSubDS1 = iterSwath.first;
1501
0
                CPLString osSubDS2;
1502
1503
0
                for (const auto &pol : iterSwath.second)
1504
0
                {
1505
0
                    if (!osSubDS2.empty())
1506
0
                        osSubDS2 += "+";
1507
0
                    osSubDS2 += pol;
1508
                    // Create single band or multiband complex SubDataset
1509
0
                    int i = 0;
1510
0
                    if (bIsSLC)
1511
0
                    {
1512
0
                        for (i = 0; i < 3; i++)
1513
0
                        {
1514
0
                            CPLString osCalibTemp = aosCalibrationValues[i];
1515
0
                            poDS->AddSubDataset(
1516
0
                                CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s_%s:%s",
1517
0
                                           osCalibTemp.c_str(),
1518
0
                                           osMDFilename.c_str(),
1519
0
                                           osSubDS1.c_str(), pol.c_str(),
1520
0
                                           aosDataUnitValues[2].c_str()),
1521
0
                                CPLSPrintf("Single band with %s swath and %s "
1522
0
                                           "polarization and %s calibration",
1523
0
                                           osSubDS1.c_str(), pol.c_str(),
1524
0
                                           osCalibTemp.c_str()));
1525
0
                        }
1526
1527
0
                        CPLString osCalibTemp = aosCalibrationValues[i];
1528
0
                        poDS->AddSubDataset(
1529
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s_%s:%s",
1530
0
                                       osCalibTemp.c_str(),
1531
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1532
0
                                       pol.c_str(),
1533
0
                                       aosDataUnitValues[1].c_str()),
1534
0
                            CPLSPrintf("Single band with %s swath and %s "
1535
0
                                       "polarization and %s calibration",
1536
0
                                       osSubDS1.c_str(), pol.c_str(),
1537
0
                                       osCalibTemp.c_str()));
1538
0
                        poDS->AddSubDataset(
1539
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s_%s:%s",
1540
0
                                       osCalibTemp.c_str(),
1541
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1542
0
                                       pol.c_str(),
1543
0
                                       aosDataUnitValues[2].c_str()),
1544
0
                            CPLSPrintf("Single band with %s swath and %s "
1545
0
                                       "polarization and %s calibration",
1546
0
                                       osSubDS1.c_str(), pol.c_str(),
1547
0
                                       osCalibTemp.c_str()));
1548
0
                    }
1549
0
                    else
1550
0
                    {
1551
0
                        i = 3;
1552
0
                        CPLString osCalibTemp = aosCalibrationValues[i];
1553
1554
0
                        poDS->AddSubDataset(
1555
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s_%s:%s",
1556
0
                                       osCalibTemp.c_str(),
1557
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1558
0
                                       pol.c_str(),
1559
0
                                       aosDataUnitValues[0].c_str()),
1560
0
                            CPLSPrintf("Single band with %s swath and %s "
1561
0
                                       "polarization and %s calibration",
1562
0
                                       osSubDS1.c_str(), pol.c_str(),
1563
0
                                       osCalibTemp.c_str()));
1564
0
                    }
1565
0
                }
1566
1567
0
                if (iterSwath.second.size() > 1)
1568
0
                {
1569
                    // Create single band subdataset with all polarizations
1570
0
                    int i = 0;
1571
0
                    if (bIsSLC)
1572
0
                    {
1573
0
                        for (i = 0; i < 3; i++)
1574
0
                        {
1575
0
                            CPLString osCalibTemp = aosCalibrationValues[i];
1576
1577
0
                            poDS->AddSubDataset(
1578
0
                                CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s:%s",
1579
0
                                           osCalibTemp.c_str(),
1580
0
                                           osMDFilename.c_str(),
1581
0
                                           osSubDS1.c_str(),
1582
0
                                           aosDataUnitValues[2].c_str()),
1583
0
                                CPLSPrintf(
1584
0
                                    "%s swath with all polarizations (%s) as "
1585
0
                                    "bands and %s calibration",
1586
0
                                    osSubDS1.c_str(), osSubDS2.c_str(),
1587
0
                                    osCalibTemp.c_str()));
1588
0
                        }
1589
1590
0
                        CPLString osCalibTemp = aosCalibrationValues[i];
1591
0
                        poDS->AddSubDataset(
1592
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s:%s",
1593
0
                                       osCalibTemp.c_str(),
1594
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1595
0
                                       aosDataUnitValues[1].c_str()),
1596
0
                            CPLSPrintf(
1597
0
                                "%s swath with all polarizations (%s) as "
1598
0
                                "bands and %s calibration",
1599
0
                                osSubDS1.c_str(), osSubDS2.c_str(),
1600
0
                                osCalibTemp.c_str()));
1601
0
                        poDS->AddSubDataset(
1602
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s:%s",
1603
0
                                       osCalibTemp.c_str(),
1604
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1605
0
                                       aosDataUnitValues[2].c_str()),
1606
0
                            CPLSPrintf(
1607
0
                                "%s swath with all polarizations (%s) as "
1608
0
                                "bands and %s calibration",
1609
0
                                osSubDS1.c_str(), osSubDS2.c_str(),
1610
0
                                osCalibTemp.c_str()));
1611
0
                    }
1612
0
                    else
1613
0
                    {
1614
0
                        i = 3;
1615
0
                        CPLString osCalibTemp = aosCalibrationValues[i];
1616
0
                        poDS->AddSubDataset(
1617
0
                            CPLSPrintf("SENTINEL1_CALIB:%s:%s:%s:%s",
1618
0
                                       osCalibTemp.c_str(),
1619
0
                                       osMDFilename.c_str(), osSubDS1.c_str(),
1620
0
                                       aosDataUnitValues[0].c_str()),
1621
0
                            CPLSPrintf(
1622
0
                                "%s swath with all polarizations (%s) as "
1623
0
                                "bands and %s calibration",
1624
0
                                osSubDS1.c_str(), osSubDS2.c_str(),
1625
0
                                osCalibTemp.c_str()));
1626
0
                    }
1627
0
                }
1628
0
            }
1629
0
        }
1630
0
        else
1631
0
        {
1632
0
            for (const CPLString &osImgSwPol : oImageNumberSwPol)
1633
0
            {
1634
0
                CPLString osImgSwPolTmp = osImgSwPol;
1635
0
                const char *pszImage = strchr(osImgSwPolTmp.c_str(), ' ');
1636
0
                CPLString osImage, osSwath, osPolarization;
1637
0
                if (pszImage != nullptr)
1638
0
                {
1639
0
                    osImage = osImgSwPolTmp;
1640
0
                    osImage.resize(pszImage - osImgSwPolTmp.c_str());
1641
0
                    osImgSwPolTmp = pszImage + strlen(" ");
1642
0
                    const char *pszSwath = strchr(osImgSwPolTmp.c_str(), ' ');
1643
0
                    if (pszSwath != nullptr)
1644
0
                    {
1645
0
                        osSwath = osImgSwPolTmp;
1646
0
                        osSwath.resize(pszSwath - osImgSwPolTmp.c_str());
1647
0
                        osPolarization = pszSwath + strlen(" ");
1648
0
                        int i = 0;
1649
1650
0
                        if (bIsSLC)
1651
0
                        {
1652
0
                            for (i = 0; i < 3; i++)
1653
0
                            {
1654
0
                                CPLString osCalibTemp = aosCalibrationValues[i];
1655
1656
0
                                poDS->AddSubDataset(
1657
0
                                    CPLSPrintf(
1658
0
                                        "SENTINEL1_CALIB:%s:%s:%s_%s_%s:%s",
1659
0
                                        osCalibTemp.c_str(),
1660
0
                                        osMDFilename.c_str(), osSwath.c_str(),
1661
0
                                        osPolarization.c_str(), osImage.c_str(),
1662
0
                                        aosDataUnitValues[2].c_str()),
1663
0
                                    CPLSPrintf(
1664
0
                                        "Single band with %s swath and %s "
1665
0
                                        "polarization and %s calibration",
1666
0
                                        osSwath.c_str(), osPolarization.c_str(),
1667
0
                                        osCalibTemp.c_str()));
1668
0
                            }
1669
1670
0
                            CPLString osCalibTemp = aosCalibrationValues[i];
1671
1672
0
                            poDS->AddSubDataset(
1673
0
                                CPLSPrintf(
1674
0
                                    "SENTINEL1_CALIB:%s:%s:%s_%s_%s:%s",
1675
0
                                    osCalibTemp.c_str(), osMDFilename.c_str(),
1676
0
                                    osSwath.c_str(), osPolarization.c_str(),
1677
0
                                    osImage.c_str(),
1678
0
                                    aosDataUnitValues[1].c_str()),
1679
0
                                CPLSPrintf("Single band with %s swath and %s "
1680
0
                                           "polarization and %s calibration",
1681
0
                                           osSwath.c_str(),
1682
0
                                           osPolarization.c_str(),
1683
0
                                           osCalibTemp.c_str()));
1684
1685
0
                            poDS->AddSubDataset(
1686
0
                                CPLSPrintf(
1687
0
                                    "SENTINEL1_CALIB:%s:%s:%s_%s_%s:%s",
1688
0
                                    osCalibTemp.c_str(), osMDFilename.c_str(),
1689
0
                                    osSwath.c_str(), osPolarization.c_str(),
1690
0
                                    osImage.c_str(),
1691
0
                                    aosDataUnitValues[2].c_str()),
1692
0
                                CPLSPrintf("Single band with %s swath and %s "
1693
0
                                           "polarization and %s calibration",
1694
0
                                           osSwath.c_str(),
1695
0
                                           osPolarization.c_str(),
1696
0
                                           osCalibTemp.c_str()));
1697
0
                        }
1698
0
                        else
1699
0
                        {
1700
0
                            i = 3;
1701
0
                            CPLString osCalibTemp = aosCalibrationValues[i];
1702
1703
0
                            poDS->AddSubDataset(
1704
0
                                CPLSPrintf(
1705
0
                                    "SENTINEL1_CALIB:%s:%s:%s_%s_%s:%s",
1706
0
                                    osCalibTemp.c_str(), osMDFilename.c_str(),
1707
0
                                    osSwath.c_str(), osPolarization.c_str(),
1708
0
                                    osImage.c_str(),
1709
0
                                    aosDataUnitValues[0].c_str()),
1710
0
                                CPLSPrintf("Single band with %s swath and %s "
1711
0
                                           "polarization and %s calibration",
1712
0
                                           osSwath.c_str(),
1713
0
                                           osPolarization.c_str(),
1714
0
                                           osCalibTemp.c_str()));
1715
0
                        }
1716
0
                    }
1717
0
                }
1718
0
            }
1719
0
        }
1720
0
    }
1721
0
    if (poDS->GetRasterCount() == 0)
1722
0
    {
1723
0
        CPLError(CE_Failure, CPLE_OpenFailed, "Measurement bands not found.");
1724
0
        return nullptr;
1725
0
    }
1726
1727
    /* -------------------------------------------------------------------- */
1728
    /*      Collect more metadata elements                                  */
1729
    /* -------------------------------------------------------------------- */
1730
1731
    /* -------------------------------------------------------------------- */
1732
    /*      Platform information                                            */
1733
    /* -------------------------------------------------------------------- */
1734
0
    const CPLXMLNode *psPlatformAttrs =
1735
0
        SAFEDataset::GetMetaDataObject(psMetaDataObjects, "platform");
1736
1737
0
    if (psPlatformAttrs != nullptr)
1738
0
    {
1739
0
        const char *pszItem =
1740
0
            CPLGetXMLValue(psPlatformAttrs,
1741
0
                           "metadataWrap.xmlData.safe:platform"
1742
0
                           ".safe:familyName",
1743
0
                           "");
1744
0
        poDS->SetMetadataItem("SATELLITE_IDENTIFIER", pszItem);
1745
1746
0
        pszItem =
1747
0
            CPLGetXMLValue(psPlatformAttrs,
1748
0
                           "metadataWrap.xmlData.safe:platform"
1749
0
                           ".safe:instrument.safe:familyName.abbreviation",
1750
0
                           "");
1751
0
        poDS->SetMetadataItem("SENSOR_IDENTIFIER", pszItem);
1752
1753
0
        pszItem = CPLGetXMLValue(psPlatformAttrs,
1754
0
                                 "metadataWrap.xmlData.safe:platform"
1755
0
                                 ".safe:instrument.safe:extension"
1756
0
                                 ".s1sarl1:instrumentMode.s1sarl1:mode",
1757
0
                                 "UNK");
1758
0
        poDS->SetMetadataItem("BEAM_MODE", pszItem);
1759
1760
0
        pszItem = CPLGetXMLValue(psPlatformAttrs,
1761
0
                                 "metadataWrap.xmlData.safe:platform"
1762
0
                                 ".safe:instrument.safe:extension"
1763
0
                                 ".s1sarl1:instrumentMode.s1sarl1:swath",
1764
0
                                 "UNK");
1765
0
        poDS->SetMetadataItem("BEAM_SWATH", pszItem);
1766
0
    }
1767
1768
    /* -------------------------------------------------------------------- */
1769
    /*      Acquisition Period information                                  */
1770
    /* -------------------------------------------------------------------- */
1771
0
    const CPLXMLNode *psAcquisitionAttrs =
1772
0
        SAFEDataset::GetMetaDataObject(psMetaDataObjects, "acquisitionPeriod");
1773
1774
0
    if (psAcquisitionAttrs != nullptr)
1775
0
    {
1776
0
        const char *pszItem =
1777
0
            CPLGetXMLValue(psAcquisitionAttrs,
1778
0
                           "metadataWrap.xmlData.safe:acquisitionPeriod"
1779
0
                           ".safe:startTime",
1780
0
                           "UNK");
1781
0
        poDS->SetMetadataItem("ACQUISITION_START_TIME", pszItem);
1782
0
        pszItem = CPLGetXMLValue(psAcquisitionAttrs,
1783
0
                                 "metadataWrap.xmlData.safe:acquisitionPeriod"
1784
0
                                 ".safe:stopTime",
1785
0
                                 "UNK");
1786
0
        poDS->SetMetadataItem("ACQUISITION_STOP_TIME", pszItem);
1787
0
    }
1788
1789
    /* -------------------------------------------------------------------- */
1790
    /*      Processing information                                          */
1791
    /* -------------------------------------------------------------------- */
1792
0
    const CPLXMLNode *psProcessingAttrs =
1793
0
        SAFEDataset::GetMetaDataObject(psMetaDataObjects, "processing");
1794
1795
0
    if (psProcessingAttrs != nullptr)
1796
0
    {
1797
0
        const char *pszItem = CPLGetXMLValue(
1798
0
            psProcessingAttrs,
1799
0
            "metadataWrap.xmlData.safe:processing.safe:facility.name", "UNK");
1800
0
        poDS->SetMetadataItem("FACILITY_IDENTIFIER", pszItem);
1801
0
    }
1802
1803
    /* -------------------------------------------------------------------- */
1804
    /*      Measurement Orbit Reference information                         */
1805
    /* -------------------------------------------------------------------- */
1806
0
    const CPLXMLNode *psOrbitAttrs = SAFEDataset::GetMetaDataObject(
1807
0
        psMetaDataObjects, "measurementOrbitReference");
1808
1809
0
    if (psOrbitAttrs != nullptr)
1810
0
    {
1811
0
        const char *pszItem =
1812
0
            CPLGetXMLValue(psOrbitAttrs,
1813
0
                           "metadataWrap.xmlData.safe:orbitReference"
1814
0
                           ".safe:orbitNumber",
1815
0
                           "UNK");
1816
0
        poDS->SetMetadataItem("ORBIT_NUMBER", pszItem);
1817
0
        pszItem = CPLGetXMLValue(psOrbitAttrs,
1818
0
                                 "metadataWrap.xmlData.safe:orbitReference"
1819
0
                                 ".safe:extension.s1:orbitProperties.s1:pass",
1820
0
                                 "UNK");
1821
0
        poDS->SetMetadataItem("ORBIT_DIRECTION", pszItem);
1822
0
    }
1823
1824
    /* -------------------------------------------------------------------- */
1825
    /*      Footprint                                                       */
1826
    /* -------------------------------------------------------------------- */
1827
0
    const CPLXMLNode *psFrameSet = SAFEDataset::GetMetaDataObject(
1828
0
        psMetaDataObjects, "measurementFrameSet");
1829
1830
0
    if (psFrameSet)
1831
0
    {
1832
0
        const auto psFootPrint = CPLGetXMLNode(psFrameSet, "metadataWrap."
1833
0
                                                           "xmlData."
1834
0
                                                           "safe:frameSet."
1835
0
                                                           "safe:frame."
1836
0
                                                           "safe:footPrint");
1837
0
        if (psFootPrint)
1838
0
        {
1839
0
            const char *pszSRSName =
1840
0
                CPLGetXMLValue(psFootPrint, "srsName", nullptr);
1841
0
            const char *pszCoordinates =
1842
0
                CPLGetXMLValue(psFootPrint, "gml:coordinates", nullptr);
1843
0
            if (pszSRSName &&
1844
0
                EQUAL(pszSRSName,
1845
0
                      "http://www.opengis.net/gml/srs/epsg.xml#4326") &&
1846
0
                pszCoordinates)
1847
0
            {
1848
0
                const CPLStringList aosValues(
1849
0
                    CSLTokenizeString2(pszCoordinates, " ,", 0));
1850
0
                if (aosValues.size() == 8)
1851
0
                {
1852
0
                    poDS->SetMetadataItem(
1853
0
                        "FOOTPRINT",
1854
0
                        CPLSPrintf("POLYGON((%s %s,%s %s,%s %s,%s %s, %s %s))",
1855
0
                                   aosValues[1], aosValues[0], aosValues[3],
1856
0
                                   aosValues[2], aosValues[5], aosValues[4],
1857
0
                                   aosValues[7], aosValues[6], aosValues[1],
1858
0
                                   aosValues[0]));
1859
0
                }
1860
0
            }
1861
0
        }
1862
0
    }
1863
1864
    /* -------------------------------------------------------------------- */
1865
    /*      Initialize any PAM information.                                 */
1866
    /* -------------------------------------------------------------------- */
1867
0
    const CPLString osDescription = osMDFilename;
1868
1869
    /* -------------------------------------------------------------------- */
1870
    /*      Initialize any PAM information.                                 */
1871
    /* -------------------------------------------------------------------- */
1872
0
    poDS->SetDescription(osDescription);
1873
1874
0
    poDS->SetPhysicalFilename(osMDFilename);
1875
0
    if (!osSubdatasetName.empty())
1876
0
    {
1877
0
        poDS->SetDescription(poOpenInfo->pszFilename);
1878
0
        poDS->SetSubdatasetName(osSubdatasetName);
1879
0
    }
1880
1881
0
    poDS->TryLoadXML();
1882
1883
    /* -------------------------------------------------------------------- */
1884
    /*      Check for overviews.                                            */
1885
    /* -------------------------------------------------------------------- */
1886
0
    poDS->oOvManager.Initialize(poDS.get(), ":::VIRTUAL:::");
1887
1888
0
    return poDS.release();
1889
0
}
1890
1891
/************************************************************************/
1892
/*                           AddSubDataset()                            */
1893
/************************************************************************/
1894
void SAFEDataset::AddSubDataset(const CPLString &osName,
1895
                                const CPLString &osDesc)
1896
0
{
1897
0
    ++m_nSubDSNum;
1898
0
    papszSubDatasets = CSLAddNameValue(
1899
0
        papszSubDatasets, CPLSPrintf("SUBDATASET_%d_NAME", m_nSubDSNum),
1900
0
        osName.c_str());
1901
0
    papszSubDatasets = CSLAddNameValue(
1902
0
        papszSubDatasets, CPLSPrintf("SUBDATASET_%d_DESC", m_nSubDSNum),
1903
0
        osDesc.c_str());
1904
0
}
1905
1906
/************************************************************************/
1907
/*                            GetGCPCount()                             */
1908
/************************************************************************/
1909
1910
int SAFEDataset::GetGCPCount()
1911
0
{
1912
0
    return nGCPCount;
1913
0
}
1914
1915
/************************************************************************/
1916
/*                          GetGCPSpatialRef()                          */
1917
/************************************************************************/
1918
1919
const OGRSpatialReference *SAFEDataset::GetGCPSpatialRef() const
1920
1921
0
{
1922
0
    return m_oGCPSRS.IsEmpty() ? nullptr : &m_oGCPSRS;
1923
0
}
1924
1925
/************************************************************************/
1926
/*                              GetGCPs()                               */
1927
/************************************************************************/
1928
1929
const GDAL_GCP *SAFEDataset::GetGCPs()
1930
1931
0
{
1932
0
    return pasGCPList;
1933
0
}
1934
1935
/************************************************************************/
1936
/*                       GetMetadataDomainList()                        */
1937
/************************************************************************/
1938
1939
char **SAFEDataset::GetMetadataDomainList()
1940
0
{
1941
0
    return BuildMetadataDomainList(GDALDataset::GetMetadataDomainList(), TRUE,
1942
0
                                   "SUBDATASETS", nullptr);
1943
0
}
1944
1945
/************************************************************************/
1946
/*                            GetMetadata()                             */
1947
/************************************************************************/
1948
1949
CSLConstList SAFEDataset::GetMetadata(const char *pszDomain)
1950
0
{
1951
0
    if (pszDomain != nullptr && STARTS_WITH_CI(pszDomain, "SUBDATASETS"))
1952
0
        return papszSubDatasets;
1953
1954
0
    return GDALDataset::GetMetadata(pszDomain);
1955
0
}
1956
1957
/************************************************************************/
1958
/*                         GDALRegister_SAFE()                          */
1959
/************************************************************************/
1960
1961
void GDALRegister_SAFE()
1962
22
{
1963
22
    if (GDALGetDriverByName("SAFE") != nullptr)
1964
0
        return;
1965
1966
22
    GDALDriver *poDriver = new GDALDriver();
1967
1968
22
    poDriver->SetDescription("SAFE");
1969
22
    poDriver->SetMetadataItem(GDAL_DCAP_RASTER, "YES");
1970
22
    poDriver->SetMetadataItem(GDAL_DCAP_VIRTUALIO, "YES");
1971
22
    poDriver->SetMetadataItem(GDAL_DMD_LONGNAME, "Sentinel-1 SAR SAFE Product");
1972
22
    poDriver->SetMetadataItem(GDAL_DMD_HELPTOPIC, "drivers/raster/safe.html");
1973
1974
22
    poDriver->pfnOpen = SAFEDataset::Open;
1975
22
    poDriver->pfnIdentify = SAFEDataset::Identify;
1976
1977
22
    GetGDALDriverManager()->RegisterDriver(poDriver);
1978
22
}
1979
1980
/************************************************************************/
1981
/*                   Azimuth time handling functions                    */
1982
/************************************************************************/
1983
1984
SAFECalibratedRasterBand::TimePoint
1985
SAFECalibratedRasterBand::getTimePoint(const char *pszTime)
1986
0
{
1987
0
    int nYear, nMonth, nDay, nHour, nMinute, nSeconds;
1988
0
    long nMicroSeconds;
1989
0
    sscanf(pszTime, "%d-%d-%dT%d:%d:%d.%ld", &nYear, &nMonth, &nDay, &nHour,
1990
0
           &nMinute, &nSeconds, &nMicroSeconds);
1991
1992
0
    struct tm oTm;
1993
0
    oTm.tm_sec = nSeconds;
1994
0
    oTm.tm_min = nMinute;
1995
0
    oTm.tm_hour = nHour;
1996
0
    oTm.tm_mday = nDay;
1997
0
    oTm.tm_mon = nMonth - 1;
1998
0
    oTm.tm_year = nYear - 1900;
1999
0
    oTm.tm_isdst = -1;
2000
2001
0
    std::time_t oTt = static_cast<std::time_t>(CPLYMDHMSToUnixTime(&oTm));
2002
2003
0
    TimePoint oTp = std::chrono::system_clock::from_time_t(oTt);
2004
0
    TimePoint oTp1 = oTp + std::chrono::microseconds(nMicroSeconds);
2005
0
    return oTp1;
2006
0
}
2007
2008
double SAFECalibratedRasterBand::getTimeDiff(TimePoint oT1, TimePoint oT2)
2009
0
{
2010
0
    std::chrono::duration<double> oResult = (oT2 - oT1);
2011
0
    return oResult.count();
2012
0
}
2013
2014
SAFECalibratedRasterBand::TimePoint
2015
SAFECalibratedRasterBand::getazTime(TimePoint oStart, TimePoint oStop,
2016
                                    long nNumOfLines, int nOffset)
2017
0
{
2018
0
    double dfTemp = getTimeDiff(oStart, oStop);
2019
0
    dfTemp /= (nNumOfLines - 1);
2020
0
    unsigned long nTimeDiffMicro = static_cast<unsigned long>(dfTemp * 1000000);
2021
0
    TimePoint oResult =
2022
0
        oStart + (nOffset * std::chrono::microseconds(nTimeDiffMicro));
2023
0
    return oResult;
2024
0
}
2025
2026
/************************************************************************/
2027
/*               Utility functions used in interpolation                */
2028
/************************************************************************/
2029
2030
int SAFECalibratedRasterBand::getCalibrationVectorIndex(int nLineNo)
2031
0
{
2032
0
    for (size_t i = 1; i < m_anLineLUT.size(); i++)
2033
0
    {
2034
0
        if (nLineNo < m_anLineLUT[i])
2035
0
            return static_cast<int>(i - 1);
2036
0
    }
2037
0
    return 0;
2038
0
}
2039
2040
int SAFECalibratedRasterBand::getPixelIndex(int nPixelNo)
2041
0
{
2042
0
    for (int i = 1; i < m_nNumPixels; i++)
2043
0
    {
2044
0
        if (nPixelNo < m_anPixelLUT[i])
2045
0
            return i - 1;
2046
0
    }
2047
0
    return 0;
2048
0
}