Line data Source code
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 : * Permission is hereby granted, free of charge, to any person obtaining a
12 : * copy of this software and associated documentation files (the "Software"),
13 : * to deal in the Software without restriction, including without limitation
14 : * the rights to use, copy, modify, merge, publish, distribute, sublicense,
15 : * and/or sell copies of the Software, and to permit persons to whom the
16 : * Software is furnished to do so, subject to the following conditions:
17 : *
18 : * The above copyright notice and this permission notice shall be included
19 : * in all copies or substantial portions of the Software.
20 : *
21 : * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
22 : * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
23 : * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
24 : * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
25 : * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
26 : * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
27 : * DEALINGS IN THE SOFTWARE.
28 : ****************************************************************************/
29 :
30 : #include "cpl_port.h"
31 : #include "gdal_alg.h"
32 :
33 : #include <cmath>
34 : #include <cstddef>
35 : #include <algorithm>
36 :
37 : #include "cpl_conv.h"
38 : #include "cpl_error.h"
39 : #include "cpl_vsi.h"
40 : #include "gdal.h"
41 : #include "gdal_priv.h"
42 :
43 : /************************************************************************/
44 : /* GDALChecksumImage() */
45 : /************************************************************************/
46 :
47 : /**
48 : * Compute checksum for image region.
49 : *
50 : * Computes a 16bit (0-65535) checksum from a region of raster data on a GDAL
51 : * supported band. Floating point data is converted to 32bit integer
52 : * so decimal portions of such raster data will not affect the checksum.
53 : * Real and Imaginary components of complex bands influence the result.
54 : *
55 : * @param hBand the raster band to read from.
56 : * @param nXOff pixel offset of window to read.
57 : * @param nYOff line offset of window to read.
58 : * @param nXSize pixel size of window to read.
59 : * @param nYSize line size of window to read.
60 : *
61 : * @return Checksum value, or -1 in case of error (starting with GDAL 3.6)
62 : */
63 :
64 8991 : int CPL_STDCALL GDALChecksumImage(GDALRasterBandH hBand, int nXOff, int nYOff,
65 : int nXSize, int nYSize)
66 :
67 : {
68 8991 : VALIDATE_POINTER1(hBand, "GDALChecksumImage", 0);
69 :
70 : const static int anPrimes[11] = {7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43};
71 :
72 8991 : int nChecksum = 0;
73 8991 : int iPrime = 0;
74 8991 : const GDALDataType eDataType = GDALGetRasterDataType(hBand);
75 8991 : const bool bComplex = CPL_TO_BOOL(GDALDataTypeIsComplex(eDataType));
76 8991 : const bool bIsFloatingPoint =
77 8404 : (eDataType == GDT_Float32 || eDataType == GDT_Float64 ||
78 17395 : eDataType == GDT_CFloat32 || eDataType == GDT_CFloat64);
79 :
80 75271800 : const auto IntFromDouble = [](double dfVal)
81 : {
82 : int nVal;
83 75271800 : if (CPLIsNan(dfVal) || CPLIsInf(dfVal))
84 : {
85 : // Most compilers seem to cast NaN or Inf to 0x80000000.
86 : // but VC7 is an exception. So we force the result
87 : // of such a cast.
88 18029 : nVal = 0x80000000;
89 : }
90 : else
91 : {
92 : // Standard behavior of GDALCopyWords when converting
93 : // from floating point to Int32.
94 75253800 : dfVal += 0.5;
95 :
96 75253800 : if (dfVal < -2147483647.0)
97 1215 : nVal = -2147483647;
98 75252600 : else if (dfVal > 2147483647)
99 2690 : nVal = 2147483647;
100 : else
101 75249900 : nVal = static_cast<GInt32>(floor(dfVal));
102 : }
103 75271800 : return nVal;
104 : };
105 :
106 8991 : if (bIsFloatingPoint && nXOff == 0 && nYOff == 0)
107 : {
108 1064 : const GDALDataType eDstDataType = bComplex ? GDT_CFloat64 : GDT_Float64;
109 1064 : int nBlockXSize = 0;
110 1064 : int nBlockYSize = 0;
111 1064 : GDALGetBlockSize(hBand, &nBlockXSize, &nBlockYSize);
112 1064 : const int nDstDataTypeSize = GDALGetDataTypeSizeBytes(eDstDataType);
113 1064 : int nChunkXSize = nBlockXSize;
114 1064 : const int nChunkYSize = nBlockYSize;
115 1064 : if (nBlockXSize < nXSize)
116 : {
117 : const GIntBig nMaxChunkSize =
118 96 : std::max(static_cast<GIntBig>(10 * 1000 * 1000),
119 48 : GDALGetCacheMax64() / 10);
120 48 : if (nDstDataTypeSize > 0 &&
121 48 : static_cast<GIntBig>(nXSize) * nChunkYSize <
122 48 : nMaxChunkSize / nDstDataTypeSize)
123 : {
124 : // A full line of height nChunkYSize can fit in the maximum
125 : // allowed memory
126 48 : nChunkXSize = nXSize;
127 : }
128 : else
129 : {
130 : // Otherwise compute a size that is a multiple of nBlockXSize
131 0 : nChunkXSize = static_cast<int>(std::min(
132 0 : static_cast<GIntBig>(nXSize),
133 0 : nBlockXSize *
134 0 : std::max(static_cast<GIntBig>(1),
135 0 : nMaxChunkSize /
136 0 : (static_cast<GIntBig>(nBlockXSize) *
137 0 : nChunkYSize * nDstDataTypeSize))));
138 : }
139 : }
140 :
141 : double *padfLineData = static_cast<double *>(
142 1064 : VSI_MALLOC3_VERBOSE(nChunkXSize, nChunkYSize, nDstDataTypeSize));
143 1064 : if (padfLineData == nullptr)
144 : {
145 0 : return -1;
146 : }
147 1064 : const int nValsPerIter = bComplex ? 2 : 1;
148 :
149 1064 : const int nYBlocks = DIV_ROUND_UP(nYSize, nChunkYSize);
150 1064 : const int nXBlocks = DIV_ROUND_UP(nXSize, nChunkXSize);
151 26744 : for (int iYBlock = 0; iYBlock < nYBlocks; ++iYBlock)
152 : {
153 25680 : const int iYStart = iYBlock * nChunkYSize;
154 25680 : const int iYEnd =
155 25680 : iYBlock == nYBlocks - 1 ? nYSize : iYStart + nChunkYSize;
156 25680 : const int nChunkActualHeight = iYEnd - iYStart;
157 51359 : for (int iXBlock = 0; iXBlock < nXBlocks; ++iXBlock)
158 : {
159 25680 : const int iXStart = iXBlock * nChunkXSize;
160 25680 : const int iXEnd =
161 25680 : iXBlock == nXBlocks - 1 ? nXSize : iXStart + nChunkXSize;
162 25680 : const int nChunkActualXSize = iXEnd - iXStart;
163 25680 : if (GDALRasterIO(
164 : hBand, GF_Read, iXStart, iYStart, nChunkActualXSize,
165 : nChunkActualHeight, padfLineData, nChunkActualXSize,
166 25680 : nChunkActualHeight, eDstDataType, 0, 0) != CE_None)
167 : {
168 1 : CPLError(CE_Failure, CPLE_FileIO,
169 : "Checksum value could not be computed due to I/O "
170 : "read error.");
171 1 : nChecksum = -1;
172 1 : iYBlock = nYBlocks;
173 1 : break;
174 : }
175 25679 : const size_t xIters =
176 25679 : static_cast<size_t>(nValsPerIter) * nChunkActualXSize;
177 77556 : for (int iY = iYStart; iY < iYEnd; ++iY)
178 : {
179 : // Initialize iPrime so that it is consistent with a
180 : // per full line iteration strategy
181 51877 : iPrime = (nValsPerIter *
182 51877 : (static_cast<int64_t>(iY) * nXSize + iXStart)) %
183 : 11;
184 51877 : const size_t nOffset = nValsPerIter *
185 51877 : static_cast<size_t>(iY - iYStart) *
186 51877 : nChunkActualXSize;
187 75323300 : for (size_t i = 0; i < xIters; ++i)
188 : {
189 75271400 : const double dfVal = padfLineData[nOffset + i];
190 75271400 : nChecksum += IntFromDouble(dfVal) % anPrimes[iPrime++];
191 75271400 : if (iPrime > 10)
192 6842400 : iPrime = 0;
193 : }
194 51877 : nChecksum &= 0xffff;
195 : }
196 : }
197 : }
198 :
199 1064 : CPLFree(padfLineData);
200 : }
201 7927 : else if (bIsFloatingPoint)
202 : {
203 3 : const GDALDataType eDstDataType = bComplex ? GDT_CFloat64 : GDT_Float64;
204 :
205 3 : double *padfLineData = static_cast<double *>(VSI_MALLOC2_VERBOSE(
206 : nXSize, GDALGetDataTypeSizeBytes(eDstDataType)));
207 3 : if (padfLineData == nullptr)
208 : {
209 0 : return -1;
210 : }
211 :
212 33 : for (int iLine = nYOff; iLine < nYOff + nYSize; iLine++)
213 : {
214 30 : if (GDALRasterIO(hBand, GF_Read, nXOff, iLine, nXSize, 1,
215 : padfLineData, nXSize, 1, eDstDataType, 0,
216 30 : 0) != CE_None)
217 : {
218 0 : CPLError(CE_Failure, CPLE_FileIO,
219 : "Checksum value couldn't be computed due to "
220 : "I/O read error.");
221 0 : nChecksum = -1;
222 0 : break;
223 : }
224 30 : const size_t nCount = bComplex ? static_cast<size_t>(nXSize) * 2
225 : : static_cast<size_t>(nXSize);
226 :
227 480 : for (size_t i = 0; i < nCount; i++)
228 : {
229 450 : const double dfVal = padfLineData[i];
230 450 : nChecksum += IntFromDouble(dfVal) % anPrimes[iPrime++];
231 450 : if (iPrime > 10)
232 40 : iPrime = 0;
233 :
234 450 : nChecksum &= 0xffff;
235 : }
236 : }
237 :
238 3 : CPLFree(padfLineData);
239 : }
240 7924 : else if (nXOff == 0 && nYOff == 0)
241 : {
242 7861 : const GDALDataType eDstDataType = bComplex ? GDT_CInt32 : GDT_Int32;
243 7861 : int nBlockXSize = 0;
244 7861 : int nBlockYSize = 0;
245 7861 : GDALGetBlockSize(hBand, &nBlockXSize, &nBlockYSize);
246 7861 : const int nDstDataTypeSize = GDALGetDataTypeSizeBytes(eDstDataType);
247 7861 : int nChunkXSize = nBlockXSize;
248 7861 : const int nChunkYSize = nBlockYSize;
249 7861 : if (nBlockXSize < nXSize)
250 : {
251 : const GIntBig nMaxChunkSize =
252 1170 : std::max(static_cast<GIntBig>(10 * 1000 * 1000),
253 585 : GDALGetCacheMax64() / 10);
254 585 : if (nDstDataTypeSize > 0 &&
255 585 : static_cast<GIntBig>(nXSize) * nChunkYSize <
256 585 : nMaxChunkSize / nDstDataTypeSize)
257 : {
258 : // A full line of height nChunkYSize can fit in the maximum
259 : // allowed memory
260 585 : nChunkXSize = nXSize;
261 : }
262 : else
263 : {
264 : // Otherwise compute a size that is a multiple of nBlockXSize
265 0 : nChunkXSize = static_cast<int>(std::min(
266 0 : static_cast<GIntBig>(nXSize),
267 0 : nBlockXSize *
268 0 : std::max(static_cast<GIntBig>(1),
269 0 : nMaxChunkSize /
270 0 : (static_cast<GIntBig>(nBlockXSize) *
271 0 : nChunkYSize * nDstDataTypeSize))));
272 : }
273 : }
274 :
275 : int *panChunkData = static_cast<GInt32 *>(
276 7861 : VSI_MALLOC3_VERBOSE(nChunkXSize, nChunkYSize, nDstDataTypeSize));
277 7861 : if (panChunkData == nullptr)
278 : {
279 0 : return -1;
280 : }
281 7861 : const int nValsPerIter = bComplex ? 2 : 1;
282 :
283 7861 : const int nYBlocks = DIV_ROUND_UP(nYSize, nChunkYSize);
284 7861 : const int nXBlocks = DIV_ROUND_UP(nXSize, nChunkXSize);
285 2341180 : for (int iYBlock = 0; iYBlock < nYBlocks; ++iYBlock)
286 : {
287 2333320 : const int iYStart = iYBlock * nChunkYSize;
288 2333320 : const int iYEnd =
289 2333320 : iYBlock == nYBlocks - 1 ? nYSize : iYStart + nChunkYSize;
290 2333320 : const int nChunkActualHeight = iYEnd - iYStart;
291 4666550 : for (int iXBlock = 0; iXBlock < nXBlocks; ++iXBlock)
292 : {
293 2333320 : const int iXStart = iXBlock * nChunkXSize;
294 2333320 : const int iXEnd =
295 2333320 : iXBlock == nXBlocks - 1 ? nXSize : iXStart + nChunkXSize;
296 2333320 : const int nChunkActualXSize = iXEnd - iXStart;
297 2333320 : if (GDALRasterIO(
298 : hBand, GF_Read, iXStart, iYStart, nChunkActualXSize,
299 : nChunkActualHeight, panChunkData, nChunkActualXSize,
300 2333320 : nChunkActualHeight, eDstDataType, 0, 0) != CE_None)
301 : {
302 81 : CPLError(CE_Failure, CPLE_FileIO,
303 : "Checksum value could not be computed due to I/O "
304 : "read error.");
305 82 : nChecksum = -1;
306 82 : iYBlock = nYBlocks;
307 82 : break;
308 : }
309 2333240 : const size_t xIters =
310 2333240 : static_cast<size_t>(nValsPerIter) * nChunkActualXSize;
311 5337460 : for (int iY = iYStart; iY < iYEnd; ++iY)
312 : {
313 : // Initialize iPrime so that it is consistent with a
314 : // per full line iteration strategy
315 3004230 : iPrime = (nValsPerIter *
316 3004230 : (static_cast<int64_t>(iY) * nXSize + iXStart)) %
317 : 11;
318 3004230 : const size_t nOffset = nValsPerIter *
319 3004230 : static_cast<size_t>(iY - iYStart) *
320 3004230 : nChunkActualXSize;
321 855132000 : for (size_t i = 0; i < xIters; ++i)
322 : {
323 852128000 : nChecksum +=
324 852128000 : panChunkData[nOffset + i] % anPrimes[iPrime++];
325 852128000 : if (iPrime > 10)
326 80772400 : iPrime = 0;
327 : }
328 3004230 : nChecksum &= 0xffff;
329 : }
330 : }
331 : }
332 :
333 7861 : CPLFree(panChunkData);
334 : }
335 : else
336 : {
337 63 : const GDALDataType eDstDataType = bComplex ? GDT_CInt32 : GDT_Int32;
338 :
339 63 : int *panLineData = static_cast<GInt32 *>(VSI_MALLOC2_VERBOSE(
340 : nXSize, GDALGetDataTypeSizeBytes(eDstDataType)));
341 63 : if (panLineData == nullptr)
342 : {
343 0 : return -1;
344 : }
345 :
346 153 : for (int iLine = nYOff; iLine < nYOff + nYSize; iLine++)
347 : {
348 90 : if (GDALRasterIO(hBand, GF_Read, nXOff, iLine, nXSize, 1,
349 : panLineData, nXSize, 1, eDstDataType, 0,
350 90 : 0) != CE_None)
351 : {
352 0 : CPLError(CE_Failure, CPLE_FileIO,
353 : "Checksum value could not be computed due to I/O "
354 : "read error.");
355 0 : nChecksum = -1;
356 0 : break;
357 : }
358 90 : const size_t nCount = bComplex ? static_cast<size_t>(nXSize) * 2
359 : : static_cast<size_t>(nXSize);
360 :
361 600 : for (size_t i = 0; i < nCount; i++)
362 : {
363 510 : nChecksum += panLineData[i] % anPrimes[iPrime++];
364 510 : if (iPrime > 10)
365 40 : iPrime = 0;
366 :
367 510 : nChecksum &= 0xffff;
368 : }
369 : }
370 :
371 63 : CPLFree(panLineData);
372 : }
373 :
374 8990 : return nChecksum;
375 : }
|