Coverage Report

Created: 2025-06-09 08:44

/src/gdal/alg/gdalchecksum.cpp
Line
Count
Source (jump to first uncovered line)
1
/******************************************************************************
2
 *
3
 * Project:  GDAL
4
 * Purpose:  Compute simple checksum for a region of image data.
5
 * Author:   Frank Warmerdam, warmerdam@pobox.com
6
 *
7
 ******************************************************************************
8
 * Copyright (c) 2003, Frank Warmerdam
9
 * Copyright (c) 2007-2008, Even Rouault <even dot rouault at spatialys.com>
10
 *
11
 * SPDX-License-Identifier: MIT
12
 ****************************************************************************/
13
14
#include "cpl_port.h"
15
#include "gdal_alg.h"
16
17
#include <cmath>
18
#include <cstddef>
19
#include <algorithm>
20
21
#include "cpl_conv.h"
22
#include "cpl_error.h"
23
#include "cpl_vsi.h"
24
#include "gdal.h"
25
#include "gdal_priv.h"
26
27
/************************************************************************/
28
/*                         GDALChecksumImage()                          */
29
/************************************************************************/
30
31
/**
32
 * Compute checksum for image region.
33
 *
34
 * Computes a 16bit (0-65535) checksum from a region of raster data on a GDAL
35
 * supported band.   Floating point data is converted to 32bit integer
36
 * so decimal portions of such raster data will not affect the checksum.
37
 * Real and Imaginary components of complex bands influence the result.
38
 *
39
 * @param hBand the raster band to read from.
40
 * @param nXOff pixel offset of window to read.
41
 * @param nYOff line offset of window to read.
42
 * @param nXSize pixel size of window to read.
43
 * @param nYSize line size of window to read.
44
 *
45
 * @return Checksum value, or -1 in case of error (starting with GDAL 3.6)
46
 */
47
48
int CPL_STDCALL GDALChecksumImage(GDALRasterBandH hBand, int nXOff, int nYOff,
49
                                  int nXSize, int nYSize)
50
51
457k
{
52
457k
    VALIDATE_POINTER1(hBand, "GDALChecksumImage", 0);
53
54
457k
    const static int anPrimes[11] = {7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43};
55
56
457k
    int nChecksum = 0;
57
457k
    int iPrime = 0;
58
457k
    const GDALDataType eDataType = GDALGetRasterDataType(hBand);
59
457k
    const bool bComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eDataType));
60
457k
    const bool bIsFloatingPoint =
61
457k
        (eDataType == GDT_Float16 || eDataType == GDT_Float32 ||
62
457k
         eDataType == GDT_Float64 || eDataType == GDT_CFloat16 ||
63
457k
         eDataType == GDT_CFloat32 || eDataType == GDT_CFloat64);
64
65
457k
    const auto IntFromDouble = [](double dfVal)
66
693M
    {
67
693M
        int nVal;
68
693M
        if (!std::isfinite(dfVal))
69
81.9M
        {
70
81.9M
            nVal = INT_MIN;
71
81.9M
        }
72
611M
        else
73
611M
        {
74
            // Standard behavior of GDALCopyWords when converting
75
            // from floating point to Int32.
76
611M
            dfVal += 0.5;
77
78
611M
            if (dfVal < -2147483647.0)
79
157M
                nVal = -2147483647;
80
454M
            else if (dfVal > 2147483647)
81
4.78M
                nVal = 2147483647;
82
449M
            else
83
449M
                nVal = static_cast<GInt32>(floor(dfVal));
84
611M
        }
85
693M
        return nVal;
86
693M
    };
87
88
457k
    if (bIsFloatingPoint && nXOff == 0 && nYOff == 0)
89
35.4k
    {
90
35.4k
        const GDALDataType eDstDataType = bComplex ? GDT_CFloat64 : GDT_Float64;
91
35.4k
        int nBlockXSize = 0;
92
35.4k
        int nBlockYSize = 0;
93
35.4k
        GDALGetBlockSize(hBand, &nBlockXSize, &nBlockYSize);
94
35.4k
        const int nDstDataTypeSize = GDALGetDataTypeSizeBytes(eDstDataType);
95
35.4k
        int nChunkXSize = nBlockXSize;
96
35.4k
        const int nChunkYSize = nBlockYSize;
97
35.4k
        if (nBlockXSize < nXSize)
98
2.98k
        {
99
2.98k
            const GIntBig nMaxChunkSize =
100
2.98k
                std::max(static_cast<GIntBig>(10 * 1000 * 1000),
101
2.98k
                         GDALGetCacheMax64() / 10);
102
2.98k
            if (nDstDataTypeSize > 0 &&
103
2.98k
                static_cast<GIntBig>(nXSize) * nChunkYSize <
104
2.98k
                    nMaxChunkSize / nDstDataTypeSize)
105
2.98k
            {
106
                // A full line of height nChunkYSize can fit in the maximum
107
                // allowed memory
108
2.98k
                nChunkXSize = nXSize;
109
2.98k
            }
110
0
            else if (nDstDataTypeSize > 0)
111
0
            {
112
                // Otherwise compute a size that is a multiple of nBlockXSize
113
0
                nChunkXSize = static_cast<int>(std::min(
114
0
                    static_cast<GIntBig>(nXSize),
115
0
                    nBlockXSize *
116
0
                        std::max(static_cast<GIntBig>(1),
117
0
                                 nMaxChunkSize /
118
0
                                     (static_cast<GIntBig>(nBlockXSize) *
119
0
                                      nChunkYSize * nDstDataTypeSize))));
120
0
            }
121
2.98k
        }
122
123
35.4k
        double *padfLineData = static_cast<double *>(
124
35.4k
            VSI_MALLOC3_VERBOSE(nChunkXSize, nChunkYSize, nDstDataTypeSize));
125
35.4k
        if (padfLineData == nullptr)
126
0
        {
127
0
            return -1;
128
0
        }
129
35.4k
        const int nValsPerIter = bComplex ? 2 : 1;
130
131
35.4k
        const int nYBlocks = DIV_ROUND_UP(nYSize, nChunkYSize);
132
35.4k
        const int nXBlocks = DIV_ROUND_UP(nXSize, nChunkXSize);
133
4.12M
        for (int iYBlock = 0; iYBlock < nYBlocks; ++iYBlock)
134
4.10M
        {
135
4.10M
            const int iYStart = iYBlock * nChunkYSize;
136
4.10M
            const int iYEnd =
137
4.10M
                iYBlock == nYBlocks - 1 ? nYSize : iYStart + nChunkYSize;
138
4.10M
            const int nChunkActualHeight = iYEnd - iYStart;
139
8.19M
            for (int iXBlock = 0; iXBlock < nXBlocks; ++iXBlock)
140
4.10M
            {
141
4.10M
                const int iXStart = iXBlock * nChunkXSize;
142
4.10M
                const int iXEnd =
143
4.10M
                    iXBlock == nXBlocks - 1 ? nXSize : iXStart + nChunkXSize;
144
4.10M
                const int nChunkActualXSize = iXEnd - iXStart;
145
4.10M
                if (GDALRasterIO(
146
4.10M
                        hBand, GF_Read, iXStart, iYStart, nChunkActualXSize,
147
4.10M
                        nChunkActualHeight, padfLineData, nChunkActualXSize,
148
4.10M
                        nChunkActualHeight, eDstDataType, 0, 0) != CE_None)
149
13.5k
                {
150
13.5k
                    CPLError(CE_Failure, CPLE_FileIO,
151
13.5k
                             "Checksum value could not be computed due to I/O "
152
13.5k
                             "read error.");
153
13.5k
                    nChecksum = -1;
154
13.5k
                    break;
155
13.5k
                }
156
4.08M
                const size_t xIters =
157
4.08M
                    static_cast<size_t>(nValsPerIter) * nChunkActualXSize;
158
10.0M
                for (int iY = iYStart; iY < iYEnd; ++iY)
159
5.92M
                {
160
                    // Initialize iPrime so that it is consistent with a
161
                    // per full line iteration strategy
162
5.92M
                    iPrime = (nValsPerIter *
163
5.92M
                              (static_cast<int64_t>(iY) * nXSize + iXStart)) %
164
5.92M
                             11;
165
5.92M
                    const size_t nOffset = nValsPerIter *
166
5.92M
                                           static_cast<size_t>(iY - iYStart) *
167
5.92M
                                           nChunkActualXSize;
168
699M
                    for (size_t i = 0; i < xIters; ++i)
169
693M
                    {
170
693M
                        const double dfVal = padfLineData[nOffset + i];
171
693M
                        nChecksum += IntFromDouble(dfVal) % anPrimes[iPrime++];
172
693M
                        if (iPrime > 10)
173
63.0M
                            iPrime = 0;
174
693M
                    }
175
5.92M
                    nChecksum &= 0xffff;
176
5.92M
                }
177
4.08M
            }
178
179
4.10M
            if (nChecksum < 0)
180
13.5k
                break;
181
4.10M
        }
182
183
35.4k
        CPLFree(padfLineData);
184
35.4k
    }
185
422k
    else if (bIsFloatingPoint)
186
0
    {
187
0
        const GDALDataType eDstDataType = bComplex ? GDT_CFloat64 : GDT_Float64;
188
189
0
        double *padfLineData = static_cast<double *>(VSI_MALLOC2_VERBOSE(
190
0
            nXSize, GDALGetDataTypeSizeBytes(eDstDataType)));
191
0
        if (padfLineData == nullptr)
192
0
        {
193
0
            return -1;
194
0
        }
195
196
0
        for (int iLine = nYOff; iLine < nYOff + nYSize; iLine++)
197
0
        {
198
0
            if (GDALRasterIO(hBand, GF_Read, nXOff, iLine, nXSize, 1,
199
0
                             padfLineData, nXSize, 1, eDstDataType, 0,
200
0
                             0) != CE_None)
201
0
            {
202
0
                CPLError(CE_Failure, CPLE_FileIO,
203
0
                         "Checksum value couldn't be computed due to "
204
0
                         "I/O read error.");
205
0
                nChecksum = -1;
206
0
                break;
207
0
            }
208
0
            const size_t nCount = bComplex ? static_cast<size_t>(nXSize) * 2
209
0
                                           : static_cast<size_t>(nXSize);
210
211
0
            for (size_t i = 0; i < nCount; i++)
212
0
            {
213
0
                const double dfVal = padfLineData[i];
214
0
                nChecksum += IntFromDouble(dfVal) % anPrimes[iPrime++];
215
0
                if (iPrime > 10)
216
0
                    iPrime = 0;
217
218
0
                nChecksum &= 0xffff;
219
0
            }
220
0
        }
221
222
0
        CPLFree(padfLineData);
223
0
    }
224
422k
    else if (nXOff == 0 && nYOff == 0)
225
422k
    {
226
422k
        const GDALDataType eDstDataType = bComplex ? GDT_CInt32 : GDT_Int32;
227
422k
        int nBlockXSize = 0;
228
422k
        int nBlockYSize = 0;
229
422k
        GDALGetBlockSize(hBand, &nBlockXSize, &nBlockYSize);
230
422k
        const int nDstDataTypeSize = GDALGetDataTypeSizeBytes(eDstDataType);
231
422k
        int nChunkXSize = nBlockXSize;
232
422k
        const int nChunkYSize = nBlockYSize;
233
422k
        if (nBlockXSize < nXSize)
234
45.6k
        {
235
45.6k
            const GIntBig nMaxChunkSize =
236
45.6k
                std::max(static_cast<GIntBig>(10 * 1000 * 1000),
237
45.6k
                         GDALGetCacheMax64() / 10);
238
45.6k
            if (nDstDataTypeSize > 0 &&
239
45.6k
                static_cast<GIntBig>(nXSize) * nChunkYSize <
240
45.6k
                    nMaxChunkSize / nDstDataTypeSize)
241
45.6k
            {
242
                // A full line of height nChunkYSize can fit in the maximum
243
                // allowed memory
244
45.6k
                nChunkXSize = nXSize;
245
45.6k
            }
246
1
            else if (nDstDataTypeSize > 0)
247
1
            {
248
                // Otherwise compute a size that is a multiple of nBlockXSize
249
1
                nChunkXSize = static_cast<int>(std::min(
250
1
                    static_cast<GIntBig>(nXSize),
251
1
                    nBlockXSize *
252
1
                        std::max(static_cast<GIntBig>(1),
253
1
                                 nMaxChunkSize /
254
1
                                     (static_cast<GIntBig>(nBlockXSize) *
255
1
                                      nChunkYSize * nDstDataTypeSize))));
256
1
            }
257
45.6k
        }
258
259
422k
        int *panChunkData = static_cast<GInt32 *>(
260
422k
            VSI_MALLOC3_VERBOSE(nChunkXSize, nChunkYSize, nDstDataTypeSize));
261
422k
        if (panChunkData == nullptr)
262
0
        {
263
0
            return -1;
264
0
        }
265
422k
        const int nValsPerIter = bComplex ? 2 : 1;
266
267
422k
        const int nYBlocks = DIV_ROUND_UP(nYSize, nChunkYSize);
268
422k
        const int nXBlocks = DIV_ROUND_UP(nXSize, nChunkXSize);
269
24.3M
        for (int iYBlock = 0; iYBlock < nYBlocks; ++iYBlock)
270
24.2M
        {
271
24.2M
            const int iYStart = iYBlock * nChunkYSize;
272
24.2M
            const int iYEnd =
273
24.2M
                iYBlock == nYBlocks - 1 ? nYSize : iYStart + nChunkYSize;
274
24.2M
            const int nChunkActualHeight = iYEnd - iYStart;
275
48.1M
            for (int iXBlock = 0; iXBlock < nXBlocks; ++iXBlock)
276
24.2M
            {
277
24.2M
                const int iXStart = iXBlock * nChunkXSize;
278
24.2M
                const int iXEnd =
279
24.2M
                    iXBlock == nXBlocks - 1 ? nXSize : iXStart + nChunkXSize;
280
24.2M
                const int nChunkActualXSize = iXEnd - iXStart;
281
24.2M
                if (GDALRasterIO(
282
24.2M
                        hBand, GF_Read, iXStart, iYStart, nChunkActualXSize,
283
24.2M
                        nChunkActualHeight, panChunkData, nChunkActualXSize,
284
24.2M
                        nChunkActualHeight, eDstDataType, 0, 0) != CE_None)
285
244k
                {
286
244k
                    CPLError(CE_Failure, CPLE_FileIO,
287
244k
                             "Checksum value could not be computed due to I/O "
288
244k
                             "read error.");
289
244k
                    nChecksum = -1;
290
244k
                    break;
291
244k
                }
292
23.9M
                const size_t xIters =
293
23.9M
                    static_cast<size_t>(nValsPerIter) * nChunkActualXSize;
294
54.5M
                for (int iY = iYStart; iY < iYEnd; ++iY)
295
30.6M
                {
296
                    // Initialize iPrime so that it is consistent with a
297
                    // per full line iteration strategy
298
30.6M
                    iPrime = (nValsPerIter *
299
30.6M
                              (static_cast<int64_t>(iY) * nXSize + iXStart)) %
300
30.6M
                             11;
301
30.6M
                    const size_t nOffset = nValsPerIter *
302
30.6M
                                           static_cast<size_t>(iY - iYStart) *
303
30.6M
                                           nChunkActualXSize;
304
8.53G
                    for (size_t i = 0; i < xIters; ++i)
305
8.50G
                    {
306
8.50G
                        nChecksum +=
307
8.50G
                            panChunkData[nOffset + i] % anPrimes[iPrime++];
308
8.50G
                        if (iPrime > 10)
309
773M
                            iPrime = 0;
310
8.50G
                    }
311
30.6M
                    nChecksum &= 0xffff;
312
30.6M
                }
313
23.9M
            }
314
315
24.2M
            if (nChecksum < 0)
316
244k
                break;
317
24.2M
        }
318
319
422k
        CPLFree(panChunkData);
320
422k
    }
321
0
    else
322
0
    {
323
0
        const GDALDataType eDstDataType = bComplex ? GDT_CInt32 : GDT_Int32;
324
325
0
        int *panLineData = static_cast<GInt32 *>(VSI_MALLOC2_VERBOSE(
326
0
            nXSize, GDALGetDataTypeSizeBytes(eDstDataType)));
327
0
        if (panLineData == nullptr)
328
0
        {
329
0
            return -1;
330
0
        }
331
332
0
        for (int iLine = nYOff; iLine < nYOff + nYSize; iLine++)
333
0
        {
334
0
            if (GDALRasterIO(hBand, GF_Read, nXOff, iLine, nXSize, 1,
335
0
                             panLineData, nXSize, 1, eDstDataType, 0,
336
0
                             0) != CE_None)
337
0
            {
338
0
                CPLError(CE_Failure, CPLE_FileIO,
339
0
                         "Checksum value could not be computed due to I/O "
340
0
                         "read error.");
341
0
                nChecksum = -1;
342
0
                break;
343
0
            }
344
0
            const size_t nCount = bComplex ? static_cast<size_t>(nXSize) * 2
345
0
                                           : static_cast<size_t>(nXSize);
346
347
0
            for (size_t i = 0; i < nCount; i++)
348
0
            {
349
0
                nChecksum += panLineData[i] % anPrimes[iPrime++];
350
0
                if (iPrime > 10)
351
0
                    iPrime = 0;
352
353
0
                nChecksum &= 0xffff;
354
0
            }
355
0
        }
356
357
0
        CPLFree(panLineData);
358
0
    }
359
360
    // coverity[return_overflow]
361
457k
    return nChecksum;
362
457k
}