Line data Source code
1 : /******************************************************************************
2 : *
3 : * Name: gdalmultidim_array_maths.cpp
4 : * Project: GDAL Core
5 : * Purpose: Mathematic operations on GDALMDArray
6 : * Author: Even Rouault <even.rouault at spatialys.com>
7 : *
8 : ******************************************************************************
9 : * Copyright (c) 2019, Even Rouault <even.rouault at spatialys.com>
10 : *
11 : * SPDX-License-Identifier: MIT
12 : ****************************************************************************/
13 :
14 : #include "gdal_multidim.h"
15 : #include "gdal_pam_multidim.h"
16 : #include "ogr_spatialref.h"
17 :
18 : #include <algorithm>
19 : #include <cmath>
20 : #include <functional>
21 : #include <limits>
22 :
23 : /************************************************************************/
24 : /* GDALMDArray::HasSameShapeAs() */
25 : /************************************************************************/
26 :
27 : /** Returns true if both arrays have the same shape.
28 : *
29 : * That is to say the same number of dimensions and the size of corresponding
30 : * dimensions is the same.
31 : *
32 : * @since 3.14
33 : */
34 165 : bool GDALMDArray::HasSameShapeAs(const GDALMDArray &other) const
35 : {
36 165 : bool bRet = (GetDimensionCount() == other.GetDimensionCount());
37 165 : if (bRet)
38 : {
39 165 : const auto &apoThisDims = GetDimensions();
40 165 : const auto &apoOtherDims = other.GetDimensions();
41 388 : for (size_t i = 0; bRet && i < apoThisDims.size(); ++i)
42 : {
43 223 : bRet = (apoThisDims[i]->GetSize() == apoOtherDims[i]->GetSize());
44 : }
45 : }
46 165 : return bRet;
47 : }
48 :
49 : /************************************************************************/
50 : /* GDALMathOperationMDArray */
51 : /************************************************************************/
52 :
53 : //! @cond Doxygen_Suppress
54 :
55 : class GDALMathOperationMDArray final : public GDALPamMDArray
56 : {
57 : public:
58 : enum class Operation
59 : {
60 : OP_ADD,
61 : OP_SUBTRACT,
62 : OP_MULTIPLY,
63 : OP_DIVIDE,
64 : };
65 :
66 298 : static const char *OperationToString(Operation op)
67 : {
68 298 : switch (op)
69 : {
70 96 : case Operation::OP_ADD:
71 96 : break;
72 84 : case Operation::OP_SUBTRACT:
73 84 : return "-";
74 67 : case Operation::OP_MULTIPLY:
75 67 : return "*";
76 51 : case Operation::OP_DIVIDE:
77 51 : return "/";
78 : }
79 96 : return "+";
80 : }
81 :
82 : static std::shared_ptr<GDALMathOperationMDArray>
83 : Create(const std::shared_ptr<GDALMDArray> &arrayLeft,
84 : const std::shared_ptr<GDALMDArray> &arrayRight, Operation op);
85 :
86 : protected:
87 292 : static std::string GetName(const std::shared_ptr<GDALMDArray> &arrayLeft,
88 : const std::shared_ptr<GDALMDArray> &arrayRight,
89 : Operation op)
90 : {
91 584 : return std::string(arrayLeft->GetFullName())
92 292 : .append(" ")
93 292 : .append(OperationToString(op))
94 292 : .append(" ")
95 584 : .append(arrayRight->GetFullName());
96 : }
97 :
98 146 : GDALMathOperationMDArray(const std::shared_ptr<GDALMDArray> &arrayLeft,
99 : const std::shared_ptr<GDALMDArray> &arrayRight,
100 : Operation op)
101 292 : : GDALAbstractMDArray(std::string(),
102 146 : GetName(arrayLeft, arrayRight, op)),
103 292 : GDALPamMDArray(std::string(), GetName(arrayLeft, arrayRight, op),
104 : // Arbitrary linking to arrayLeft for PAM purposes
105 292 : GDALPamMultiDim::GetPAM(arrayLeft),
106 : arrayLeft->GetContext()),
107 : m_arrayLeft(arrayLeft), m_arrayRight(arrayRight), m_op(op),
108 : m_dt(GDALExtendedDataType::Create(GDT_Float64)),
109 : m_dfNoDataLeft(
110 146 : m_arrayLeft->GetNoDataValueAsDouble(&m_bHasNoDataLeft)),
111 : m_dfNoDataRight(
112 876 : m_arrayRight->GetNoDataValueAsDouble(&m_bHasNoDataRight))
113 : {
114 146 : if (m_bHasNoDataLeft || m_bHasNoDataRight)
115 : {
116 50 : if (m_bHasNoDataLeft && m_bHasNoDataRight)
117 : {
118 37 : if (m_dfNoDataLeft == m_dfNoDataRight)
119 32 : m_dfNoData = m_dfNoDataLeft;
120 : }
121 13 : else if (m_bHasNoDataLeft)
122 7 : m_dfNoData = m_dfNoDataLeft;
123 : else
124 6 : m_dfNoData = m_dfNoDataRight;
125 50 : m_abyRawNoDataValue.resize(sizeof(double));
126 50 : memcpy(m_abyRawNoDataValue.data(), &m_dfNoData, sizeof(m_dfNoData));
127 : }
128 :
129 146 : const auto &leftUnit = m_arrayLeft->GetUnit();
130 146 : const auto &rightUnit = m_arrayRight->GetUnit();
131 146 : switch (m_op)
132 : {
133 90 : case Operation::OP_ADD:
134 : case Operation::OP_SUBTRACT:
135 : {
136 90 : if (leftUnit == rightUnit)
137 86 : m_osUnit = leftUnit;
138 90 : break;
139 : }
140 56 : case Operation::OP_MULTIPLY:
141 : case Operation::OP_DIVIDE:
142 : {
143 56 : if (!leftUnit.empty() && !rightUnit.empty())
144 : {
145 6 : m_osUnit = leftUnit;
146 6 : m_osUnit += ' ';
147 6 : m_osUnit += OperationToString(m_op);
148 6 : m_osUnit += ' ';
149 6 : m_osUnit += rightUnit;
150 : }
151 56 : break;
152 : }
153 : }
154 146 : }
155 :
156 1 : bool IsWritable() const override
157 : {
158 1 : return false;
159 : }
160 :
161 127 : const std::string &GetFilename() const override
162 : {
163 127 : const auto &filename1 = m_arrayLeft->GetFilename();
164 127 : const auto &filename2 = m_arrayRight->GetFilename();
165 127 : if (filename1 == filename2)
166 86 : return filename1;
167 41 : static std::string emptyString;
168 41 : return emptyString;
169 : }
170 :
171 : const std::vector<std::shared_ptr<GDALDimension>> &
172 891 : GetDimensions() const override
173 : {
174 891 : return m_arrayLeft->GetDimensions();
175 : }
176 :
177 : std::vector<std::shared_ptr<GDALMDArray>>
178 17 : GetCoordinateVariables() const override
179 : {
180 34 : auto left = m_arrayLeft->GetCoordinateVariables();
181 34 : auto right = m_arrayRight->GetCoordinateVariables();
182 17 : if (left == right)
183 17 : return left;
184 0 : return GDALMDArray::GetCoordinateVariables();
185 : }
186 :
187 422 : const GDALExtendedDataType &GetDataType() const override
188 : {
189 422 : return m_dt;
190 : }
191 :
192 71 : const void *GetRawNoDataValue() const override
193 : {
194 71 : return m_abyRawNoDataValue.empty() ? nullptr
195 71 : : m_abyRawNoDataValue.data();
196 : }
197 :
198 40 : const std::string &GetUnit() const override
199 : {
200 40 : return m_osUnit;
201 : }
202 :
203 20 : std::shared_ptr<OGRSpatialReference> GetSpatialRef() const override
204 : {
205 40 : auto leftSRS = m_arrayLeft->GetSpatialRef();
206 40 : auto rightSRS = m_arrayRight->GetSpatialRef();
207 20 : if (!leftSRS || !rightSRS || !leftSRS->IsSame(rightSRS.get()))
208 19 : return nullptr;
209 1 : return leftSRS;
210 : }
211 :
212 20 : double GetOffset(bool *pbHasOffset,
213 : GDALDataType *peStorageType) const override
214 : {
215 20 : bool bHasLeftOffset = false;
216 20 : GDALDataType leftStorageType = GDT_Unknown;
217 : const double dfLeftOffset =
218 20 : m_arrayLeft->GetOffset(&bHasLeftOffset, &leftStorageType);
219 :
220 20 : bool bHasRightOffset = false;
221 20 : GDALDataType rightStorageType = GDT_Unknown;
222 : const double dfRightOffset =
223 20 : m_arrayRight->GetOffset(&bHasRightOffset, &rightStorageType);
224 :
225 20 : bool bHasLeftScale = false;
226 20 : const double dfLeftScale = m_arrayLeft->GetScale(&bHasLeftScale);
227 :
228 20 : bool bHasRightScale = false;
229 20 : const double dfRightScale = m_arrayRight->GetScale(&bHasRightScale);
230 :
231 20 : bool bRet = false;
232 20 : double dfRes = 0.0;
233 :
234 20 : switch (m_op)
235 : {
236 7 : case Operation::OP_ADD:
237 : {
238 10 : bRet = bHasLeftOffset && bHasRightOffset &&
239 3 : (bHasLeftScale == bHasRightScale &&
240 3 : (!bHasLeftScale || dfLeftScale == dfRightScale));
241 7 : if (bRet)
242 2 : dfRes = dfLeftOffset + dfRightOffset;
243 7 : break;
244 : }
245 :
246 7 : case Operation::OP_SUBTRACT:
247 : {
248 10 : bRet = bHasLeftOffset && bHasRightOffset &&
249 3 : (bHasLeftScale == bHasRightScale &&
250 3 : (!bHasLeftScale || dfLeftScale == dfRightScale));
251 7 : if (bRet)
252 2 : dfRes = dfLeftOffset - dfRightOffset;
253 7 : break;
254 : }
255 :
256 6 : case Operation::OP_MULTIPLY:
257 : case Operation::OP_DIVIDE:
258 6 : break;
259 : }
260 :
261 20 : if (pbHasOffset)
262 20 : *pbHasOffset = bRet;
263 20 : if (bRet && peStorageType)
264 : {
265 0 : *peStorageType =
266 0 : GDALDataTypeUnion(leftStorageType, rightStorageType);
267 : }
268 :
269 20 : return dfRes;
270 : }
271 :
272 8 : double GetScale(bool *pbHasScale,
273 : GDALDataType *peStorageType) const override
274 : {
275 8 : bool bHasLeftOffset = false;
276 8 : const double dfLeftOffset = m_arrayLeft->GetOffset(&bHasLeftOffset);
277 :
278 8 : bool bHasRightOffset = false;
279 8 : const double dfRightOffset = m_arrayRight->GetOffset(&bHasRightOffset);
280 :
281 8 : bool bHasLeftScale = false;
282 8 : GDALDataType leftStorageType = GDT_Unknown;
283 : const double dfLeftScale =
284 8 : m_arrayLeft->GetScale(&bHasLeftScale, &leftStorageType);
285 :
286 8 : bool bHasRightScale = false;
287 8 : GDALDataType rightStorageType = GDT_Unknown;
288 : const double dfRightScale =
289 8 : m_arrayRight->GetScale(&bHasRightScale, &rightStorageType);
290 :
291 8 : bool bRet = false;
292 8 : double dfRes = 1.0;
293 :
294 8 : const bool bZeroOffset =
295 16 : (bHasLeftOffset == bHasRightOffset &&
296 8 : (!bHasLeftOffset ||
297 4 : (dfLeftOffset == dfRightOffset && dfLeftOffset == 0)));
298 :
299 8 : switch (m_op)
300 : {
301 4 : case Operation::OP_ADD:
302 : case Operation::OP_SUBTRACT:
303 : {
304 6 : bRet = bZeroOffset &&
305 2 : (bHasLeftScale == bHasRightScale &&
306 2 : (!bHasLeftScale || dfLeftScale == dfRightScale));
307 4 : if (bRet)
308 2 : dfRes = dfLeftScale;
309 4 : break;
310 : }
311 :
312 2 : case Operation::OP_MULTIPLY:
313 : {
314 2 : bRet = bZeroOffset && bHasLeftScale && bHasRightScale;
315 2 : if (bRet)
316 1 : dfRes = dfLeftScale * dfRightScale;
317 2 : break;
318 : }
319 :
320 2 : case Operation::OP_DIVIDE:
321 : {
322 2 : bRet = bZeroOffset && bHasLeftScale && bHasRightScale;
323 2 : if (bRet)
324 1 : dfRes = dfLeftScale / dfRightScale;
325 2 : break;
326 : }
327 : }
328 :
329 8 : if (pbHasScale)
330 8 : *pbHasScale = bRet;
331 8 : if (bRet && peStorageType)
332 : {
333 0 : *peStorageType =
334 0 : GDALDataTypeUnion(leftStorageType, rightStorageType);
335 : }
336 :
337 8 : return dfRes;
338 : }
339 :
340 32 : std::vector<GUInt64> GetBlockSize() const override
341 : {
342 64 : const auto leftBlockSize = m_arrayLeft->GetBlockSize();
343 64 : const auto rightBlockSize = m_arrayRight->GetBlockSize();
344 32 : if (leftBlockSize != rightBlockSize)
345 2 : return GDALMDArray::GetBlockSize();
346 30 : return leftBlockSize;
347 : }
348 :
349 : bool IRead(const GUInt64 *arrayStartIdx, const size_t *count,
350 : const GInt64 *arrayStep, const GPtrDiff_t *bufferStride,
351 : const GDALExtendedDataType &bufferDataType,
352 : void *pDstBuffer) const override;
353 :
354 : private:
355 : const std::shared_ptr<GDALMDArray> m_arrayLeft;
356 : const std::shared_ptr<GDALMDArray> m_arrayRight;
357 : const Operation m_op;
358 : const GDALExtendedDataType m_dt;
359 : bool m_bHasNoDataLeft = false;
360 : bool m_bHasNoDataRight = false;
361 : const double m_dfNoDataLeft;
362 : const double m_dfNoDataRight;
363 : double m_dfNoData = std::numeric_limits<double>::quiet_NaN();
364 : std::vector<GByte> m_abyRawNoDataValue{};
365 : std::string m_osUnit{};
366 : mutable std::vector<double> m_leftValues{};
367 : mutable std::vector<double> m_rightValues{};
368 :
369 896 : inline bool IsInvalid(double dfLeft) const
370 : {
371 1792 : return std::isnan(dfLeft) ||
372 1792 : (m_bHasNoDataLeft && m_dfNoDataLeft == dfLeft);
373 : }
374 :
375 899 : inline bool IsInvalidTuple(double dfLeft, double dfRight) const
376 : {
377 899 : return std::isnan(dfLeft) ||
378 1785 : (m_bHasNoDataLeft && m_dfNoDataLeft == dfLeft) ||
379 2684 : std::isnan(dfRight) ||
380 1785 : (m_bHasNoDataRight && m_dfNoDataRight == dfRight);
381 : }
382 :
383 : template <class Op>
384 : void PerformOperation(double *leftValues, const double *rightValues,
385 : size_t nElts, bool bIntegerAndAllValid) const;
386 : };
387 :
388 : /************************************************************************/
389 : /* GDALMathOperationMDArray::Create() */
390 : /************************************************************************/
391 :
392 : /* static */ std::shared_ptr<GDALMathOperationMDArray>
393 152 : GDALMathOperationMDArray::Create(const std::shared_ptr<GDALMDArray> &arrayLeft,
394 : const std::shared_ptr<GDALMDArray> &arrayRight,
395 : Operation op)
396 : {
397 152 : if (!arrayLeft)
398 : {
399 0 : CPLError(CE_Failure, CPLE_AppDefined, "arrayLeft is null");
400 0 : return nullptr;
401 : }
402 :
403 152 : if (!arrayRight)
404 : {
405 0 : CPLError(CE_Failure, CPLE_AppDefined, "arrayRight is null");
406 0 : return nullptr;
407 : }
408 :
409 152 : if (!arrayLeft->HasSameShapeAs(*arrayRight))
410 : {
411 4 : CPLError(CE_Failure, CPLE_AppDefined, "%s",
412 8 : ("Arrays " + arrayLeft->GetFullName() + " and " +
413 12 : arrayRight->GetFullName() + " do not have the same shape")
414 : .c_str());
415 4 : return nullptr;
416 : }
417 :
418 737 : for (const auto &array : {arrayLeft, arrayRight})
419 : {
420 295 : if (array->GetDataType().GetClass() != GEDTC_NUMERIC)
421 : {
422 2 : CPLError(
423 : CE_Failure, CPLE_AppDefined, "%s",
424 4 : ("Array " + array->GetFullName() + " is not numeric").c_str());
425 2 : return nullptr;
426 : }
427 : }
428 :
429 : auto newAr(std::shared_ptr<GDALMathOperationMDArray>(
430 292 : new GDALMathOperationMDArray(arrayLeft, arrayRight, op)));
431 146 : newAr->SetSelf(newAr);
432 146 : return newAr;
433 : }
434 :
435 : /************************************************************************/
436 : /* GDALMathOperationMDArray::PerformOperation() */
437 : /************************************************************************/
438 :
439 : template <class Op>
440 : #if defined(__GNUC__) && !defined(__clang__)
441 : __attribute__((optimize("tree-vectorize")))
442 : #endif
443 111 : void GDALMathOperationMDArray::PerformOperation(double *leftValues,
444 : const double *rightValues,
445 : size_t nElts,
446 : bool bIntegerAndAllValid) const
447 : {
448 111 : if (bIntegerAndAllValid)
449 : {
450 28 : if (!rightValues)
451 : {
452 21 : for (size_t i = 0; i < nElts; ++i)
453 : {
454 18 : leftValues[i] = Op()(leftValues[i], leftValues[i]);
455 : }
456 : }
457 : else
458 : {
459 4467 : for (size_t i = 0; i < nElts; ++i)
460 : {
461 4442 : leftValues[i] = Op()(leftValues[i], rightValues[i]);
462 : }
463 : }
464 : }
465 : else
466 : {
467 83 : if (!rightValues)
468 : {
469 54 : if ((!m_bHasNoDataLeft || std::isnan(m_dfNoDataLeft)) &&
470 18 : std::isnan(m_dfNoData))
471 : {
472 1012 : for (size_t i = 0; i < nElts; ++i)
473 : {
474 994 : leftValues[i] = Op()(leftValues[i], leftValues[i]);
475 : }
476 : }
477 : else
478 : {
479 914 : for (size_t i = 0; i < nElts; ++i)
480 : {
481 896 : if (IsInvalid(leftValues[i]))
482 8 : leftValues[i] = m_dfNoData;
483 : else
484 888 : leftValues[i] = Op()(leftValues[i], leftValues[i]);
485 : }
486 : }
487 : }
488 : else
489 : {
490 23 : if ((!m_bHasNoDataLeft || std::isnan(m_dfNoDataLeft)) &&
491 90 : (!m_bHasNoDataRight || std::isnan(m_dfNoDataRight)) &&
492 20 : std::isnan(m_dfNoData))
493 : {
494 364 : for (size_t i = 0; i < nElts; ++i)
495 : {
496 344 : leftValues[i] = Op()(leftValues[i], rightValues[i]);
497 : }
498 : }
499 : else
500 : {
501 926 : for (size_t i = 0; i < nElts; ++i)
502 : {
503 899 : if (IsInvalidTuple(leftValues[i], rightValues[i]))
504 26 : leftValues[i] = m_dfNoData;
505 : else
506 873 : leftValues[i] = Op()(leftValues[i], rightValues[i]);
507 : }
508 : }
509 : }
510 : }
511 111 : }
512 :
513 : /************************************************************************/
514 : /* GDALMathOperationMDArray::IRead() */
515 : /************************************************************************/
516 :
517 114 : bool GDALMathOperationMDArray::IRead(const GUInt64 *arrayStartIdx,
518 : const size_t *count,
519 : const GInt64 *arrayStep,
520 : const GPtrDiff_t *bufferStride,
521 : const GDALExtendedDataType &bufferDataType,
522 : void *pDstBuffer) const
523 : {
524 114 : if (bufferDataType.GetClass() != GEDTC_NUMERIC)
525 : {
526 0 : CPLError(CE_Failure, CPLE_AppDefined,
527 : "GDALMathOperationMDArray::IRead(): not supported with "
528 : "non-numeric buffer data type");
529 0 : return false;
530 : }
531 :
532 114 : const auto &apoDims = GetDimensions();
533 114 : const size_t nDims = apoDims.size();
534 :
535 : const bool bIntegerAndAllValid =
536 114 : m_abyRawNoDataValue.empty() &&
537 69 : GDALDataTypeIsInteger(
538 183 : m_arrayLeft->GetDataType().GetNumericDataType()) &&
539 31 : GDALDataTypeIsInteger(m_arrayRight->GetDataType().GetNumericDataType());
540 :
541 131 : if (bIntegerAndAllValid && m_op == Operation::OP_SUBTRACT &&
542 17 : m_arrayLeft == m_arrayRight)
543 : {
544 3 : CopyContiguousBufferToBuffer(nDims, count, nullptr,
545 6 : GDALExtendedDataType::Create(GDT_Unknown),
546 : pDstBuffer, bufferDataType, bufferStride);
547 3 : return true;
548 : }
549 :
550 : // Check if the output buffer is in contiguous row major order
551 : // If this is true and its type is also the type of the array, we can
552 : // directly read the left array into the output buffer and update it in
553 : // place.
554 111 : bool bRowMajorBufferStride = true;
555 111 : GPtrDiff_t nExpectedStrideValue = 1;
556 277 : for (size_t i = nDims; i > 0; /* decremented in loop */)
557 : {
558 169 : --i;
559 169 : if (bufferStride[i] != nExpectedStrideValue)
560 : {
561 3 : bRowMajorBufferStride = false;
562 3 : break;
563 : }
564 : // The expected stride derives from the requested count, not
565 : // from the size of the array dimensions, which are only equal when
566 : // the whole array is requested.
567 166 : nExpectedStrideValue =
568 166 : static_cast<GPtrDiff_t>(nExpectedStrideValue * count[i]);
569 : }
570 :
571 111 : size_t nElts = 1;
572 111 : bool bFullArrayRequested = true;
573 280 : for (size_t i = 0; i < nDims; ++i)
574 : {
575 169 : nElts *= count[i];
576 169 : if (bFullArrayRequested)
577 169 : bFullArrayRequested = (count[i] == apoDims[i]->GetSize());
578 : }
579 :
580 : const bool bLeftValuesAllocNeeded =
581 111 : !bRowMajorBufferStride || bufferDataType != m_dt;
582 :
583 : // We can skip both I/O and memory allocations if operating on the same
584 : // array.
585 111 : const bool bRightValuesAllocNeeded = m_arrayLeft != m_arrayRight;
586 :
587 : try
588 : {
589 111 : if (bLeftValuesAllocNeeded)
590 : {
591 5 : if (nElts > m_leftValues.size())
592 5 : m_leftValues.resize(nElts);
593 : }
594 111 : if (bRightValuesAllocNeeded)
595 : {
596 72 : if (nElts > m_rightValues.size())
597 72 : m_rightValues.resize(nElts);
598 : }
599 : }
600 0 : catch (const std::exception &)
601 : {
602 0 : CPLError(CE_Failure, CPLE_OutOfMemory,
603 : "GDALMathOperationMDArray::IRead(): out of memory");
604 0 : return false;
605 : }
606 :
607 : double *leftValues = bLeftValuesAllocNeeded
608 111 : ? m_leftValues.data()
609 111 : : static_cast<double *>(pDstBuffer);
610 :
611 111 : const bool ret = m_arrayLeft->Read(arrayStartIdx, count, arrayStep, nullptr,
612 222 : m_dt, leftValues) &&
613 111 : (!bRightValuesAllocNeeded ||
614 144 : m_arrayRight->Read(arrayStartIdx, count, arrayStep,
615 72 : nullptr, m_dt, m_rightValues.data()));
616 :
617 : const double *rightValues =
618 111 : bRightValuesAllocNeeded ? m_rightValues.data() : nullptr;
619 :
620 111 : if (ret)
621 : {
622 111 : switch (m_op)
623 : {
624 21 : case Operation::OP_ADD:
625 : {
626 21 : PerformOperation<std::plus<double>>(leftValues, rightValues,
627 : nElts, bIntegerAndAllValid);
628 21 : break;
629 : }
630 44 : case Operation::OP_SUBTRACT:
631 : {
632 44 : PerformOperation<std::minus<double>>(
633 : leftValues, rightValues, nElts, bIntegerAndAllValid);
634 44 : break;
635 : }
636 29 : case Operation::OP_MULTIPLY:
637 : {
638 29 : PerformOperation<std::multiplies<double>>(
639 : leftValues, rightValues, nElts, bIntegerAndAllValid);
640 29 : break;
641 : }
642 17 : case Operation::OP_DIVIDE:
643 : {
644 17 : PerformOperation<std::divides<double>>(
645 : leftValues, rightValues, nElts, bIntegerAndAllValid);
646 17 : break;
647 : }
648 : }
649 :
650 111 : if (bLeftValuesAllocNeeded)
651 : {
652 5 : CopyContiguousBufferToBuffer(nDims, count, leftValues, m_dt,
653 : pDstBuffer, bufferDataType,
654 : bufferStride);
655 : }
656 : }
657 :
658 111 : if (bFullArrayRequested)
659 : {
660 108 : m_leftValues.clear();
661 108 : m_rightValues.clear();
662 : }
663 :
664 111 : return ret;
665 : }
666 :
667 : //! @endcond
668 :
669 : /************************************************************************/
670 : /* operator+() */
671 : /************************************************************************/
672 :
673 : /** Add this array with another one of the same shape.
674 : *
675 : * The resulting array is lazy evaluated.
676 : *
677 : * The resulting array type is Float64.
678 : *
679 : * The operation is nodata-aware.
680 : *
681 : * This is the same as C function GDALMDArrayBinaryOperation().
682 : *
683 : * @since 3.14
684 : * @return a new array, or nullptr in case of error
685 : */
686 : std::shared_ptr<GDALMDArray>
687 51 : GDALMDArray::operator+(const std::shared_ptr<GDALMDArray> &other) const
688 : {
689 102 : return GDALMathOperationMDArray::Create(
690 153 : GetSelf(), other, GDALMathOperationMDArray::Operation::OP_ADD);
691 : }
692 :
693 : /************************************************************************/
694 : /* operator-() */
695 : /************************************************************************/
696 :
697 : /** Subtract this array with another one of the same shape.
698 : *
699 : * The resulting array is lazy evaluated.
700 : *
701 : * The resulting array type is Float64.
702 : *
703 : * The operation is nodata-aware.
704 : *
705 : * This is the same as C function GDALMDArrayBinaryOperation().
706 : *
707 : * @since 3.14
708 : * @return a new array, or nullptr in case of error
709 : */
710 : std::shared_ptr<GDALMDArray>
711 43 : GDALMDArray::operator-(const std::shared_ptr<GDALMDArray> &other) const
712 : {
713 86 : return GDALMathOperationMDArray::Create(
714 129 : GetSelf(), other, GDALMathOperationMDArray::Operation::OP_SUBTRACT);
715 : }
716 :
717 : /************************************************************************/
718 : /* operator*() */
719 : /************************************************************************/
720 :
721 : /** Multiply this array with another one of the same shape.
722 : *
723 : * The resulting array is lazy evaluated.
724 : *
725 : * The resulting array type is Float64.
726 : *
727 : * The operation is nodata-aware.
728 : *
729 : * This is the same as C function GDALMDArrayBinaryOperation().
730 : *
731 : * @since 3.14
732 : * @return a new array, or nullptr in case of error
733 : */
734 : std::shared_ptr<GDALMDArray>
735 32 : GDALMDArray::operator*(const std::shared_ptr<GDALMDArray> &other) const
736 : {
737 64 : return GDALMathOperationMDArray::Create(
738 96 : GetSelf(), other, GDALMathOperationMDArray::Operation::OP_MULTIPLY);
739 : }
740 :
741 : /************************************************************************/
742 : /* operator/() */
743 : /************************************************************************/
744 :
745 : /** Divide this array by another one of the same shape.
746 : *
747 : * The resulting array is lazy evaluated.
748 : *
749 : * The resulting array type is Float64.
750 : *
751 : * The operation is nodata-aware.
752 : *
753 : * This is the same as C function GDALMDArrayBinaryOperation().
754 : *
755 : * @since 3.14
756 : * @return a new array, or nullptr in case of error
757 : */
758 : std::shared_ptr<GDALMDArray>
759 26 : GDALMDArray::operator/(const std::shared_ptr<GDALMDArray> &other) const
760 : {
761 52 : return GDALMathOperationMDArray::Create(
762 78 : GetSelf(), other, GDALMathOperationMDArray::Operation::OP_DIVIDE);
763 : }
|
