Line data Source code
1 : #include "DataContainer.h"
2 : #include "DataContainerFormatter.hpp"
3 : #include "FormatConfig.h"
4 : #include "DataHandlerCPU.h"
5 : #include "DataHandlerMapCPU.h"
6 : #include "BlockDescriptor.h"
7 : #include "RandomBlocksDescriptor.h"
8 : #include "PartitionDescriptor.h"
9 : #include "Error.h"
10 : #include "TypeCasts.hpp"
11 : #include "Assertions.h"
12 :
13 : #include <utility>
14 : #include <algorithm>
15 :
16 : namespace elsa
17 : {
18 :
19 : template <typename data_t>
20 : DataContainer<data_t>::DataContainer(const DataDescriptor& dataDescriptor,
21 : DataHandlerType handlerType)
22 : : _dataDescriptor{dataDescriptor.clone()},
23 : _dataHandler{createDataHandler(handlerType, _dataDescriptor->getNumberOfCoefficients())},
24 : _dataHandlerType{handlerType}
25 63780 : {
26 63780 : }
27 :
28 : template <typename data_t>
29 : DataContainer<data_t>::DataContainer(const DataDescriptor& dataDescriptor,
30 : const Eigen::Matrix<data_t, Eigen::Dynamic, 1>& data,
31 : DataHandlerType handlerType)
32 : : _dataDescriptor{dataDescriptor.clone()},
33 : _dataHandler{createDataHandler(handlerType, _dataDescriptor->getNumberOfCoefficients())},
34 : _dataHandlerType{handlerType}
35 1527 : {
36 1527 : if (_dataHandler->getSize() != data.size())
37 0 : throw InvalidArgumentError("DataContainer: initialization vector has invalid size");
38 :
39 3602655 : for (index_t i = 0; i < _dataHandler->getSize(); ++i)
40 3601128 : (*_dataHandler)[i] = data[i];
41 1527 : }
42 :
43 : template <typename data_t>
44 : DataContainer<data_t>::DataContainer(const DataContainer<data_t>& other)
45 : : _dataDescriptor{other._dataDescriptor->clone()},
46 : _dataHandler{other._dataHandler->clone()},
47 : _dataHandlerType{other._dataHandlerType}
48 15751 : {
49 15751 : }
50 :
51 : template <typename data_t>
52 : DataContainer<data_t>& DataContainer<data_t>::operator=(const DataContainer<data_t>& other)
53 18330 : {
54 18330 : if (this != &other) {
55 18330 : _dataDescriptor = other._dataDescriptor->clone();
56 :
57 18330 : if (_dataHandler && canAssign(other._dataHandlerType)) {
58 18320 : *_dataHandler = *other._dataHandler;
59 18320 : } else {
60 10 : _dataHandler = other._dataHandler->clone();
61 10 : _dataHandlerType = other._dataHandlerType;
62 10 : }
63 18330 : }
64 :
65 18330 : return *this;
66 18330 : }
67 :
68 : template <typename data_t>
69 : DataContainer<data_t>::DataContainer(DataContainer<data_t>&& other) noexcept
70 : : _dataDescriptor{std::move(other._dataDescriptor)},
71 : _dataHandler{std::move(other._dataHandler)},
72 : _dataHandlerType{std::move(other._dataHandlerType)}
73 77 : {
74 : // leave other in a valid state
75 77 : other._dataDescriptor = nullptr;
76 77 : other._dataHandler = nullptr;
77 77 : }
78 :
79 : template <typename data_t>
80 : DataContainer<data_t>& DataContainer<data_t>::operator=(DataContainer<data_t>&& other)
81 8218 : {
82 8218 : _dataDescriptor = std::move(other._dataDescriptor);
83 :
84 8218 : if (_dataHandler && canAssign(other._dataHandlerType)) {
85 8218 : *_dataHandler = std::move(*other._dataHandler);
86 8218 : } else {
87 0 : _dataHandler = std::move(other._dataHandler);
88 0 : _dataHandlerType = std::move(other._dataHandlerType);
89 0 : }
90 :
91 : // leave other in a valid state
92 8218 : other._dataDescriptor = nullptr;
93 8218 : other._dataHandler = nullptr;
94 :
95 8218 : return *this;
96 8218 : }
97 :
98 : template <typename data_t>
99 : const DataDescriptor& DataContainer<data_t>::getDataDescriptor() const
100 1357875 : {
101 1357875 : return *_dataDescriptor;
102 1357875 : }
103 :
104 : template <typename data_t>
105 : index_t DataContainer<data_t>::getSize() const
106 18555240 : {
107 18555240 : return _dataHandler->getSize();
108 18555240 : }
109 :
110 : template <typename data_t>
111 : data_t& DataContainer<data_t>::operator[](index_t index)
112 9643048 : {
113 9643048 : ELSA_VERIFY(index >= 0);
114 9643048 : ELSA_VERIFY(index < getSize());
115 :
116 9643048 : return (*_dataHandler)[index];
117 9643048 : }
118 :
119 : template <typename data_t>
120 : const data_t& DataContainer<data_t>::operator[](index_t index) const
121 4131438 : {
122 4131438 : ELSA_VERIFY(index >= 0);
123 4131438 : ELSA_VERIFY(index < getSize());
124 :
125 4131438 : return static_cast<const DataHandler<data_t>&>(*_dataHandler)[index];
126 4131438 : }
127 :
128 : template <typename data_t>
129 : data_t DataContainer<data_t>::at(const IndexVector_t& coordinate) const
130 0 : {
131 0 : const auto arr = coordinate.array();
132 0 : if ((arr < 0).any()
133 0 : || (arr >= _dataDescriptor->getNumberOfCoefficientsPerDimension().array()).any()) {
134 0 : return 0;
135 0 : }
136 :
137 0 : return (*this)[_dataDescriptor->getIndexFromCoordinate(coordinate)];
138 0 : }
139 :
140 : template <typename data_t>
141 : data_t& DataContainer<data_t>::operator()(const IndexVector_t& coordinate)
142 1299528 : {
143 : // const auto arr = coordinate.array();
144 : // const auto shape = _dataDescriptor->getNumberOfCoefficientsPerDimension().array();
145 : // ELSA_VERIFY((arr >= 0).all());
146 : // ELSA_VERIFY((arr < shape).all());
147 :
148 1299528 : return (*this)[_dataDescriptor->getIndexFromCoordinate(coordinate)];
149 1299528 : }
150 :
151 : template <typename data_t>
152 : const data_t& DataContainer<data_t>::operator()(const IndexVector_t& coordinate) const
153 203495 : {
154 : // const auto arr = coordinate.array();
155 : // const auto shape = _dataDescriptor->getNumberOfCoefficientsPerDimension().array();
156 : // ELSA_VERIFY((arr >= 0).all());
157 : // ELSA_VERIFY((arr < shape).all());
158 :
159 203495 : return (*this)[_dataDescriptor->getIndexFromCoordinate(coordinate)];
160 203495 : }
161 :
162 : template <typename data_t>
163 : data_t DataContainer<data_t>::dot(const DataContainer<data_t>& other) const
164 632 : {
165 632 : return _dataHandler->dot(*other._dataHandler);
166 632 : }
167 :
168 : template <typename data_t>
169 : GetFloatingPointType_t<data_t> DataContainer<data_t>::squaredL2Norm() const
170 21528 : {
171 21528 : return _dataHandler->squaredL2Norm();
172 21528 : }
173 :
174 : template <typename data_t>
175 : GetFloatingPointType_t<data_t> DataContainer<data_t>::l2Norm() const
176 2729 : {
177 2729 : return _dataHandler->l2Norm();
178 2729 : }
179 :
180 : template <typename data_t>
181 : index_t DataContainer<data_t>::l0PseudoNorm() const
182 9 : {
183 9 : return _dataHandler->l0PseudoNorm();
184 9 : }
185 :
186 : template <typename data_t>
187 : GetFloatingPointType_t<data_t> DataContainer<data_t>::l1Norm() const
188 13 : {
189 13 : return _dataHandler->l1Norm();
190 13 : }
191 :
192 : template <typename data_t>
193 : GetFloatingPointType_t<data_t> DataContainer<data_t>::lInfNorm() const
194 13 : {
195 13 : return _dataHandler->lInfNorm();
196 13 : }
197 :
198 : template <typename data_t>
199 : data_t DataContainer<data_t>::sum() const
200 5 : {
201 5 : return _dataHandler->sum();
202 5 : }
203 :
204 : template <typename data_t>
205 : data_t DataContainer<data_t>::minElement() const
206 36 : {
207 36 : return _dataHandler->minElement();
208 36 : }
209 :
210 : template <typename data_t>
211 : data_t DataContainer<data_t>::maxElement() const
212 14 : {
213 14 : return _dataHandler->maxElement();
214 14 : }
215 :
216 : template <typename data_t>
217 : void DataContainer<data_t>::fft(FFTNorm norm) const
218 146 : {
219 146 : this->_dataHandler->fft(*this->_dataDescriptor, norm);
220 146 : }
221 :
222 : template <typename data_t>
223 : void DataContainer<data_t>::ifft(FFTNorm norm) const
224 250 : {
225 250 : this->_dataHandler->ifft(*this->_dataDescriptor, norm);
226 250 : }
227 :
228 : template <typename data_t>
229 : DataContainer<data_t>& DataContainer<data_t>::operator+=(const DataContainer<data_t>& dc)
230 8327 : {
231 8327 : *_dataHandler += *dc._dataHandler;
232 8327 : return *this;
233 8327 : }
234 :
235 : template <typename data_t>
236 : DataContainer<data_t>& DataContainer<data_t>::operator-=(const DataContainer<data_t>& dc)
237 11740 : {
238 11740 : *_dataHandler -= *dc._dataHandler;
239 11740 : return *this;
240 11740 : }
241 :
242 : template <typename data_t>
243 : DataContainer<data_t>& DataContainer<data_t>::operator*=(const DataContainer<data_t>& dc)
244 5 : {
245 5 : *_dataHandler *= *dc._dataHandler;
246 5 : return *this;
247 5 : }
248 :
249 : template <typename data_t>
250 : DataContainer<data_t>& DataContainer<data_t>::operator/=(const DataContainer<data_t>& dc)
251 5 : {
252 5 : *_dataHandler /= *dc._dataHandler;
253 5 : return *this;
254 5 : }
255 :
256 : template <typename data_t>
257 : DataContainer<data_t>& DataContainer<data_t>::operator+=(data_t scalar)
258 5 : {
259 5 : *_dataHandler += scalar;
260 5 : return *this;
261 5 : }
262 :
263 : template <typename data_t>
264 : DataContainer<data_t>& DataContainer<data_t>::operator-=(data_t scalar)
265 5 : {
266 5 : *_dataHandler -= scalar;
267 5 : return *this;
268 5 : }
269 :
270 : template <typename data_t>
271 : DataContainer<data_t>& DataContainer<data_t>::operator*=(data_t scalar)
272 504 : {
273 504 : *_dataHandler *= scalar;
274 504 : return *this;
275 504 : }
276 :
277 : template <typename data_t>
278 : DataContainer<data_t>& DataContainer<data_t>::operator/=(data_t scalar)
279 388 : {
280 388 : *_dataHandler /= scalar;
281 388 : return *this;
282 388 : }
283 :
284 : template <typename data_t>
285 : DataContainer<data_t>& DataContainer<data_t>::operator=(data_t scalar)
286 3585 : {
287 3585 : *_dataHandler = scalar;
288 3585 : return *this;
289 3585 : }
290 :
291 : template <typename data_t>
292 : template <typename... Args>
293 : std::unique_ptr<DataHandler<data_t>>
294 : DataContainer<data_t>::createDataHandler(DataHandlerType handlerType, Args&&... args)
295 65307 : {
296 65307 : switch (handlerType) {
297 64447 : case DataHandlerType::CPU:
298 64447 : return std::make_unique<DataHandlerCPU<data_t>>(std::forward<Args>(args)...);
299 860 : case DataHandlerType::MAP_CPU:
300 860 : return std::make_unique<DataHandlerCPU<data_t>>(std::forward<Args>(args)...);
301 : #ifdef ELSA_CUDA_VECTOR
302 : case DataHandlerType::GPU:
303 : return std::make_unique<DataHandlerGPU<data_t>>(std::forward<Args>(args)...);
304 : case DataHandlerType::MAP_GPU:
305 : return std::make_unique<DataHandlerGPU<data_t>>(std::forward<Args>(args)...);
306 : #endif
307 0 : default:
308 0 : throw InvalidArgumentError("DataContainer: unknown handler type");
309 65307 : }
310 65307 : }
311 :
312 : template <typename data_t>
313 : DataContainer<data_t>::DataContainer(const DataDescriptor& dataDescriptor,
314 : std::unique_ptr<DataHandler<data_t>> dataHandler,
315 : DataHandlerType dataType)
316 : : _dataDescriptor{dataDescriptor.clone()},
317 : _dataHandler{std::move(dataHandler)},
318 : _dataHandlerType{dataType}
319 4727 : {
320 4727 : }
321 :
322 : template <typename data_t>
323 : bool DataContainer<data_t>::operator==(const DataContainer<data_t>& other) const
324 1729 : {
325 1729 : if (*_dataDescriptor != *other._dataDescriptor)
326 0 : return false;
327 :
328 1729 : if (*_dataHandler != *other._dataHandler)
329 0 : return false;
330 :
331 1729 : return true;
332 1729 : }
333 :
334 : template <typename data_t>
335 : bool DataContainer<data_t>::operator!=(const DataContainer<data_t>& other) const
336 312 : {
337 312 : return !(*this == other);
338 312 : }
339 :
340 : template <typename data_t>
341 : DataContainer<data_t> DataContainer<data_t>::getBlock(index_t i)
342 2011 : {
343 2011 : const auto blockDesc = downcast_safe<BlockDescriptor>(_dataDescriptor.get());
344 2011 : if (!blockDesc)
345 5 : throw LogicError("DataContainer: cannot get block from not-blocked container");
346 :
347 2006 : if (i >= blockDesc->getNumberOfBlocks() || i < 0)
348 0 : throw InvalidArgumentError("DataContainer: block index out of bounds");
349 :
350 2006 : index_t startIndex = blockDesc->getOffsetOfBlock(i);
351 2006 : const auto& ithDesc = blockDesc->getDescriptorOfBlock(i);
352 2006 : index_t blockSize = ithDesc.getNumberOfCoefficients();
353 :
354 2006 : DataHandlerType newHandlerType = (_dataHandlerType == DataHandlerType::CPU
355 2006 : || _dataHandlerType == DataHandlerType::MAP_CPU)
356 2006 : ? DataHandlerType::MAP_CPU
357 2006 : : DataHandlerType::MAP_GPU;
358 :
359 2006 : return DataContainer<data_t>{ithDesc, _dataHandler->getBlock(startIndex, blockSize),
360 2006 : newHandlerType};
361 2006 : }
362 :
363 : template <typename data_t>
364 : const DataContainer<data_t> DataContainer<data_t>::getBlock(index_t i) const
365 832 : {
366 832 : const auto blockDesc = downcast_safe<BlockDescriptor>(_dataDescriptor.get());
367 832 : if (!blockDesc)
368 5 : throw LogicError("DataContainer: cannot get block from not-blocked container");
369 :
370 827 : if (i >= blockDesc->getNumberOfBlocks() || i < 0)
371 0 : throw InvalidArgumentError("DataContainer: block index out of bounds");
372 :
373 827 : index_t startIndex = blockDesc->getOffsetOfBlock(i);
374 827 : const auto& ithDesc = blockDesc->getDescriptorOfBlock(i);
375 827 : index_t blockSize = ithDesc.getNumberOfCoefficients();
376 :
377 827 : DataHandlerType newHandlerType = (_dataHandlerType == DataHandlerType::CPU
378 827 : || _dataHandlerType == DataHandlerType::MAP_CPU)
379 827 : ? DataHandlerType::MAP_CPU
380 827 : : DataHandlerType::MAP_GPU;
381 :
382 : // getBlock() returns a pointer to non-const DH, but that's fine as it gets wrapped in a
383 : // constant container
384 827 : return DataContainer<data_t>{ithDesc, _dataHandler->getBlock(startIndex, blockSize),
385 827 : newHandlerType};
386 827 : }
387 :
388 : template <typename data_t>
389 : DataContainer<data_t> DataContainer<data_t>::viewAs(const DataDescriptor& dataDescriptor)
390 1154 : {
391 1154 : if (dataDescriptor.getNumberOfCoefficients() != getSize())
392 0 : throw InvalidArgumentError("DataContainer: view must have same size as container");
393 :
394 1154 : DataHandlerType newHandlerType = (_dataHandlerType == DataHandlerType::CPU
395 1154 : || _dataHandlerType == DataHandlerType::MAP_CPU)
396 1154 : ? DataHandlerType::MAP_CPU
397 1154 : : DataHandlerType::MAP_GPU;
398 :
399 1154 : return DataContainer<data_t>{dataDescriptor, _dataHandler->getBlock(0, getSize()),
400 1154 : newHandlerType};
401 1154 : }
402 :
403 : template <typename data_t>
404 : const DataContainer<data_t>
405 : DataContainer<data_t>::viewAs(const DataDescriptor& dataDescriptor) const
406 740 : {
407 740 : if (dataDescriptor.getNumberOfCoefficients() != getSize())
408 0 : throw InvalidArgumentError("DataContainer: view must have same size as container");
409 :
410 740 : DataHandlerType newHandlerType = (_dataHandlerType == DataHandlerType::CPU
411 740 : || _dataHandlerType == DataHandlerType::MAP_CPU)
412 740 : ? DataHandlerType::MAP_CPU
413 740 : : DataHandlerType::MAP_GPU;
414 :
415 : // getBlock() returns a pointer to non-const DH, but that's fine as it gets wrapped in a
416 : // constant container
417 740 : return DataContainer<data_t>{dataDescriptor, _dataHandler->getBlock(0, getSize()),
418 740 : newHandlerType};
419 740 : }
420 :
421 : template <typename data_t>
422 : const DataContainer<data_t> DataContainer<data_t>::slice(index_t i) const
423 614 : {
424 614 : auto& desc = getDataDescriptor();
425 614 : auto dim = desc.getNumberOfDimensions();
426 614 : auto sizeOfLastDim = desc.getNumberOfCoefficientsPerDimension()[dim - 1];
427 :
428 614 : if (i >= sizeOfLastDim) {
429 8 : throw LogicError("Trying to access out of bound slice");
430 8 : }
431 :
432 606 : if (sizeOfLastDim == 1) {
433 4 : return *this;
434 4 : }
435 :
436 602 : auto sliceDesc = PartitionDescriptor(desc, sizeOfLastDim);
437 :
438 : // Now set the slice
439 602 : return viewAs(sliceDesc).getBlock(i);
440 602 : }
441 :
442 : template <typename data_t>
443 : DataContainer<data_t> DataContainer<data_t>::slice(index_t i)
444 890 : {
445 890 : auto& desc = getDataDescriptor();
446 890 : auto dim = desc.getNumberOfDimensions();
447 890 : auto sizeOfLastDim = desc.getNumberOfCoefficientsPerDimension()[dim - 1];
448 :
449 890 : if (i >= sizeOfLastDim) {
450 4 : throw LogicError("Trying to access out of bound slice");
451 4 : }
452 :
453 886 : if (sizeOfLastDim == 1) {
454 4 : return *this;
455 4 : }
456 :
457 882 : auto sliceDesc = PartitionDescriptor(desc, sizeOfLastDim);
458 :
459 : // Now set the slice
460 882 : return viewAs(sliceDesc).getBlock(i);
461 882 : }
462 :
463 : template <typename data_t>
464 : typename DataContainer<data_t>::iterator DataContainer<data_t>::begin()
465 904 : {
466 904 : return iterator(&(*this)[0]);
467 904 : }
468 :
469 : template <typename data_t>
470 : typename DataContainer<data_t>::const_iterator DataContainer<data_t>::begin() const
471 866 : {
472 866 : return cbegin();
473 866 : }
474 :
475 : template <typename data_t>
476 : typename DataContainer<data_t>::const_iterator DataContainer<data_t>::cbegin() const
477 9296 : {
478 9296 : return const_iterator(&(*this)[0]);
479 9296 : }
480 :
481 : template <typename data_t>
482 : typename DataContainer<data_t>::iterator DataContainer<data_t>::end()
483 519228 : {
484 519228 : return iterator(&(*this)[0] + getSize());
485 519228 : }
486 :
487 : template <typename data_t>
488 : typename DataContainer<data_t>::const_iterator DataContainer<data_t>::end() const
489 520056 : {
490 520056 : return cend();
491 520056 : }
492 :
493 : template <typename data_t>
494 : typename DataContainer<data_t>::const_iterator DataContainer<data_t>::cend() const
495 520086 : {
496 520086 : return const_iterator(&(*this)[0] + getSize());
497 520086 : }
498 :
499 : template <typename data_t>
500 : typename DataContainer<data_t>::reverse_iterator DataContainer<data_t>::rbegin()
501 0 : {
502 0 : return reverse_iterator(end());
503 0 : }
504 :
505 : template <typename data_t>
506 : typename DataContainer<data_t>::const_reverse_iterator DataContainer<data_t>::rbegin() const
507 0 : {
508 0 : return crbegin();
509 0 : }
510 :
511 : template <typename data_t>
512 : typename DataContainer<data_t>::const_reverse_iterator DataContainer<data_t>::crbegin() const
513 3 : {
514 3 : return const_reverse_iterator(cend());
515 3 : }
516 :
517 : template <typename data_t>
518 : typename DataContainer<data_t>::reverse_iterator DataContainer<data_t>::rend()
519 0 : {
520 0 : return reverse_iterator(begin());
521 0 : }
522 :
523 : template <typename data_t>
524 : typename DataContainer<data_t>::const_reverse_iterator DataContainer<data_t>::rend() const
525 0 : {
526 0 : return crend();
527 0 : }
528 :
529 : template <typename data_t>
530 : typename DataContainer<data_t>::const_reverse_iterator DataContainer<data_t>::crend() const
531 8423 : {
532 8423 : return const_reverse_iterator(cbegin());
533 8423 : }
534 :
535 : template <typename data_t>
536 : DataHandlerType DataContainer<data_t>::getDataHandlerType() const
537 146098 : {
538 146098 : return _dataHandlerType;
539 146098 : }
540 :
541 : template <typename data_t>
542 : void DataContainer<data_t>::format(std::ostream& os, format_config cfg) const
543 0 : {
544 0 : DataContainerFormatter<data_t> fmt{cfg};
545 0 : fmt.format(os, *this);
546 0 : }
547 :
548 : template <typename data_t>
549 : DataContainer<data_t> DataContainer<data_t>::loadToCPU()
550 0 : {
551 0 : if (_dataHandlerType == DataHandlerType::CPU
552 0 : || _dataHandlerType == DataHandlerType::MAP_CPU) {
553 0 : throw LogicError(
554 0 : "DataContainer: cannot load data to CPU with already CPU based container");
555 0 : }
556 :
557 0 : DataContainer<data_t> dcCPU(*_dataDescriptor, DataHandlerType::CPU);
558 :
559 0 : for (index_t i = 0; i < getSize(); i++) {
560 0 : dcCPU[i] = this->operator[](i);
561 0 : }
562 :
563 0 : return dcCPU;
564 0 : }
565 :
566 : template <typename data_t>
567 : DataContainer<data_t> DataContainer<data_t>::loadToGPU()
568 0 : {
569 0 : if (_dataHandlerType == DataHandlerType::GPU
570 0 : || _dataHandlerType == DataHandlerType::MAP_GPU) {
571 0 : throw LogicError(
572 0 : "DataContainer: cannot load data to GPU with already GPU based container");
573 0 : }
574 :
575 0 : DataContainer<data_t> dcGPU(*_dataDescriptor, DataHandlerType::GPU);
576 :
577 0 : for (index_t i = 0; i < getSize(); i++) {
578 0 : dcGPU[i] = this->operator[](i);
579 0 : }
580 :
581 0 : return dcGPU;
582 0 : }
583 :
584 : template <typename data_t>
585 : bool DataContainer<data_t>::canAssign(DataHandlerType handlerType)
586 26538 : {
587 26538 : if (_dataHandlerType == DataHandlerType::CPU
588 26538 : || _dataHandlerType == DataHandlerType::MAP_CPU) {
589 26538 : switch (handlerType) {
590 25915 : case DataHandlerType::CPU:
591 25915 : return true;
592 0 : break;
593 623 : case DataHandlerType::MAP_CPU:
594 623 : return true;
595 0 : break;
596 0 : default:
597 0 : return false;
598 0 : }
599 0 : } else {
600 0 : switch (handlerType) {
601 0 : case DataHandlerType::GPU:
602 0 : return true;
603 0 : break;
604 0 : case DataHandlerType::MAP_GPU:
605 0 : return true;
606 0 : break;
607 0 : default:
608 0 : return false;
609 0 : }
610 0 : }
611 26538 : }
612 :
613 : template <typename data_t>
614 : DataContainer<data_t> clip(DataContainer<data_t> dc, data_t min, data_t max)
615 20 : {
616 136 : std::transform(dc.begin(), dc.end(), dc.begin(), [&](auto x) {
617 136 : if (x < min) {
618 28 : return min;
619 108 : } else if (x > max) {
620 28 : return max;
621 80 : } else {
622 80 : return x;
623 80 : }
624 136 : });
625 :
626 20 : return dc;
627 20 : }
628 :
629 : template <typename data_t>
630 : DataContainer<data_t> concatenate(const DataContainer<data_t>& dc1,
631 : const DataContainer<data_t>& dc2)
632 28 : {
633 28 : auto desc1 = dc1.getDataDescriptor().clone();
634 28 : auto desc2 = dc2.getDataDescriptor().clone();
635 :
636 28 : if (desc1->getNumberOfDimensions() != desc2->getNumberOfDimensions()) {
637 4 : throw LogicError("Can't concatenate two DataContainers with different dimension");
638 4 : }
639 :
640 24 : std::vector<std::unique_ptr<DataDescriptor>> descriptors;
641 24 : descriptors.reserve(2);
642 24 : descriptors.push_back(std::move(desc1));
643 24 : descriptors.push_back(std::move(desc2));
644 :
645 24 : auto blockDesc = RandomBlocksDescriptor(descriptors);
646 24 : auto concatenated = DataContainer<data_t>(blockDesc);
647 :
648 24 : concatenated.getBlock(0) = dc1;
649 24 : concatenated.getBlock(1) = dc2;
650 24 : return concatenated;
651 24 : }
652 :
653 : template <typename data_t>
654 : DataContainer<data_t> fftShift2D(DataContainer<data_t> dc)
655 24 : {
656 24 : assert(dc.getDataDescriptor().getNumberOfDimensions() == 2
657 24 : && "DataContainer::fftShift2D: currently only supporting 2D signals");
658 :
659 24 : const DataDescriptor& dataDescriptor = dc.getDataDescriptor();
660 24 : IndexVector_t numOfCoeffsPerDim = dataDescriptor.getNumberOfCoefficientsPerDimension();
661 24 : index_t m = numOfCoeffsPerDim[0];
662 24 : index_t n = numOfCoeffsPerDim[1];
663 :
664 24 : index_t firstShift = m / 2;
665 24 : index_t secondShift = n / 2;
666 :
667 24 : DataContainer<data_t> copyDC(dataDescriptor);
668 :
669 96 : for (index_t i = 0; i < m; ++i) {
670 352 : for (index_t j = 0; j < n; ++j) {
671 280 : copyDC((i + firstShift) % m, (j + secondShift) % n) = dc(i, j);
672 280 : }
673 72 : }
674 :
675 24 : return copyDC;
676 24 : }
677 :
678 : template <typename data_t>
679 : DataContainer<data_t> ifftShift2D(DataContainer<data_t> dc)
680 24 : {
681 24 : assert(dc.getDataDescriptor().getNumberOfDimensions() == 2
682 24 : && "DataContainer::ifftShift2D: currently only supporting 2D signals");
683 :
684 24 : const DataDescriptor& dataDescriptor = dc.getDataDescriptor();
685 24 : IndexVector_t numOfCoeffsPerDim = dataDescriptor.getNumberOfCoefficientsPerDimension();
686 24 : index_t m = numOfCoeffsPerDim[0];
687 24 : index_t n = numOfCoeffsPerDim[1];
688 :
689 24 : index_t firstShift = -m / 2;
690 24 : index_t secondShift = -n / 2;
691 :
692 24 : DataContainer<data_t> copyDC(dataDescriptor);
693 :
694 96 : for (index_t i = 0; i < m; ++i) {
695 352 : for (index_t j = 0; j < n; ++j) {
696 280 : index_t leftIndex = (((i + firstShift) % m) + m) % m;
697 280 : index_t rightIndex = (((j + secondShift) % n) + n) % n;
698 280 : copyDC(leftIndex, rightIndex) = dc(i, j);
699 280 : }
700 72 : }
701 :
702 24 : return copyDC;
703 24 : }
704 :
705 : // ------------------------------------------
706 : // explicit template instantiation
707 : template class DataContainer<float>;
708 : template class DataContainer<complex<float>>;
709 : template class DataContainer<double>;
710 : template class DataContainer<complex<double>>;
711 : template class DataContainer<index_t>;
712 :
713 : template DataContainer<float> clip<float>(DataContainer<float> dc, float min, float max);
714 : template DataContainer<double> clip<double>(DataContainer<double> dc, double min, double max);
715 :
716 : template DataContainer<float> concatenate<float>(const DataContainer<float>&,
717 : const DataContainer<float>&);
718 : template DataContainer<double> concatenate<double>(const DataContainer<double>&,
719 : const DataContainer<double>&);
720 : template DataContainer<complex<float>>
721 : concatenate<complex<float>>(const DataContainer<complex<float>>&,
722 : const DataContainer<complex<float>>&);
723 : template DataContainer<complex<double>>
724 : concatenate<complex<double>>(const DataContainer<complex<double>>&,
725 : const DataContainer<complex<double>>&);
726 :
727 : template DataContainer<float> fftShift2D<float>(DataContainer<float>);
728 : template DataContainer<complex<float>>
729 : fftShift2D<complex<float>>(DataContainer<complex<float>>);
730 : template DataContainer<double> fftShift2D<double>(DataContainer<double>);
731 : template DataContainer<complex<double>>
732 : fftShift2D<complex<double>>(DataContainer<complex<double>>);
733 :
734 : template DataContainer<float> ifftShift2D<float>(DataContainer<float>);
735 : template DataContainer<complex<float>>
736 : ifftShift2D<complex<float>>(DataContainer<complex<float>>);
737 : template DataContainer<double> ifftShift2D<double>(DataContainer<double>);
738 : template DataContainer<complex<double>>
739 : ifftShift2D<complex<double>>(DataContainer<complex<double>>);
740 :
741 : } // namespace elsa
|
