Coverage Report

Created: 2025-11-15 08:43

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/gdal/frmts/mrf/mrf_band.cpp
Line
Count
Source
1
/*
2
 * Copyright (c) 2002-2012, California Institute of Technology.
3
 * All rights reserved.  Based on Government Sponsored Research under contracts
4
 * NAS7-1407 and/or NAS7-03001.
5
 *
6
 * Redistribution and use in source and binary forms, with or without
7
 * modification, are permitted provided that the following conditions are met:
8
 *   1. Redistributions of source code must retain the above copyright notice,
9
 * this list of conditions and the following disclaimer.
10
 *   2. Redistributions in binary form must reproduce the above copyright
11
 * notice, this list of conditions and the following disclaimer in the
12
 * documentation and/or other materials provided with the distribution.
13
 *   3. Neither the name of the California Institute of Technology (Caltech),
14
 * its operating division the Jet Propulsion Laboratory (JPL), the National
15
 * Aeronautics and Space Administration (NASA), nor the names of its
16
 * contributors may be used to endorse or promote products derived from this
17
 * software without specific prior written permission.
18
 *
19
 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
20
 * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
21
 * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
22
 * ARE DISCLAIMED. IN NO EVENT SHALL THE CALIFORNIA INSTITUTE OF TECHNOLOGY BE
23
 * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
24
 * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
25
 * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
26
 * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
27
 * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
28
 * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
29
 * POSSIBILITY OF SUCH DAMAGE.
30
 *
31
 * Copyright 2014-2021 Esri
32
 *
33
 * Licensed under the Apache License, Version 2.0 (the "License");
34
 * you may not use this file except in compliance with the License.
35
 * You may obtain a copy of the License at
36
 *
37
 * http://www.apache.org/licenses/LICENSE-2.0
38
 *
39
 * Unless required by applicable law or agreed to in writing, software
40
 * distributed under the License is distributed on an "AS IS" BASIS,
41
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
42
 * See the License for the specific language governing permissions and
43
 * limitations under the License.
44
 */
45
46
/******************************************************************************
47
 *
48
 * Project:  Meta Raster File Format Driver Implementation, RasterBand
49
 * Purpose:  Implementation of MRF band
50
 *
51
 * Author:   Lucian Plesea, Lucian.Plesea jpl.nasa.gov, lplesea esri.com
52
 *
53
 ****************************************************************************/
54
55
#include "marfa.h"
56
#include "gdal_priv.h"
57
#include "ogr_srs_api.h"
58
#include "ogr_spatialref.h"
59
60
#include <vector>
61
#include <algorithm>
62
#include <cassert>
63
#include <zlib.h>
64
#if defined(ZSTD_SUPPORT)
65
#include <zstd.h>
66
#endif
67
68
using namespace std::chrono;
69
70
NAMESPACE_MRF_START
71
72
// packs a block of a given type, with a stride
73
// Count is the number of items that need to be copied
74
// These are separate to allow for optimization
75
76
template <typename T>
77
static void cpy_stride_in(void *dst, void *src, int c, int stride)
78
535k
{
79
535k
    T *s = reinterpret_cast<T *>(src);
80
535k
    T *d = reinterpret_cast<T *>(dst);
81
82
198M
    while (c--)
83
197M
    {
84
197M
        *d++ = *s;
85
197M
        s += stride;
86
197M
    }
87
535k
}
mrf_band.cpp:void GDAL_MRF::cpy_stride_in<unsigned char>(void*, void*, int, int)
Line
Count
Source
78
141k
{
79
141k
    T *s = reinterpret_cast<T *>(src);
80
141k
    T *d = reinterpret_cast<T *>(dst);
81
82
148M
    while (c--)
83
147M
    {
84
147M
        *d++ = *s;
85
147M
        s += stride;
86
147M
    }
87
141k
}
mrf_band.cpp:void GDAL_MRF::cpy_stride_in<short>(void*, void*, int, int)
Line
Count
Source
78
149k
{
79
149k
    T *s = reinterpret_cast<T *>(src);
80
149k
    T *d = reinterpret_cast<T *>(dst);
81
82
28.5M
    while (c--)
83
28.4M
    {
84
28.4M
        *d++ = *s;
85
28.4M
        s += stride;
86
28.4M
    }
87
149k
}
mrf_band.cpp:void GDAL_MRF::cpy_stride_in<int>(void*, void*, int, int)
Line
Count
Source
78
152k
{
79
152k
    T *s = reinterpret_cast<T *>(src);
80
152k
    T *d = reinterpret_cast<T *>(dst);
81
82
14.1M
    while (c--)
83
14.0M
    {
84
14.0M
        *d++ = *s;
85
14.0M
        s += stride;
86
14.0M
    }
87
152k
}
mrf_band.cpp:void GDAL_MRF::cpy_stride_in<long long>(void*, void*, int, int)
Line
Count
Source
78
91.7k
{
79
91.7k
    T *s = reinterpret_cast<T *>(src);
80
91.7k
    T *d = reinterpret_cast<T *>(dst);
81
82
7.38M
    while (c--)
83
7.29M
    {
84
7.29M
        *d++ = *s;
85
7.29M
        s += stride;
86
7.29M
    }
87
91.7k
}
88
89
template <typename T>
90
static void cpy_stride_out(void *dst, void *src, int c, int stride)
91
0
{
92
0
    T *s = reinterpret_cast<T *>(src);
93
0
    T *d = reinterpret_cast<T *>(dst);
94
95
0
    while (c--)
96
0
    {
97
0
        *d = *s++;
98
0
        d += stride;
99
0
    }
100
0
}
Unexecuted instantiation: mrf_band.cpp:void GDAL_MRF::cpy_stride_out<unsigned char>(void*, void*, int, int)
Unexecuted instantiation: mrf_band.cpp:void GDAL_MRF::cpy_stride_out<short>(void*, void*, int, int)
Unexecuted instantiation: mrf_band.cpp:void GDAL_MRF::cpy_stride_out<int>(void*, void*, int, int)
Unexecuted instantiation: mrf_band.cpp:void GDAL_MRF::cpy_stride_out<long long>(void*, void*, int, int)
101
102
// Does every value in the buffer have the same value, using strict comparison
103
template <typename T>
104
inline int isAllVal(const T *b, size_t bytecount, double ndv)
105
0
{
106
0
    T val = static_cast<T>(ndv);
107
0
    size_t count = bytecount / sizeof(T);
108
0
    for (; count; --count)
109
0
    {
110
0
        if (*(b++) != val)
111
0
        {
112
0
            return FALSE;
113
0
        }
114
0
    }
115
0
    return TRUE;
116
0
}
Unexecuted instantiation: int GDAL_MRF::isAllVal<unsigned char>(unsigned char const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<signed char>(signed char const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<unsigned short>(unsigned short const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<short>(short const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<unsigned int>(unsigned int const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<int>(int const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<unsigned long long>(unsigned long long const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<long long>(long long const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<float>(float const*, unsigned long, double)
Unexecuted instantiation: int GDAL_MRF::isAllVal<double>(double const*, unsigned long, double)
117
118
// Dispatcher based on gdal data type
119
static int isAllVal(GDALDataType gt, void *b, size_t bytecount, double ndv)
120
0
{
121
    // Test to see if it has data
122
0
    int isempty = false;
123
124
    // A case branch in a temporary macro, conversion from gdal enum to type
125
0
#define TEST_T(GType, T)                                                       \
126
0
    case GType:                                                                \
127
0
        isempty = isAllVal(reinterpret_cast<T *>(b), bytecount, ndv);          \
128
0
        break
129
130
0
    switch (gt)
131
0
    {
132
0
        TEST_T(GDT_Byte, GByte);
133
0
        TEST_T(GDT_Int8, GInt8);
134
0
        TEST_T(GDT_UInt16, GUInt16);
135
0
        TEST_T(GDT_Int16, GInt16);
136
0
        TEST_T(GDT_UInt32, GUInt32);
137
0
        TEST_T(GDT_Int32, GInt32);
138
0
        TEST_T(GDT_UInt64, GUInt64);
139
0
        TEST_T(GDT_Int64, GInt64);
140
0
        TEST_T(GDT_Float32, float);
141
0
        TEST_T(GDT_Float64, double);
142
0
        default:
143
0
            break;
144
0
    }
145
0
#undef TEST_T
146
147
0
    return isempty;
148
0
}
149
150
// Swap bytes in place, unconditional
151
// cppcheck-suppress constParameterReference
152
static void swab_buff(buf_mgr &src, const ILImage &img)
153
0
{
154
0
    size_t i;
155
0
    switch (GDALGetDataTypeSizeBytes(img.dt))
156
0
    {
157
0
        case 2:
158
0
        {
159
0
            uint16_t *b = reinterpret_cast<uint16_t *>(src.buffer);
160
0
            for (i = src.size / 2; i; b++, i--)
161
0
                *b = swab16(*b);
162
0
            break;
163
0
        }
164
0
        case 4:
165
0
        {
166
0
            uint32_t *b = reinterpret_cast<uint32_t *>(src.buffer);
167
0
            for (i = src.size / 4; i; b++, i--)
168
0
                *b = swab32(*b);
169
0
            break;
170
0
        }
171
0
        case 8:
172
0
        {
173
0
            uint64_t *b = reinterpret_cast<uint64_t *>(src.buffer);
174
0
            for (i = src.size / 8; i; b++, i--)
175
0
                *b = swab64(*b);
176
0
            break;
177
0
        }
178
0
    }
179
0
}
180
181
// Similar to compress2() but with flags to control zlib features
182
// Returns true if it worked
183
static int ZPack(const buf_mgr &src, buf_mgr &dst, int flags)
184
0
{
185
0
    z_stream stream;
186
0
    int err;
187
188
0
    memset(&stream, 0, sizeof(stream));
189
0
    stream.next_in = reinterpret_cast<Bytef *>(src.buffer);
190
0
    stream.avail_in = (uInt)src.size;
191
0
    stream.next_out = reinterpret_cast<Bytef *>(dst.buffer);
192
0
    stream.avail_out = (uInt)dst.size;
193
194
0
    int level = std::clamp(flags & ZFLAG_LMASK, 1, 9);
195
0
    int wb = MAX_WBITS;
196
    // if gz flag is set, ignore raw request
197
0
    if (flags & ZFLAG_GZ)
198
0
        wb += 16;
199
0
    else if (flags & ZFLAG_RAW)
200
0
        wb = -wb;
201
0
    int memlevel = 8;  // Good compromise
202
0
    int strategy = (flags & ZFLAG_SMASK) >> 6;
203
0
    if (strategy > 4)
204
0
        strategy = 0;
205
206
0
    err = deflateInit2(&stream, level, Z_DEFLATED, wb, memlevel, strategy);
207
0
    if (err != Z_OK)
208
0
    {
209
0
        deflateEnd(&stream);
210
0
        return err;
211
0
    }
212
213
0
    err = deflate(&stream, Z_FINISH);
214
0
    if (err != Z_STREAM_END)
215
0
    {
216
0
        deflateEnd(&stream);
217
0
        return false;
218
0
    }
219
0
    dst.size = stream.total_out;
220
0
    err = deflateEnd(&stream);
221
0
    return err == Z_OK;
222
0
}
223
224
// Similar to uncompress() from zlib, accepts the ZFLAG_RAW
225
// Return true if it worked
226
static int ZUnPack(const buf_mgr &src, buf_mgr &dst, int flags)
227
7
{
228
229
7
    z_stream stream;
230
7
    int err;
231
232
7
    memset(&stream, 0, sizeof(stream));
233
7
    stream.next_in = reinterpret_cast<Bytef *>(src.buffer);
234
7
    stream.avail_in = (uInt)src.size;
235
7
    stream.next_out = reinterpret_cast<Bytef *>(dst.buffer);
236
7
    stream.avail_out = (uInt)dst.size;
237
238
    // 32 means autodetec gzip or zlib header, negative 15 is for raw
239
7
    int wb = (ZFLAG_RAW & flags) ? -MAX_WBITS : 32 + MAX_WBITS;
240
7
    err = inflateInit2(&stream, wb);
241
7
    if (err != Z_OK)
242
0
        return false;
243
244
7
    err = inflate(&stream, Z_FINISH);
245
7
    if (err != Z_STREAM_END)
246
6
    {
247
6
        inflateEnd(&stream);
248
6
        return false;
249
6
    }
250
1
    dst.size = stream.total_out;
251
1
    err = inflateEnd(&stream);
252
1
    return err == Z_OK;
253
7
}
254
255
/*
256
 * Deflates a buffer, extrasize is the available size in the buffer past the
257
 * input If the output fits past the data, it uses that area, otherwise it uses
258
 * a temporary buffer and copies the data over the input on return, returning a
259
 * pointer to it. The output size is returned in src.size Returns nullptr when
260
 * compression failed
261
 */
262
static void *DeflateBlock(buf_mgr &src, size_t extrasize, int flags)
263
0
{
264
    // The one we might need to allocate
265
0
    void *dbuff = nullptr;
266
0
    buf_mgr dst = {src.buffer + src.size, extrasize};
267
268
    // Allocate a temp buffer if there is not sufficient space,
269
    // We need to have a bit more than half the buffer available
270
0
    if (extrasize < (src.size + 64))
271
0
    {
272
0
        dst.size = src.size + 64;
273
0
        dbuff = VSIMalloc(dst.size);
274
0
        dst.buffer = (char *)dbuff;
275
0
        if (!dst.buffer)
276
0
            return nullptr;
277
0
    }
278
279
0
    if (!ZPack(src, dst, flags))
280
0
    {
281
0
        CPLFree(dbuff);  // Safe to call with NULL
282
0
        return nullptr;
283
0
    }
284
285
0
    if (src.size + extrasize < dst.size)
286
0
    {
287
0
        CPLError(CE_Failure, CPLE_AppDefined,
288
0
                 "DeflateBlock(): too small buffer");
289
0
        CPLFree(dbuff);
290
0
        return nullptr;
291
0
    }
292
293
    // source size is used to hold the output size
294
0
    src.size = dst.size;
295
    // If we didn't allocate a buffer, the receiver can use it already
296
0
    if (!dbuff)
297
0
        return dst.buffer;
298
299
    // If we allocated a buffer, we need to copy the data to the input buffer.
300
0
    memcpy(src.buffer, dbuff, src.size);
301
0
    CPLFree(dbuff);
302
0
    return src.buffer;
303
0
}
304
305
#if defined(ZSTD_SUPPORT)
306
307
CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static void rankfilter(buf_mgr &src,
308
                                                            size_t factor)
309
0
{
310
    // Arange bytes by rank
311
0
    if (factor > 1)
312
0
    {
313
0
        std::vector<char> tempb(src.size);
314
0
        char *d = tempb.data();
315
0
        for (size_t j = 0; j < factor; j++)
316
0
            for (size_t i = j; i < src.size; i += factor)
317
0
                *d++ = src.buffer[i];
318
0
        memcpy(src.buffer, tempb.data(), src.size);
319
0
    }
320
    // byte delta
321
0
    auto p = reinterpret_cast<GByte *>(src.buffer);
322
0
    auto guard = p + src.size;
323
0
    GByte b(0);
324
0
    while (p < guard)
325
0
    {
326
0
        GByte temp = *p;
327
0
        *p -= b;
328
0
        b = temp;
329
0
        p++;
330
0
    }
331
0
}
332
333
CPL_NOSANITIZE_UNSIGNED_INT_OVERFLOW static void derank(buf_mgr &src,
334
                                                        size_t factor)
335
0
{
336
    // undo delta
337
0
    auto p = reinterpret_cast<GByte *>(src.buffer);
338
0
    auto guard = p + src.size;
339
0
    GByte b(0);
340
0
    while (p < guard)
341
0
    {
342
0
        b += *p;
343
0
        *p = b;
344
0
        p++;
345
0
    }
346
0
    if (factor > 1)
347
0
    {  // undo rank separation
348
0
        std::vector<char> tempb(src.size);
349
0
        char *d = tempb.data();
350
0
        size_t chunk = src.size / factor;
351
0
        for (size_t i = 0; i < chunk; i++)
352
0
            for (size_t j = 0; j < factor; j++)
353
0
                *d++ = src.buffer[chunk * j + i];
354
0
        memcpy(src.buffer, tempb.data(), src.size);
355
0
    }
356
0
}
357
358
/*
359
 * Compress a buffer using zstd, extrasize is the available size in the buffer
360
 * past the input If ranks > 0, apply the rank filter If the output fits past
361
 * the data, it uses that area, otherwise it uses a temporary buffer and copies
362
 * the data over the input on return, returning a pointer to it. The output size
363
 * is returned in src.size Returns nullptr when compression failed
364
 */
365
static void *ZstdCompBlock(buf_mgr &src, size_t extrasize, int c_level,
366
                           ZSTD_CCtx *cctx, size_t ranks)
367
0
{
368
0
    if (!cctx)
369
0
        return nullptr;
370
0
    if (ranks && (src.size % ranks) == 0)
371
0
        rankfilter(src, ranks);
372
373
    // might need a buffer for the zstd output
374
0
    std::vector<char> dbuff;
375
0
    void *dst = src.buffer + src.size;
376
0
    size_t size = extrasize;
377
    // Allocate a temp buffer if there is not sufficient space.
378
    // Zstd bound is about (size * 1.004 + 64)
379
0
    if (size < ZSTD_compressBound(src.size))
380
0
    {
381
0
        size = ZSTD_compressBound(src.size);
382
0
        dbuff.resize(size);
383
0
        dst = dbuff.data();
384
0
    }
385
386
    // Use the streaming interface, it's faster and better
387
    // See discussion at https://github.com/facebook/zstd/issues/3729
388
0
    ZSTD_outBuffer output = {dst, size, 0};
389
0
    ZSTD_inBuffer input = {src.buffer, src.size, 0};
390
    // Set level
391
0
    ZSTD_CCtx_setParameter(cctx, ZSTD_c_compressionLevel, c_level);
392
    // First, pass a continue flag, otherwise it will compress in one go
393
0
    size_t val = ZSTD_compressStream2(cctx, &output, &input, ZSTD_e_continue);
394
    // If it worked, pass the end flag to flush the buffer
395
0
    if (val == 0)
396
0
        val = ZSTD_compressStream2(cctx, &output, &input, ZSTD_e_end);
397
0
    if (ZSTD_isError(val))
398
0
        return nullptr;
399
0
    val = output.pos;
400
401
    // If we didn't need the buffer, packed data is already in the user buffer
402
0
    if (dbuff.empty())
403
0
    {
404
0
        src.size = val;
405
0
        return dst;
406
0
    }
407
408
0
    if (val > (src.size + extrasize))
409
0
    {  // Doesn't fit in user buffer
410
0
        CPLError(CE_Failure, CPLE_AssertionFailed,
411
0
                 "MRF: ZSTD compression buffer too small");
412
0
        return nullptr;  // Error
413
0
    }
414
415
0
    memcpy(src.buffer, dbuff.data(), val);
416
0
    src.size = val;
417
0
    return src.buffer;
418
0
}
419
#endif
420
421
//
422
// The deflate_flags are available in all bands even if the DEFLATE option
423
// itself is not set.  This allows for PNG features to be controlled, as well
424
// as any other bands that use zlib by itself
425
//
426
MRFRasterBand::MRFRasterBand(MRFDataset *parent_dataset, const ILImage &image,
427
                             int band, int ov)
428
1.56M
    : poMRFDS(parent_dataset),
429
1.56M
      dodeflate(GetOptlist().FetchBoolean("DEFLATE", FALSE)),
430
      // Bring the quality to 0 to 9
431
1.56M
      deflate_flags(image.quality / 10),
432
1.56M
      dozstd(GetOptlist().FetchBoolean("ZSTD", FALSE)), zstd_level(9), m_l(ov),
433
1.56M
      img(image)
434
1.56M
{
435
1.56M
    nBand = band;
436
1.56M
    eDataType = parent_dataset->current.dt;
437
1.56M
    nRasterXSize = img.size.x;
438
1.56M
    nRasterYSize = img.size.y;
439
1.56M
    nBlockXSize = img.pagesize.x;
440
1.56M
    nBlockYSize = img.pagesize.y;
441
1.56M
    nBlocksPerRow = img.pagecount.x;
442
1.56M
    nBlocksPerColumn = img.pagecount.y;
443
1.56M
    img.NoDataValue = MRFRasterBand::GetNoDataValue(&img.hasNoData);
444
445
    // Pick up the twists, aka GZ, RAWZ headers
446
1.56M
    if (GetOptlist().FetchBoolean("GZ", FALSE))
447
0
        deflate_flags |= ZFLAG_GZ;
448
1.56M
    else if (GetOptlist().FetchBoolean("RAWZ", FALSE))
449
0
        deflate_flags |= ZFLAG_RAW;
450
    // And Pick up the ZLIB strategy, if any
451
1.56M
    const char *zstrategy = GetOptlist().FetchNameValueDef("Z_STRATEGY", "");
452
1.56M
    int zv = Z_DEFAULT_STRATEGY;
453
1.56M
    if (EQUAL(zstrategy, "Z_HUFFMAN_ONLY"))
454
0
        zv = Z_HUFFMAN_ONLY;
455
1.56M
    else if (EQUAL(zstrategy, "Z_RLE"))
456
0
        zv = Z_RLE;
457
1.56M
    else if (EQUAL(zstrategy, "Z_FILTERED"))
458
0
        zv = Z_FILTERED;
459
1.56M
    else if (EQUAL(zstrategy, "Z_FIXED"))
460
0
        zv = Z_FIXED;
461
1.56M
    deflate_flags |= (zv << 6);
462
1.56M
    if (image.quality < 23 && image.quality > 0)
463
0
        zstd_level = image.quality;
464
465
#if !defined(ZSTD_SUPPORT)
466
    if (dozstd)
467
    {  // signal error condition to caller
468
        CPLError(CE_Failure, CPLE_AssertionFailed,
469
                 "MRF: ZSTD support is not available");
470
        dozstd = FALSE;
471
    }
472
#endif
473
    // Chose zstd over deflate if both are enabled and available
474
1.56M
    if (dozstd && dodeflate)
475
0
        dodeflate = FALSE;
476
1.56M
}
477
478
// Clean up the overviews if they exist
479
MRFRasterBand::~MRFRasterBand()
480
1.56M
{
481
1.56M
    while (!overviews.empty())
482
0
    {
483
0
        delete overviews.back();
484
0
        overviews.pop_back();
485
0
    }
486
1.56M
}
487
488
// Look for a string from the dataset options or from the environment
489
const char *MRFRasterBand::GetOptionValue(const char *opt,
490
                                          const char *def) const
491
1.54M
{
492
1.54M
    const char *optValue = poMRFDS->optlist.FetchNameValue(opt);
493
1.54M
    if (optValue)
494
0
        return optValue;
495
1.54M
    return CPLGetConfigOption(opt, def);
496
1.54M
}
497
498
// Utility function, returns a value from a vector corresponding to the band
499
// index or the first entry
500
static double getBandValue(const std::vector<double> &v, int idx)
501
0
{
502
0
    return (static_cast<int>(v.size()) > idx) ? v[idx] : v[0];
503
0
}
504
505
// Maybe we should check against the type range?
506
// It is not keeping track of how many values have been set,
507
// so the application should set none or all the bands
508
// This call is only valid during Create
509
CPLErr MRFRasterBand::SetNoDataValue(double val)
510
0
{
511
0
    if (poMRFDS->bCrystalized)
512
0
    {
513
0
        CPLError(CE_Failure, CPLE_AssertionFailed,
514
0
                 "MRF: NoData can be set only during file create");
515
0
        return CE_Failure;
516
0
    }
517
0
    if (GInt32(poMRFDS->vNoData.size()) < nBand)
518
0
        poMRFDS->vNoData.resize(nBand);
519
0
    poMRFDS->vNoData[nBand - 1] = val;
520
    // We also need to set it for this band
521
0
    img.NoDataValue = val;
522
0
    img.hasNoData = true;
523
0
    return CE_None;
524
0
}
525
526
double MRFRasterBand::GetNoDataValue(int *pbSuccess)
527
1.72M
{
528
1.72M
    const std::vector<double> &v = poMRFDS->vNoData;
529
1.72M
    if (v.empty())
530
1.72M
        return GDALPamRasterBand::GetNoDataValue(pbSuccess);
531
0
    if (pbSuccess)
532
0
        *pbSuccess = TRUE;
533
0
    return getBandValue(v, nBand - 1);
534
1.72M
}
535
536
double MRFRasterBand::GetMinimum(int *pbSuccess)
537
17.9k
{
538
17.9k
    const std::vector<double> &v = poMRFDS->vMin;
539
17.9k
    if (v.empty())
540
17.9k
        return GDALPamRasterBand::GetMinimum(pbSuccess);
541
0
    if (pbSuccess)
542
0
        *pbSuccess = TRUE;
543
0
    return getBandValue(v, nBand - 1);
544
17.9k
}
545
546
double MRFRasterBand::GetMaximum(int *pbSuccess)
547
17.9k
{
548
17.9k
    const std::vector<double> &v = poMRFDS->vMax;
549
17.9k
    if (v.empty())
550
17.9k
        return GDALPamRasterBand::GetMaximum(pbSuccess);
551
0
    if (pbSuccess)
552
0
        *pbSuccess = TRUE;
553
0
    return getBandValue(v, nBand - 1);
554
17.9k
}
555
556
// Fill with typed ndv, count is always in bytes
557
template <typename T>
558
static CPLErr buff_fill(void *b, size_t count, const T ndv)
559
0
{
560
0
    T *buffer = static_cast<T *>(b);
561
0
    count /= sizeof(T);
562
0
    while (count--)
563
0
        *buffer++ = ndv;
564
0
    return CE_None;
565
0
}
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<unsigned short>(void*, unsigned long, unsigned short)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<short>(void*, unsigned long, short)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<unsigned int>(void*, unsigned long, unsigned int)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<int>(void*, unsigned long, int)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<unsigned long long>(void*, unsigned long, unsigned long long)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<long long>(void*, unsigned long, long long)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<float>(void*, unsigned long, float)
Unexecuted instantiation: mrf_band.cpp:CPLErr GDAL_MRF::buff_fill<double>(void*, unsigned long, double)
566
567
/**
568
 *\brief Fills a buffer with no data
569
 *
570
 */
571
CPLErr MRFRasterBand::FillBlock(void *buffer)
572
2.46k
{
573
2.46k
    int success;
574
2.46k
    double ndv = GetNoDataValue(&success);
575
2.46k
    if (!success)
576
2.46k
        ndv = 0.0;
577
2.46k
    size_t bsb = blockSizeBytes();
578
579
    // use memset for speed for bytes, or if nodata is zeros
580
2.46k
    if (0.0 == ndv || eDataType == GDT_Byte || eDataType == GDT_Int8)
581
2.46k
    {
582
2.46k
        memset(buffer, int(ndv), bsb);
583
2.46k
        return CE_None;
584
2.46k
    }
585
586
0
#define bf(T) buff_fill<T>(buffer, bsb, T(ndv));
587
0
    switch (eDataType)
588
0
    {
589
0
        case GDT_UInt16:
590
0
            return bf(GUInt16);
591
0
        case GDT_Int16:
592
0
            return bf(GInt16);
593
0
        case GDT_UInt32:
594
0
            return bf(GUInt32);
595
0
        case GDT_Int32:
596
0
            return bf(GInt32);
597
0
        case GDT_UInt64:
598
0
            return bf(GUInt64);
599
0
        case GDT_Int64:
600
0
            return bf(GInt64);
601
0
        case GDT_Float32:
602
0
            return bf(float);
603
0
        case GDT_Float64:
604
0
            return bf(double);
605
0
        default:
606
0
            break;
607
0
    }
608
0
#undef bf
609
    // Should exit before
610
0
    return CE_Failure;
611
0
}
612
613
/*\brief Interleave block fill
614
 *
615
 *  Acquire space for all the other bands, fill each one then drop the locks
616
 *  The current band output goes directly into the buffer
617
 */
618
619
CPLErr MRFRasterBand::FillBlock(int xblk, int yblk, void *buffer)
620
0
{
621
0
    std::vector<GDALRasterBlock *> blocks;
622
623
0
    for (int i = 0; i < poMRFDS->nBands; i++)
624
0
    {
625
0
        GDALRasterBand *b = poMRFDS->GetRasterBand(i + 1);
626
0
        if (b->GetOverviewCount() && 0 != m_l)
627
0
            b = b->GetOverview(m_l - 1);
628
629
        // Get the other band blocks, keep them around until later
630
0
        if (b == this)
631
0
        {
632
0
            FillBlock(buffer);
633
0
        }
634
0
        else
635
0
        {
636
0
            GDALRasterBlock *poBlock = b->GetLockedBlockRef(xblk, yblk, 1);
637
0
            if (poBlock == nullptr)  // Didn't get this block
638
0
                break;
639
0
            FillBlock(poBlock->GetDataRef());
640
0
            blocks.push_back(poBlock);
641
0
        }
642
0
    }
643
644
    // Drop the locks for blocks we acquired
645
0
    for (int i = 0; i < int(blocks.size()); i++)
646
0
        blocks[i]->DropLock();
647
648
0
    return CE_None;
649
0
}
650
651
/*\brief Interleave block read
652
 *
653
 *  Acquire space for all the other bands, unpack from the dataset buffer, then
654
 * drop the locks The current band output goes directly into the buffer
655
 */
656
657
CPLErr MRFRasterBand::ReadInterleavedBlock(int xblk, int yblk, void *buffer)
658
2.86k
{
659
2.86k
    std::vector<GDALRasterBlock *> blocks;
660
661
538k
    for (int i = 0; i < poMRFDS->nBands; i++)
662
535k
    {
663
535k
        GDALRasterBand *b = poMRFDS->GetRasterBand(i + 1);
664
535k
        if (b->GetOverviewCount() && 0 != m_l)
665
0
            b = b->GetOverview(m_l - 1);
666
667
535k
        void *ob = buffer;
668
        // Get the other band blocks, keep them around until later
669
535k
        if (b != this)
670
532k
        {
671
532k
            GDALRasterBlock *poBlock = b->GetLockedBlockRef(xblk, yblk, 1);
672
532k
            if (poBlock == nullptr)
673
0
                break;
674
532k
            ob = poBlock->GetDataRef();
675
532k
            blocks.push_back(poBlock);
676
532k
        }
677
678
        // Just the right mix of templates and macros make deinterleaving tidy
679
535k
        void *pbuffer = poMRFDS->GetPBuffer();
680
535k
#define CpySI(T)                                                               \
681
535k
    cpy_stride_in<T>(ob, reinterpret_cast<T *>(pbuffer) + i,                   \
682
535k
                     blockSizeBytes() / sizeof(T), img.pagesize.c)
683
684
        // Page is already in poMRFDS->pbuffer, not empty
685
        // There are only four cases, since only the data size matters
686
535k
        switch (GDALGetDataTypeSizeBytes(eDataType))
687
535k
        {
688
141k
            case 1:
689
141k
                CpySI(GByte);
690
141k
                break;
691
149k
            case 2:
692
149k
                CpySI(GInt16);
693
149k
                break;
694
152k
            case 4:
695
152k
                CpySI(GInt32);
696
152k
                break;
697
91.7k
            case 8:
698
91.7k
                CpySI(GIntBig);
699
91.7k
                break;
700
535k
        }
701
535k
    }
702
703
2.86k
#undef CpySI
704
705
    // Drop the locks we acquired
706
535k
    for (int i = 0; i < int(blocks.size()); i++)
707
532k
        blocks[i]->DropLock();
708
709
2.86k
    return CE_None;
710
2.86k
}
711
712
/**
713
 *\brief Fetch a block from the backing store dataset and keep a copy in the
714
 *cache
715
 *
716
 * @param xblk The X block number, zero based
717
 * @param yblk The Y block number, zero based
718
 * @param buffer buffer
719
 *
720
 */
721
CPLErr MRFRasterBand::FetchBlock(int xblk, int yblk, void *buffer)
722
0
{
723
0
    assert(!poMRFDS->source.empty());
724
0
    CPLDebug("MRF_IB", "FetchBlock %d,%d,0,%d, level  %d\n", xblk, yblk, nBand,
725
0
             m_l);
726
727
0
    if (poMRFDS->clonedSource)  // This is a clone
728
0
        return FetchClonedBlock(xblk, yblk, buffer);
729
730
0
    const GInt32 cstride = img.pagesize.c;  // 1 if band separate
731
0
    ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
732
0
    GUIntBig infooffset = IdxOffset(req, img);
733
734
0
    GDALDataset *poSrcDS = nullptr;
735
0
    if (nullptr == (poSrcDS = poMRFDS->GetSrcDS()))
736
0
    {
737
0
        CPLError(CE_Failure, CPLE_AppDefined, "MRF: Can't open source file %s",
738
0
                 poMRFDS->source.c_str());
739
0
        return CE_Failure;
740
0
    }
741
742
    // Scale to base resolution
743
0
    double scl = pow(poMRFDS->scale, m_l);
744
0
    if (0 == m_l)
745
0
        scl = 1;  // To allow for precision issues
746
747
    // Prepare parameters for RasterIO, they might be different from a full page
748
0
    const GSpacing vsz = GDALGetDataTypeSizeBytes(eDataType);
749
0
    int Xoff = int(xblk * img.pagesize.x * scl + 0.5);
750
0
    int Yoff = int(yblk * img.pagesize.y * scl + 0.5);
751
0
    int readszx = int(img.pagesize.x * scl + 0.5);
752
0
    int readszy = int(img.pagesize.y * scl + 0.5);
753
754
    // Compare with the full size and clip to the right and bottom if needed
755
0
    int clip = 0;
756
0
    if (Xoff + readszx > poMRFDS->full.size.x)
757
0
    {
758
0
        clip |= 1;
759
0
        readszx = poMRFDS->full.size.x - Xoff;
760
0
    }
761
0
    if (Yoff + readszy > poMRFDS->full.size.y)
762
0
    {
763
0
        clip |= 1;
764
0
        readszy = poMRFDS->full.size.y - Yoff;
765
0
    }
766
767
    // This is where the whole page fits
768
0
    void *ob = buffer;
769
0
    if (cstride != 1)
770
0
        ob = poMRFDS->GetPBuffer();
771
772
    // Fill buffer with NoData if clipping
773
0
    if (clip)
774
0
        FillBlock(ob);
775
776
    // Use the dataset RasterIO to read one or all bands if interleaved
777
0
    CPLErr ret = poSrcDS->RasterIO(
778
0
        GF_Read, Xoff, Yoff, readszx, readszy, ob, pcount(readszx, int(scl)),
779
0
        pcount(readszy, int(scl)), eDataType, cstride,
780
0
        (1 == cstride) ? &nBand : nullptr, vsz * cstride,
781
0
        vsz * cstride * img.pagesize.x,
782
0
        (cstride != 1) ? vsz : vsz * img.pagesize.x * img.pagesize.y, nullptr);
783
784
0
    if (ret != CE_None)
785
0
        return ret;
786
787
    // Might have the block in the pbuffer, mark it anyhow
788
0
    poMRFDS->tile = req;
789
0
    buf_mgr filesrc;
790
0
    filesrc.buffer = (char *)ob;
791
0
    filesrc.size = static_cast<size_t>(img.pageSizeBytes);
792
793
0
    if (poMRFDS->bypass_cache)
794
0
    {  // No local caching, just return the data
795
0
        if (1 == cstride)
796
0
            return CE_None;
797
0
        return ReadInterleavedBlock(xblk, yblk, buffer);
798
0
    }
799
800
    // Test to see if it needs to be written, or just marked as checked
801
0
    int success;
802
0
    double val = GetNoDataValue(&success);
803
0
    if (!success)
804
0
        val = 0.0;
805
806
    // TODO: test band by band if data is interleaved
807
0
    if (isAllVal(eDataType, ob, img.pageSizeBytes, val))
808
0
    {
809
        // Mark it empty and checked, ignore the possible write error
810
0
        poMRFDS->WriteTile((void *)1, infooffset, 0);
811
0
        if (1 == cstride)
812
0
            return CE_None;
813
0
        return ReadInterleavedBlock(xblk, yblk, buffer);
814
0
    }
815
816
    // Write the page in the local cache
817
818
    // Have to use a separate buffer for compression output.
819
0
    void *outbuff = VSIMalloc(poMRFDS->pbsize);
820
0
    if (nullptr == outbuff)
821
0
    {
822
0
        CPLError(CE_Failure, CPLE_AppDefined,
823
0
                 "Can't get buffer for writing page");
824
        // This is not really an error for a cache, the data is fine
825
0
        return CE_Failure;
826
0
    }
827
828
0
    buf_mgr filedst = {static_cast<char *>(outbuff), poMRFDS->pbsize};
829
0
    auto start_time = steady_clock::now();
830
0
    if (Compress(filedst, filesrc) != CE_None)
831
0
    {
832
0
        return CE_Failure;
833
0
    }
834
835
    // Where the output is, in case we deflate
836
0
    void *usebuff = outbuff;
837
0
    if (dodeflate)
838
0
    {
839
0
        CPLAssert(poMRFDS->pbsize <= filedst.size);
840
0
        usebuff = DeflateBlock(filedst, poMRFDS->pbsize - filedst.size,
841
0
                               deflate_flags);
842
0
        if (!usebuff)
843
0
        {
844
0
            CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
845
0
            return CE_Failure;
846
0
        }
847
0
    }
848
849
0
#if defined(ZSTD_SUPPORT)
850
0
    else if (dozstd)
851
0
    {
852
0
        size_t ranks = 0;  // Assume no need for byte rank sort
853
0
        if (img.comp == IL_NONE || img.comp == IL_ZSTD)
854
0
            ranks =
855
0
                static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) * cstride;
856
0
        usebuff = ZstdCompBlock(filedst, poMRFDS->pbsize - filedst.size,
857
0
                                zstd_level, poMRFDS->getzsc(), ranks);
858
0
        if (!usebuff)
859
0
        {
860
0
            CPLError(CE_Failure, CPLE_AppDefined,
861
0
                     "MRF: ZSTD compression error");
862
0
            return CE_Failure;
863
0
        }
864
0
    }
865
0
#endif
866
867
0
    poMRFDS->write_timer +=
868
0
        duration_cast<nanoseconds>(steady_clock::now() - start_time);
869
870
    // Write and update the tile index
871
0
    ret = poMRFDS->WriteTile(usebuff, infooffset, filedst.size);
872
0
    CPLFree(outbuff);
873
874
    // If we hit an error or if unpaking is not needed
875
0
    if (ret != CE_None || cstride == 1)
876
0
        return ret;
877
878
    // data is already in DS buffer, deinterlace it in pixel blocks
879
0
    return ReadInterleavedBlock(xblk, yblk, buffer);
880
0
}
881
882
/**
883
 *\brief Fetch for a cloned MRF
884
 *
885
 * @param xblk The X block number, zero based
886
 * @param yblk The Y block number, zero based
887
 * @param buffer buffer
888
 *
889
 */
890
891
CPLErr MRFRasterBand::FetchClonedBlock(int xblk, int yblk, void *buffer)
892
0
{
893
0
    CPLDebug("MRF_IB", "FetchClonedBlock %d,%d,0,%d, level  %d\n", xblk, yblk,
894
0
             nBand, m_l);
895
896
    // Paranoid check
897
0
    assert(poMRFDS->clonedSource);
898
0
    MRFDataset *poSrc = static_cast<MRFDataset *>(poMRFDS->GetSrcDS());
899
0
    if (nullptr == poSrc)
900
0
    {
901
0
        CPLError(CE_Failure, CPLE_AppDefined, "MRF: Can't open source file %s",
902
0
                 poMRFDS->source.c_str());
903
0
        return CE_Failure;
904
0
    }
905
906
0
    if (poMRFDS->bypass_cache || GF_Read == DataMode())
907
0
    {
908
        // Can't store, so just fetch from source, which is an MRF with
909
        // identical structure
910
0
        MRFRasterBand *b =
911
0
            static_cast<MRFRasterBand *>(poSrc->GetRasterBand(nBand));
912
0
        if (b->GetOverviewCount() && m_l)
913
0
            b = static_cast<MRFRasterBand *>(b->GetOverview(m_l - 1));
914
0
        if (b == nullptr)
915
0
            return CE_Failure;
916
0
        return b->IReadBlock(xblk, yblk, buffer);
917
0
    }
918
919
0
    ILSize req(xblk, yblk, 0, (nBand - 1) / img.pagesize.c, m_l);
920
0
    ILIdx tinfo;
921
922
    // Get the cloned source tile info
923
    // The cloned source index is after the current one
924
0
    if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img, poMRFDS->idxSize))
925
0
    {
926
0
        CPLError(CE_Failure, CPLE_AppDefined,
927
0
                 "MRF: Unable to read cloned index entry");
928
0
        return CE_Failure;
929
0
    }
930
931
0
    GUIntBig infooffset = IdxOffset(req, img);
932
0
    CPLErr err;
933
934
    // Does the source have this tile?
935
0
    if (tinfo.size == 0)
936
0
    {  // Nope, mark it empty and return fill
937
0
        err = poMRFDS->WriteTile((void *)1, infooffset, 0);
938
0
        if (CE_None != err)
939
0
            return err;
940
0
        return FillBlock(buffer);
941
0
    }
942
943
0
    VSILFILE *srcfd = poSrc->DataFP();
944
0
    if (nullptr == srcfd)
945
0
    {
946
0
        CPLError(CE_Failure, CPLE_AppDefined,
947
0
                 "MRF: Can't open source data file %s",
948
0
                 poMRFDS->source.c_str());
949
0
        return CE_Failure;
950
0
    }
951
952
    // Need to read the tile from the source
953
0
    if (tinfo.size <= 0 || tinfo.size > INT_MAX)
954
0
    {
955
0
        CPLError(CE_Failure, CPLE_OutOfMemory,
956
0
                 "Invalid tile size " CPL_FRMT_GIB, tinfo.size);
957
0
        return CE_Failure;
958
0
    }
959
0
    char *buf = static_cast<char *>(VSIMalloc(static_cast<size_t>(tinfo.size)));
960
0
    if (buf == nullptr)
961
0
    {
962
0
        CPLError(CE_Failure, CPLE_OutOfMemory,
963
0
                 "Cannot allocate " CPL_FRMT_GIB " bytes", tinfo.size);
964
0
        return CE_Failure;
965
0
    }
966
967
0
    VSIFSeekL(srcfd, tinfo.offset, SEEK_SET);
968
0
    if (tinfo.size !=
969
0
        GIntBig(VSIFReadL(buf, 1, static_cast<size_t>(tinfo.size), srcfd)))
970
0
    {
971
0
        CPLFree(buf);
972
0
        CPLError(CE_Failure, CPLE_AppDefined,
973
0
                 "MRF: Can't read data from source %s",
974
0
                 poSrc->current.datfname.c_str());
975
0
        return CE_Failure;
976
0
    }
977
978
    // Write it then reissue the read
979
0
    err = poMRFDS->WriteTile(buf, infooffset, tinfo.size);
980
0
    CPLFree(buf);
981
0
    if (CE_None != err)
982
0
        return err;
983
    // Reissue read, it will work from the cloned data
984
0
    return IReadBlock(xblk, yblk, buffer);
985
0
}
986
987
/**
988
 *\brief read a block in the provided buffer
989
 *
990
 *  For separate band model, the DS buffer is not used, the read is direct in
991
 *the buffer For pixel interleaved model, the DS buffer holds the temp copy and
992
 *all the other bands are force read
993
 *
994
 */
995
996
CPLErr MRFRasterBand::IReadBlock(int xblk, int yblk, void *buffer)
997
49.1k
{
998
49.1k
    GInt32 cstride = img.pagesize.c;
999
49.1k
    ILIdx tinfo;
1000
49.1k
    ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
1001
49.1k
    CPLDebug("MRF_IB",
1002
49.1k
             "IReadBlock %d,%d,0,%d, level %d, idxoffset " CPL_FRMT_GIB "\n",
1003
49.1k
             xblk, yblk, nBand - 1, m_l, IdxOffset(req, img));
1004
1005
    // If this is a caching file and bypass is on, just do the fetch
1006
49.1k
    if (poMRFDS->bypass_cache && !poMRFDS->source.empty())
1007
0
        return FetchBlock(xblk, yblk, buffer);
1008
1009
49.1k
    tinfo.size = 0;  // Just in case it is missing
1010
49.1k
    if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img))
1011
7.14k
    {
1012
7.14k
        if (!poMRFDS->no_errors)
1013
7.14k
        {
1014
7.14k
            CPLError(CE_Failure, CPLE_AppDefined,
1015
7.14k
                     "MRF: Unable to read index at offset " CPL_FRMT_GIB,
1016
7.14k
                     IdxOffset(req, img));
1017
7.14k
            return CE_Failure;
1018
7.14k
        }
1019
0
        return FillBlock(buffer);
1020
7.14k
    }
1021
1022
42.0k
    if (0 == tinfo.size)
1023
2.46k
    {  // Could be missing or it could be caching
1024
        // Offset != 0 means no data, Update mode is for local MRFs only
1025
        // if caching index mode is RO don't try to fetch
1026
        // Also, caching MRFs can't be opened in update mode
1027
2.46k
        if (0 != tinfo.offset || GA_Update == poMRFDS->eAccess ||
1028
2.30k
            poMRFDS->source.empty() || IdxMode() == GF_Read)
1029
2.46k
            return FillBlock(buffer);
1030
1031
        // caching MRF, need to fetch a block
1032
0
        return FetchBlock(xblk, yblk, buffer);
1033
2.46k
    }
1034
1035
39.5k
    CPLDebug("MRF_IB", "Tinfo offset " CPL_FRMT_GIB ", size " CPL_FRMT_GIB "\n",
1036
39.5k
             tinfo.offset, tinfo.size);
1037
    // If we have a tile, read it
1038
1039
    // Should use a permanent buffer, like the pbuffer mechanism
1040
    // Get a large buffer, in case we need to unzip
1041
1042
    // We add a padding of 3 bytes since in LERC1 decompression, we can
1043
    // dereference a unsigned int at the end of the buffer, that can be
1044
    // partially out of the buffer.
1045
    // Can be reproduced with :
1046
    // gdal_translate ../autotest/gcore/data/byte.tif out.mrf -of MRF -co
1047
    // COMPRESS=LERC -co OPTIONS=V1:YES -ot Float32 valgrind gdalinfo -checksum
1048
    // out.mrf Invalid read of size 4 at BitStuffer::read(unsigned char**,
1049
    // std::vector<unsigned int, std::allocator<unsigned int> >&) const
1050
    // (BitStuffer.cpp:153)
1051
1052
    // No stored tile should be larger than twice the raw size.
1053
39.5k
    if (tinfo.size <= 0 || tinfo.size > poMRFDS->pbsize * 2)
1054
447
    {
1055
447
        if (!poMRFDS->no_errors)
1056
447
        {
1057
447
            CPLError(CE_Failure, CPLE_OutOfMemory,
1058
447
                     "Stored tile is too large: " CPL_FRMT_GIB, tinfo.size);
1059
447
            return CE_Failure;
1060
447
        }
1061
0
        return FillBlock(buffer);
1062
447
    }
1063
1064
39.1k
    VSILFILE *dfp = DataFP();
1065
1066
    // No data file to read from
1067
39.1k
    if (dfp == nullptr)
1068
78
        return CE_Failure;
1069
1070
39.0k
    void *data = VSIMalloc(static_cast<size_t>(tinfo.size + PADDING_BYTES));
1071
39.0k
    if (data == nullptr)
1072
0
    {
1073
0
        CPLError(CE_Failure, CPLE_OutOfMemory,
1074
0
                 "Could not allocate memory for tile size: " CPL_FRMT_GIB,
1075
0
                 tinfo.size);
1076
0
        return CE_Failure;
1077
0
    }
1078
1079
    // This part is not thread safe, but it is what GDAL expects
1080
39.0k
    VSIFSeekL(dfp, tinfo.offset, SEEK_SET);
1081
39.0k
    if (1 != VSIFReadL(data, static_cast<size_t>(tinfo.size), 1, dfp))
1082
363
    {
1083
363
        CPLFree(data);
1084
363
        if (poMRFDS->no_errors)
1085
0
            return FillBlock(buffer);
1086
363
        CPLError(CE_Failure, CPLE_AppDefined, "Unable to read data page, %d@%x",
1087
363
                 static_cast<int>(tinfo.size), static_cast<int>(tinfo.offset));
1088
363
        return CE_Failure;
1089
363
    }
1090
1091
    /* initialize padding bytes */
1092
38.6k
    memset(((char *)data) + static_cast<size_t>(tinfo.size), 0, PADDING_BYTES);
1093
38.6k
    buf_mgr src = {(char *)data, static_cast<size_t>(tinfo.size)};
1094
38.6k
    buf_mgr dst;
1095
1096
38.6k
    auto start_time = steady_clock::now();
1097
1098
    // We got the data, do we need to decompress it before decoding?
1099
38.6k
    if (dodeflate)
1100
7
    {
1101
7
        if (img.pageSizeBytes > INT_MAX - 1440)
1102
0
        {
1103
0
            CPLFree(data);
1104
0
            CPLError(CE_Failure, CPLE_AppDefined, "Page size is too big at %d",
1105
0
                     img.pageSizeBytes);
1106
0
            return CE_Failure;
1107
0
        }
1108
7
        dst.size =
1109
7
            img.pageSizeBytes +
1110
7
            1440;  // in case the packed page is a bit larger than the raw one
1111
7
        dst.buffer = static_cast<char *>(VSIMalloc(dst.size));
1112
7
        if (nullptr == dst.buffer)
1113
0
        {
1114
0
            CPLFree(data);
1115
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "Cannot allocate %d bytes",
1116
0
                     static_cast<int>(dst.size));
1117
0
            return CE_Failure;
1118
0
        }
1119
1120
7
        if (ZUnPack(src, dst, deflate_flags))
1121
1
        {  // Got it unpacked, update the pointers
1122
1
            CPLFree(data);
1123
1
            data = dst.buffer;
1124
1
            tinfo.size = dst.size;
1125
1
        }
1126
6
        else
1127
6
        {  // assume the page was not gzipped, warn only
1128
6
            CPLFree(dst.buffer);
1129
6
            if (!poMRFDS->no_errors)
1130
6
                CPLError(CE_Warning, CPLE_AppDefined, "Can't inflate page!");
1131
6
        }
1132
7
    }
1133
1134
38.6k
#if defined(ZSTD_SUPPORT)
1135
    // undo ZSTD
1136
38.6k
    else if (dozstd)
1137
0
    {
1138
0
        auto ctx = poMRFDS->getzsd();
1139
0
        if (!ctx)
1140
0
        {
1141
0
            CPLFree(data);
1142
0
            CPLError(CE_Failure, CPLE_AppDefined, "Can't acquire ZSTD context");
1143
0
            return CE_Failure;
1144
0
        }
1145
0
        if (img.pageSizeBytes > INT_MAX - 1440)
1146
0
        {
1147
0
            CPLFree(data);
1148
0
            CPLError(CE_Failure, CPLE_AppDefined, "Page is too large at %d",
1149
0
                     img.pageSizeBytes);
1150
0
            return CE_Failure;
1151
0
        }
1152
0
        dst.size =
1153
0
            img.pageSizeBytes +
1154
0
            1440;  // Allow for a slight increase from previous compressions
1155
0
        dst.buffer = static_cast<char *>(VSIMalloc(dst.size));
1156
0
        if (nullptr == dst.buffer)
1157
0
        {
1158
0
            CPLFree(data);
1159
0
            CPLError(CE_Failure, CPLE_OutOfMemory, "Cannot allocate %d bytes",
1160
0
                     static_cast<int>(dst.size));
1161
0
            return CE_Failure;
1162
0
        }
1163
1164
0
        auto raw_size = ZSTD_decompressDCtx(ctx, dst.buffer, dst.size,
1165
0
                                            src.buffer, src.size);
1166
0
        if (ZSTD_isError(raw_size))
1167
0
        {  // assume page was not packed, warn only
1168
0
            CPLFree(dst.buffer);
1169
0
            if (!poMRFDS->no_errors)
1170
0
                CPLError(CE_Warning, CPLE_AppDefined,
1171
0
                         "Can't unpack ZSTD page!");
1172
0
        }
1173
0
        else
1174
0
        {
1175
0
            CPLFree(data);  // The compressed data
1176
0
            data = dst.buffer;
1177
0
            tinfo.size = raw_size;
1178
            // Might need to undo the rank sort
1179
0
            size_t ranks = 0;
1180
0
            if (img.comp == IL_NONE || img.comp == IL_ZSTD)
1181
0
                ranks = static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) *
1182
0
                        img.pagesize.c;
1183
0
            if (ranks)
1184
0
            {
1185
0
                src.buffer = static_cast<char *>(data);
1186
0
                src.size = static_cast<size_t>(tinfo.size);
1187
0
                derank(src, ranks);
1188
0
            }
1189
0
        }
1190
0
    }
1191
38.6k
#endif
1192
1193
38.6k
    src.buffer = static_cast<char *>(data);
1194
38.6k
    src.size = static_cast<size_t>(tinfo.size);
1195
1196
    // After unpacking, the size has to be pageSizeBytes
1197
    // If pages are interleaved, use the dataset page buffer instead
1198
38.6k
    dst.buffer = reinterpret_cast<char *>(
1199
38.6k
        (1 == cstride) ? buffer : poMRFDS->GetPBuffer());
1200
38.6k
    dst.size = img.pageSizeBytes;
1201
1202
38.6k
    if (poMRFDS->no_errors)
1203
0
        CPLPushErrorHandler(CPLQuietErrorHandler);
1204
38.6k
    CPLErr ret = Decompress(dst, src);
1205
1206
38.6k
    poMRFDS->read_timer +=
1207
38.6k
        duration_cast<nanoseconds>(steady_clock::now() - start_time);
1208
1209
38.6k
    dst.size =
1210
38.6k
        img.pageSizeBytes;  // In case the decompress failed, force it back
1211
1212
    // Swap whatever we decompressed if we need to
1213
38.6k
    if (is_Endianness_Dependent(img.dt, img.comp) && (img.nbo != NET_ORDER))
1214
0
        swab_buff(dst, img);
1215
1216
38.6k
    CPLFree(data);
1217
38.6k
    if (poMRFDS->no_errors)
1218
0
    {
1219
0
        CPLPopErrorHandler();
1220
0
        if (ret != CE_None)  // Set each page buffer to the correct no data
1221
                             // value, then proceed
1222
0
            return (1 == cstride) ? FillBlock(buffer)
1223
0
                                  : FillBlock(xblk, yblk, buffer);
1224
0
    }
1225
1226
    // If pages are separate or we had errors, we're done
1227
38.6k
    if (1 == cstride || CE_None != ret)
1228
35.7k
        return ret;
1229
1230
    // De-interleave page from dataset buffer and return
1231
2.86k
    return ReadInterleavedBlock(xblk, yblk, buffer);
1232
38.6k
}
1233
1234
/**
1235
 *\brief Write a block from the provided buffer
1236
 *
1237
 * Same trick as read, use a temporary tile buffer for pixel interleave
1238
 * For band separate, use a
1239
 * Write the block once it has all the bands, report
1240
 * if a new block is started before the old one was completed
1241
 *
1242
 */
1243
1244
CPLErr MRFRasterBand::IWriteBlock(int xblk, int yblk, void *buffer)
1245
0
{
1246
0
    GInt32 cstride = img.pagesize.c;
1247
0
    ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
1248
0
    GUIntBig infooffset = IdxOffset(req, img);
1249
1250
0
    CPLDebug("MRF_IB", "IWriteBlock %d,%d,0,%d, level %d, stride %d\n", xblk,
1251
0
             yblk, nBand, m_l, cstride);
1252
1253
    // Finish the Create call
1254
0
    if (!poMRFDS->bCrystalized && !poMRFDS->Crystalize())
1255
0
    {
1256
0
        CPLError(CE_Failure, CPLE_AppDefined, "MRF: Error creating files");
1257
0
        return CE_Failure;
1258
0
    }
1259
1260
0
    if (1 == cstride)
1261
0
    {  // Separate bands, we can write it as is
1262
        // Empty page skip
1263
0
        int success;
1264
0
        double val = GetNoDataValue(&success);
1265
0
        if (!success)
1266
0
            val = 0.0;
1267
0
        if (isAllVal(eDataType, buffer, img.pageSizeBytes, val))
1268
0
            return poMRFDS->WriteTile(nullptr, infooffset, 0);
1269
1270
        // Use the pbuffer to hold the compressed page before writing it
1271
0
        poMRFDS->tile = ILSize();  // Mark it corrupt
1272
1273
0
        buf_mgr src;
1274
0
        src.buffer = (char *)buffer;
1275
0
        src.size = static_cast<size_t>(img.pageSizeBytes);
1276
0
        buf_mgr dst = {(char *)poMRFDS->GetPBuffer(),
1277
0
                       poMRFDS->GetPBufferSize()};
1278
1279
        // Swab the source before encoding if we need to
1280
0
        if (is_Endianness_Dependent(img.dt, img.comp) && (img.nbo != NET_ORDER))
1281
0
            swab_buff(src, img);
1282
1283
0
        auto start_time = steady_clock::now();
1284
1285
        // Compress functions need to return the compressed size in
1286
        // the bytes in buffer field
1287
0
        if (Compress(dst, src) != CE_None)
1288
0
            return CE_Failure;
1289
1290
0
        void *usebuff = dst.buffer;
1291
0
        if (dodeflate)
1292
0
        {
1293
0
            CPLAssert(dst.size <= poMRFDS->pbsize);
1294
0
            usebuff =
1295
0
                DeflateBlock(dst, poMRFDS->pbsize - dst.size, deflate_flags);
1296
0
            if (!usebuff)
1297
0
            {
1298
0
                CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
1299
0
                return CE_Failure;
1300
0
            }
1301
0
        }
1302
1303
0
#if defined(ZSTD_SUPPORT)
1304
0
        else if (dozstd)
1305
0
        {
1306
0
            size_t ranks = 0;  // Assume no need for byte rank sort
1307
0
            if (img.comp == IL_NONE || img.comp == IL_ZSTD)
1308
0
                ranks = static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt));
1309
0
            usebuff = ZstdCompBlock(dst, poMRFDS->pbsize - dst.size, zstd_level,
1310
0
                                    poMRFDS->getzsc(), ranks);
1311
0
            if (!usebuff)
1312
0
            {
1313
0
                CPLError(CE_Failure, CPLE_AppDefined,
1314
0
                         "MRF: Zstd Compression error");
1315
0
                return CE_Failure;
1316
0
            }
1317
0
        }
1318
0
#endif
1319
0
        poMRFDS->write_timer +=
1320
0
            duration_cast<nanoseconds>(steady_clock::now() - start_time);
1321
0
        return poMRFDS->WriteTile(usebuff, infooffset, dst.size);
1322
0
    }
1323
1324
    // Multiple bands per page, use a temporary to assemble the page
1325
    // Temporary is large because we use it to hold both the uncompressed and
1326
    // the compressed
1327
0
    poMRFDS->tile = req;
1328
0
    poMRFDS->bdirty = 0;
1329
1330
    // Keep track of what bands are empty
1331
0
    GUIntBig empties = 0;
1332
1333
0
    void *tbuffer = VSIMalloc(img.pageSizeBytes + poMRFDS->pbsize);
1334
1335
0
    if (!tbuffer)
1336
0
    {
1337
0
        CPLError(CE_Failure, CPLE_AppDefined,
1338
0
                 "MRF: Can't allocate write buffer");
1339
0
        return CE_Failure;
1340
0
    }
1341
1342
    // Get the other bands from the block cache
1343
0
    for (int iBand = 0; iBand < poMRFDS->nBands; iBand++)
1344
0
    {
1345
0
        char *pabyThisImage = nullptr;
1346
0
        GDALRasterBlock *poBlock = nullptr;
1347
1348
0
        if (iBand == nBand - 1)
1349
0
        {
1350
0
            pabyThisImage = reinterpret_cast<char *>(buffer);
1351
0
            poMRFDS->bdirty |= bandbit();
1352
0
        }
1353
0
        else
1354
0
        {
1355
0
            GDALRasterBand *band = poMRFDS->GetRasterBand(iBand + 1);
1356
            // Pick the right overview
1357
0
            if (m_l)
1358
0
                band = band->GetOverview(m_l - 1);
1359
0
            poBlock = (reinterpret_cast<MRFRasterBand *>(band))
1360
0
                          ->TryGetLockedBlockRef(xblk, yblk);
1361
0
            if (nullptr == poBlock)
1362
0
                continue;
1363
            // This is where the image data is for this band
1364
1365
0
            pabyThisImage = reinterpret_cast<char *>(poBlock->GetDataRef());
1366
0
            poMRFDS->bdirty |= bandbit(iBand);
1367
0
        }
1368
1369
        // Keep track of empty bands, but encode them anyhow, in case some are
1370
        // not empty
1371
0
        int success;
1372
0
        double val = GetNoDataValue(&success);
1373
0
        if (!success)
1374
0
            val = 0.0;
1375
0
        if (isAllVal(eDataType, pabyThisImage, blockSizeBytes(), val))
1376
0
            empties |= bandbit(iBand);
1377
1378
            // Copy the data into the dataset buffer here
1379
            // Just the right mix of templates and macros make this real tidy
1380
0
#define CpySO(T)                                                               \
1381
0
    cpy_stride_out<T>((reinterpret_cast<T *>(tbuffer)) + iBand, pabyThisImage, \
1382
0
                      blockSizeBytes() / sizeof(T), cstride)
1383
1384
        // Build the page in tbuffer
1385
0
        switch (GDALGetDataTypeSizeBytes(eDataType))
1386
0
        {
1387
0
            case 1:
1388
0
                CpySO(GByte);
1389
0
                break;
1390
0
            case 2:
1391
0
                CpySO(GInt16);
1392
0
                break;
1393
0
            case 4:
1394
0
                CpySO(GInt32);
1395
0
                break;
1396
0
            case 8:
1397
0
                CpySO(GIntBig);
1398
0
                break;
1399
0
            default:
1400
0
            {
1401
0
                CPLError(CE_Failure, CPLE_AppDefined,
1402
0
                         "MRF: Write datatype of %d bytes "
1403
0
                         "not implemented",
1404
0
                         GDALGetDataTypeSizeBytes(eDataType));
1405
0
                if (poBlock != nullptr)
1406
0
                {
1407
0
                    poBlock->MarkClean();
1408
0
                    poBlock->DropLock();
1409
0
                }
1410
0
                CPLFree(tbuffer);
1411
0
                return CE_Failure;
1412
0
            }
1413
0
        }
1414
1415
0
        if (poBlock != nullptr)
1416
0
        {
1417
0
            poBlock->MarkClean();
1418
0
            poBlock->DropLock();
1419
0
        }
1420
0
    }
1421
1422
    // Should keep track of the individual band buffers and only mix them if
1423
    // this is not an empty page ( move the Copy with Stride Out from above
1424
    // below this test This way works fine, but it does work extra for empty
1425
    // pages
1426
1427
0
    if (GIntBig(empties) == AllBandMask())
1428
0
    {
1429
0
        CPLFree(tbuffer);
1430
0
        return poMRFDS->WriteTile(nullptr, infooffset, 0);
1431
0
    }
1432
1433
0
    if (poMRFDS->bdirty != AllBandMask())
1434
0
        CPLError(CE_Warning, CPLE_AppDefined,
1435
0
                 "MRF: IWrite, band dirty mask is " CPL_FRMT_GIB
1436
0
                 " instead of " CPL_FRMT_GIB,
1437
0
                 poMRFDS->bdirty, AllBandMask());
1438
1439
0
    buf_mgr src;
1440
0
    src.buffer = (char *)tbuffer;
1441
0
    src.size = static_cast<size_t>(img.pageSizeBytes);
1442
1443
    // Use the space after pagesizebytes for compressed output, it is of pbsize
1444
0
    char *outbuff = (char *)tbuffer + img.pageSizeBytes;
1445
1446
0
    buf_mgr dst = {outbuff, poMRFDS->pbsize};
1447
0
    CPLErr ret;
1448
1449
0
    auto start_time = steady_clock::now();
1450
1451
0
    ret = Compress(dst, src);
1452
0
    if (ret != CE_None)
1453
0
    {
1454
        // Compress failed, write it as an empty tile
1455
0
        CPLFree(tbuffer);
1456
0
        poMRFDS->WriteTile(nullptr, infooffset, 0);
1457
0
        return CE_None;  // Should report the error, but it triggers partial
1458
                         // band attempts
1459
0
    }
1460
1461
    // Where the output is, in case we deflate
1462
0
    void *usebuff = outbuff;
1463
0
    if (dodeflate)
1464
0
    {
1465
        // Move the packed part at the start of tbuffer, to make more space
1466
        // available
1467
0
        memcpy(tbuffer, outbuff, dst.size);
1468
0
        dst.buffer = static_cast<char *>(tbuffer);
1469
0
        usebuff = DeflateBlock(dst,
1470
0
                               static_cast<size_t>(img.pageSizeBytes) +
1471
0
                                   poMRFDS->pbsize - dst.size,
1472
0
                               deflate_flags);
1473
0
        if (!usebuff)
1474
0
            CPLError(CE_Failure, CPLE_AppDefined, "MRF: Deflate error");
1475
0
    }
1476
1477
0
#if defined(ZSTD_SUPPORT)
1478
0
    else if (dozstd)
1479
0
    {
1480
0
        memcpy(tbuffer, outbuff, dst.size);
1481
0
        dst.buffer = static_cast<char *>(tbuffer);
1482
0
        size_t ranks = 0;  // Assume no need for byte rank sort
1483
0
        if (img.comp == IL_NONE || img.comp == IL_ZSTD)
1484
0
            ranks =
1485
0
                static_cast<size_t>(GDALGetDataTypeSizeBytes(img.dt)) * cstride;
1486
0
        usebuff = ZstdCompBlock(dst,
1487
0
                                static_cast<size_t>(img.pageSizeBytes) +
1488
0
                                    poMRFDS->pbsize - dst.size,
1489
0
                                zstd_level, poMRFDS->getzsc(), ranks);
1490
0
        if (!usebuff)
1491
0
            CPLError(CE_Failure, CPLE_AppDefined,
1492
0
                     "MRF: ZStd compression error");
1493
0
    }
1494
0
#endif
1495
1496
0
    poMRFDS->write_timer +=
1497
0
        duration_cast<nanoseconds>(steady_clock::now() - start_time);
1498
1499
0
    if (!usebuff)
1500
0
    {  // Error was signaled
1501
0
        CPLFree(tbuffer);
1502
0
        poMRFDS->WriteTile(nullptr, infooffset, 0);
1503
0
        poMRFDS->bdirty = 0;
1504
0
        return CE_Failure;
1505
0
    }
1506
1507
0
    ret = poMRFDS->WriteTile(usebuff, infooffset, dst.size);
1508
0
    CPLFree(tbuffer);
1509
1510
0
    poMRFDS->bdirty = 0;
1511
0
    return ret;
1512
0
}
1513
1514
//
1515
// Tests if a given block exists without reading it
1516
// returns false only when it is definitely not existing
1517
//
1518
bool MRFRasterBand::TestBlock(int xblk, int yblk)
1519
0
{
1520
    // When bypassing the cache, assume all blocks are valid
1521
0
    if (poMRFDS->bypass_cache && !poMRFDS->source.empty())
1522
0
        return true;
1523
1524
    // Blocks outside of image have no data by default
1525
0
    if (xblk < 0 || yblk < 0 || xblk >= img.pagecount.x ||
1526
0
        yblk >= img.pagecount.y)
1527
0
        return false;
1528
1529
0
    ILIdx tinfo;
1530
0
    GInt32 cstride = img.pagesize.c;
1531
0
    ILSize req(xblk, yblk, 0, (nBand - 1) / cstride, m_l);
1532
1533
0
    if (CE_None != poMRFDS->ReadTileIdx(tinfo, req, img))
1534
        // Got an error reading the tile index
1535
0
        return !poMRFDS->no_errors;
1536
1537
    // Got an index, if the size is readable, the block does exist
1538
0
    if (0 < tinfo.size && tinfo.size < poMRFDS->pbsize * 2)
1539
0
        return true;
1540
1541
    // We are caching, but the tile has not been checked, so it could exist
1542
0
    return (!poMRFDS->source.empty() && 0 == tinfo.offset);
1543
0
}
1544
1545
int MRFRasterBand::GetOverviewCount()
1546
684k
{
1547
    // First try internal overviews
1548
684k
    int nInternalOverviewCount = static_cast<int>(overviews.size());
1549
684k
    if (nInternalOverviewCount > 0)
1550
0
        return nInternalOverviewCount;
1551
684k
    return GDALPamRasterBand::GetOverviewCount();
1552
684k
}
1553
1554
GDALRasterBand *MRFRasterBand::GetOverview(int n)
1555
1.22k
{
1556
    // First try internal overviews
1557
1.22k
    if (n >= 0 && n < (int)overviews.size())
1558
0
        return overviews[n];
1559
1.22k
    return GDALPamRasterBand::GetOverview(n);
1560
1.22k
}
1561
1562
CPLErr Raw_Band::Decompress(buf_mgr &dst, buf_mgr &src)
1563
10.4k
{
1564
10.4k
    if (src.size > dst.size)
1565
1
        return CE_Failure;
1566
10.4k
    memcpy(dst.buffer, src.buffer, src.size);
1567
10.4k
    dst.size = src.size;
1568
10.4k
    return CE_None;
1569
10.4k
}
1570
1571
CPLErr MRFLRasterBand::IReadBlock(int xblk, int yblk, void *buffer)
1572
0
{
1573
0
    return pBand->IReadBlock(xblk, yblk, buffer);
1574
0
}
1575
1576
NAMESPACE_MRF_END