Line data Source code
1 : /******************************************************************************
2 : *
3 : * Project: GDAL Core
4 : * Purpose: Inline C++ templates
5 : * Author: Phil Vachon, <philippe at cowpig.ca>
6 : *
7 : ******************************************************************************
8 : * Copyright (c) 2009, Phil Vachon, <philippe at cowpig.ca>
9 : * Copyright (c) 2025, Even Rouault, <even.rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #ifndef GDAL_PRIV_TEMPLATES_HPP_INCLUDED
15 : #define GDAL_PRIV_TEMPLATES_HPP_INCLUDED
16 :
17 : #include "cpl_port.h"
18 :
19 : #include <algorithm>
20 : #include <cmath>
21 : #include <cstdint>
22 : #include <limits>
23 : #include <type_traits>
24 :
25 : #include "cpl_float.h"
26 :
27 : /************************************************************************/
28 : /* GDALGetDataLimits() */
29 : /************************************************************************/
30 : /**
31 : * Compute the limits of values that can be placed in Tout in terms of
32 : * Tin. Usually used for output clamping, when the output data type's
33 : * limits are stable relative to the input type (i.e. no roundoff error).
34 : *
35 : * @param tMaxValue the returned maximum value
36 : * @param tMinValue the returned minimum value
37 : */
38 :
39 : template <class Tin, class Tout>
40 371555664 : inline void GDALGetDataLimits(Tin &tMaxValue, Tin &tMinValue)
41 : {
42 371555664 : tMaxValue = cpl::NumericLimits<Tin>::max();
43 372000344 : tMinValue = cpl::NumericLimits<Tin>::lowest();
44 :
45 : // Compute the actual minimum value of Tout in terms of Tin.
46 : if constexpr (cpl::NumericLimits<Tout>::is_signed &&
47 : cpl::NumericLimits<Tout>::is_integer)
48 : {
49 : // the minimum value is less than zero
50 : // cppcheck-suppress knownConditionTrueFalse
51 : if constexpr (cpl::NumericLimits<Tout>::digits <
52 : cpl::NumericLimits<Tin>::digits ||
53 : !cpl::NumericLimits<Tin>::is_integer)
54 : {
55 : // Tout is smaller than Tin, so we need to clamp values in input
56 : // to the range of Tout's min/max values
57 : if constexpr (cpl::NumericLimits<Tin>::is_signed)
58 : {
59 78397869 : tMinValue =
60 78397919 : static_cast<Tin>(cpl::NumericLimits<Tout>::lowest());
61 : }
62 78912347 : tMaxValue = static_cast<Tin>(cpl::NumericLimits<Tout>::max());
63 : }
64 : }
65 : else if constexpr (cpl::NumericLimits<Tout>::is_integer)
66 : {
67 : // the output is unsigned, so we just need to determine the max
68 : /* coverity[same_on_both_sides] */
69 : if constexpr (cpl::NumericLimits<Tout>::digits <=
70 : cpl::NumericLimits<Tin>::digits)
71 : {
72 : // Tout is smaller than Tin, so we need to clamp the input values
73 : // to the range of Tout's max
74 141533086 : tMaxValue = static_cast<Tin>(cpl::NumericLimits<Tout>::max());
75 : }
76 143409476 : tMinValue = 0;
77 : }
78 371898374 : }
79 :
80 : /************************************************************************/
81 : /* GDALClampValue() */
82 : /************************************************************************/
83 : /**
84 : * Clamp values of type T to a specified range
85 : *
86 : * @param tValue the value
87 : * @param tMax the max value
88 : * @param tMin the min value
89 : */
90 : template <class T>
91 370850980 : inline T GDALClampValue(const T tValue, const T tMax, const T tMin)
92 : {
93 370850980 : return tValue > tMax ? tMax : tValue < tMin ? tMin : tValue;
94 : }
95 :
96 : /************************************************************************/
97 : /* GDALClampDoubleValue() */
98 : /************************************************************************/
99 : /**
100 : * Clamp double values to a specified range, this uses the same
101 : * argument ordering as std::clamp, returns TRUE if the value was clamped.
102 : *
103 : * @param tValue the value
104 : * @param tMin the min value
105 : * @param tMax the max value
106 : *
107 : */
108 : template <class T2, class T3>
109 235 : inline bool GDALClampDoubleValue(double &tValue, const T2 tMin, const T3 tMax)
110 : {
111 235 : const double tMin2{static_cast<double>(tMin)};
112 235 : const double tMax2{static_cast<double>(tMax)};
113 235 : if (tValue > tMax2 || tValue < tMin2)
114 : {
115 26 : tValue = tValue > tMax2 ? tMax2 : tValue < tMin2 ? tMin2 : tValue;
116 26 : return true;
117 : }
118 : else
119 : {
120 209 : return false;
121 : }
122 : }
123 :
124 : /************************************************************************/
125 : /* GDALIsValueInRange() */
126 : /************************************************************************/
127 : /**
128 : * Returns whether a value is in the type range.
129 : * NaN is considered not to be in type range.
130 : *
131 : * @param dfValue the value
132 : * @return whether the value is in the type range.
133 : */
134 145437 : template <class T> inline bool GDALIsValueInRange(double dfValue)
135 : {
136 290826 : return dfValue >= static_cast<double>(cpl::NumericLimits<T>::lowest()) &&
137 290826 : dfValue <= static_cast<double>(cpl::NumericLimits<T>::max());
138 : }
139 :
140 26 : template <> inline bool GDALIsValueInRange<double>(double dfValue)
141 : {
142 26 : return !CPLIsNan(dfValue);
143 : }
144 :
145 42088 : template <> inline bool GDALIsValueInRange<float>(double dfValue)
146 : {
147 84098 : return CPLIsInf(dfValue) ||
148 42018 : (dfValue >= -std::numeric_limits<float>::max() &&
149 84088 : dfValue <= std::numeric_limits<float>::max());
150 : }
151 :
152 0 : template <> inline bool GDALIsValueInRange<GFloat16>(double dfValue)
153 : {
154 0 : return CPLIsInf(dfValue) ||
155 0 : (dfValue >= -cpl::NumericLimits<GFloat16>::max() &&
156 0 : dfValue <= cpl::NumericLimits<GFloat16>::max());
157 : }
158 :
159 14808 : template <> inline bool GDALIsValueInRange<int64_t>(double dfValue)
160 : {
161 : // Values in the range [INT64_MAX - 1023, INT64_MAX - 1]
162 : // get converted to a double that once cast to int64_t is
163 : // INT64_MAX + 1, hence the < strict comparison.
164 : return dfValue >=
165 29615 : static_cast<double>(cpl::NumericLimits<int64_t>::lowest()) &&
166 29615 : dfValue < static_cast<double>(cpl::NumericLimits<int64_t>::max());
167 : }
168 :
169 8717 : template <> inline bool GDALIsValueInRange<uint64_t>(double dfValue)
170 : {
171 : // Values in the range [UINT64_MAX - 2047, UINT64_MAX - 1]
172 : // get converted to a double that once cast to uint64_t is
173 : // UINT64_MAX + 1, hence the < strict comparison.
174 17431 : return dfValue >= 0 &&
175 17431 : dfValue < static_cast<double>(cpl::NumericLimits<uint64_t>::max());
176 : }
177 :
178 : /************************************************************************/
179 : /* GDALIsValueExactAs() */
180 : /************************************************************************/
181 : /**
182 : * Returns whether a value can be exactly represented on type T.
183 : *
184 : * That is static_cast\<double\>(static_cast\<T\>(dfValue)) is legal and is
185 : * equal to dfValue.
186 : *
187 : * Note: for T=float or double, a NaN input leads to true
188 : *
189 : * @param dfValue the value
190 : * @return whether the value can be exactly represented on type T.
191 : */
192 567 : template <class T> inline bool GDALIsValueExactAs(double dfValue)
193 : {
194 1087 : return GDALIsValueInRange<T>(dfValue) &&
195 1087 : static_cast<double>(static_cast<T>(dfValue)) == dfValue;
196 : }
197 :
198 102 : template <> inline bool GDALIsValueExactAs<float>(double dfValue)
199 : {
200 199 : return CPLIsNan(dfValue) ||
201 97 : (GDALIsValueInRange<float>(dfValue) &&
202 195 : static_cast<double>(static_cast<float>(dfValue)) == dfValue);
203 : }
204 :
205 0 : template <> inline bool GDALIsValueExactAs<GFloat16>(double dfValue)
206 : {
207 0 : return CPLIsNan(dfValue) ||
208 0 : (GDALIsValueInRange<GFloat16>(dfValue) &&
209 0 : static_cast<double>(static_cast<GFloat16>(dfValue)) == dfValue);
210 : }
211 :
212 16 : template <> inline bool GDALIsValueExactAs<double>(double)
213 : {
214 16 : return true;
215 : }
216 :
217 : /************************************************************************/
218 : /* GDALCopyWord() */
219 : /************************************************************************/
220 :
221 : // Integer input and output: clamp the input
222 :
223 : template <class Tin, class Tout> struct sGDALCopyWord
224 : {
225 208211774 : static inline void f(const Tin tValueIn, Tout &tValueOut)
226 : {
227 : Tin tMaxVal, tMinVal;
228 208211774 : GDALGetDataLimits<Tin, Tout>(tMaxVal, tMinVal);
229 209172674 : tValueOut =
230 209059274 : static_cast<Tout>(GDALClampValue(tValueIn, tMaxVal, tMinVal));
231 209172674 : }
232 : };
233 :
234 : // Integer input and floating point output: simply convert
235 :
236 : template <class Tin> struct sGDALCopyWord<Tin, GFloat16>
237 : {
238 945 : static inline void f(const Tin tValueIn, GFloat16 &hfValueOut)
239 : {
240 945 : hfValueOut = static_cast<GFloat16>(tValueIn);
241 945 : }
242 : };
243 :
244 : template <class Tin> struct sGDALCopyWord<Tin, float>
245 : {
246 10886549 : static inline void f(const Tin tValueIn, float &fValueOut)
247 : {
248 10886549 : fValueOut = static_cast<float>(tValueIn);
249 10886549 : }
250 : };
251 :
252 : template <class Tin> struct sGDALCopyWord<Tin, double>
253 : {
254 73489140 : static inline void f(const Tin tValueIn, double &dfValueOut)
255 : {
256 73489140 : dfValueOut = static_cast<double>(tValueIn);
257 73489140 : }
258 : };
259 :
260 : // Floating point input and output, converting between identical types: simply copy
261 :
262 : template <> struct sGDALCopyWord<GFloat16, GFloat16>
263 : {
264 : static inline void f(const GFloat16 hfValueIn, GFloat16 &hfValueOut)
265 : {
266 : hfValueOut = hfValueIn;
267 : }
268 : };
269 :
270 : template <> struct sGDALCopyWord<float, float>
271 : {
272 : static inline void f(const float fValueIn, float &fValueOut)
273 : {
274 : fValueOut = fValueIn;
275 : }
276 : };
277 :
278 : template <> struct sGDALCopyWord<double, double>
279 : {
280 : static inline void f(const double dfValueIn, double &dfValueOut)
281 : {
282 : dfValueOut = dfValueIn;
283 : }
284 : };
285 :
286 : // Floating point input and output, converting to a larger type: use implicit conversion
287 :
288 : template <> struct sGDALCopyWord<GFloat16, float>
289 : {
290 225 : static inline void f(const GFloat16 hfValueIn, float &dfValueOut)
291 : {
292 225 : dfValueOut = hfValueIn;
293 225 : }
294 : };
295 :
296 : template <> struct sGDALCopyWord<GFloat16, double>
297 : {
298 322 : static inline void f(const GFloat16 hfValueIn, double &dfValueOut)
299 : {
300 322 : dfValueOut = hfValueIn;
301 322 : }
302 : };
303 :
304 : template <> struct sGDALCopyWord<float, double>
305 : {
306 46469200 : static inline void f(const float fValueIn, double &dfValueOut)
307 : {
308 46469200 : dfValueOut = fValueIn;
309 46469200 : }
310 : };
311 :
312 : // Floating point input and out, converting to a smaller type: ensure overflow results in infinity
313 :
314 : template <> struct sGDALCopyWord<float, GFloat16>
315 : {
316 165 : static inline void f(const float fValueIn, GFloat16 &hfValueOut)
317 : {
318 165 : if (fValueIn > cpl::NumericLimits<GFloat16>::max())
319 : {
320 0 : hfValueOut = cpl::NumericLimits<GFloat16>::infinity();
321 0 : return;
322 : }
323 165 : if (fValueIn < -cpl::NumericLimits<GFloat16>::max())
324 : {
325 0 : hfValueOut = -cpl::NumericLimits<GFloat16>::infinity();
326 0 : return;
327 : }
328 :
329 165 : hfValueOut = static_cast<GFloat16>(fValueIn);
330 : }
331 : };
332 :
333 : template <> struct sGDALCopyWord<double, GFloat16>
334 : {
335 167 : static inline void f(const double dfValueIn, GFloat16 &hfValueOut)
336 : {
337 167 : if (dfValueIn > cpl::NumericLimits<GFloat16>::max())
338 : {
339 0 : hfValueOut = cpl::NumericLimits<GFloat16>::infinity();
340 0 : return;
341 : }
342 167 : if (dfValueIn < -cpl::NumericLimits<GFloat16>::max())
343 : {
344 0 : hfValueOut = -cpl::NumericLimits<GFloat16>::infinity();
345 0 : return;
346 : }
347 :
348 167 : hfValueOut = static_cast<GFloat16>(dfValueIn);
349 : }
350 : };
351 :
352 : template <> struct sGDALCopyWord<double, float>
353 : {
354 2279650 : static inline void f(const double dfValueIn, float &fValueOut)
355 : {
356 2279650 : if (dfValueIn > std::numeric_limits<float>::max())
357 : {
358 48 : fValueOut = std::numeric_limits<float>::infinity();
359 48 : return;
360 : }
361 2279600 : if (dfValueIn < -std::numeric_limits<float>::max())
362 : {
363 54 : fValueOut = -std::numeric_limits<float>::infinity();
364 54 : return;
365 : }
366 :
367 2279550 : fValueOut = static_cast<float>(dfValueIn);
368 : }
369 : };
370 :
371 : // Floating point input to a small unsigned integer type: nan becomes zero, otherwise round and clamp
372 :
373 : template <class Tout> struct sGDALCopyWord<GFloat16, Tout>
374 : {
375 425 : static inline void f(const GFloat16 hfValueIn, Tout &tValueOut)
376 : {
377 425 : if (CPLIsNan(hfValueIn))
378 : {
379 0 : tValueOut = 0;
380 0 : return;
381 : }
382 : GFloat16 hfMaxVal, hfMinVal;
383 425 : GDALGetDataLimits<GFloat16, Tout>(hfMaxVal, hfMinVal);
384 425 : tValueOut = static_cast<Tout>(
385 850 : GDALClampValue(hfValueIn + GFloat16(0.5f), hfMaxVal, hfMinVal));
386 : }
387 : };
388 :
389 : template <class Tout> struct sGDALCopyWord<float, Tout>
390 : {
391 3979060 : static inline void f(const float fValueIn, Tout &tValueOut)
392 : {
393 3979060 : if (CPLIsNan(fValueIn))
394 : {
395 0 : tValueOut = 0;
396 0 : return;
397 : }
398 : float fMaxVal, fMinVal;
399 3979100 : GDALGetDataLimits<float, Tout>(fMaxVal, fMinVal);
400 3979010 : tValueOut = static_cast<Tout>(
401 3979010 : GDALClampValue(fValueIn + 0.5f, fMaxVal, fMinVal));
402 : }
403 : };
404 :
405 : template <class Tout> struct sGDALCopyWord<double, Tout>
406 : {
407 82492520 : static inline void f(const double dfValueIn, Tout &tValueOut)
408 : {
409 82492520 : if (CPLIsNan(dfValueIn))
410 : {
411 0 : tValueOut = 0;
412 0 : return;
413 : }
414 : double dfMaxVal, dfMinVal;
415 81290920 : GDALGetDataLimits<double, Tout>(dfMaxVal, dfMinVal);
416 80102720 : tValueOut = static_cast<Tout>(
417 80216320 : GDALClampValue(dfValueIn + 0.5, dfMaxVal, dfMinVal));
418 : }
419 : };
420 :
421 : // Floating point input to a large unsigned integer type: nan becomes zero, otherwise round and clamp.
422 : // Avoid roundoff while clamping.
423 :
424 : template <> struct sGDALCopyWord<GFloat16, std::uint64_t>
425 : {
426 168 : static inline void f(const GFloat16 hfValueIn, std::uint64_t &nValueOut)
427 : {
428 168 : if (!(hfValueIn > 0))
429 : {
430 3 : nValueOut = 0;
431 : }
432 165 : else if (CPLIsInf(hfValueIn))
433 : {
434 1 : nValueOut = cpl::NumericLimits<std::uint64_t>::max();
435 : }
436 : else
437 : {
438 164 : nValueOut = static_cast<std::uint64_t>(hfValueIn + GFloat16(0.5f));
439 : }
440 168 : }
441 : };
442 :
443 : template <> struct sGDALCopyWord<float, unsigned int>
444 : {
445 203 : static inline void f(const float fValueIn, unsigned int &nValueOut)
446 : {
447 203 : if (!(fValueIn > 0))
448 : {
449 20 : nValueOut = 0;
450 : }
451 183 : else if (fValueIn >=
452 183 : static_cast<float>(cpl::NumericLimits<unsigned int>::max()))
453 : {
454 20 : nValueOut = cpl::NumericLimits<unsigned int>::max();
455 : }
456 : else
457 : {
458 163 : nValueOut = static_cast<unsigned int>(fValueIn + 0.5f);
459 : }
460 203 : }
461 : };
462 :
463 : template <> struct sGDALCopyWord<float, std::uint64_t>
464 : {
465 168 : static inline void f(const float fValueIn, std::uint64_t &nValueOut)
466 : {
467 168 : if (!(fValueIn > 0))
468 : {
469 3 : nValueOut = 0;
470 : }
471 165 : else if (fValueIn >=
472 165 : static_cast<float>(cpl::NumericLimits<std::uint64_t>::max()))
473 : {
474 2 : nValueOut = cpl::NumericLimits<std::uint64_t>::max();
475 : }
476 : else
477 : {
478 163 : nValueOut = static_cast<std::uint64_t>(fValueIn + 0.5f);
479 : }
480 168 : }
481 : };
482 :
483 : template <> struct sGDALCopyWord<double, std::uint64_t>
484 : {
485 604 : static inline void f(const double dfValueIn, std::uint64_t &nValueOut)
486 : {
487 604 : if (!(dfValueIn > 0))
488 : {
489 164 : nValueOut = 0;
490 : }
491 440 : else if (dfValueIn >
492 440 : static_cast<double>(cpl::NumericLimits<uint64_t>::max()))
493 : {
494 4 : nValueOut = cpl::NumericLimits<uint64_t>::max();
495 : }
496 : else
497 : {
498 436 : nValueOut = static_cast<std::uint64_t>(dfValueIn + 0.5);
499 : }
500 604 : }
501 : };
502 :
503 : // Floating point input to a very large unsigned integer type: nan becomes zero, otherwise round and clamp.
504 : // Avoid infinity while clamping when the maximum integer is too large for the floating-point type.
505 : // Avoid roundoff while clamping.
506 :
507 : template <> struct sGDALCopyWord<GFloat16, unsigned short>
508 : {
509 213 : static inline void f(const GFloat16 hfValueIn, unsigned short &nValueOut)
510 : {
511 213 : if (!(hfValueIn > 0))
512 : {
513 30 : nValueOut = 0;
514 : }
515 183 : else if (CPLIsInf(hfValueIn))
516 : {
517 10 : nValueOut = cpl::NumericLimits<unsigned short>::max();
518 : }
519 : else
520 : {
521 173 : nValueOut = static_cast<unsigned short>(hfValueIn + GFloat16(0.5f));
522 : }
523 213 : }
524 : };
525 :
526 : template <> struct sGDALCopyWord<GFloat16, unsigned int>
527 : {
528 203 : static inline void f(const GFloat16 hfValueIn, unsigned int &nValueOut)
529 : {
530 203 : if (!(hfValueIn > 0))
531 : {
532 20 : nValueOut = 0;
533 : }
534 183 : else if (CPLIsInf(hfValueIn))
535 : {
536 0 : nValueOut = cpl::NumericLimits<unsigned int>::max();
537 : }
538 : else
539 : {
540 183 : nValueOut = static_cast<unsigned int>(hfValueIn + GFloat16(0.5f));
541 : }
542 203 : }
543 : };
544 :
545 : // Floating point input to a small signed integer type: nan becomes zero, otherwise round and clamp.
546 : // Rounding for signed integers is different than for the unsigned integers above.
547 :
548 : template <> struct sGDALCopyWord<GFloat16, signed char>
549 : {
550 233 : static inline void f(const GFloat16 hfValueIn, signed char &nValueOut)
551 : {
552 233 : if (CPLIsNan(hfValueIn))
553 : {
554 0 : nValueOut = 0;
555 0 : return;
556 : }
557 : GFloat16 hfMaxVal, hfMinVal;
558 233 : GDALGetDataLimits<GFloat16, signed char>(hfMaxVal, hfMinVal);
559 233 : GFloat16 hfValue = hfValueIn >= GFloat16(0.0f)
560 163 : ? hfValueIn + GFloat16(0.5f)
561 233 : : hfValueIn - GFloat16(0.5f);
562 233 : nValueOut = static_cast<signed char>(
563 466 : GDALClampValue(hfValue, hfMaxVal, hfMinVal));
564 : }
565 : };
566 :
567 : template <> struct sGDALCopyWord<float, signed char>
568 : {
569 297 : static inline void f(const float fValueIn, signed char &nValueOut)
570 : {
571 297 : if (CPLIsNan(fValueIn))
572 : {
573 0 : nValueOut = 0;
574 0 : return;
575 : }
576 : float fMaxVal, fMinVal;
577 297 : GDALGetDataLimits<float, signed char>(fMaxVal, fMinVal);
578 297 : float fValue = fValueIn >= 0.0f ? fValueIn + 0.5f : fValueIn - 0.5f;
579 297 : nValueOut =
580 297 : static_cast<signed char>(GDALClampValue(fValue, fMaxVal, fMinVal));
581 : }
582 : };
583 :
584 : template <> struct sGDALCopyWord<float, short>
585 : {
586 2928940 : static inline void f(const float fValueIn, short &nValueOut)
587 : {
588 2928940 : if (CPLIsNan(fValueIn))
589 : {
590 0 : nValueOut = 0;
591 0 : return;
592 : }
593 : float fMaxVal, fMinVal;
594 2928940 : GDALGetDataLimits<float, short>(fMaxVal, fMinVal);
595 2928940 : float fValue = fValueIn >= 0.0f ? fValueIn + 0.5f : fValueIn - 0.5f;
596 2928940 : nValueOut =
597 2928940 : static_cast<short>(GDALClampValue(fValue, fMaxVal, fMinVal));
598 : }
599 : };
600 :
601 : template <> struct sGDALCopyWord<double, signed char>
602 : {
603 465 : static inline void f(const double dfValueIn, signed char &nValueOut)
604 : {
605 465 : if (CPLIsNan(dfValueIn))
606 : {
607 0 : nValueOut = 0;
608 0 : return;
609 : }
610 : double dfMaxVal, dfMinVal;
611 465 : GDALGetDataLimits<double, signed char>(dfMaxVal, dfMinVal);
612 465 : double dfValue = dfValueIn > 0.0 ? dfValueIn + 0.5 : dfValueIn - 0.5;
613 465 : nValueOut = static_cast<signed char>(
614 465 : GDALClampValue(dfValue, dfMaxVal, dfMinVal));
615 : }
616 : };
617 :
618 : template <> struct sGDALCopyWord<double, short>
619 : {
620 5101410 : static inline void f(const double dfValueIn, short &nValueOut)
621 : {
622 5101410 : if (CPLIsNan(dfValueIn))
623 : {
624 0 : nValueOut = 0;
625 0 : return;
626 : }
627 : double dfMaxVal, dfMinVal;
628 5101390 : GDALGetDataLimits<double, short>(dfMaxVal, dfMinVal);
629 5101320 : double dfValue = dfValueIn > 0.0 ? dfValueIn + 0.5 : dfValueIn - 0.5;
630 5101310 : nValueOut =
631 5101320 : static_cast<short>(GDALClampValue(dfValue, dfMaxVal, dfMinVal));
632 : }
633 : };
634 :
635 : template <> struct sGDALCopyWord<double, int>
636 : {
637 70358300 : static inline void f(const double dfValueIn, int &nValueOut)
638 : {
639 70358300 : if (CPLIsNan(dfValueIn))
640 : {
641 0 : nValueOut = 0;
642 0 : return;
643 : }
644 : double dfMaxVal, dfMinVal;
645 70358300 : GDALGetDataLimits<double, int>(dfMaxVal, dfMinVal);
646 70358300 : double dfValue = dfValueIn >= 0.0 ? dfValueIn + 0.5 : dfValueIn - 0.5;
647 70358300 : nValueOut =
648 70358300 : static_cast<int>(GDALClampValue(dfValue, dfMaxVal, dfMinVal));
649 : }
650 : };
651 :
652 : // Floating point input to a large signed integer type: nan becomes zero, otherwise round and clamp.
653 : // Rounding for signed integers is different than for the unsigned integers above.
654 : // Avoid roundoff while clamping.
655 :
656 : template <> struct sGDALCopyWord<GFloat16, short>
657 : {
658 397 : static inline void f(const GFloat16 hfValueIn, short &nValueOut)
659 : {
660 397 : if (CPLIsNan(hfValueIn))
661 : {
662 0 : nValueOut = 0;
663 : }
664 397 : else if (hfValueIn >=
665 397 : static_cast<GFloat16>(cpl::NumericLimits<short>::max()))
666 : {
667 32 : nValueOut = cpl::NumericLimits<short>::max();
668 : }
669 365 : else if (hfValueIn <=
670 365 : static_cast<GFloat16>(cpl::NumericLimits<short>::lowest()))
671 : {
672 42 : nValueOut = cpl::NumericLimits<short>::lowest();
673 : }
674 : else
675 : {
676 646 : nValueOut = static_cast<short>(hfValueIn > GFloat16(0.0f)
677 323 : ? hfValueIn + GFloat16(0.5f)
678 431 : : hfValueIn - GFloat16(0.5f));
679 : }
680 397 : }
681 : };
682 :
683 : template <> struct sGDALCopyWord<float, int>
684 : {
685 155543000 : static inline void f(const float fValueIn, int &nValueOut)
686 : {
687 155543000 : if (CPLIsNan(fValueIn))
688 : {
689 0 : nValueOut = 0;
690 : }
691 158013000 : else if (fValueIn >= static_cast<float>(cpl::NumericLimits<int>::max()))
692 : {
693 160 : nValueOut = cpl::NumericLimits<int>::max();
694 : }
695 155374000 : else if (fValueIn <=
696 157778000 : static_cast<float>(cpl::NumericLimits<int>::lowest()))
697 : {
698 0 : nValueOut = cpl::NumericLimits<int>::lowest();
699 : }
700 : else
701 : {
702 158585000 : nValueOut = static_cast<int>(fValueIn > 0.0f ? fValueIn + 0.5f
703 120973 : : fValueIn - 0.5f);
704 : }
705 158464000 : }
706 : };
707 :
708 : template <> struct sGDALCopyWord<float, std::int64_t>
709 : {
710 238 : static inline void f(const float fValueIn, std::int64_t &nValueOut)
711 : {
712 238 : if (CPLIsNan(fValueIn))
713 : {
714 1 : nValueOut = 0;
715 : }
716 237 : else if (fValueIn >=
717 237 : static_cast<float>(cpl::NumericLimits<std::int64_t>::max()))
718 : {
719 2 : nValueOut = cpl::NumericLimits<std::int64_t>::max();
720 : }
721 235 : else if (fValueIn <=
722 235 : static_cast<float>(cpl::NumericLimits<std::int64_t>::lowest()))
723 : {
724 2 : nValueOut = cpl::NumericLimits<std::int64_t>::lowest();
725 : }
726 : else
727 : {
728 466 : nValueOut = static_cast<std::int64_t>(
729 233 : fValueIn > 0.0f ? fValueIn + 0.5f : fValueIn - 0.5f);
730 : }
731 238 : }
732 : };
733 :
734 : template <> struct sGDALCopyWord<double, std::int64_t>
735 : {
736 688 : static inline void f(const double dfValueIn, std::int64_t &nValueOut)
737 : {
738 688 : if (CPLIsNan(dfValueIn))
739 : {
740 1 : nValueOut = 0;
741 : }
742 687 : else if (dfValueIn >=
743 687 : static_cast<double>(cpl::NumericLimits<std::int64_t>::max()))
744 : {
745 6 : nValueOut = cpl::NumericLimits<std::int64_t>::max();
746 : }
747 681 : else if (dfValueIn <=
748 681 : static_cast<double>(cpl::NumericLimits<std::int64_t>::min()))
749 : {
750 4 : nValueOut = cpl::NumericLimits<std::int64_t>::min();
751 : }
752 : else
753 : {
754 1354 : nValueOut = static_cast<std::int64_t>(
755 677 : dfValueIn > 0.0 ? dfValueIn + 0.5 : dfValueIn - 0.5);
756 : }
757 688 : }
758 : };
759 :
760 : // Floating point input to a very large signed integer type: nan becomes zero, otherwise round and clamp.
761 : // Rounding for signed integers is different than for the unsigned integers above.
762 : // Avoid infinity while clamping when the maximum integer is too large for the floating-point type.
763 : // Avoid roundoff while clamping.
764 :
765 : template <> struct sGDALCopyWord<GFloat16, int>
766 : {
767 387 : static inline void f(const GFloat16 hfValueIn, int &nValueOut)
768 : {
769 387 : if (CPLIsNan(hfValueIn))
770 : {
771 0 : nValueOut = 0;
772 : }
773 387 : else if (CPLIsInf(hfValueIn))
774 : {
775 0 : nValueOut = hfValueIn > GFloat16(0.0f)
776 0 : ? cpl::NumericLimits<int>::max()
777 0 : : cpl::NumericLimits<int>::lowest();
778 : }
779 : else
780 : {
781 774 : nValueOut = static_cast<int>(hfValueIn > GFloat16(0.0f)
782 387 : ? hfValueIn + GFloat16(0.5f)
783 527 : : hfValueIn - GFloat16(0.5f));
784 : }
785 387 : }
786 : };
787 :
788 : template <> struct sGDALCopyWord<GFloat16, std::int64_t>
789 : {
790 238 : static inline void f(const GFloat16 hfValueIn, std::int64_t &nValueOut)
791 : {
792 238 : if (CPLIsNan(hfValueIn))
793 : {
794 1 : nValueOut = 0;
795 : }
796 237 : else if (CPLIsInf(hfValueIn))
797 : {
798 4 : nValueOut = hfValueIn > GFloat16(0.0f)
799 2 : ? cpl::NumericLimits<std::int64_t>::max()
800 1 : : cpl::NumericLimits<std::int64_t>::lowest();
801 : }
802 : else
803 : {
804 235 : nValueOut = static_cast<std::int64_t>(
805 235 : hfValueIn > GFloat16(0.0f) ? hfValueIn + GFloat16(0.5f)
806 306 : : hfValueIn - GFloat16(0.5f));
807 : }
808 238 : }
809 : };
810 :
811 : /**
812 : * Copy a single word, optionally rounding if appropriate (i.e. going
813 : * from the float to the integer case). Note that this is the function
814 : * you should specialize if you're adding a new data type.
815 : *
816 : * @param tValueIn value of type Tin; the input value to be converted
817 : * @param tValueOut value of type Tout; the output value
818 : */
819 :
820 : template <class Tin, class Tout>
821 697344401 : inline void GDALCopyWord(const Tin tValueIn, Tout &tValueOut)
822 : {
823 : if constexpr (std::is_same<Tin, Tout>::value)
824 32426697 : tValueOut = tValueIn;
825 : else
826 664917704 : sGDALCopyWord<Tin, Tout>::f(tValueIn, tValueOut);
827 695280591 : }
828 :
829 : /************************************************************************/
830 : /* GDALCopy4Words() */
831 : /************************************************************************/
832 : /**
833 : * Copy 4 packed words to 4 packed words, optionally rounding if appropriate
834 : * (i.e. going from the float to the integer case).
835 : *
836 : * @param pValueIn pointer to 4 input values of type Tin.
837 : * @param pValueOut pointer to 4 output values of type Tout.
838 : */
839 :
840 : template <class Tin, class Tout>
841 16 : inline void GDALCopy4Words(const Tin *pValueIn, Tout *const pValueOut)
842 : {
843 16 : GDALCopyWord(pValueIn[0], pValueOut[0]);
844 16 : GDALCopyWord(pValueIn[1], pValueOut[1]);
845 16 : GDALCopyWord(pValueIn[2], pValueOut[2]);
846 16 : GDALCopyWord(pValueIn[3], pValueOut[3]);
847 16 : }
848 :
849 : /************************************************************************/
850 : /* GDALCopy8Words() */
851 : /************************************************************************/
852 : /**
853 : * Copy 8 packed words to 8 packed words, optionally rounding if appropriate
854 : * (i.e. going from the float to the integer case).
855 : *
856 : * @param pValueIn pointer to 8 input values of type Tin.
857 : * @param pValueOut pointer to 8 output values of type Tout.
858 : */
859 :
860 : template <class Tin, class Tout>
861 14782979 : inline void GDALCopy8Words(const Tin *pValueIn, Tout *const pValueOut)
862 : {
863 14782979 : GDALCopy4Words(pValueIn, pValueOut);
864 14784379 : GDALCopy4Words(pValueIn + 4, pValueOut + 4);
865 14790379 : }
866 :
867 : // Needs SSE2
868 : #if defined(__x86_64) || defined(_M_X64) || defined(USE_SSE2) || \
869 : defined(USE_NEON_OPTIMIZATIONS)
870 :
871 : #ifdef USE_NEON_OPTIMIZATIONS
872 : #include "include_sse2neon.h"
873 : #else
874 : #include <emmintrin.h>
875 : #endif
876 :
877 32365335 : static inline void GDALCopyXMMToInt32(const __m128i xmm, void *pDest)
878 : {
879 32365335 : int n32 = _mm_cvtsi128_si32(xmm); // Extract lower 32 bit word
880 32365335 : memcpy(pDest, &n32, sizeof(n32));
881 32365335 : }
882 :
883 77195273 : static inline void GDALCopyXMMToInt64(const __m128i xmm, void *pDest)
884 : {
885 : _mm_storel_epi64(reinterpret_cast<__m128i *>(pDest), xmm);
886 77195273 : }
887 :
888 : #if __SSSE3__
889 : #include <tmmintrin.h>
890 : #endif
891 :
892 : #if defined(__SSE4_1__) || defined(__AVX__)
893 : #include <smmintrin.h>
894 : #endif
895 :
896 : template <>
897 26198802 : inline void GDALCopy4Words(const float *pValueIn, GByte *const pValueOut)
898 : {
899 26198802 : __m128 xmm = _mm_loadu_ps(pValueIn);
900 :
901 : // The following clamping would be useless due to the final saturating
902 : // packing if we could guarantee the input range in [INT_MIN,INT_MAX]
903 26198802 : const __m128 p0d5 = _mm_set1_ps(0.5f);
904 26198802 : const __m128 xmm_max = _mm_set1_ps(255);
905 26198802 : xmm = _mm_add_ps(xmm, p0d5);
906 52442904 : xmm = _mm_min_ps(_mm_max_ps(xmm, p0d5), xmm_max);
907 :
908 26223602 : __m128i xmm_i = _mm_cvttps_epi32(xmm);
909 :
910 : #if defined(__SSSE3__) || defined(USE_NEON_OPTIMIZATIONS)
911 : xmm_i = _mm_shuffle_epi8(
912 : xmm_i, _mm_cvtsi32_si128(0 | (4 << 8) | (8 << 16) | (12 << 24)));
913 : #else
914 26225102 : xmm_i = _mm_packs_epi32(xmm_i, xmm_i); // Pack int32 to int16
915 26225902 : xmm_i = _mm_packus_epi16(xmm_i, xmm_i); // Pack int16 to uint8
916 : #endif
917 26225902 : GDALCopyXMMToInt32(xmm_i, pValueOut);
918 26217102 : }
919 :
920 : template <>
921 3355730 : inline void GDALCopy4Words(const float *pValueIn, GInt16 *const pValueOut)
922 : {
923 3355730 : __m128 xmm = _mm_loadu_ps(pValueIn);
924 :
925 3355730 : const __m128 xmm_min = _mm_set1_ps(-32768);
926 3355730 : const __m128 xmm_max = _mm_set1_ps(32767);
927 6711460 : xmm = _mm_min_ps(_mm_max_ps(xmm, xmm_min), xmm_max);
928 :
929 3355730 : const __m128 p0d5 = _mm_set1_ps(0.5f);
930 3355730 : const __m128 m0d5 = _mm_set1_ps(-0.5f);
931 3355730 : const __m128 mask = _mm_cmpge_ps(xmm, p0d5);
932 : // f >= 0.5f ? f + 0.5f : f - 0.5f
933 13422900 : xmm = _mm_add_ps(
934 : xmm, _mm_or_ps(_mm_and_ps(mask, p0d5), _mm_andnot_ps(mask, m0d5)));
935 :
936 3355730 : __m128i xmm_i = _mm_cvttps_epi32(xmm);
937 :
938 3355730 : xmm_i = _mm_packs_epi32(xmm_i, xmm_i); // Pack int32 to int16
939 3355730 : GDALCopyXMMToInt64(xmm_i, pValueOut);
940 3355730 : }
941 :
942 : template <>
943 1 : inline void GDALCopy4Words(const float *pValueIn, GUInt16 *const pValueOut)
944 : {
945 1 : __m128 xmm = _mm_loadu_ps(pValueIn);
946 :
947 1 : const __m128 p0d5 = _mm_set1_ps(0.5f);
948 1 : const __m128 xmm_max = _mm_set1_ps(65535);
949 1 : xmm = _mm_add_ps(xmm, p0d5);
950 2 : xmm = _mm_min_ps(_mm_max_ps(xmm, p0d5), xmm_max);
951 :
952 1 : __m128i xmm_i = _mm_cvttps_epi32(xmm);
953 :
954 : #if defined(__SSE4_1__) || defined(__AVX__) || defined(USE_NEON_OPTIMIZATIONS)
955 : xmm_i = _mm_packus_epi32(xmm_i, xmm_i); // Pack int32 to uint16
956 : #else
957 : // Translate to int16 range because _mm_packus_epi32 is SSE4.1 only
958 2 : xmm_i = _mm_add_epi32(xmm_i, _mm_set1_epi32(-32768));
959 1 : xmm_i = _mm_packs_epi32(xmm_i, xmm_i); // Pack int32 to int16
960 : // Translate back to uint16 range (actually -32768==32768 in int16)
961 1 : xmm_i = _mm_add_epi16(xmm_i, _mm_set1_epi16(-32768));
962 : #endif
963 1 : GDALCopyXMMToInt64(xmm_i, pValueOut);
964 1 : }
965 :
966 : #ifdef __AVX2__
967 :
968 : #include <immintrin.h>
969 :
970 : template <>
971 : inline void GDALCopy8Words(const float *pValueIn, GByte *const pValueOut)
972 : {
973 : __m256 ymm = _mm256_loadu_ps(pValueIn);
974 :
975 : const __m256 p0d5 = _mm256_set1_ps(0.5f);
976 : const __m256 ymm_max = _mm256_set1_ps(255);
977 : ymm = _mm256_add_ps(ymm, p0d5);
978 : ymm = _mm256_min_ps(_mm256_max_ps(ymm, p0d5), ymm_max);
979 :
980 : __m256i ymm_i = _mm256_cvttps_epi32(ymm);
981 :
982 : ymm_i = _mm256_packus_epi32(ymm_i, ymm_i); // Pack int32 to uint16
983 : ymm_i = _mm256_permute4x64_epi64(ymm_i, 0 | (2 << 2)); // AVX2
984 :
985 : __m128i xmm_i = _mm256_castsi256_si128(ymm_i);
986 : xmm_i = _mm_packus_epi16(xmm_i, xmm_i);
987 : GDALCopyXMMToInt64(xmm_i, pValueOut);
988 : }
989 :
990 : template <>
991 : inline void GDALCopy8Words(const float *pValueIn, GUInt16 *const pValueOut)
992 : {
993 : __m256 ymm = _mm256_loadu_ps(pValueIn);
994 :
995 : const __m256 p0d5 = _mm256_set1_ps(0.5f);
996 : const __m256 ymm_max = _mm256_set1_ps(65535);
997 : ymm = _mm256_add_ps(ymm, p0d5);
998 : ymm = _mm256_min_ps(_mm256_max_ps(ymm, p0d5), ymm_max);
999 :
1000 : __m256i ymm_i = _mm256_cvttps_epi32(ymm);
1001 :
1002 : ymm_i = _mm256_packus_epi32(ymm_i, ymm_i); // Pack int32 to uint16
1003 : ymm_i = _mm256_permute4x64_epi64(ymm_i, 0 | (2 << 2)); // AVX2
1004 :
1005 : _mm_storeu_si128(reinterpret_cast<__m128i *>(pValueOut),
1006 : _mm256_castsi256_si128(ymm_i));
1007 : }
1008 : #else
1009 : template <>
1010 7763521 : inline void GDALCopy8Words(const float *pValueIn, GUInt16 *const pValueOut)
1011 : {
1012 7763521 : __m128 xmm = _mm_loadu_ps(pValueIn);
1013 15527002 : __m128 xmm1 = _mm_loadu_ps(pValueIn + 4);
1014 :
1015 7763521 : const __m128 p0d5 = _mm_set1_ps(0.5f);
1016 7763521 : const __m128 xmm_max = _mm_set1_ps(65535);
1017 7763521 : xmm = _mm_add_ps(xmm, p0d5);
1018 7763521 : xmm1 = _mm_add_ps(xmm1, p0d5);
1019 15525202 : xmm = _mm_min_ps(_mm_max_ps(xmm, p0d5), xmm_max);
1020 15518402 : xmm1 = _mm_min_ps(_mm_max_ps(xmm1, p0d5), xmm_max);
1021 :
1022 7756551 : __m128i xmm_i = _mm_cvttps_epi32(xmm);
1023 7764761 : __m128i xmm1_i = _mm_cvttps_epi32(xmm1);
1024 :
1025 : #if defined(__SSE4_1__) || defined(__AVX__) || defined(USE_NEON_OPTIMIZATIONS)
1026 : xmm_i = _mm_packus_epi32(xmm_i, xmm1_i); // Pack int32 to uint16
1027 : #else
1028 : // Translate to int16 range because _mm_packus_epi32 is SSE4.1 only
1029 15529502 : xmm_i = _mm_add_epi32(xmm_i, _mm_set1_epi32(-32768));
1030 15529502 : xmm1_i = _mm_add_epi32(xmm1_i, _mm_set1_epi32(-32768));
1031 7765191 : xmm_i = _mm_packs_epi32(xmm_i, xmm1_i); // Pack int32 to int16
1032 : // Translate back to uint16 range (actually -32768==32768 in int16)
1033 15530402 : xmm_i = _mm_add_epi16(xmm_i, _mm_set1_epi16(-32768));
1034 : #endif
1035 : _mm_storeu_si128(reinterpret_cast<__m128i *>(pValueOut), xmm_i);
1036 7765191 : }
1037 : #endif
1038 :
1039 : #ifdef notdef_because_slightly_slower_than_default_implementation
1040 : template <>
1041 : inline void GDALCopy4Words(const double *pValueIn, float *const pValueOut)
1042 : {
1043 : __m128d float_posmax = _mm_set1_pd(std::numeric_limits<float>::max());
1044 : __m128d float_negmax = _mm_set1_pd(-std::numeric_limits<float>::max());
1045 : __m128d float_posinf = _mm_set1_pd(std::numeric_limits<float>::infinity());
1046 : __m128d float_neginf = _mm_set1_pd(-std::numeric_limits<float>::infinity());
1047 : __m128d val01 = _mm_loadu_pd(pValueIn);
1048 : __m128d val23 = _mm_loadu_pd(pValueIn + 2);
1049 : __m128d mask_max = _mm_cmpge_pd(val01, float_posmax);
1050 : __m128d mask_max23 = _mm_cmpge_pd(val23, float_posmax);
1051 : val01 = _mm_or_pd(_mm_and_pd(mask_max, float_posinf),
1052 : _mm_andnot_pd(mask_max, val01));
1053 : val23 = _mm_or_pd(_mm_and_pd(mask_max23, float_posinf),
1054 : _mm_andnot_pd(mask_max23, val23));
1055 : __m128d mask_min = _mm_cmple_pd(val01, float_negmax);
1056 : __m128d mask_min23 = _mm_cmple_pd(val23, float_negmax);
1057 : val01 = _mm_or_pd(_mm_and_pd(mask_min, float_neginf),
1058 : _mm_andnot_pd(mask_min, val01));
1059 : val23 = _mm_or_pd(_mm_and_pd(mask_min23, float_neginf),
1060 : _mm_andnot_pd(mask_min23, val23));
1061 : __m128 val01_s = _mm_cvtpd_ps(val01);
1062 : __m128 val23_s = _mm_cvtpd_ps(val23);
1063 : __m128i val01_i = _mm_castps_si128(val01_s);
1064 : __m128i val23_i = _mm_castps_si128(val23_s);
1065 : GDALCopyXMMToInt64(val01_i, pValueOut);
1066 : GDALCopyXMMToInt64(val23_i, pValueOut + 2);
1067 : }
1068 : #endif
1069 :
1070 : #endif // defined(__x86_64) || defined(_M_X64)
1071 :
1072 : #endif // GDAL_PRIV_TEMPLATES_HPP_INCLUDED
|
