diff --git a/src/buildblock/CMakeLists.txt b/src/buildblock/CMakeLists.txt index b7f0882614..1b17ee9c4b 100644 --- a/src/buildblock/CMakeLists.txt +++ b/src/buildblock/CMakeLists.txt @@ -15,7 +15,6 @@ set(${dir_LIB_SOURCES} ParsingObject.cxx num_threads.cxx Array.cxx - IndexRange.cxx PatientPosition.cxx TimeFrameDefinitions.cxx ExamInfo.cxx diff --git a/src/buildblock/IndexRange.cxx b/src/buildblock/IndexRange.cxx deleted file mode 100644 index c47b5a9188..0000000000 --- a/src/buildblock/IndexRange.cxx +++ /dev/null @@ -1,88 +0,0 @@ -// -// -/*! - \file - \ingroup Array - \brief implementations for the stir::IndexRange class - - \author Kris Thielemans - \author PARAPET project - - -*/ -/* - Copyright (C) 2000 PARAPET partners - Copyright (C) 2000- 2009, Hammersmith Imanet Ltd - This file is part of STIR. - - SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license - - See STIR/LICENSE.txt for details -*/ -#include "stir/IndexRange.h" -#include - -START_NAMESPACE_STIR - -template -bool -IndexRange::get_regular_range(BasicCoordinate& min, - BasicCoordinate& max) const -{ - // check if empty range - if (base_type::begin() == base_type::end()) - { - std::fill(min.begin(), min.end(), 0); - std::fill(max.begin(), max.end(), -1); - return true; - } - - // if not a regular range, exit - if (is_regular_range == regular_false) - return false; - - typename base_type::const_iterator iter = base_type::begin(); - - BasicCoordinate lower_dim_min; - BasicCoordinate lower_dim_max; - if (!iter->get_regular_range(lower_dim_min, lower_dim_max)) - return false; - - if (is_regular_range == regular_to_do) - { - // check if all lower dimensional ranges have same regular range - BasicCoordinate lower_dim_min_try; - BasicCoordinate lower_dim_max_try; - - for (++iter; iter != base_type::end(); ++iter) - { - if (!iter->get_regular_range(lower_dim_min_try, lower_dim_max_try)) - { - is_regular_range = regular_false; - return false; - } - if (!std::equal(lower_dim_min.begin(), lower_dim_min.end(), lower_dim_min_try.begin()) - || !std::equal(lower_dim_max.begin(), lower_dim_max.end(), lower_dim_max_try.begin())) - { - is_regular_range = regular_false; - return false; - } - } - // yes, they do - is_regular_range = regular_true; - } - - min = join(base_type::get_min_index(), lower_dim_min); - max = join(base_type::get_max_index(), lower_dim_max); - - return true; -} - -/*************************************************** - instantiations - ***************************************************/ - -template class IndexRange<2>; -template class IndexRange<3>; -template class IndexRange<4>; -END_NAMESPACE_STIR diff --git a/src/buildblock/ProjDataInMemory.cxx b/src/buildblock/ProjDataInMemory.cxx index 71b0eaaa36..28d8d5ca28 100644 --- a/src/buildblock/ProjDataInMemory.cxx +++ b/src/buildblock/ProjDataInMemory.cxx @@ -82,7 +82,7 @@ namespace detail { template void -copy_data_from_buffer(const Array<1, float>& buffer, Array& array, std::streamoff offset) +copy_data_from_buffer(const Array<1, float, long long>& buffer, Array& array, std::streamoff offset) { #ifdef STIR_OPENMP # pragma omp critical(PROJDATAINMEMORYCOPY) @@ -96,7 +96,7 @@ copy_data_from_buffer(const Array<1, float>& buffer, Array void -copy_data_to_buffer(Array<1, float>& buffer, const Array& array, std::streamoff offset) +copy_data_to_buffer(Array<1, float, long long>& buffer, const Array& array, std::streamoff offset) { #ifdef STIR_OPENMP # pragma omp critical(PROJDATAINMEMORYCOPY) diff --git a/src/buildblock/overlap_interpolate.cxx b/src/buildblock/overlap_interpolate.cxx index 5978f49a04..c0ec69fc35 100644 --- a/src/buildblock/overlap_interpolate.cxx +++ b/src/buildblock/overlap_interpolate.cxx @@ -97,10 +97,10 @@ START_NAMESPACE_STIR */ -template +template void -overlap_interpolate(VectorWithOffset& out_data, - const VectorWithOffset& in_data, +overlap_interpolate(VectorWithOffset& out_data, + const VectorWithOffset& in_data, const float zoom, const float offset, const bool assign_rest_with_zeroes) @@ -124,8 +124,8 @@ overlap_interpolate(VectorWithOffset& out_data, // the position of the left edge of the first 'out' bin: // left_edge = (out_data.get_min_index() - .5)/zoom + offset // x1 -.5 <= left_edge < x1+.5 - int x2 = out_data.get_min_index(); - int x1 = (int)floor((x2 - .5) / zoom + offset + .5); + auto x2 = out_data.get_min_index(); + auto x1 = (indexT)floor((x2 - .5) / zoom + offset + .5); // the next variable holds the difference between the coordinates // of the right edges of the 'in' bin and the 'out' bin, computed @@ -217,9 +217,9 @@ overlap_interpolate(VectorWithOffset& out_data, const float inverse_zoom = 1 / zoom; // current coordinate in out_data - int x2 = out_data.get_min_index(); + auto x2 = out_data.get_min_index(); // current coordinate in in_data - int x1 = (int)floor((x2 - .5) * inverse_zoom + offset + .5); + auto x1 = (indexT)floor((x2 - .5) * inverse_zoom + offset + .5); // the next variable holds the difference between the coordinates // of the right edges of the 'in' bin and the 'out' bin, computed @@ -329,26 +329,7 @@ template void overlap_interpolate<>(VectorWithOffset& out_data, const bool assign_rest_with_zeroes); END_NAMESPACE_STIR -// TODO remove -#if defined(OLDDESIGN) - -# include "stir/Tensor2D.h" - -template void overlap_interpolate<>(VectorWithOffset>& out_data, - const VectorWithOffset>& in_data, - const float zoom, - const float offset, - const bool assign_rest_with_zeroes); - -template void overlap_interpolate<>(VectorWithOffset>& out_data, - const VectorWithOffset>& in_data, - const float zoom, - const float offset, - const bool assign_rest_with_zeroes); - -#else - -# include "stir/Array.h" +#include "stir/Array.h" START_NAMESPACE_STIR @@ -363,6 +344,5 @@ template void overlap_interpolate<>(VectorWithOffset>& out_data, const float zoom, const float offset, const bool assign_rest_with_zeroes); -#endif END_NAMESPACE_STIR diff --git a/src/include/stir/Array.h b/src/include/stir/Array.h index c40fe348da..6a84ec27d8 100644 --- a/src/include/stir/Array.h +++ b/src/include/stir/Array.h @@ -3,7 +3,7 @@ Copyright (C) 2000 PARAPET partners Copyright (C) 2000 - 2011-10-14, Hammersmith Imanet Ltd Copyright (C) 2011-07-01 - 2012, Kris Thielemans - Copyright (C) 2023 - 2025, University College London + Copyright (C) 2023 - 2026, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -73,16 +73,16 @@ the object. However, as grow() is a virtual function, Array::grow is called, which initialises new elements first to 0. */ -template -class Array : public NumericVectorWithOffset, elemT> +template +class Array : public NumericVectorWithOffset, elemT, indexT> { #ifdef STIR_COMPILING_SWIG_WRAPPER // work-around swig problem. It gets confused when using a private (or protected) // typedef in a definition of a public typedef/member public: #endif - typedef Array self; - typedef NumericVectorWithOffset, elemT> base_type; + typedef Array self; + typedef NumericVectorWithOffset, elemT, indexT> base_type; public: //@{ @@ -110,7 +110,7 @@ class Array : public NumericVectorWithOffset, e # ifndef ARRAY_FULL2 //! This defines an iterator type that iterates through all elements. typedef FullArrayIterator::full_iterator, + typename Array::full_iterator, elemT, full_reference, full_pointer> @@ -118,7 +118,7 @@ class Array : public NumericVectorWithOffset, e //! As full_iterator, but for const objects. typedef FullArrayIterator::const_full_iterator, + typename Array::const_full_iterator, elemT, const_full_reference, const_full_pointer> @@ -136,7 +136,7 @@ class Array : public NumericVectorWithOffset, e inline Array(); //! Construct an Array of given range of indices, elements are initialised to 0 - inline explicit Array(const IndexRange&); + inline explicit Array(const IndexRange&); //! Construct an Array pointing to existing contiguous data /*! @@ -148,7 +148,7 @@ class Array : public NumericVectorWithOffset, e The C-array \a data_ptr will be accessed with the last dimension running fastest ("row-major" order). */ - inline Array(const IndexRange& range, shared_ptr data_sptr); + inline Array(const IndexRange& range, shared_ptr data_sptr); #ifndef SWIG // swig 2.0.4 gets confused by base_type (due to numeric template arguments) @@ -200,7 +200,9 @@ class Array : public NumericVectorWithOffset, e inline const_full_iterator end_all_const() const; //@} - inline IndexRange get_index_range() const; + inline IndexRange get_index_range() const; + + inline bool empty() const override; //! return the total number of elements in this array inline size_t size_all() const; @@ -215,10 +217,10 @@ class Array : public NumericVectorWithOffset, e \warning In most cases, calling resize() will result in the array using non-contiguous memory. */ - inline virtual void resize(const IndexRange& range); + inline virtual void resize(const IndexRange& range); //! alias for resize() - virtual inline void grow(const IndexRange& range); + virtual inline void grow(const IndexRange& range); //! return sum of all elements inline elemT sum() const; @@ -254,30 +256,30 @@ class Array : public NumericVectorWithOffset, e //! find regular range, returns \c false if the range is not regular /*! \see class IndexRange for a definition of (ir)regular ranges */ - bool get_regular_range(BasicCoordinate& min, BasicCoordinate& max) const; + bool get_regular_range(BasicCoordinate& min, BasicCoordinate& max) const; //! allow array-style access, read/write - inline Array& operator[](int i); + inline Array& operator[](indexT i); //! array access, read-only - inline const Array& operator[](int i) const; + inline const Array& operator[](indexT i) const; //! allow array-style access given a BasicCoordinate to specify the indices, read/write - inline elemT& operator[](const BasicCoordinate& c); + inline elemT& operator[](const BasicCoordinate& c); //! array access given a BasicCoordinate to specify the indices, read-only // TODO alternative return value: elemT - inline const elemT& operator[](const BasicCoordinate& c) const; + inline const elemT& operator[](const BasicCoordinate& c) const; //! \name indexed access with range checking (throw std:out_of_range) //@{ - inline Array& at(int i); + inline Array& at(indexT i); - inline const Array& at(int i) const; + inline const Array& at(indexT i) const; - inline elemT& at(const BasicCoordinate& c); + inline elemT& at(const BasicCoordinate& c); - inline const elemT& at(const BasicCoordinate& c) const; + inline const elemT& at(const BasicCoordinate& c) const; //@} //! \name Numerical operations @@ -413,7 +415,7 @@ class Array : public NumericVectorWithOffset, e The C-array \a data_ptr will be accessed with the last dimension running fastest ("row-major" order). */ - inline void init_with_copy(const IndexRange& range, elemT const* const data_ptr); + inline void init_with_copy(const IndexRange& range, elemT const* const data_ptr); //! Set the array to a range of indices, and point to/copy from \c data_ptr /*! \arg data_ptr should point to a contiguous block of correct size @@ -424,9 +426,9 @@ class Array : public NumericVectorWithOffset, e \warning This function should only be called from within a constructor. It will ignore any existing content and therefore would cause memory leaks. */ - inline void init(const IndexRange& range, elemT* const data_ptr, bool copy_data); + inline void init(const IndexRange& range, elemT* const data_ptr, bool copy_data); // Make sure that we can access init() recursively - template + template friend class Array; using base_type::grow; @@ -438,16 +440,16 @@ class Array : public NumericVectorWithOffset, e **************************************************/ //! The 1-dimensional (partial) specialisation of Array. -template -class Array<1, elemT> : public NumericVectorWithOffset +template +class Array<1, elemT, indexT> : public NumericVectorWithOffset { #ifdef STIR_COMPILING_SWIG_WRAPPER // work-around swig problem. It gets confused when using a private (or protected) // typedef in a definition of a public typedef/member public: #endif - typedef NumericVectorWithOffset base_type; - typedef Array<1, elemT> self; + typedef NumericVectorWithOffset base_type; + typedef Array<1, elemT, indexT> self; public: //! typedefs such that we do not need to have \c typename wherever we use these types defined in the base class @@ -474,29 +476,29 @@ class Array<1, elemT> : public NumericVectorWithOffset //! default constructor: array of length 0 inline Array(); - //! constructor given an IndexRange<1>, initialising elements to 0 - inline explicit Array(const IndexRange<1>& range); + //! constructor given an IndexRange<1, indexT>, initialising elements to 0 + inline explicit Array(const IndexRange<1, indexT>& range); //! constructor given first and last indices, initialising elements to 0 - inline Array(const int min_index, const int max_index); + inline Array(const indexT min_index, const indexT max_index); - //! constructor given an IndexRange<1>, pointing to existing contiguous data + //! constructor given an IndexRange<1, indexT>, pointing to existing contiguous data /*! \arg data_ptr should point to a contiguous block of correct size. The constructed Array will essentially be a "view" of the \c data_sptr block. Therefore, any modifications to the array will modify the data at \a data_sptr. This will be the case until the Array is resized. */ - inline Array(const IndexRange<1>& range, shared_ptr data_sptr); + inline Array(const IndexRange<1, indexT>& range, shared_ptr data_sptr); - //! constructor given an IndexRange<1> from existing contiguous data (will copy) + //! constructor given an IndexRange<1, indexT> from existing contiguous data (will copy) /*! \arg data_ptr should point to a contiguous block of correct size. */ - inline Array(const IndexRange<1>& range, const elemT* const data_ptr); + inline Array(const IndexRange<1, indexT>& range, const elemT* const data_ptr); //! constructor from basetype - inline Array(const NumericVectorWithOffset& il); + inline Array(const NumericVectorWithOffset& il); //! Copy constructor // implementation needed as the above doesn't replace the normal copy-constructor @@ -537,23 +539,23 @@ class Array<1, elemT> : public NumericVectorWithOffset //@} //! return the range of indices used - inline IndexRange<1> get_index_range() const; + inline IndexRange<1, indexT> get_index_range() const; //! return the total number of elements in this array inline size_t size_all() const; //! Array::grow initialises new elements to 0 - inline virtual void grow(const IndexRange<1>& range); + inline virtual void grow(const IndexRange<1, indexT>& range); // Array::grow initialises new elements to 0 - inline void grow(const int min_index, const int max_index) override; + inline void grow(const indexT min_index, const indexT max_index) override; //! Array::resize initialises new elements to 0 - inline virtual void resize(const IndexRange<1>& range); + inline virtual void resize(const IndexRange<1, indexT>& range); - inline void resize(const int min_index, const int max_index, bool initialise_with_0); + inline void resize(const indexT min_index, const indexT max_index, bool initialise_with_0); //! resize, initialising new elements to 0 - inline void resize(const int min_index, const int max_index) override; + inline void resize(const indexT min_index, const indexT max_index) override; //! \name access to the data via a pointer //@{ @@ -603,7 +605,7 @@ class Array<1, elemT> : public NumericVectorWithOffset inline bool is_regular() const; //! find regular range, returns \c false if the range is not regular - bool get_regular_range(BasicCoordinate<1, int>& min, BasicCoordinate<1, int>& max) const; + bool get_regular_range(BasicCoordinate<1, indexT>& min, BasicCoordinate<1, indexT>& max) const; /* Add numerical operators with correct return value, as opposed to those from the base class */ @@ -632,42 +634,42 @@ class Array<1, elemT> : public NumericVectorWithOffset inline self operator/(const elemT a) const; //! allow array-style access, read/write - inline elemT& operator[](int i); + inline elemT& operator[](indexT i); //! array access, read-only - inline const elemT& operator[](int i) const; + inline const elemT& operator[](indexT i) const; //! allow array-style access giving its BasicCoordinate, read/write - inline const elemT& operator[](const BasicCoordinate<1, int>& c) const; + inline const elemT& operator[](const BasicCoordinate<1, indexT>& c) const; //! array access giving its BasicCoordinate, read-only - inline elemT& operator[](const BasicCoordinate<1, int>& c); + inline elemT& operator[](const BasicCoordinate<1, indexT>& c); //! \name indexed access with range checking (throw std:out_of_range) //@{ - inline elemT& at(int i); + inline elemT& at(indexT i); - inline const elemT& at(int i) const; + inline const elemT& at(indexT i) const; - inline elemT& at(const BasicCoordinate<1, int>& c); + inline elemT& at(const BasicCoordinate<1, indexT>& c); - inline const elemT& at(const BasicCoordinate<1, int>& c) const; + inline const elemT& at(const BasicCoordinate<1, indexT>& c) const; //@} private: // Make sure we can call init() recursively. - template + template friend class Array; //! change vector with new index range and copy data from \c data_ptr /*! \arg data_ptr should start to a contiguous block of correct size */ - inline void init_with_copy(const IndexRange<1>& range, elemT const* const data_ptr); + inline void init_with_copy(const IndexRange<1, indexT>& range, elemT const* const data_ptr); //! change vector with new index range and point to \c data_ptr /*! \arg data_ptr should start to a contiguous block of correct size */ - inline void init(const IndexRange<1>& range, elemT* const data_ptr, bool copy_data); + inline void init(const IndexRange<1, indexT>& range, elemT* const data_ptr, bool copy_data); }; END_NAMESPACE_STIR diff --git a/src/include/stir/Array.inl b/src/include/stir/Array.inl index 0dd5c21566..24505b89b4 100644 --- a/src/include/stir/Array.inl +++ b/src/include/stir/Array.inl @@ -4,7 +4,7 @@ Copyright (C) 2000 PARAPET partners Copyright (C) 2000 - 2011-01-11, Hammersmith Imanet Ltd Copyright (C) 2011-07-01 - 2012, Kris Thielemans - Copyright (C) 2023 - 2025, University College London + Copyright (C) 2023 - 2026, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -37,11 +37,11 @@ START_NAMESPACE_STIR /********************************************** - inlines for Array + inlines for Array **********************************************/ -template +template bool -Array::is_contiguous() const +Array::is_contiguous() const { auto mem = &(*this->begin_all()); for (auto i = this->get_min_index(); i <= this->get_max_index(); ++i) @@ -57,21 +57,26 @@ Array::is_contiguous() const return true; } -template +template void -Array::resize(const IndexRange& range) +Array::resize(const IndexRange& range) { base_type::resize(range.get_min_index(), range.get_max_index()); typename base_type::iterator iter = this->begin(); - typename IndexRange::const_iterator range_iter = range.begin(); + typename IndexRange::const_iterator range_iter = range.begin(); for (; iter != this->end(); ++iter, ++range_iter) (*iter).resize(*range_iter); } -template +template void -Array::init(const IndexRange& range, elemT* const data_ptr, bool copy_data) +Array::init(const IndexRange& range, elemT* const data_ptr, bool copy_data) { + if (range.empty()) + { + this->recycle(); + return; + } base_type::resize(range.get_min_index(), range.get_max_index()); auto iter = this->begin(); auto range_iter = range.begin(); @@ -83,10 +88,15 @@ Array::init(const IndexRange& range, elem } } -template +template void -Array::init_with_copy(const IndexRange& range, elemT const* const data_ptr) +Array::init_with_copy(const IndexRange& range, elemT const* const data_ptr) { + if (range.empty()) + { + this->recycle(); + return; + } base_type::resize(range.get_min_index(), range.get_max_index()); auto iter = this->begin(); auto range_iter = range.begin(); @@ -98,21 +108,21 @@ Array::init_with_copy(const IndexRange& r } } -template +template void -Array::grow(const IndexRange& range) +Array::grow(const IndexRange& range) { resize(range); } -template -Array::Array() +template +Array::Array() : base_type(), _allocated_full_data_ptr(nullptr) {} -template -Array::Array(const IndexRange& range) +template +Array::Array(const IndexRange& range) : base_type(), _allocated_full_data_ptr(new elemT[range.size_all()]) { @@ -125,15 +135,15 @@ Array::Array(const IndexRange& range) this->init(range, this->_allocated_full_data_ptr.get(), false); } -template -Array::Array(const IndexRange& range, shared_ptr data_sptr) +template +Array::Array(const IndexRange& range, shared_ptr data_sptr) { this->_allocated_full_data_ptr = data_sptr; this->init(range, this->_allocated_full_data_ptr.get(), false); } -template -Array::Array(const self& t) +template +Array::Array(const self& t) : base_type(), _allocated_full_data_ptr(new elemT[t.size_all()]) { @@ -144,8 +154,8 @@ Array::Array(const self& t) #ifndef SWIG // swig cannot parse this ATM, but we don't need it anyway in the wrappers -template -Array::Array(const base_type& t) +template +Array::Array(const base_type& t) : base_type(t), _allocated_full_data_ptr(nullptr) { @@ -153,8 +163,8 @@ Array::Array(const base_type& t) } #endif -template -Array::~Array() +template +Array::~Array() { if (this->_allocated_full_data_ptr) { @@ -163,103 +173,97 @@ Array::~Array() } } -template -Array::Array(Array&& other) noexcept +template +Array::Array(Array&& other) noexcept : Array() { swap(*this, other); DEBINFO("move constructor " + std::to_string(num_dimensions) + "copy of size " + std::to_string(this->size_all())); } -template -Array& -Array::operator=(Array other) +template +Array& +Array::operator=(Array other) { swap(*this, other); DEBINFO("Array= " + std::to_string(num_dimensions) + "copy of size " + std::to_string(this->size_all())); return *this; } -template -typename Array::full_iterator -Array::end_all() +template +typename Array::full_iterator +Array::end_all() { // note this value is fixed by the current convention in full_iterator::operator++() return full_iterator(this->end(), this->end(), - typename Array::full_iterator(0), - typename Array::full_iterator(0)); + typename Array::full_iterator(0), + typename Array::full_iterator(0)); } -template -typename Array::const_full_iterator -Array::end_all_const() const +template +typename Array::const_full_iterator +Array::end_all_const() const { return const_full_iterator(this->end(), this->end(), - typename Array::const_full_iterator(0), - typename Array::const_full_iterator(0)); + typename Array::const_full_iterator(0), + typename Array::const_full_iterator(0)); } -template -typename Array::const_full_iterator -Array::end_all() const +template +typename Array::const_full_iterator +Array::end_all() const { return this->end_all_const(); } -template -typename Array::full_iterator -Array::begin_all() +template +typename Array::full_iterator +Array::begin_all() { - if (this->begin() == this->end()) - { - // empty array - return end_all(); - } + if (this->empty()) + return end_all(); else return full_iterator(this->begin(), this->end(), this->begin()->begin_all(), this->begin()->end_all()); } -template -typename Array::const_full_iterator -Array::begin_all_const() const +template +typename Array::const_full_iterator +Array::begin_all_const() const { - if (this->begin() == this->end()) - { - // empty array - return end_all(); - } + if (this->empty()) + return end_all(); else return const_full_iterator(this->begin(), this->end(), this->begin()->begin_all_const(), this->begin()->end_all_const()); } -template -typename Array::const_full_iterator -Array::begin_all() const +template +typename Array::const_full_iterator +Array::begin_all() const { return begin_all_const(); } -template -IndexRange -Array::get_index_range() const +template +IndexRange +Array::get_index_range() const { - VectorWithOffset> range(this->get_min_index(), this->get_max_index()); + VectorWithOffset, indexT> range(this->get_min_index(), this->get_max_index()); - typename VectorWithOffset>::iterator range_iter = range.begin(); - const_iterator array_iter = this->begin(); + auto range_iter = range.begin(); + auto array_iter = this->begin(); for (; range_iter != range.end(); range_iter++, array_iter++) { *range_iter = (*array_iter).get_index_range(); } - return IndexRange(range); + return IndexRange(range); } -template +template size_t -Array::size_all() const +Array::size_all() const { this->check_state(); size_t acc = 0; @@ -268,11 +272,25 @@ Array::size_all() const # pragma omp parallel for reduction(+ : acc) # endif #endif - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) acc += this->num[i].size_all(); return acc; } +template +bool +Array::empty() const +{ + this->check_state(); + if (base_type::empty()) + return true; + // else + for (auto i : *this) + if (i.empty()) + return true; + return false; +} + /*! If is_contiguous() is \c false, calls error(). Otherwise, return a \c elemT* to the first element of the array. @@ -287,9 +305,9 @@ Array::size_all() const get_const_full_data_ptr() and release_const_full_data_ptr(). (This is checked with assert() in DEBUG mode.) */ -template +template elemT* -Array::get_full_data_ptr() +Array::get_full_data_ptr() { this->_full_pointer_access = true; if (!this->is_contiguous()) @@ -306,9 +324,9 @@ Array::get_full_data_ptr() \see get_full_data_ptr() */ -template +template const elemT* -Array::get_const_full_data_ptr() const +Array::get_const_full_data_ptr() const { this->_full_pointer_access = true; if (!this->is_contiguous()) @@ -324,9 +342,9 @@ Array::get_const_full_data_ptr() const \see get_full_data_ptr() */ -template +template void -Array::release_full_data_ptr() +Array::release_full_data_ptr() { assert(this->_full_pointer_access); @@ -340,17 +358,17 @@ Array::release_full_data_ptr() \see get_const_full_data_ptr() */ -template +template void -Array::release_const_full_data_ptr() const +Array::release_const_full_data_ptr() const { assert(this->_full_pointer_access); this->_full_pointer_access = false; } -template +template elemT -Array::sum() const +Array::sum() const { this->check_state(); typename HigherPrecision::type acc; @@ -360,14 +378,14 @@ Array::sum() const # pragma omp parallel for reduction(+ : acc) # endif #endif - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) acc += this->num[i].sum(); return static_cast(acc); } -template +template elemT -Array::sum_positive() const +Array::sum_positive() const { this->check_state(); typename HigherPrecision::type acc; @@ -377,14 +395,14 @@ Array::sum_positive() const # pragma omp parallel for reduction(+ : acc) # endif #endif - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) acc += this->num[i].sum_positive(); return static_cast(acc); } -template +template elemT -Array::find_max() const +Array::find_max() const { this->check_state(); if (this->size() > 0) @@ -395,7 +413,7 @@ Array::find_max() const # pragma omp parallel for reduction(max : maxval) # endif #endif - for (int i = this->get_min_index() + 1; i <= this->get_max_index(); i++) + for (auto i = this->get_min_index() + 1; i <= this->get_max_index(); i++) { maxval = std::max(this->num[i].find_max(), maxval); } @@ -408,9 +426,9 @@ Array::find_max() const } } -template +template elemT -Array::find_min() const +Array::find_min() const { this->check_state(); if (this->size() > 0) @@ -421,7 +439,7 @@ Array::find_min() const # pragma omp parallel for reduction(min : minval) # endif #endif - for (int i = this->get_min_index() + 1; i <= this->get_max_index(); i++) + for (auto i = this->get_min_index() + 1; i <= this->get_max_index(); i++) { minval = std::min(this->num[i].find_min(), minval); } @@ -434,116 +452,116 @@ Array::find_min() const } } -template +template void -Array::fill(const elemT& n) +Array::fill(const elemT& n) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i].fill(n); this->check_state(); } -template +template void -Array::apply_lower_threshold(const elemT& l) +Array::apply_lower_threshold(const elemT& l) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i].apply_lower_threshold(l); this->check_state(); } -template +template void -Array::apply_upper_threshold(const elemT& u) +Array::apply_upper_threshold(const elemT& u) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i].apply_upper_threshold(u); this->check_state(); } -template +template bool -Array::is_regular() const +Array::is_regular() const { return get_index_range().is_regular(); } // TODO terribly inefficient at the moment -template +template bool -Array::get_regular_range(BasicCoordinate& min, - BasicCoordinate& max) const +Array::get_regular_range(BasicCoordinate& min, + BasicCoordinate& max) const { - const IndexRange range = get_index_range(); + const auto range = get_index_range(); return range.get_regular_range(min, max); } -template -Array& -Array::operator[](int i) +template +Array& +Array::operator[](indexT i) { return base_type::operator[](i); } -template -const Array& -Array::operator[](int i) const +template +const Array& +Array::operator[](indexT i) const { return base_type::operator[](i); } -template +template elemT& -Array::operator[](const BasicCoordinate& c) +Array::operator[](const BasicCoordinate& c) { return (*this)[c[1]][cut_first_dimension(c)]; } -template +template const elemT& -Array::operator[](const BasicCoordinate& c) const +Array::operator[](const BasicCoordinate& c) const { return (*this)[c[1]][cut_first_dimension(c)]; } -template -Array& -Array::at(int i) +template +Array& +Array::at(indexT i) { return base_type::at(i); } -template -const Array& -Array::at(int i) const +template +const Array& +Array::at(indexT i) const { return base_type::at(i); } -template +template elemT& -Array::at(const BasicCoordinate& c) +Array::at(const BasicCoordinate& c) { return (*this).at(c[1]).at(cut_first_dimension(c)); } -template +template const elemT& -Array::at(const BasicCoordinate& c) const +Array::at(const BasicCoordinate& c) const { return (*this).at(c[1]).at(cut_first_dimension(c)); } -template +template template void -Array::axpby(const elemT2 a, const Array& x, const elemT2 b, const Array& y) +Array::axpby(const elemT2 a, const Array& x, const elemT2 b, const Array& y) { - Array::xapyb(x, a, y, b); + Array::xapyb(x, a, y, b); } -template +template void -Array::xapyb(const Array& x, const elemT a, const Array& y, const elemT b) +Array::xapyb(const Array& x, const elemT a, const Array& y, const elemT b) { this->check_state(); if ((this->get_index_range() != x.get_index_range()) || (this->get_index_range() != y.get_index_range())) @@ -558,9 +576,9 @@ Array::xapyb(const Array& x, const elemT a, const Array& } } -template +template void -Array::xapyb(const Array& x, const Array& a, const Array& y, const Array& b) +Array::xapyb(const Array& x, const Array& a, const Array& y, const Array& b) { this->check_state(); if ((this->get_index_range() != x.get_index_range()) || (this->get_index_range() != y.get_index_range()) @@ -579,42 +597,42 @@ Array::xapyb(const Array& x, const Array& a, const Array& } } -template +template template void -Array::sapyb(const T& a, const Array& y, const T& b) +Array::sapyb(const T& a, const Array& y, const T& b) { this->xapyb(*this, a, y, b); } /********************************************** - inlines for Array<1, elemT> + inlines for Array<1, elemT, indexT> **********************************************/ -template +template void -Array<1, elemT>::init_with_copy(const IndexRange<1>& range, elemT const* const data_ptr) +Array<1, elemT, indexT>::init_with_copy(const IndexRange<1, indexT>& range, elemT const* const data_ptr) { base_type::init_with_copy(range.get_min_index(), range.get_max_index(), data_ptr); } -template +template void -Array<1, elemT>::init(const IndexRange<1>& range, elemT* const data_ptr, bool copy_data) +Array<1, elemT, indexT>::init(const IndexRange<1, indexT>& range, elemT* const data_ptr, bool copy_data) { base_type::init(range.get_min_index(), range.get_max_index(), data_ptr, copy_data); } -template +template void -Array<1, elemT>::resize(const int min_index, const int max_index, bool initialise_with_0) +Array<1, elemT, indexT>::resize(const indexT min_index, const indexT max_index, bool initialise_with_0) { this->check_state(); - const int oldstart = this->get_min_index(); + const indexT oldstart = this->get_min_index(); const size_type oldlength = this->size(); base_type::resize(min_index, max_index); - if (!initialise_with_0) + if (!initialise_with_0 || this->size() == 0) { this->check_state(); return; @@ -622,165 +640,165 @@ Array<1, elemT>::resize(const int min_index, const int max_index, bool initialis if (oldlength == 0) { - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) - assign(this->num[i], 0); + for (auto& i : *this) + assign(i, 0); } else { - for (int i = this->get_min_index(); i < oldstart && i <= this->get_max_index(); ++i) + for (auto i = this->get_min_index(); i < oldstart && i <= this->get_max_index(); ++i) assign(this->num[i], 0); - for (int i = std::max(static_cast(oldstart + oldlength), this->get_min_index()); i <= this->get_max_index(); ++i) + for (auto i = std::max(static_cast(oldstart + oldlength), this->get_min_index()); i <= this->get_max_index(); ++i) assign(this->num[i], 0); } this->check_state(); } -template +template void -Array<1, elemT>::resize(const int min_index, const int max_index) +Array<1, elemT, indexT>::resize(const indexT min_index, const indexT max_index) { resize(min_index, max_index, true); } -template +template void -Array<1, elemT>::resize(const IndexRange<1>& range) +Array<1, elemT, indexT>::resize(const IndexRange<1, indexT>& range) { resize(range.get_min_index(), range.get_max_index()); } -template +template void -Array<1, elemT>::grow(const int min_index, const int max_index) +Array<1, elemT, indexT>::grow(const indexT min_index, const indexT max_index) { resize(min_index, max_index); } -template +template void -Array<1, elemT>::grow(const IndexRange<1>& range) +Array<1, elemT, indexT>::grow(const IndexRange<1, indexT>& range) { grow(range.get_min_index(), range.get_max_index()); } -template -Array<1, elemT>::Array() +template +Array<1, elemT, indexT>::Array() : base_type() {} -template -Array<1, elemT>::Array(const IndexRange<1>& range) +template +Array<1, elemT, indexT>::Array(const IndexRange<1, indexT>& range) : base_type() { grow(range); } -template -Array<1, elemT>::Array(const int min_index, const int max_index) +template +Array<1, elemT, indexT>::Array(const indexT min_index, const indexT max_index) : base_type() { grow(min_index, max_index); } -template -Array<1, elemT>::Array(const IndexRange<1>& range, shared_ptr data_sptr) +template +Array<1, elemT, indexT>::Array(const IndexRange<1, indexT>& range, shared_ptr data_sptr) : base_type(range.get_min_index(), range.get_max_index(), data_sptr) {} -template -Array<1, elemT>::Array(const IndexRange<1>& range, const elemT* const data_ptr) +template +Array<1, elemT, indexT>::Array(const IndexRange<1, indexT>& range, const elemT* const data_ptr) : base_type(range.get_min_index(), range.get_max_index(), data_ptr) {} -template -Array<1, elemT>::Array(const base_type& il) +template +Array<1, elemT, indexT>::Array(const base_type& il) : base_type(il) {} -template -Array<1, elemT>::Array(const Array<1, elemT>& other) +template +Array<1, elemT, indexT>::Array(const Array<1, elemT, indexT>& other) : base_type(other) {} -template -Array<1, elemT>::~Array() +template +Array<1, elemT, indexT>::~Array() {} -template -Array<1, elemT>::Array(Array<1, elemT>&& other) noexcept +template +Array<1, elemT, indexT>::Array(Array<1, elemT, indexT>&& other) noexcept : Array() { swap(*this, other); } -template -Array<1, elemT>& -Array<1, elemT>::operator=(const Array<1, elemT>& other) +template +Array<1, elemT, indexT>& +Array<1, elemT, indexT>::operator=(const Array<1, elemT, indexT>& other) { // use the base_type assignment, as this tries to avoid reallocating memory base_type::operator=(other); return *this; } -template -typename Array<1, elemT>::full_iterator -Array<1, elemT>::begin_all() +template +typename Array<1, elemT, indexT>::full_iterator +Array<1, elemT, indexT>::begin_all() { return this->begin(); } -template -typename Array<1, elemT>::const_full_iterator -Array<1, elemT>::begin_all() const +template +typename Array<1, elemT, indexT>::const_full_iterator +Array<1, elemT, indexT>::begin_all() const { return this->begin(); } -template -typename Array<1, elemT>::full_iterator -Array<1, elemT>::end_all() +template +typename Array<1, elemT, indexT>::full_iterator +Array<1, elemT, indexT>::end_all() { return this->end(); } -template -typename Array<1, elemT>::const_full_iterator -Array<1, elemT>::end_all() const +template +typename Array<1, elemT, indexT>::const_full_iterator +Array<1, elemT, indexT>::end_all() const { return this->end(); } -template -typename Array<1, elemT>::const_full_iterator -Array<1, elemT>::begin_all_const() const +template +typename Array<1, elemT, indexT>::const_full_iterator +Array<1, elemT, indexT>::begin_all_const() const { return this->begin(); } -template -typename Array<1, elemT>::const_full_iterator -Array<1, elemT>::end_all_const() const +template +typename Array<1, elemT, indexT>::const_full_iterator +Array<1, elemT, indexT>::end_all_const() const { return this->end(); } -template -IndexRange<1> -Array<1, elemT>::get_index_range() const +template +IndexRange<1, indexT> +Array<1, elemT, indexT>::get_index_range() const { - return IndexRange<1>(this->get_min_index(), this->get_max_index()); + return IndexRange<1, indexT>(this->get_min_index(), this->get_max_index()); } -template +template size_t -Array<1, elemT>::size_all() const +Array<1, elemT, indexT>::size_all() const { return this->size(); } -template +template elemT -Array<1, elemT>::sum() const +Array<1, elemT, indexT>::sum() const { this->check_state(); typename HigherPrecision::type acc; @@ -790,14 +808,14 @@ Array<1, elemT>::sum() const # pragma omp parallel for reduction(+ : acc) # endif #endif - for (int i = this->get_min_index(); i <= this->get_max_index(); ++i) + for (auto i = this->get_min_index(); i <= this->get_max_index(); ++i) acc += this->num[i]; return static_cast(acc); }; -template +template elemT -Array<1, elemT>::sum_positive() const +Array<1, elemT, indexT>::sum_positive() const { this->check_state(); typename HigherPrecision::type acc; @@ -807,7 +825,7 @@ Array<1, elemT>::sum_positive() const # pragma omp parallel for reduction(+ : acc) # endif #endif - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) { if (this->num[i] > 0) acc += this->num[i]; @@ -815,9 +833,9 @@ Array<1, elemT>::sum_positive() const return static_cast(acc); }; -template +template elemT -Array<1, elemT>::find_max() const +Array<1, elemT, indexT>::find_max() const { this->check_state(); if (this->size() > 0) @@ -832,9 +850,9 @@ Array<1, elemT>::find_max() const this->check_state(); }; -template +template elemT -Array<1, elemT>::find_min() const +Array<1, elemT, indexT>::find_min() const { this->check_state(); if (this->size() > 0) @@ -849,18 +867,18 @@ Array<1, elemT>::find_min() const this->check_state(); }; -template +template bool -Array<1, elemT>::is_regular() const +Array<1, elemT, indexT>::is_regular() const { return true; } -template +template bool -Array<1, elemT>::get_regular_range(BasicCoordinate<1, int>& min, BasicCoordinate<1, int>& max) const +Array<1, elemT, indexT>::get_regular_range(BasicCoordinate<1, indexT>& min, BasicCoordinate<1, indexT>& max) const { - const IndexRange<1> range = get_index_range(); + const IndexRange<1, indexT> range = get_index_range(); return range.get_regular_range(min, max); } @@ -873,129 +891,129 @@ its base_type (which happens if these function are not repeated in this class). Complicated... */ -template -Array<1, elemT> -Array<1, elemT>::operator+(const base_type& iv) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator+(const base_type& iv) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return retval += iv; }; -template -Array<1, elemT> -Array<1, elemT>::operator-(const base_type& iv) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator-(const base_type& iv) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return retval -= iv; } -template -Array<1, elemT> -Array<1, elemT>::operator*(const base_type& iv) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator*(const base_type& iv) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return retval *= iv; } -template -Array<1, elemT> -Array<1, elemT>::operator/(const base_type& iv) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator/(const base_type& iv) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return retval /= iv; } -template -Array<1, elemT> -Array<1, elemT>::operator+(const elemT a) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator+(const elemT a) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return (retval += a); }; -template -Array<1, elemT> -Array<1, elemT>::operator-(const elemT a) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator-(const elemT a) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return (retval -= a); }; -template -Array<1, elemT> -Array<1, elemT>::operator*(const elemT a) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator*(const elemT a) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return (retval *= a); }; -template -Array<1, elemT> -Array<1, elemT>::operator/(const elemT a) const +template +Array<1, elemT, indexT> +Array<1, elemT, indexT>::operator/(const elemT a) const { this->check_state(); - Array<1, elemT> retval(*this); + Array<1, elemT, indexT> retval(*this); return (retval /= a); }; -template +template const elemT& -Array<1, elemT>::operator[](int i) const +Array<1, elemT, indexT>::operator[](indexT i) const { return base_type::operator[](i); }; -template +template elemT& -Array<1, elemT>::operator[](int i) +Array<1, elemT, indexT>::operator[](indexT i) { return base_type::operator[](i); }; -template +template const elemT& -Array<1, elemT>::operator[](const BasicCoordinate<1, int>& c) const +Array<1, elemT, indexT>::operator[](const BasicCoordinate<1, indexT>& c) const { return (*this)[c[1]]; }; -template +template elemT& -Array<1, elemT>::operator[](const BasicCoordinate<1, int>& c) +Array<1, elemT, indexT>::operator[](const BasicCoordinate<1, indexT>& c) { return (*this)[c[1]]; }; -template +template const elemT& -Array<1, elemT>::at(int i) const +Array<1, elemT, indexT>::at(indexT i) const { return base_type::at(i); }; -template +template elemT& -Array<1, elemT>::at(int i) +Array<1, elemT, indexT>::at(indexT i) { return base_type::at(i); }; -template +template const elemT& -Array<1, elemT>::at(const BasicCoordinate<1, int>& c) const +Array<1, elemT, indexT>::at(const BasicCoordinate<1, indexT>& c) const { return (*this).at(c[1]); }; -template +template elemT& -Array<1, elemT>::at(const BasicCoordinate<1, int>& c) +Array<1, elemT, indexT>::at(const BasicCoordinate<1, indexT>& c) { return (*this).at(c[1]); }; diff --git a/src/include/stir/ArrayFilter1DUsingConvolution.h b/src/include/stir/ArrayFilter1DUsingConvolution.h index 32a4d7849a..a326c93fcf 100644 --- a/src/include/stir/ArrayFilter1DUsingConvolution.h +++ b/src/include/stir/ArrayFilter1DUsingConvolution.h @@ -27,9 +27,6 @@ START_NAMESPACE_STIR -template -class VectorWithOffset; - /*! \ingroup Array \brief This class implements convolution of a 1D array with an diff --git a/src/include/stir/ArrayFilter1DUsingConvolutionSymmetricKernel.h b/src/include/stir/ArrayFilter1DUsingConvolutionSymmetricKernel.h index fd490c1f2f..1dab910a7a 100644 --- a/src/include/stir/ArrayFilter1DUsingConvolutionSymmetricKernel.h +++ b/src/include/stir/ArrayFilter1DUsingConvolutionSymmetricKernel.h @@ -26,9 +26,6 @@ START_NAMESPACE_STIR -template -class VectorWithOffset; - /*! \ingroup Array \brief This class implements convolution of a 1D array with a symmetric kernel. diff --git a/src/include/stir/ArrayFilter2DUsingConvolution.h b/src/include/stir/ArrayFilter2DUsingConvolution.h index 4c8b72c2b0..88c2f1ccb2 100644 --- a/src/include/stir/ArrayFilter2DUsingConvolution.h +++ b/src/include/stir/ArrayFilter2DUsingConvolution.h @@ -19,13 +19,9 @@ #define __stir_ArrayFilter2DUsingConvolution_H__ #include "stir/ArrayFunctionObject_2ArgumentImplementation.h" -//#include "stir/VectorWithOffset.h" START_NAMESPACE_STIR -template -class VectorWithOffset; - template class ArrayFilter2DUsingConvolution : public ArrayFunctionObject_2ArgumentImplementation<2, elemT> { diff --git a/src/include/stir/ArrayFilter3DUsingConvolution.h b/src/include/stir/ArrayFilter3DUsingConvolution.h index a32d1bc1d1..fbeb30e054 100644 --- a/src/include/stir/ArrayFilter3DUsingConvolution.h +++ b/src/include/stir/ArrayFilter3DUsingConvolution.h @@ -23,9 +23,6 @@ START_NAMESPACE_STIR -template -class VectorWithOffset; - template class ArrayFilter3DUsingConvolution : public ArrayFunctionObject_2ArgumentImplementation<3, elemT> { diff --git a/src/include/stir/ArrayFunction.h b/src/include/stir/ArrayFunction.h index 2d7b20f172..d0308fd1e6 100644 --- a/src/include/stir/ArrayFunction.h +++ b/src/include/stir/ArrayFunction.h @@ -68,37 +68,37 @@ START_NAMESPACE_STIR //! Replace elements by their logarithm, 1D version /*! \ingroup Array */ -template -inline Array<1, elemT>& in_place_log(Array<1, elemT>& v); +template +inline Array<1, elemT, indexT>& in_place_log(Array<1, elemT, indexT>& v); //! apply log to each element of the multi-dimensional array /*! \ingroup Array */ -template -inline Array& in_place_log(Array& v); +template +inline Array& in_place_log(Array& v); //! Replace elements by their exponentiation, 1D version /*! \ingroup Array */ -template -inline Array<1, elemT>& in_place_exp(Array<1, elemT>& v); +template +inline Array<1, elemT, indexT>& in_place_exp(Array<1, elemT, indexT>& v); //! apply exp to each element of the multi-dimensional array /*! \ingroup Array */ -template -inline Array& in_place_exp(Array& v); +template +inline Array& in_place_exp(Array& v); //! Replace elements by their absolute value, 1D version /*! \ingroup Array \warning The implementation does not work with complex numbers. */ -template -inline Array<1, elemT>& in_place_abs(Array<1, elemT>& v); +template +inline Array<1, elemT, indexT>& in_place_abs(Array<1, elemT, indexT>& v); //! store absolute value of each element of the multi-dimensional array /*! \ingroup Array \warning The implementation does not work with complex numbers. */ -template -inline Array& in_place_abs(Array& v); +template +inline Array& in_place_abs(Array& v); //! apply any function(object) to each element of the multi-dimensional array /*! \ingroup Array @@ -141,8 +141,8 @@ inline T& in_place_apply_function(T& v, FUNCTION f); \todo Add a specialisation such that this function would handle function objects and (smart) pointers to function objects. At the moment, it's only the latter. */ -template -inline void in_place_apply_array_function_on_1st_index(Array& array, FunctionObjectPtr f); +template +inline void in_place_apply_array_function_on_1st_index(Array& array, FunctionObjectPtr f); //! apply any function(object) to each element of the multi-dimensional array, storing results in a different array /*! \ingroup Array @@ -153,9 +153,10 @@ inline void in_place_apply_array_function_on_1st_index(Array& ar \endcode \a f should not modify the index range of the output argument. */ -template -inline void -apply_array_function_on_1st_index(Array& out_array, const Array& in_array, FunctionObjectPtr f); +template +inline void apply_array_function_on_1st_index(Array& out_array, + const Array& in_array, + FunctionObjectPtr f); /* local #define used for 2 purposes: - in partial template specialisation that uses ArrayFunctionObject types @@ -171,7 +172,7 @@ apply_array_function_on_1st_index(Array& out_array, const Array< Ideally, the code should be rewritten to work with any kind of (smart) ptr. TODO */ #if !defined(__GNUC__) && !defined(_MSC_VER) -# define ActualFunctionObjectPtrIter VectorWithOffset>>::const_iterator +# define ActualFunctionObjectPtrIter VectorWithOffset>>::const_iterator #else /* Puzzlingly, although the code is actually called with iterators of the type above, @@ -180,7 +181,7 @@ apply_array_function_on_1st_index(Array& out_array, const Array< VC also refuses to compile it. A work-around is to use the following type */ -# define ActualFunctionObjectPtrIter shared_ptr> const* +# define ActualFunctionObjectPtrIter shared_ptr> const* #endif //! Apply a sequence of 1d array-function objects on every dimension of the input array @@ -198,15 +199,16 @@ apply_array_function_on_1st_index(Array& out_array, const Array< objects and (smart) pointers to function objects. At the moment, it's only the latter. */ // TODO add specialisation that uses ArrayFunctionObject::is_trivial -template -inline void in_place_apply_array_functions_on_each_index(Array& array, +template +inline void in_place_apply_array_functions_on_each_index(Array& array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop); //! 1d specialisation of the above. -template -inline void -in_place_apply_array_functions_on_each_index(Array<1, elemT>& array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop); +template +inline void in_place_apply_array_functions_on_each_index(Array<1, elemT, indexT>& array, + FunctionObjectPtrIter start, + FunctionObjectPtrIter stop); //! Apply a sequence of 1d array-function objects on every dimension of the input array, store in output array /*! \ingroup Array @@ -223,9 +225,9 @@ in_place_apply_array_functions_on_each_index(Array<1, elemT>& array, FunctionObj \todo Add a specialisation such that this function would handle iterators of function objects and (smart) pointers to function objects. At the moment, it's only the latter. */ -template -inline void apply_array_functions_on_each_index(Array& out_array, - const Array& in_array, +template +inline void apply_array_functions_on_each_index(Array& out_array, + const Array& in_array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop); @@ -237,9 +239,9 @@ inline void apply_array_functions_on_each_index(Array& out_array \todo Modify such that this function would handle function objects and (smart) pointers to ArrayFunctionObject objects. At the moment, it's only the latter. */ -template -inline void apply_array_functions_on_each_index(Array& out_array, - const Array& in_array, +template +inline void apply_array_functions_on_each_index(Array& out_array, + const Array& in_array, ActualFunctionObjectPtrIter start, ActualFunctionObjectPtrIter stop); @@ -247,25 +249,27 @@ inline void apply_array_functions_on_each_index(Array& out_array /*! \ingroup Array */ // has to be here to get general 1D specialisation to compile -template -inline void apply_array_functions_on_each_index(Array<1, elemT>& out_array, - const Array<1, elemT>& in_array, +template +inline void apply_array_functions_on_each_index(Array<1, elemT, indexT>& out_array, + const Array<1, elemT, indexT>& in_array, ActualFunctionObjectPtrIter start, ActualFunctionObjectPtrIter stop); -template //! 1d specialisation for general function objects /*! \ingroup Array */ -inline void apply_array_functions_on_each_index(Array<1, elemT>& out_array, - const Array<1, elemT>& in_array, +template +inline void apply_array_functions_on_each_index(Array<1, elemT, indexT>& out_array, + const Array<1, elemT, indexT>& in_array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop); -template -inline void transform_array_to_periodic_indices(Array& out_array, const Array& in_array); -template -inline void transform_array_from_periodic_indices(Array& out_array, const Array& in_array); +template +inline void transform_array_to_periodic_indices(Array& out_array, + const Array& in_array); +template +inline void transform_array_from_periodic_indices(Array& out_array, + const Array& in_array); END_NAMESPACE_STIR diff --git a/src/include/stir/ArrayFunction.inl b/src/include/stir/ArrayFunction.inl index 612b7d8c03..c7452a4d8d 100644 --- a/src/include/stir/ArrayFunction.inl +++ b/src/include/stir/ArrayFunction.inl @@ -27,13 +27,6 @@ #include #include -#ifdef BOOST_NO_STDC_NAMESPACE -namespace std -{ -using ::log; -using ::exp; -} // namespace std -#endif START_NAMESPACE_STIR @@ -41,57 +34,57 @@ START_NAMESPACE_STIR // element wise and in place numeric functions //---------------------------------------------------------------------- -template -inline Array<1, elemT>& -in_place_log(Array<1, elemT>& v) +template +inline Array<1, elemT, indexT>& +in_place_log(Array<1, elemT, indexT>& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) v[i] = std::log(v[i]); return v; } -template -inline Array& -in_place_log(Array& v) +template +inline Array& +in_place_log(Array& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) in_place_log(v[i]); return v; } -template -inline Array<1, elemT>& -in_place_exp(Array<1, elemT>& v) +template +inline Array<1, elemT, indexT>& +in_place_exp(Array<1, elemT, indexT>& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) v[i] = std::exp(v[i]); return v; } -template -inline Array& -in_place_exp(Array& v) +template +inline Array& +in_place_exp(Array& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) in_place_exp(v[i]); return v; } -template -inline Array<1, elemT>& -in_place_abs(Array<1, elemT>& v) +template +inline Array<1, elemT, indexT>& +in_place_abs(Array<1, elemT, indexT>& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) if (v[i] < 0) v[i] = -v[i]; return v; } -template -inline Array& -in_place_abs(Array& v) +template +inline Array& +in_place_abs(Array& v) { - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) in_place_abs(v[i]); return v; } @@ -100,8 +93,8 @@ template inline T& in_place_apply_function(T& v, FUNCTION f) { - typename T::full_iterator iter = v.begin_all(); - const typename T::full_iterator end_iter = v.end_all(); + auto iter = v.begin_all(); + const auto end_iter = v.end_all(); while (iter != end_iter) { *iter = f(*iter); @@ -110,31 +103,31 @@ in_place_apply_function(T& v, FUNCTION f) return v; } -template +template inline void -in_place_apply_array_function_on_1st_index(Array& array, FunctionObjectPtr f) +in_place_apply_array_function_on_1st_index(Array& array, FunctionObjectPtr f) { assert(array.is_regular()); - const int outer_min_index = array.get_min_index(); - const int outer_max_index = array.get_max_index(); + const auto outer_min_index = array.get_min_index(); + const auto outer_max_index = array.get_max_index(); // construct a vector with a full_iterator for every array[i] VectorWithOffset< #ifndef _MSC_VER typename #endif - Array::full_iterator> + Array::full_iterator> full_iterators(outer_min_index, outer_max_index); - for (int i = outer_min_index; i <= outer_max_index; ++i) + for (auto i = outer_min_index; i <= outer_max_index; ++i) full_iterators[i] = array[i].begin_all(); // allocate 1d array - Array<1, elemT> array1d(outer_min_index, outer_max_index); + Array<1, elemT, indexT> array1d(outer_min_index, outer_max_index); while (full_iterators[outer_min_index] != array[outer_min_index].end_all()) { // copy elements into 1d array - for (int i = outer_min_index; i <= outer_max_index; ++i) + for (auto i = outer_min_index; i <= outer_max_index; ++i) array1d[i] = *full_iterators[i]; // apply function @@ -142,41 +135,43 @@ in_place_apply_array_function_on_1st_index(Array& array, Functio // put results back // and increment full_iterators to do next index - for (int i = outer_min_index; i <= outer_max_index; ++i) + for (auto i = outer_min_index; i <= outer_max_index; ++i) *full_iterators[i]++ = array1d[i]; } } -template +template inline void -apply_array_function_on_1st_index(Array& out_array, const Array& in_array, FunctionObjectPtr f) +apply_array_function_on_1st_index(Array& out_array, + const Array& in_array, + FunctionObjectPtr f) { assert(in_array.is_regular()); assert(out_array.is_regular()); - const int in_min_index = in_array.get_min_index(); - const int in_max_index = in_array.get_max_index(); - const int out_min_index = out_array.get_min_index(); - const int out_max_index = out_array.get_max_index(); + const auto in_min_index = in_array.get_min_index(); + const auto in_max_index = in_array.get_max_index(); + const auto out_min_index = out_array.get_min_index(); + const auto out_max_index = out_array.get_max_index(); // construct a vector with a full_iterator for every in_array[i] - VectorWithOffset::const_full_iterator> in_full_iterators(in_min_index, in_max_index); - for (int i = in_min_index; i <= in_max_index; ++i) + VectorWithOffset::const_full_iterator> in_full_iterators(in_min_index, in_max_index); + for (auto i = in_min_index; i <= in_max_index; ++i) in_full_iterators[i] = in_array[i].begin_all(); // same for out_array[i] - VectorWithOffset::full_iterator> out_full_iterators(out_min_index, out_max_index); - for (int i = out_min_index; i <= out_max_index; ++i) + VectorWithOffset::full_iterator> out_full_iterators(out_min_index, out_max_index); + for (auto i = out_min_index; i <= out_max_index; ++i) out_full_iterators[i] = out_array[i].begin_all(); // allocate 1d array - Array<1, elemT> in_array1d(in_min_index, in_max_index); - Array<1, elemT> out_array1d(out_min_index, out_max_index); + Array<1, elemT, indexT> in_array1d(in_min_index, in_max_index); + Array<1, elemT, indexT> out_array1d(out_min_index, out_max_index); while (in_full_iterators[in_min_index] != in_array[in_min_index].end_all()) { assert(out_full_iterators[out_min_index] != out_array[out_min_index].end_all()); // copy elements into 1d array // increment in_full_iterators for next index - for (int i = in_min_index; i <= in_max_index; ++i) + for (auto i = in_min_index; i <= in_max_index; ++i) in_array1d[i] = *(in_full_iterators[i]++); // apply function @@ -186,15 +181,15 @@ apply_array_function_on_1st_index(Array& out_array, const Array< // put results back // increment out_full_iterators for next index - for (int i = out_min_index; i <= out_max_index; ++i) + for (auto i = out_min_index; i <= out_max_index; ++i) *(out_full_iterators[i]++) = out_array1d[i]; } assert(out_full_iterators[out_min_index] == out_array[out_min_index].end_all()); } -template +template inline void -in_place_apply_array_functions_on_each_index(Array& array, +in_place_apply_array_functions_on_each_index(Array& array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop) { @@ -203,22 +198,24 @@ in_place_apply_array_functions_on_each_index(Array& array, in_place_apply_array_function_on_1st_index(array, *start); ++start; - for (typename Array::iterator restiter = array.begin(); restiter != array.end(); ++restiter) + for (auto restiter = array.begin(); restiter != array.end(); ++restiter) in_place_apply_array_functions_on_each_index(*restiter, start, stop); } -template +template inline void -in_place_apply_array_functions_on_each_index(Array<1, elemT>& array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop) +in_place_apply_array_functions_on_each_index(Array<1, elemT, indexT>& array, + FunctionObjectPtrIter start, + FunctionObjectPtrIter stop) { assert(start + 1 == stop); (**start)(array); } -template +template inline void -apply_array_functions_on_each_index(Array& out_array, - const Array& in_array, +apply_array_functions_on_each_index(Array& out_array, + const Array& in_array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop) { @@ -226,23 +223,23 @@ apply_array_functions_on_each_index(Array& out_array, assert(num_dim > 1); assert(out_array.is_regular()); - BasicCoordinate tmp_out_min_indices, tmp_out_max_indices; + BasicCoordinate tmp_out_min_indices, tmp_out_max_indices; out_array.get_regular_range(tmp_out_min_indices, tmp_out_max_indices); tmp_out_min_indices[1] = in_array.get_min_index(); tmp_out_max_indices[1] = in_array.get_max_index(); - Array tmp_out_array(IndexRange(tmp_out_min_indices, tmp_out_max_indices)); + Array tmp_out_array(IndexRange(tmp_out_min_indices, tmp_out_max_indices)); - for (int i = in_array.get_min_index(); i <= in_array.get_max_index(); ++i) + for (auto i = in_array.get_min_index(); i <= in_array.get_max_index(); ++i) apply_array_functions_on_each_index(tmp_out_array[i], in_array[i], start + 1, stop); apply_array_function_on_1st_index(out_array, tmp_out_array, *start); } // specialisation that uses ArrayFunctionObject::is_trivial etc -template +template inline void -apply_array_functions_on_each_index(Array& out_array, - const Array& in_array, +apply_array_functions_on_each_index(Array& out_array, + const Array& in_array, ActualFunctionObjectPtrIter start, ActualFunctionObjectPtrIter stop) { @@ -252,7 +249,7 @@ apply_array_functions_on_each_index(Array& out_array, // cerr << "apply_array_functions_on_each_index dim " << num_dim << std::endl; if ((**start).is_trivial()) { - for (int i = max(out_array.get_min_index(), in_array.get_min_index()); + for (auto i = max(out_array.get_min_index(), in_array.get_min_index()); i <= min(out_array.get_max_index(), in_array.get_max_index()); ++i) apply_array_functions_on_each_index(out_array[i], in_array[i], start + 1, stop); @@ -261,23 +258,23 @@ apply_array_functions_on_each_index(Array& out_array, { assert(out_array.is_regular()); - IndexRange<1> influencing_indices; + IndexRange<1, indexT> influencing_indices; if ((**start).get_influencing_indices(influencing_indices, - IndexRange<1>(out_array.get_min_index(), out_array.get_max_index())) + IndexRange<1, indexT>(out_array.get_min_index(), out_array.get_max_index())) == Succeeded::no) - influencing_indices = IndexRange<1>(influencing_indices.get_min_index(), in_array.get_max_index()); + influencing_indices = IndexRange<1, indexT>(influencing_indices.get_min_index(), in_array.get_max_index()); else { - influencing_indices = IndexRange<1>(max(influencing_indices.get_min_index(), in_array.get_min_index()), - min(influencing_indices.get_max_index(), in_array.get_max_index())); + influencing_indices = IndexRange<1, indexT>(max(influencing_indices.get_min_index(), in_array.get_min_index()), + min(influencing_indices.get_max_index(), in_array.get_max_index())); } - BasicCoordinate tmp_out_min_indices, tmp_out_max_indices; + BasicCoordinate tmp_out_min_indices, tmp_out_max_indices; out_array.get_regular_range(tmp_out_min_indices, tmp_out_max_indices); tmp_out_min_indices[1] = influencing_indices.get_min_index(); tmp_out_max_indices[1] = influencing_indices.get_max_index(); - Array tmp_out_array(IndexRange(tmp_out_min_indices, tmp_out_max_indices)); + Array tmp_out_array(IndexRange(tmp_out_min_indices, tmp_out_max_indices)); - for (int i = influencing_indices.get_min_index(); i <= influencing_indices.get_max_index(); ++i) + for (auto i = influencing_indices.get_min_index(); i <= influencing_indices.get_max_index(); ++i) apply_array_functions_on_each_index(tmp_out_array[i], in_array[i], start + 1, stop); apply_array_function_on_1st_index(out_array, tmp_out_array, *start); @@ -285,10 +282,10 @@ apply_array_functions_on_each_index(Array& out_array, } // has to be here to get general 1D specialisation to compile -template +template inline void -apply_array_functions_on_each_index(Array<1, elemT>& out_array, - const Array<1, elemT>& in_array, +apply_array_functions_on_each_index(Array<1, elemT, indexT>& out_array, + const Array<1, elemT, indexT>& in_array, ActualFunctionObjectPtrIter start, ActualFunctionObjectPtrIter stop) { @@ -296,10 +293,10 @@ apply_array_functions_on_each_index(Array<1, elemT>& out_array, (**start)(out_array, in_array); } -template +template inline void -apply_array_functions_on_each_index(Array<1, elemT>& out_array, - const Array<1, elemT>& in_array, +apply_array_functions_on_each_index(Array<1, elemT, indexT>& out_array, + const Array<1, elemT, indexT>& in_array, FunctionObjectPtrIter start, FunctionObjectPtrIter stop) { @@ -309,19 +306,20 @@ apply_array_functions_on_each_index(Array<1, elemT>& out_array, /******************* functions that copy arrays using transformed indices *****/ -template +template inline void -transform_array_to_periodic_indices(Array& out_array, const Array& in_array) +transform_array_to_periodic_indices(Array& out_array, + const Array& in_array) { assert(norm(get_min_indices(out_array)) < .01); // check out_array is 0-based - BasicCoordinate min_indices, max_indices; + BasicCoordinate min_indices, max_indices; assert(out_array.get_regular_range(min_indices, max_indices)); out_array.get_regular_range(min_indices, max_indices); - const BasicCoordinate out_sizes = max_indices - min_indices + 1; + const auto out_sizes = max_indices - min_indices + 1; { - BasicCoordinate index = get_min_indices(in_array); + BasicCoordinate index = get_min_indices(in_array); do { out_array[modulo(index, out_sizes)] = in_array[index]; @@ -329,19 +327,20 @@ transform_array_to_periodic_indices(Array& out_array, con } } -template +template inline void -transform_array_from_periodic_indices(Array& out_array, const Array& in_array) +transform_array_from_periodic_indices(Array& out_array, + const Array& in_array) { assert(norm(get_min_indices(in_array)) < .01); // check in_array is 0-based - BasicCoordinate min_indices, max_indices; + BasicCoordinate min_indices, max_indices; assert(in_array.get_regular_range(min_indices, max_indices)); in_array.get_regular_range(min_indices, max_indices); - const BasicCoordinate in_sizes = max_indices - min_indices + 1; + const BasicCoordinate in_sizes = max_indices - min_indices + 1; - BasicCoordinate index = get_min_indices(out_array); + auto index = get_min_indices(out_array); do { out_array[index] = in_array[modulo(index, in_sizes)]; diff --git a/src/include/stir/ArrayFunctionObject.h b/src/include/stir/ArrayFunctionObject.h index 4cc6694629..55b19bb635 100644 --- a/src/include/stir/ArrayFunctionObject.h +++ b/src/include/stir/ArrayFunctionObject.h @@ -27,13 +27,11 @@ START_NAMESPACE_STIR -template -class IndexRange; /*! \ingroup Array \brief A class for operations on n-dimensional Arrays */ -template +template class ArrayFunctionObject { public: @@ -42,12 +40,13 @@ class ArrayFunctionObject /*! \warning Not all derived classes will be able to handle arbitrary index ranges for \a in_array. */ - virtual void operator()(Array& array) const = 0; + virtual void operator()(Array& array) const = 0; //! result stored in another array /*! \warning Not all derived classes will be able to handle arbitrary index ranges in \a out_array and \a in_array. */ - virtual void operator()(Array& out_array, const Array& in_array) const = 0; + virtual void operator()(Array& out_array, + const Array& in_array) const = 0; //! Should return true when the operations won't modify the object at all /*! For the 2 argument version, elements in \a out_array will be set to corresponding elements in \a in_array. Elements in \a out_array that do not @@ -63,8 +62,8 @@ class ArrayFunctionObject going to affect the \a output_indices (independent of the size of the input array) or of it is too difficult for the derived class to return a sensible index range. */ - virtual Succeeded get_influencing_indices(IndexRange& influencing_indices, - const IndexRange& output_indices) const + virtual Succeeded get_influencing_indices(IndexRange& influencing_indices, + const IndexRange& output_indices) const { return Succeeded::no; } @@ -76,8 +75,8 @@ class ArrayFunctionObject going to be affected by the \a input_indices (independent of the size of the output array) or of it is too difficult for the derived class to return a sensible index range. */ - virtual Succeeded get_influenced_indices(IndexRange& influenced_indices, - const IndexRange& input_indices) const + virtual Succeeded get_influenced_indices(IndexRange& influenced_indices, + const IndexRange& input_indices) const { return Succeeded::no; } diff --git a/src/include/stir/ArrayFwd.h b/src/include/stir/ArrayFwd.h index c3a09a5359..39cebb76d2 100644 --- a/src/include/stir/ArrayFwd.h +++ b/src/include/stir/ArrayFwd.h @@ -8,19 +8,26 @@ See STIR/LICENSE.txt for details */ +#ifndef ___stir_ArrayFwd_H___ +#define ___stir_ArrayFwd_H___ /*! \file \ingroup Array \brief forward declaration of stir::Array class for multi-dimensional (numeric) arrays */ -#include "stir/common.h" +#include "stir/VectorWithOffsetFwd.h" +#include "stir/IndexRangeFwd.h" +#include "stir/BasicCoordinateFwd.h" namespace stir { -template +template class Array; //! type alias for future-proofing for "large" rectangular arrays -template -using ArrayType = Array; +template +using ArrayType = Array; + } // namespace stir + +#endif diff --git a/src/include/stir/BasicCoordinate.h b/src/include/stir/BasicCoordinate.h index 80e2cf4513..38fec89b25 100644 --- a/src/include/stir/BasicCoordinate.h +++ b/src/include/stir/BasicCoordinate.h @@ -30,7 +30,7 @@ */ -#include "stir/common.h" +#include "stir/BasicCoordinateFwd.h" #include #include diff --git a/src/include/stir/BasicCoordinateFwd.h b/src/include/stir/BasicCoordinateFwd.h new file mode 100644 index 0000000000..b2c56e2c89 --- /dev/null +++ b/src/include/stir/BasicCoordinateFwd.h @@ -0,0 +1,26 @@ + +/* + Copyright (C) 2025, University College London + This file is part of STIR. + + SPDX-License-Identifier: Apache-2.0 + + See STIR/LICENSE.txt for details +*/ + +#ifndef ___stir_BasicCoordinateFwd_H___ +#define ___stir_BasicCoordinateFwd_H___ +/*! + \file + \ingroup Array + \brief forward declaration of stir::BasicCoordinate class +*/ +#include "stir/common.h" + +namespace stir +{ +template +class BasicCoordinate; +} // namespace stir + +#endif diff --git a/src/include/stir/IO/interfile.h b/src/include/stir/IO/interfile.h index 37108781b9..a235191f20 100644 --- a/src/include/stir/IO/interfile.h +++ b/src/include/stir/IO/interfile.h @@ -37,13 +37,9 @@ START_NAMESPACE_STIR -template -class IndexRange; template class DiscretisedDensity; template -class VectorWithOffset; -template class CartesianCoordinate3D; template class Coordinate3D; @@ -53,8 +49,6 @@ class ProjDataFromStream; class DynamicDiscretisedDensity; template class ParametricDiscretisedDensity; -template -class VoxelsOnCartesianGrid; template class KineticParameters; diff --git a/src/include/stir/IO/read_data.h b/src/include/stir/IO/read_data.h index 9af4015a82..9c32922ab3 100644 --- a/src/include/stir/IO/read_data.h +++ b/src/include/stir/IO/read_data.h @@ -40,8 +40,9 @@ class NumericInfo; \warning When an error occurs, the function immediately returns. However, the data might have been partially read from \a s. */ -template -inline Succeeded read_data(IStreamT& s, ArrayType& data, const ByteOrder byte_order = ByteOrder::native); +template +inline Succeeded +read_data(IStreamT& s, ArrayType& data, const ByteOrder byte_order = ByteOrder::native); /*! \ingroup Array_IO \brief Read the data of an Array from file as a different type. @@ -55,9 +56,9 @@ inline Succeeded read_data(IStreamT& s, ArrayType& data, \see find_scale_factor() for the meaning of \a scale_factor. */ -template +template inline Succeeded read_data(IStreamT& s, - ArrayType& data, + ArrayType& data, NumericInfo input_type, ScaleT& scale_factor, const ByteOrder byte_order = ByteOrder::native); @@ -72,9 +73,9 @@ inline Succeeded read_data(IStreamT& s, const bool) The only difference is that the input type is now specified using NumericType. */ -template +template inline Succeeded read_data(IStreamT& s, - ArrayType& data, + ArrayType& data, NumericType type, ScaleT& scale, const ByteOrder byte_order = ByteOrder::native); diff --git a/src/include/stir/IO/read_data.inl b/src/include/stir/IO/read_data.inl index 4daca13e09..ac1c7b3032 100644 --- a/src/include/stir/IO/read_data.inl +++ b/src/include/stir/IO/read_data.inl @@ -1,6 +1,6 @@ /* Copyright (C) 2004- 2007, Hammersmith Imanet Ltd - Copyright (C) 2024, University College London + Copyright (C) 2024, 2025, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 @@ -31,15 +31,15 @@ START_NAMESPACE_STIR namespace detail { /* Generic implementation of read_data(). See test_if_1d.h for info why we do this.*/ -template +template inline Succeeded -read_data_help(is_not_1d, IStreamT& s, ArrayType& data, const ByteOrder byte_order) +read_data_help(is_not_1d, IStreamT& s, ArrayType& data, const ByteOrder byte_order) { if (data.is_contiguous()) return read_data_1d(s, data, byte_order); // otherwise, recurse - for (typename ArrayType::iterator iter = data.begin(); iter != data.end(); ++iter) + for (auto iter = data.begin(); iter != data.end(); ++iter) { if (read_data(s, *iter, byte_order) == Succeeded::no) return Succeeded::no; @@ -49,26 +49,26 @@ read_data_help(is_not_1d, IStreamT& s, ArrayType& data, c /* 1D implementation of read_data(). See test_if_1d.h for info why we do this.*/ // specialisation for 1D case -template +template inline Succeeded -read_data_help(is_1d, IStreamT& s, ArrayType<1, elemT>& data, const ByteOrder byte_order) +read_data_help(is_1d, IStreamT& s, ArrayType<1, elemT, indexT>& data, const ByteOrder byte_order) { return read_data_1d(s, data, byte_order); } } // end of namespace detail -template +template inline Succeeded -read_data(IStreamT& s, ArrayType& data, const ByteOrder byte_order) +read_data(IStreamT& s, ArrayType& data, const ByteOrder byte_order) { return detail::read_data_help(detail::test_if_1d(), s, data, byte_order); } -template +template inline Succeeded read_data(IStreamT& s, - ArrayType& data, + ArrayType& data, NumericInfo input_type, ScaleT& scale_factor, const ByteOrder byte_order) @@ -82,7 +82,7 @@ read_data(IStreamT& s, } else { - ArrayType in_data(data.get_index_range()); + ArrayType in_data(data.get_index_range()); Succeeded success = read_data(s, in_data, byte_order); if (success == Succeeded::no) return Succeeded::no; @@ -91,9 +91,10 @@ read_data(IStreamT& s, } } -template +template inline Succeeded -read_data(IStreamT& s, ArrayType& data, NumericType type, ScaleT& scale, const ByteOrder byte_order) +read_data( + IStreamT& s, ArrayType& data, NumericType type, ScaleT& scale, const ByteOrder byte_order) { switch (type.id) { diff --git a/src/include/stir/IO/read_data_1d.h b/src/include/stir/IO/read_data_1d.h index 1df5de7018..d4e5172522 100644 --- a/src/include/stir/IO/read_data_1d.h +++ b/src/include/stir/IO/read_data_1d.h @@ -17,7 +17,7 @@ See STIR/LICENSE.txt for details */ #include "stir/ArrayFwd.h" -#include +#include #include START_NAMESPACE_STIR @@ -32,15 +32,15 @@ namespace detail This function might propagate any exceptions by std::istream::read. */ -template -inline Succeeded read_data_1d(std::istream& s, ArrayType& data, const ByteOrder byte_order); +template +inline Succeeded read_data_1d(std::istream& s, ArrayType& data, const ByteOrder byte_order); /* \ingroup Array_IO_detail \brief This is the (internal) function that does the actual reading from a FILE*. \internal */ -template -inline Succeeded read_data_1d(FILE*&, ArrayType& data, const ByteOrder byte_order); +template +inline Succeeded read_data_1d(FILE*&, ArrayType& data, const ByteOrder byte_order); } // end namespace detail END_NAMESPACE_STIR diff --git a/src/include/stir/IO/read_data_1d.inl b/src/include/stir/IO/read_data_1d.inl index ed8c4e38e7..8221bd1780 100644 --- a/src/include/stir/IO/read_data_1d.inl +++ b/src/include/stir/IO/read_data_1d.inl @@ -28,9 +28,9 @@ namespace detail /***************** version for istream *******************************/ -template +template Succeeded -read_data_1d(std::istream& s, ArrayType& data, const ByteOrder byte_order) +read_data_1d(std::istream& s, ArrayType& data, const ByteOrder byte_order) { if (!s || (dynamic_cast(&s) != 0 && !dynamic_cast(&s)->is_open()) || (dynamic_cast(&s) != 0 && !dynamic_cast(&s)->is_open())) @@ -64,9 +64,9 @@ read_data_1d(std::istream& s, ArrayType& data, const Byte /***************** version for FILE *******************************/ // largely a copy of above, but with calls to stdio function -template +template Succeeded -read_data_1d(FILE*& fptr_ref, ArrayType& data, const ByteOrder byte_order) +read_data_1d(FILE*& fptr_ref, ArrayType& data, const ByteOrder byte_order) { FILE* fptr = fptr_ref; if (fptr == NULL || ferror(fptr)) diff --git a/src/include/stir/IO/write_data.h b/src/include/stir/IO/write_data.h index f6f2d6b3c5..664a72f2e4 100644 --- a/src/include/stir/IO/write_data.h +++ b/src/include/stir/IO/write_data.h @@ -48,9 +48,9 @@ class NumericInfo; is then not templated in \c num_dimensions, but explicitly defined for a few dimensions (see the source). */ -template +template inline Succeeded write_data(OStreamT& s, - const ArrayType& data, + const ArrayType& data, const ByteOrder byte_order = ByteOrder::native, const bool can_corrupt_data = false); /*! \ingroup Array_IO @@ -69,9 +69,9 @@ inline Succeeded write_data(OStreamT& s, is then not templated in \c num_dimensions, but explicitly defined for a few dimensions (see the source). */ -template +template inline Succeeded write_data(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericInfo output_type, ScaleT& scale_factor, const ByteOrder byte_order = ByteOrder::native, @@ -94,9 +94,9 @@ inline Succeeded write_data(OStreamT& s, is then not templated in \c num_dimensions, but explicitly defined for a few dimensions (see the source). */ -template +template inline Succeeded write_data_with_fixed_scale_factor(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericInfo output_type, const ScaleT scale_factor, const ByteOrder byte_order = ByteOrder::native, @@ -118,9 +118,9 @@ inline Succeeded write_data_with_fixed_scale_factor(OStreamT& s, defined for a few dimensions (see the source). */ -template +template inline Succeeded write_data(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericType type, ScaleT& scale, const ByteOrder byte_order = ByteOrder::native, diff --git a/src/include/stir/IO/write_data.inl b/src/include/stir/IO/write_data.inl index c4d71688af..a0305c0094 100644 --- a/src/include/stir/IO/write_data.inl +++ b/src/include/stir/IO/write_data.inl @@ -32,11 +32,11 @@ namespace detail /* Generic implementation of write_data_with_fixed_scale_factor(). See test_if_1d.h for info why we do this this way. */ -template +template inline Succeeded write_data_with_fixed_scale_factor_help(is_not_1d, OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericInfo output_type, const ScaleT scale_factor, const ByteOrder byte_order, @@ -51,11 +51,11 @@ write_data_with_fixed_scale_factor_help(is_not_1d, } // specialisation for 1D case -template +template inline Succeeded write_data_with_fixed_scale_factor_help(is_1d, OStreamT& s, - const ArrayType<1, elemT>& data, + const ArrayType<1, elemT, indexT>& data, NumericInfo, const ScaleT scale_factor, const ByteOrder byte_order, @@ -77,10 +77,10 @@ write_data_with_fixed_scale_factor_help(is_1d, } // end of namespace detail -template +template Succeeded write_data_with_fixed_scale_factor(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericInfo output_type, const ScaleT scale_factor, const ByteOrder byte_order, @@ -90,10 +90,10 @@ write_data_with_fixed_scale_factor(OStreamT& s, detail::test_if_1d(), s, data, output_type, scale_factor, byte_order, can_corrupt_data); } -template +template Succeeded write_data(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericInfo output_type, ScaleT& scale_factor, const ByteOrder byte_order, @@ -103,17 +103,20 @@ write_data(OStreamT& s, return write_data_with_fixed_scale_factor(s, data, output_type, scale_factor, byte_order, can_corrupt_data); } -template +template inline Succeeded -write_data(OStreamT& s, const ArrayType& data, const ByteOrder byte_order, const bool can_corrupt_data) +write_data(OStreamT& s, + const ArrayType& data, + const ByteOrder byte_order, + const bool can_corrupt_data) { return write_data_with_fixed_scale_factor(s, data, NumericInfo(), 1.F, byte_order, can_corrupt_data); } -template +template Succeeded write_data(OStreamT& s, - const ArrayType& data, + const ArrayType& data, NumericType type, ScaleT& scale, const ByteOrder byte_order, diff --git a/src/include/stir/IO/write_data_1d.h b/src/include/stir/IO/write_data_1d.h index a174f1445d..c7f58234e7 100644 --- a/src/include/stir/IO/write_data_1d.h +++ b/src/include/stir/IO/write_data_1d.h @@ -18,7 +18,7 @@ See STIR/LICENSE.txt for details */ #include "stir/ArrayFwd.h" -#include +#include #include START_NAMESPACE_STIR @@ -34,9 +34,11 @@ namespace detail This function does not throw any exceptions. Exceptions thrown by std::ostream::write are caught. */ -template -inline Succeeded -write_data_1d(std::ostream& s, const Array& data, const ByteOrder byte_order, const bool can_corrupt_data); +template +inline Succeeded write_data_1d(std::ostream& s, + const Array& data, + const ByteOrder byte_order, + const bool can_corrupt_data); /*! \ingroup Array_IO_detail \brief This is an internal function called by \c write_data(). It does the actual writing to \c FILE* using stdio functions. @@ -44,9 +46,9 @@ write_data_1d(std::ostream& s, const Array& data, const B This function does not throw any exceptions. */ -template +template inline Succeeded write_data_1d(FILE*& fptr_ref, - const ArrayType& data, + const ArrayType& data, const ByteOrder byte_order, const bool can_corrupt_data); } // namespace detail diff --git a/src/include/stir/IO/write_data_1d.inl b/src/include/stir/IO/write_data_1d.inl index b1a0a9c291..b506551b2a 100644 --- a/src/include/stir/IO/write_data_1d.inl +++ b/src/include/stir/IO/write_data_1d.inl @@ -28,9 +28,12 @@ namespace detail /***************** version for ostream *******************************/ -template +template inline Succeeded -write_data_1d(std::ostream& s, const Array& data, const ByteOrder byte_order, const bool can_corrupt_data) +write_data_1d(std::ostream& s, + const Array& data, + const ByteOrder byte_order, + const bool can_corrupt_data) { if (!s || (dynamic_cast(&s) != 0 && !dynamic_cast(&s)->is_open()) || (dynamic_cast(&s) != 0 && !dynamic_cast(&s)->is_open())) @@ -45,7 +48,7 @@ write_data_1d(std::ostream& s, const Array& data, const B /* if (!byte_order.is_native_order()) { - Array a_copy(data); + Array a_copy(data); for(int i=data.get_min_index(); i<=data.get_max_index(); i++) ByteOrder::swap_order(a_copy[i]); return write_data(s, a_copy, ByteOrder::native, true); @@ -53,7 +56,7 @@ write_data_1d(std::ostream& s, const Array& data, const B */ if (!byte_order.is_native_order()) { - Array& data_ref = const_cast&>(data); + auto& data_ref = const_cast&>(data); for (auto iter = data_ref.begin_all(); iter != data_ref.end_all(); ++iter) ByteOrder::swap_order(*iter); } @@ -76,7 +79,7 @@ write_data_1d(std::ostream& s, const Array& data, const B if (!can_corrupt_data && !byte_order.is_native_order()) { - Array& data_ref = const_cast&>(data); + auto& data_ref = const_cast&>(data); for (auto iter = data_ref.begin_all(); iter != data_ref.end_all(); ++iter) ByteOrder::swap_order(*iter); } @@ -93,9 +96,12 @@ write_data_1d(std::ostream& s, const Array& data, const B /***************** version for FILE *******************************/ // largely a copy of above, but with calls to stdio function -template +template inline Succeeded -write_data_1d(FILE*& fptr_ref, const Array& data, const ByteOrder byte_order, const bool can_corrupt_data) +write_data_1d(FILE*& fptr_ref, + const Array& data, + const ByteOrder byte_order, + const bool can_corrupt_data) { FILE* fptr = fptr_ref; if (fptr == 0 || ferror(fptr)) @@ -110,7 +116,7 @@ write_data_1d(FILE*& fptr_ref, const Array& data, const B /* if (!byte_order.is_native_order()) { - Array a_copy(data); + Array a_copy(data); for(int i=data.get_min_index(); i<=data.get_max_index(); i++) ByteOrder::swap_order(a_copy[i]); return write_data(s, a_copy, ByteOrder::native, true); @@ -118,7 +124,7 @@ write_data_1d(FILE*& fptr_ref, const Array& data, const B */ if (!byte_order.is_native_order()) { - Array& data_ref = const_cast&>(data); + auto& data_ref = const_cast&>(data); for (auto iter = data_ref.begin_all(); iter != data_ref.end_all(); ++iter) ByteOrder::swap_order(*iter); } @@ -134,7 +140,7 @@ write_data_1d(FILE*& fptr_ref, const Array& data, const B if (!can_corrupt_data && !byte_order.is_native_order()) { - Array& data_ref = const_cast&>(data); + auto& data_ref = const_cast&>(data); for (auto iter = data_ref.begin_all(); iter != data_ref.end_all(); ++iter) ByteOrder::swap_order(*iter); } diff --git a/src/include/stir/IndexRange.h b/src/include/stir/IndexRange.h index 882cfb33cc..9c2213666f 100644 --- a/src/include/stir/IndexRange.h +++ b/src/include/stir/IndexRange.h @@ -26,6 +26,7 @@ */ #include "stir/VectorWithOffset.h" #include "stir/BasicCoordinate.h" +#include "stir/IndexRangeFwd.h" START_NAMESPACE_STIR @@ -64,11 +65,11 @@ START_NAMESPACE_STIR } \endcode */ -template -class IndexRange : public VectorWithOffset> +template +class IndexRange : public VectorWithOffset, indexT> { protected: - typedef VectorWithOffset> base_type; + typedef VectorWithOffset, indexT> base_type; public: //! typedefs such that we do not need to have \a typename wherever we use iterators @@ -78,17 +79,20 @@ class IndexRange : public VectorWithOffset> //! Empty range inline IndexRange(); +#ifndef SWIG + // SWIG bug prevents using base_type here. Leads to problems with num_dimensions-1 //! Make an IndexRange from the base type - inline IndexRange(const VectorWithOffset>& range); + inline IndexRange(const base_type& range); +#endif //! Copy constructor - inline IndexRange(const IndexRange& range); + inline IndexRange(const IndexRange& range); //! Construct a regular range given by all minimum indices and all maximum indices. - inline IndexRange(const BasicCoordinate& min, const BasicCoordinate& max); + inline IndexRange(const BasicCoordinate& min, const BasicCoordinate& max); //! Construct a regular range given by sizes (minimum indices will be 0) - inline IndexRange(const BasicCoordinate& sizes); + inline IndexRange(const BasicCoordinate& sizes); //! return the total number of elements in this range inline size_t size_all() const; @@ -96,22 +100,23 @@ class IndexRange : public VectorWithOffset> // these are derived from VectorWithOffset // TODO these should be overloaded, to set regular_range as well. /* - const IndexRange& operator[](int i) const + const IndexRange& operator[](indexT i) const { return range[i]; } - IndexRange& operator[](int i) + IndexRange& operator[](indexT i) { return range[i]; } */ //! comparison operator - inline bool operator==(const IndexRange&) const; - inline bool operator!=(const IndexRange&) const; + inline bool operator==(const IndexRange&) const; + inline bool operator!=(const IndexRange&) const; //! checks if the range is 'regular' inline bool is_regular() const; + inline bool empty() const override; //! find regular range, returns false if the range is not regular - bool get_regular_range(BasicCoordinate& min, BasicCoordinate& max) const; + inline bool get_regular_range(BasicCoordinate& min, BasicCoordinate& max) const; private: //! enum to encode the current knowledge about regularity @@ -127,37 +132,41 @@ class IndexRange : public VectorWithOffset> }; //! The (simple) 1 dimensional specialisation of IndexRange -template <> -class IndexRange<1> +template +class IndexRange<1, indexT> { public: + typedef size_t size_type; inline IndexRange(); - inline IndexRange(const int min, const int max); + inline IndexRange(const indexT min, const indexT max); - inline IndexRange(const BasicCoordinate<1, int>& min, const BasicCoordinate<1, int>& max); + inline IndexRange(const BasicCoordinate<1, indexT>& min, const BasicCoordinate<1, indexT>& max); - inline IndexRange(const int length); - inline IndexRange(const BasicCoordinate<1, int>& size); + inline IndexRange(const size_type length); + inline IndexRange(const BasicCoordinate<1, indexT>& size); - inline int get_min_index() const; - inline int get_max_index() const; - inline int get_length() const; + inline indexT get_min_index() const; + inline indexT get_max_index() const; + inline size_type get_length() const; + inline size_type size() const { return get_length(); } //! return the total number of elements in this range inline size_t size_all() const; + inline bool empty() const; - inline bool operator==(const IndexRange<1>& range2) const; + inline bool operator==(const IndexRange<1, indexT>& range2) const; //! checks if the range is 'regular' (always true for the 1d case) inline bool is_regular() const; //! fills in min and max, and returns true - inline bool get_regular_range(BasicCoordinate<1, int>& min, BasicCoordinate<1, int>& max) const; + inline bool get_regular_range(BasicCoordinate<1, indexT>& min, BasicCoordinate<1, indexT>& max) const; //! resets to new index range - inline void resize(const int min_index, const int max_index); + inline void resize(const indexT min_index, const indexT max_index); private: - int min; - int max; + indexT min; + indexT max; + inline void recycle(); }; END_NAMESPACE_STIR diff --git a/src/include/stir/IndexRange.inl b/src/include/stir/IndexRange.inl index 8e02a65a26..0595d56c1c 100644 --- a/src/include/stir/IndexRange.inl +++ b/src/include/stir/IndexRange.inl @@ -3,6 +3,7 @@ /* Copyright (C) 2000 PARAPET partners Copyright (C) 2000- 2005, Hammersmith Imanet Ltd + Copyright (C) 2025-2026, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -13,7 +14,7 @@ /*! \file \ingroup Array - \brief inline definitions for the IndexRange class + \brief inline definitions for the stir::IndexRange class \author Kris Thielemans \author PARAPET project @@ -22,6 +23,7 @@ */ #include +#include START_NAMESPACE_STIR @@ -29,75 +31,146 @@ START_NAMESPACE_STIR n-D version ***************************************/ -template -IndexRange::IndexRange() +template +IndexRange::IndexRange() : base_type(), is_regular_range(regular_true) {} -template -IndexRange::IndexRange(const IndexRange& range) +template +IndexRange::IndexRange(const IndexRange& range) : base_type(range), is_regular_range(range.is_regular_range) {} -template -IndexRange::IndexRange(const base_type& range) +template +IndexRange::IndexRange(const base_type& range) : base_type(range), is_regular_range(regular_to_do) {} -template -IndexRange::IndexRange(const BasicCoordinate& min_v, - const BasicCoordinate& max_v) +template +IndexRange::IndexRange(const BasicCoordinate& min_v, + const BasicCoordinate& max_v) : base_type(min_v[1], max_v[1]), is_regular_range(regular_true) { - const IndexRange lower_dims(cut_first_dimension(min_v), cut_first_dimension(max_v)); + const IndexRange lower_dims(cut_first_dimension(min_v), cut_first_dimension(max_v)); this->fill(lower_dims); } -template -IndexRange::IndexRange(const BasicCoordinate& sizes) +template +IndexRange::IndexRange(const BasicCoordinate& sizes) : base_type(sizes[1]), is_regular_range(regular_true) { - const IndexRange lower_dims(cut_first_dimension(sizes)); + const IndexRange lower_dims(cut_first_dimension(sizes)); this->fill(lower_dims); } -template +template +bool +IndexRange::empty() const +{ + this->check_state(); + if (base_type::empty()) + return true; + if (this->is_regular_range == regular_true) + return this->begin()->empty(); + // else + for (auto i : *this) + if (i.empty()) + return true; + return false; +} +template std::size_t -IndexRange::size_all() const +IndexRange::size_all() const { this->check_state(); - if (this->is_regular_range == regular_true && this->get_length() > 0) + if (this->empty()) + return std::size_t(0); + if (this->is_regular_range == regular_true) return this->get_length() * this->begin()->size_all(); // else size_t acc = 0; - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (indexT i = this->get_min_index(); i <= this->get_max_index(); i++) acc += this->num[i].size_all(); return acc; } -template +template bool -IndexRange::operator==(const IndexRange& range2) const +IndexRange::operator==(const IndexRange& range2) const { return this->get_min_index() == range2.get_min_index() && this->get_length() == range2.get_length() && std::equal(this->begin(), this->end(), range2.begin()); } -template +template bool -IndexRange::operator!=(const IndexRange& range2) const +IndexRange::operator!=(const IndexRange& range2) const { return !(*this == range2); } -template +template +bool +IndexRange::get_regular_range(BasicCoordinate& min, + BasicCoordinate& max) const +{ + if (base_type::empty()) + { + // use empty range (suitable for unsigned indexT) + std::fill(min.begin(), min.end(), 1); + std::fill(max.begin(), max.end(), 0); + return true; + } + + // if not a regular range, exit + if (is_regular_range == regular_false) + return false; + + typename base_type::const_iterator iter = base_type::begin(); + + BasicCoordinate lower_dim_min; + BasicCoordinate lower_dim_max; + if (!iter->get_regular_range(lower_dim_min, lower_dim_max)) + return false; + + if (is_regular_range == regular_to_do) + { + // check if all lower dimensional ranges have same regular range + BasicCoordinate lower_dim_min_try; + BasicCoordinate lower_dim_max_try; + + for (++iter; iter != base_type::end(); ++iter) + { + if (!iter->get_regular_range(lower_dim_min_try, lower_dim_max_try)) + { + is_regular_range = regular_false; + return false; + } + if (!std::equal(lower_dim_min.begin(), lower_dim_min.end(), lower_dim_min_try.begin()) + || !std::equal(lower_dim_max.begin(), lower_dim_max.end(), lower_dim_max_try.begin())) + { + is_regular_range = regular_false; + return false; + } + } + // yes, they do + is_regular_range = regular_true; + } + + min = join(base_type::get_min_index(), lower_dim_min); + max = join(base_type::get_max_index(), lower_dim_max); + + return true; +} + +template bool -IndexRange::is_regular() const +IndexRange::is_regular() const { switch (is_regular_range) { @@ -106,8 +179,8 @@ IndexRange::is_regular() const case regular_false: return false; case regular_to_do: { - BasicCoordinate min; - BasicCoordinate max; + BasicCoordinate min; + BasicCoordinate max; return get_regular_range(min, max); } } @@ -121,73 +194,103 @@ IndexRange::is_regular() const 1D version ***************************************/ -IndexRange<1>::IndexRange() - : min(0), - max(0) -{} +template +bool +IndexRange<1, indexT>::empty() const +{ + return max < min; +} -IndexRange<1>::IndexRange(const int min_v, const int max_v) +template +void +IndexRange<1, indexT>::recycle() +{ + // note: max - min + 1 needs to be 0 (as we don't do a check in size()) + // note: use (1,0) as opposed to (0,-1) to avoid problem with unsigned indexT + this->min = indexT(1); + this->max = indexT(0); +} + +template +IndexRange<1, indexT>::IndexRange() +{ + this->recycle(); +} + +template +IndexRange<1, indexT>::IndexRange(const indexT min_v, const indexT max_v) : min(min_v), max(max_v) -{} +{ + if (max_v < min_v) + this->recycle(); +} -IndexRange<1>::IndexRange(const BasicCoordinate<1, int>& min_v, const BasicCoordinate<1, int>& max_v) - : min(min_v[1]), - max(max_v[1]) +template +IndexRange<1, indexT>::IndexRange(const BasicCoordinate<1, indexT>& min_v, const BasicCoordinate<1, indexT>& max_v) + : IndexRange<1, indexT>(min_v[1], max_v[1]) {} -IndexRange<1>::IndexRange(const int length) - : min(0), - max(length - 1) +template +IndexRange<1, indexT>::IndexRange(const size_type sz) + : min((std::is_signed_v || sz > 0) ? indexT(0) : indexT(1)), + max((std::is_signed_v || sz > 0) ? indexT(sz - 1) : indexT(0)) {} -IndexRange<1>::IndexRange(const BasicCoordinate<1, int>& size) - : min(0), - max(size[1] - 1) +template +IndexRange<1, indexT>::IndexRange(const BasicCoordinate<1, indexT>& size) + : IndexRange<1, indexT>(size[1]) {} -int -IndexRange<1>::get_min_index() const +template +indexT +IndexRange<1, indexT>::get_min_index() const { return min; } -int -IndexRange<1>::get_max_index() const +template +indexT +IndexRange<1, indexT>::get_max_index() const { return max; } -int -IndexRange<1>::get_length() const +template +typename IndexRange<1, indexT>::size_type +IndexRange<1, indexT>::get_length() const { - return max - min + 1; + return static_cast((max + 1) - min); // note: this order of calculation avoids problems with unsigned indexT } +template std::size_t -IndexRange<1>::size_all() const +IndexRange<1, indexT>::size_all() const { return std::size_t(this->get_length()); } +template bool -IndexRange<1>::operator==(const IndexRange<1>& range2) const +IndexRange<1, indexT>::operator==(const IndexRange<1, indexT>& range2) const { return get_min_index() == range2.get_min_index() && get_length() == range2.get_length(); } +template bool -IndexRange<1>::is_regular() const +IndexRange<1, indexT>::is_regular() const { // 1D case: always true return true; } +template bool -IndexRange<1>::get_regular_range(BasicCoordinate<1, int>& min_v, BasicCoordinate<1, int>& max_v) const +IndexRange<1, indexT>::get_regular_range(BasicCoordinate<1, indexT>& min_v, BasicCoordinate<1, indexT>& max_v) const { - // somewhat complicated as we can't assign ints to BasicCoordinate<1,int> - BasicCoordinate<1, int> tmp; + // somewhat complicated as we can't assign ints to BasicCoordinate<1,indexT> + BasicCoordinate<1, indexT> tmp; tmp[1] = min; min_v = tmp; tmp[1] = max; @@ -195,9 +298,15 @@ IndexRange<1>::get_regular_range(BasicCoordinate<1, int>& min_v, BasicCoordinate return true; } +template void -IndexRange<1>::resize(const int min_index, const int max_index) +IndexRange<1, indexT>::resize(const indexT min_index, const indexT max_index) { + if (max_index < min_index) + { + this->recycle(); + return; + } min = min_index; max = max_index; } diff --git a/src/include/stir/IndexRangeFwd.h b/src/include/stir/IndexRangeFwd.h new file mode 100644 index 0000000000..93d13f0c26 --- /dev/null +++ b/src/include/stir/IndexRangeFwd.h @@ -0,0 +1,26 @@ + +/* + Copyright (C) 2025, University College London + This file is part of STIR. + + SPDX-License-Identifier: Apache-2.0 + + See STIR/LICENSE.txt for details +*/ + +#ifndef ___stir_IndexRangeFwd_H___ +#define ___stir_IndexRangeFwd_H___ +/*! + \file + \ingroup Array + \brief forward declaration of stir::IndexRange class +*/ +#include "stir/common.h" + +namespace stir +{ +template +class IndexRange; +} // namespace stir + +#endif diff --git a/src/include/stir/NumericVectorWithOffset.h b/src/include/stir/NumericVectorWithOffset.h index e8f975de56..d8b7b6773c 100644 --- a/src/include/stir/NumericVectorWithOffset.h +++ b/src/include/stir/NumericVectorWithOffset.h @@ -43,20 +43,20 @@ START_NAMESPACE_STIR \warning It is likely that the automatic growing feature will be removed at some point. */ -template -class NumericVectorWithOffset : public VectorWithOffset +template +class NumericVectorWithOffset : public VectorWithOffset { #ifdef SWIG public: // needs to be public for SWIG to be able to parse the "using" statement below #endif - typedef VectorWithOffset base_type; + typedef VectorWithOffset base_type; public: using base_type::base_type; //! Constructor from an object of this class' base_type - inline NumericVectorWithOffset(const VectorWithOffset& t); + inline NumericVectorWithOffset(const VectorWithOffset& t); //! Constructor from an object of this class' base_type inline NumericVectorWithOffset(const NumericVectorWithOffset& t) diff --git a/src/include/stir/NumericVectorWithOffset.inl b/src/include/stir/NumericVectorWithOffset.inl index fff2770928..15e070c97f 100644 --- a/src/include/stir/NumericVectorWithOffset.inl +++ b/src/include/stir/NumericVectorWithOffset.inl @@ -28,31 +28,31 @@ START_NAMESPACE_STIR -template -NumericVectorWithOffset::NumericVectorWithOffset(const VectorWithOffset& t) +template +NumericVectorWithOffset::NumericVectorWithOffset(const VectorWithOffset& t) : base_type(t) {} -template -NumericVectorWithOffset::NumericVectorWithOffset(NumericVectorWithOffset&& other) noexcept +template +NumericVectorWithOffset::NumericVectorWithOffset(NumericVectorWithOffset&& other) noexcept : NumericVectorWithOffset() { swap(*this, other); } // assignment -template -NumericVectorWithOffset& -NumericVectorWithOffset::operator=(const NumericVectorWithOffset& other) +template +NumericVectorWithOffset& +NumericVectorWithOffset::operator=(const NumericVectorWithOffset& other) { base_type::operator=(other); return *this; } // addition -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator+(const NumericVectorWithOffset& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator+(const NumericVectorWithOffset& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -60,9 +60,9 @@ NumericVectorWithOffset::operator+(const NumericVectorWithOffset& v) } // subtraction -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator-(const NumericVectorWithOffset& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator-(const NumericVectorWithOffset& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -70,9 +70,9 @@ NumericVectorWithOffset::operator-(const NumericVectorWithOffset& v) } // elem by elem multiplication -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator*(const NumericVectorWithOffset& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator*(const NumericVectorWithOffset& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -80,9 +80,9 @@ NumericVectorWithOffset::operator*(const NumericVectorWithOffset& v) } // elem by elem division -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator/(const NumericVectorWithOffset& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator/(const NumericVectorWithOffset& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -90,9 +90,9 @@ NumericVectorWithOffset::operator/(const NumericVectorWithOffset& v) } // Add a constant to every element -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator+(const NUMBER& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator+(const NUMBER& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -100,9 +100,9 @@ NumericVectorWithOffset::operator+(const NUMBER& v) const } // Subtract a constant from every element -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator-(const NUMBER& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator-(const NUMBER& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -110,9 +110,9 @@ NumericVectorWithOffset::operator-(const NUMBER& v) const } // Multiply every element by a constant -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator*(const NUMBER& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator*(const NUMBER& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -120,9 +120,9 @@ NumericVectorWithOffset::operator*(const NUMBER& v) const } // Divide every element by a constant -template -inline NumericVectorWithOffset -NumericVectorWithOffset::operator/(const NUMBER& v) const +template +inline NumericVectorWithOffset +NumericVectorWithOffset::operator/(const NUMBER& v) const { this->check_state(); NumericVectorWithOffset retval(*this); @@ -132,9 +132,9 @@ NumericVectorWithOffset::operator/(const NUMBER& v) const /*! This will grow the vector automatically if the 2nd argument has smaller min_index and/or larger max_index. New elements are first initialised with T() before adding.*/ -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator+=(const NumericVectorWithOffset& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator+=(const NumericVectorWithOffset& v) { this->check_state(); // first check if *this is empty @@ -143,16 +143,16 @@ NumericVectorWithOffset::operator+=(const NumericVectorWithOffset& v) return *this = v; } this->grow(std::min(this->get_min_index(), v.get_min_index()), std::max(this->get_max_index(), v.get_max_index())); - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) this->num[i] += v.num[i]; this->check_state(); return *this; } /*! See operator+= (const NumericVectorWithOffset&) for growing behaviour */ -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator-=(const NumericVectorWithOffset& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator-=(const NumericVectorWithOffset& v) { this->check_state(); // first check if *this is empty @@ -166,16 +166,16 @@ NumericVectorWithOffset::operator-=(const NumericVectorWithOffset& v) error("NumericVectorWithOffset-= called with rhs that is larger in size than lhs. This is no longer supported."); } this->grow(std::min(this->get_min_index(), v.get_min_index()), std::max(this->get_max_index(), v.get_max_index())); - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) this->num[i] -= v.num[i]; this->check_state(); return *this; } /*! See operator+= (const NumericVectorWithOffset&) for growing behaviour */ -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator*=(const NumericVectorWithOffset& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator*=(const NumericVectorWithOffset& v) { this->check_state(); // first check if *this is empty @@ -190,16 +190,16 @@ NumericVectorWithOffset::operator*=(const NumericVectorWithOffset& v) return *this; } this->grow(std::min(this->get_min_index(), v.get_min_index()), std::max(this->get_max_index(), v.get_max_index())); - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) this->num[i] *= v.num[i]; this->check_state(); return *this; } /*! See operator+= (const NumericVectorWithOffset&) for growing behaviour */ -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator/=(const NumericVectorWithOffset& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator/=(const NumericVectorWithOffset& v) { this->check_state(); // first check if *this is empty @@ -214,73 +214,73 @@ NumericVectorWithOffset::operator/=(const NumericVectorWithOffset& v) return *this; } this->grow(std::min(this->get_min_index(), v.get_min_index()), std::max(this->get_max_index(), v.get_max_index())); - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) this->num[i] /= v.num[i]; this->check_state(); return *this; } -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator+=(const NUMBER& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator+=(const NUMBER& v) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i] += v; this->check_state(); return *this; } -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator-=(const NUMBER& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator-=(const NUMBER& v) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i] -= v; this->check_state(); return *this; } -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator*=(const NUMBER& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator*=(const NUMBER& v) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i] *= v; this->check_state(); return *this; } -template -inline NumericVectorWithOffset& -NumericVectorWithOffset::operator/=(const NUMBER& v) +template +inline NumericVectorWithOffset& +NumericVectorWithOffset::operator/=(const NUMBER& v) { this->check_state(); - for (int i = this->get_min_index(); i <= this->get_max_index(); i++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) this->num[i] /= v; this->check_state(); return *this; } -template +template template inline void -NumericVectorWithOffset::axpby(const NUMBER2 a, - const NumericVectorWithOffset& x, - const NUMBER2 b, - const NumericVectorWithOffset& y) +NumericVectorWithOffset::axpby(const NUMBER2 a, + const NumericVectorWithOffset& x, + const NUMBER2 b, + const NumericVectorWithOffset& y) { - NumericVectorWithOffset::xapyb(x, a, y, b); + NumericVectorWithOffset::xapyb(x, a, y, b); } -template +template inline void -NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, - const NUMBER a, - const NumericVectorWithOffset& y, - const NUMBER b) +NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, + const NUMBER a, + const NumericVectorWithOffset& y, + const NUMBER b) { this->check_state(); if ((this->get_min_index() != x.get_min_index()) || (this->get_min_index() != y.get_min_index()) @@ -296,12 +296,12 @@ NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, } } -template +template inline void -NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, - const NumericVectorWithOffset& a, - const NumericVectorWithOffset& y, - const NumericVectorWithOffset& b) +NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, + const NumericVectorWithOffset& a, + const NumericVectorWithOffset& y, + const NumericVectorWithOffset& b) { this->check_state(); if ((this->get_min_index() != x.get_min_index()) || (this->get_min_index() != y.get_min_index()) @@ -322,10 +322,10 @@ NumericVectorWithOffset::xapyb(const NumericVectorWithOffset& x, } } -template +template template inline void -NumericVectorWithOffset::sapyb(const T2& a, const NumericVectorWithOffset& y, const T2& b) +NumericVectorWithOffset::sapyb(const T2& a, const NumericVectorWithOffset& y, const T2& b) { this->xapyb(*this, a, y, b); } diff --git a/src/include/stir/ProjDataInMemory.h b/src/include/stir/ProjDataInMemory.h index 9f3ae6a852..fa080e27fa 100644 --- a/src/include/stir/ProjDataInMemory.h +++ b/src/include/stir/ProjDataInMemory.h @@ -202,10 +202,10 @@ class ProjDataInMemory : public ProjData * iterator typedefs */ ///@{ - typedef Array<1, float>::iterator iterator; - typedef Array<1, float>::const_iterator const_iterator; - typedef Array<1, float>::full_iterator full_iterator; - typedef Array<1, float>::const_full_iterator const_full_iterator; + typedef Array<1, float, long long>::iterator iterator; + typedef Array<1, float, long long>::const_iterator const_iterator; + typedef Array<1, float, long long>::full_iterator full_iterator; + typedef Array<1, float, long long>::const_full_iterator const_full_iterator; ///@} //! start value for iterating through all elements in the array, see iterator @@ -277,7 +277,7 @@ class ProjDataInMemory : public ProjData //@} private: - Array<1, float> buffer; + Array<1, float, long long> buffer; //! allocates buffer for storing the data. Has to be called by constructors void create_buffer(const bool initialise_with_0 = false); diff --git a/src/include/stir/RunTests.h b/src/include/stir/RunTests.h index b406962608..7acc2187e6 100644 --- a/src/include/stir/RunTests.h +++ b/src/include/stir/RunTests.h @@ -138,8 +138,8 @@ class RunTests return check_if_equal(a.real(), b.real(), str) && check_if_equal(a.imag(), b.imag(), str); } //! check equality by comparing ranges and calling check_if_equal on all elements - template - bool check_if_equal(const VectorWithOffset& t1, const VectorWithOffset& t2, const std::string& str = "") + template + bool check_if_equal(const VectorWithOffset& t1, const VectorWithOffset& t2, const std::string& str = "") { if (t1.get_min_index() != t2.get_min_index() || t1.get_max_index() != t2.get_max_index()) { @@ -147,7 +147,7 @@ class RunTests return everything_ok = false; } - for (int i = t1.get_min_index(); i <= t1.get_max_index(); i++) + for (auto i = t1.get_min_index(); i <= t1.get_max_index(); i++) { if (!check_if_equal(t1[i], t2[i], str)) { @@ -183,8 +183,8 @@ class RunTests bool check_if_equal(const ProjDataInMemory& t1, const ProjDataInMemory& t2, const std::string& str = ""); // VC 6.0 needs definition of template members in the class def unfortunately. - template - bool check_if_equal(const IndexRange& t1, const IndexRange& t2, const std::string& str = "") + template + bool check_if_equal(const IndexRange& t1, const IndexRange& t2, const std::string& str = "") { if (t1 != t2) { @@ -230,10 +230,10 @@ class RunTests // VC needs definition of template members in the class def unfortunately. //! use check_if_zero on all elements - template - bool check_if_zero(const VectorWithOffset& t, const std::string& str = "") + template + bool check_if_zero(const VectorWithOffset& t, const std::string& str = "") { - for (int i = t.get_min_index(); i <= t.get_max_index(); i++) + for (auto i = t.get_min_index(); i <= t.get_max_index(); i++) { if (!check_if_zero(t[i], str)) { diff --git a/src/include/stir/VectorWithOffset.h b/src/include/stir/VectorWithOffset.h index 33416ae6ee..74e1ec71db 100644 --- a/src/include/stir/VectorWithOffset.h +++ b/src/include/stir/VectorWithOffset.h @@ -4,7 +4,7 @@ Copyright (C) 2000 PARAPET partners Copyright (C) 2000 - 2007-10-08, Hammersmith Imanet Ltd Copyright (C) 2012-06-01 - 2012, Kris Thielemans - Copyright (C) 2023 - 2025, University College London + Copyright (C) 2023 - 2026, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -23,6 +23,7 @@ */ +#include "stir/VectorWithOffsetFwd.h" #include "stir/shared_ptr.h" #include "stir/deprecated.h" #include @@ -60,7 +61,7 @@ START_NAMESPACE_STIR we would have to use uninitialized_copy etc. in some places. */ -template +template class VectorWithOffset { public: @@ -83,10 +84,10 @@ class VectorWithOffset inline VectorWithOffset(); //! Construct a VectorWithOffset of given length (initialised with \c T()) - inline explicit VectorWithOffset(const int hsz); + inline explicit VectorWithOffset(const indexT hsz); //! Construct a VectorWithOffset with offset \c min_index (initialised with \c T()) - inline VectorWithOffset(const int min_index, const int max_index); + inline VectorWithOffset(const indexT min_index, const indexT max_index); #if STIR_VERSION < 070000 //! Construct a VectorWithOffset of given length pointing to existing data @@ -96,7 +97,7 @@ class VectorWithOffset \deprecated */ - STIR_DEPRECATED VectorWithOffset(const int hsz, T* const data_ptr, T* const end_of_data_ptr); + STIR_DEPRECATED VectorWithOffset(const indexT hsz, T* const data_ptr, T* const end_of_data_ptr); //! Construct a VectorWithOffset with offset \c min_index pointing to existing data /*! @@ -105,28 +106,31 @@ class VectorWithOffset \deprecated */ - STIR_DEPRECATED inline VectorWithOffset(const int min_index, const int max_index, T* const data_ptr, T* const end_of_data_ptr); + STIR_DEPRECATED inline VectorWithOffset(const indexT min_index, + const indexT max_index, + T* const data_ptr, + T* const end_of_data_ptr); #endif //! Construct a VectorWithOffset of given length from a bare pointer (copying data) - VectorWithOffset(const int hsz, const T* const data_ptr); + VectorWithOffset(const indexT hsz, const T* const data_ptr); //! Construct a VectorWithOffset with offset \c min_index from a bare pointer (copying data) - inline VectorWithOffset(const int min_index, const int max_index, const T* const data_ptr); + inline VectorWithOffset(const indexT min_index, const indexT max_index, const T* const data_ptr); //! Construct a VectorWithOffset sharing existing data /*! \warning This refers to the original memory range, so any modifications to this object will modify the original data as well. */ - inline VectorWithOffset(const int min_index, const int max_index, shared_ptr data_sptr); + inline VectorWithOffset(const indexT min_index, const indexT max_index, shared_ptr data_sptr); //! Construct a VectorWithOffset sharing existing data /*! \warning This refers to the original memory range, so any modifications to this object will modify the original data as well. */ - inline VectorWithOffset(const int sz, shared_ptr data_sptr) + inline VectorWithOffset(const indexT sz, shared_ptr data_sptr) : VectorWithOffset(0, sz - 1, data_sptr) {} @@ -169,22 +173,22 @@ class VectorWithOffset //@{ //! return number of elements in this vector /*! \deprecated Use size() instead. */ - inline int get_length() const; + inline indexT get_length() const; //! return number of elements in this vector inline size_t size() const; //! get value of first valid index - inline int get_min_index() const; + inline indexT get_min_index() const; //! get value of last valid index - inline int get_max_index() const; + inline indexT get_max_index() const; //! change value of starting index - inline void set_offset(const int min_index); + inline void set_offset(const indexT min_index); //! identical to set_offset() - inline void set_min_index(const int min_index); + inline void set_min_index(const indexT min_index); //! grow the range of the vector, new elements are set to \c T() /*! Currently, it is only checked with assert() if old range is @@ -194,10 +198,10 @@ class VectorWithOffset resize() in a derived class, it is probably safest to overload grow() as well. */ - inline virtual void grow(const int min_index, const int max_index); + inline virtual void grow(const indexT min_index, const indexT max_index); //! grow the range of the vector from 0 to new_size-1, new elements are set to \c T() - inline void grow(const unsigned int new_size); + inline void grow(const size_type new_size); //! change the range of the vector, new elements are set to \c T() /*! New memory is allocated if the range grows outside the range specified by @@ -207,16 +211,16 @@ class VectorWithOffset \todo in principle reallocation could be avoided when the new range would fit in the old one by shifting. */ - inline virtual void resize(const int min_index, const int max_index); + inline virtual void resize(const indexT min_index, const indexT max_index); //! change the range of the vector from 0 to new_size-1, new elements are set to \c T() - inline void resize(const unsigned int new_size); + inline void resize(const size_type new_size); //! make the allocated range at least from \a min_index to \a max_index - inline void reserve(const int min_index, const int max_index); + inline void reserve(const indexT min_index, const indexT max_index); //! make the allocated range at least from 0 to new_size-1 - inline void reserve(const unsigned int new_size); + inline void reserve(const size_type new_size); //! get allocated size inline size_t capacity() const; @@ -232,7 +236,7 @@ class VectorWithOffset For a vector of 0 length, this function returns 0. */ - inline int get_capacity_min_index() const; + inline indexT get_capacity_min_index() const; //! get max_index within allocated range /*! This value depends on get_min_index() and hence will change @@ -240,23 +244,23 @@ class VectorWithOffset For a vector of 0 length, this function returns capacity()-1. */ - inline int get_capacity_max_index() const; + inline indexT get_capacity_max_index() const; //@} //! allow array-style access, read/write - inline T& operator[](int i); + inline T& operator[](indexT i); //! array access, read-only - inline const T& operator[](int i) const; + inline const T& operator[](indexT i) const; //! allow array-style access, read/write, but with range checking (throws std::out_of_range) - inline T& at(int i); + inline T& at(indexT i); //! array access, read-only, but with range checking (throws std::out_of_range) - inline const T& at(int i) const; + inline const T& at(indexT i) const; //! checks if the vector is empty - inline bool empty() const; + inline virtual bool empty() const; //! \name comparison operators //@{ @@ -377,7 +381,7 @@ class VectorWithOffset calls resize() */ - inline void init_with_copy(const int min_index, const int max_index, T const* const data_ptr); + inline void init_with_copy(const indexT min_index, const indexT max_index, T const* const data_ptr); //! initialise vector to the given index range and either copy from or point to \c data_ptr /*! \arg data_ptr should start to a contiguous block of correct size @@ -387,13 +391,13 @@ class VectorWithOffset \warning This function should only be called from within a constructor. It will ignore any existing content and therefore would cause memory leaks. */ - inline void init(const int min_index, const int max_index, T* const data_ptr, bool copy_data); + inline void init(const indexT min_index, const indexT max_index, T* const data_ptr, bool copy_data); private: //! length of vector - unsigned int length; + size_type length; //! starting index - int start; + indexT start; T* begin_allocated_memory; T* end_allocated_memory; diff --git a/src/include/stir/VectorWithOffset.inl b/src/include/stir/VectorWithOffset.inl index a9ada9b8ae..19a2d7ed6a 100644 --- a/src/include/stir/VectorWithOffset.inl +++ b/src/include/stir/VectorWithOffset.inl @@ -4,7 +4,7 @@ Copyright (C) 2000 PARAPET partners Copyright (C) 2000 - 2010-07-01, Hammersmith Imanet Ltd Copyright (C) 2012-06-01 - 2012, Kris Thielemans - Copyright (C) 2023 - 2024, University College London + Copyright (C) 2023 - 2026, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -25,14 +25,15 @@ #include "stir/IndexRange.h" #include #include +#include #include "thresholding.h" #include "stir/error.h" START_NAMESPACE_STIR -template +template void -VectorWithOffset::init() +VectorWithOffset::init() { length = 0; // i.e. an empty row of zero length, start = 0; // no offsets @@ -42,18 +43,19 @@ VectorWithOffset::init() allocated_memory_sptr = nullptr; } -template +template void -VectorWithOffset::init_with_copy(const int min_index, const int max_index, T const* const data_ptr) +VectorWithOffset::init_with_copy(const indexT min_index, const indexT max_index, T const* const data_ptr) { this->pointer_access = false; this->resize(min_index, max_index); - std::copy(data_ptr, data_ptr + this->length, this->begin()); + if (this->length > 0) + std::copy(data_ptr, data_ptr + this->length, this->begin()); } -template +template void -VectorWithOffset::init(const int min_index, const int max_index, T* const data_ptr, bool copy_data) +VectorWithOffset::init(const indexT min_index, const indexT max_index, T* const data_ptr, bool copy_data) { if (copy_data) { @@ -62,7 +64,7 @@ VectorWithOffset::init(const int min_index, const int max_index, T* const dat else { this->pointer_access = false; - this->length = max_index >= min_index ? static_cast(max_index - min_index) + 1 : 0U; + this->length = max_index >= min_index ? static_cast(max_index - min_index) + 1 : 0U; this->start = min_index; this->begin_allocated_memory = data_ptr; this->end_allocated_memory = data_ptr + this->length; @@ -71,9 +73,9 @@ VectorWithOffset::init(const int min_index, const int max_index, T* const dat } } -template +template bool -VectorWithOffset::owns_memory_for_data() const +VectorWithOffset::owns_memory_for_data() const { return this->allocated_memory_sptr ? true : false; } @@ -82,9 +84,9 @@ VectorWithOffset::owns_memory_for_data() const This function (only non-empty when debugging) is used before and after any modification of the object */ -template +template void -VectorWithOffset::check_state() const +VectorWithOffset::check_state() const { // disable for normal debugging #if _DEBUG > 1 @@ -93,13 +95,13 @@ VectorWithOffset::check_state() const #endif assert(begin_allocated_memory <= num + start); assert(end_allocated_memory >= begin_allocated_memory); - assert(static_cast(end_allocated_memory - begin_allocated_memory) >= length); + assert(static_cast(end_allocated_memory - begin_allocated_memory) >= length); assert(!allocated_memory_sptr || (allocated_memory_sptr.get() == begin_allocated_memory)); } -template +template void -VectorWithOffset::_destruct_and_deallocate() +VectorWithOffset::_destruct_and_deallocate() { // check if data is being accessed via a pointer (see get_data_ptr()) assert(pointer_access == false); @@ -110,33 +112,34 @@ VectorWithOffset::_destruct_and_deallocate() this->allocated_memory_sptr = nullptr; } -template +template void -VectorWithOffset::recycle() +VectorWithOffset::recycle() { this->check_state(); this->_destruct_and_deallocate(); this->init(); } -template -int -VectorWithOffset::get_min_index() const +template +indexT +VectorWithOffset::get_min_index() const { return start; } -template -int -VectorWithOffset::get_max_index() const +template +indexT +VectorWithOffset::get_max_index() const { + assert(std::is_signed_v || (length > 0)); return start + length - 1; } /*! Out of range errors are detected using assert() */ -template +template T& -VectorWithOffset::operator[](int i) +VectorWithOffset::operator[](indexT i) { this->check_state(); assert(i >= this->get_min_index()); @@ -146,9 +149,9 @@ VectorWithOffset::operator[](int i) } /*! Out of range errors are detected using assert() */ -template +template const T& -VectorWithOffset::operator[](int i) const +VectorWithOffset::operator[](indexT i) const { this->check_state(); assert(i >= this->get_min_index()); @@ -157,9 +160,9 @@ VectorWithOffset::operator[](int i) const return num[i]; } -template +template T& -VectorWithOffset::at(int i) +VectorWithOffset::at(indexT i) { if (length == 0 || i < this->get_min_index() || i > this->get_max_index()) throw std::out_of_range("index out of range"); @@ -167,9 +170,9 @@ VectorWithOffset::at(int i) return num[i]; } -template +template const T& -VectorWithOffset::at(int i) const +VectorWithOffset::at(indexT i) const { if (length == 0 || i < this->get_min_index() || i > this->get_max_index()) throw std::out_of_range("index out of range"); @@ -178,92 +181,94 @@ VectorWithOffset::at(int i) const return num[i]; } -template +template bool -VectorWithOffset::empty() const +VectorWithOffset::empty() const { return length == 0; } -template -typename VectorWithOffset::iterator -VectorWithOffset::begin() +template +typename VectorWithOffset::iterator +VectorWithOffset::begin() { this->check_state(); - return typename VectorWithOffset::iterator(num + this->get_min_index()); + return typename VectorWithOffset::iterator(num + this->get_min_index()); } -template -typename VectorWithOffset::const_iterator -VectorWithOffset::begin() const +template +typename VectorWithOffset::const_iterator +VectorWithOffset::begin() const { this->check_state(); - return typename VectorWithOffset::const_iterator(num + this->get_min_index()); + return typename VectorWithOffset::const_iterator(num + this->get_min_index()); } -template -typename VectorWithOffset::iterator -VectorWithOffset::end() +template +typename VectorWithOffset::iterator +VectorWithOffset::end() { - this->check_state(); - return typename VectorWithOffset::iterator(num + (this->get_max_index() + 1)); + return this->begin() + this->length; } -template -typename VectorWithOffset::const_iterator -VectorWithOffset::end() const +template +typename VectorWithOffset::const_iterator +VectorWithOffset::end() const { - this->check_state(); - return typename VectorWithOffset::const_iterator(num + (this->get_max_index() + 1)); + return this->begin() + this->length; } -template -typename VectorWithOffset::reverse_iterator -VectorWithOffset::rbegin() +template +typename VectorWithOffset::reverse_iterator +VectorWithOffset::rbegin() { this->check_state(); return std::make_reverse_iterator(end()); } -template -typename VectorWithOffset::const_reverse_iterator -VectorWithOffset::rbegin() const +template +typename VectorWithOffset::const_reverse_iterator +VectorWithOffset::rbegin() const { this->check_state(); return std::make_reverse_iterator(end()); } -template -typename VectorWithOffset::reverse_iterator -VectorWithOffset::rend() +template +typename VectorWithOffset::reverse_iterator +VectorWithOffset::rend() { this->check_state(); return std::make_reverse_iterator(begin()); } -template -typename VectorWithOffset::const_reverse_iterator -VectorWithOffset::rend() const +template +typename VectorWithOffset::const_reverse_iterator +VectorWithOffset::rend() const { this->check_state(); return std::make_reverse_iterator(begin()); } -template -VectorWithOffset::VectorWithOffset() +template +VectorWithOffset::VectorWithOffset() : pointer_access(false) { this->init(); } -template -VectorWithOffset::VectorWithOffset(const int hsz) - : VectorWithOffset(0, hsz - 1) -{} +template +VectorWithOffset::VectorWithOffset(const indexT hsz) + : VectorWithOffset(0, hsz > 0 ? hsz - 1 : 0) +{ + // note: somewhat awkward implementation to avoid problems when indexT is unsigned + if (hsz <= 0) + this->recycle(); +} -template -VectorWithOffset::VectorWithOffset(const int min_index, const int max_index) - : length(static_cast(max_index - min_index) + 1), +template +VectorWithOffset::VectorWithOffset(const indexT min_index, const indexT max_index) + : length(max_index >= min_index ? static_cast(max_index - min_index + 1) : 0), start(min_index), pointer_access(false) { @@ -280,9 +285,12 @@ VectorWithOffset::VectorWithOffset(const int min_index, const int max_index) } #if STIR_VERSION < 070000 -template -VectorWithOffset::VectorWithOffset(const int min_index, const int max_index, T* const data_ptr, T* const end_of_data_ptr) - : length(static_cast(max_index - min_index) + 1), +template +VectorWithOffset::VectorWithOffset(const indexT min_index, + const indexT max_index, + T* const data_ptr, + T* const end_of_data_ptr) + : length(static_cast(max_index - min_index) + 1), start(min_index), allocated_memory_sptr(nullptr), // we don't own the data pointer_access(false) @@ -293,106 +301,111 @@ VectorWithOffset::VectorWithOffset(const int min_index, const int max_index, this->check_state(); } -template -VectorWithOffset::VectorWithOffset(const int sz, T* const data_ptr, T* const end_of_data_ptr) - : VectorWithOffset(0, sz - 1, data_ptr, end_of_data_ptr) +template +VectorWithOffset::VectorWithOffset(const indexT sz, T* const data_ptr, T* const end_of_data_ptr) + : VectorWithOffset((std::is_signed_v || sz > 0) ? 0 : 1, + (std::is_signed_v || sz > 0) ? sz - 1 : 0, + data_ptr, + end_of_data_ptr) {} #endif // STIR_VERSION < 070000 -template -VectorWithOffset::VectorWithOffset(const int min_index, const int max_index, const T* const data_ptr) +template +VectorWithOffset::VectorWithOffset(const indexT min_index, const indexT max_index, const T* const data_ptr) { // first set empty, such that resize() will work ok this->init(); this->init_with_copy(min_index, max_index, data_ptr); } -template -VectorWithOffset::VectorWithOffset(const int sz, const T* const data_ptr) - : VectorWithOffset(0, sz - 1, data_ptr) +template +VectorWithOffset::VectorWithOffset(const indexT sz, const T* const data_ptr) + : VectorWithOffset((std::is_signed_v || sz > 0) ? 0 : 1, (std::is_signed_v || sz > 0) ? sz - 1 : 0, data_ptr) {} -template -VectorWithOffset::VectorWithOffset(const int min_index, const int max_index, shared_ptr data_sptr) +template +VectorWithOffset::VectorWithOffset(const indexT min_index, const indexT max_index, shared_ptr data_sptr) { this->allocated_memory_sptr = data_sptr; this->init(min_index, max_index, data_sptr.get(), /* copy_data = */ false); } -template -VectorWithOffset::VectorWithOffset(VectorWithOffset&& other) noexcept +template +VectorWithOffset::VectorWithOffset(VectorWithOffset&& other) noexcept : VectorWithOffset() { swap(*this, other); } -template -VectorWithOffset::~VectorWithOffset() +template +VectorWithOffset::~VectorWithOffset() { // check if data is being accessed via a pointer (see get_data_ptr()) assert(pointer_access == false); _destruct_and_deallocate(); } -template +template void -VectorWithOffset::set_offset(const int min_index) +VectorWithOffset::set_offset(const indexT min_index) { this->check_state(); // only do something when non-zero length if (length == 0) return; - num += start - min_index; + // note: num += (start - min_index), but split up in 2 steps in case indexT is unsigned + num += start; + num -= min_index; start = min_index; } -template +template void -VectorWithOffset::set_min_index(const int min_index) +VectorWithOffset::set_min_index(const indexT min_index) { this->set_offset(min_index); } -template +template size_t -VectorWithOffset::capacity() const +VectorWithOffset::capacity() const { return size_t(end_allocated_memory - begin_allocated_memory); } -template -int -VectorWithOffset::get_capacity_min_index() const +template +indexT +VectorWithOffset::get_capacity_min_index() const { // the behaviour for length==0 depends on num==begin_allocated_memory assert(length > 0 || num == begin_allocated_memory); - return static_cast(begin_allocated_memory - num); + return static_cast(begin_allocated_memory - num); } -template -int -VectorWithOffset::get_capacity_max_index() const +template +indexT +VectorWithOffset::get_capacity_max_index() const { // the behaviour for length==0 depends on num==begin_allocated_memory assert(length > 0 || num == begin_allocated_memory); - return static_cast(end_allocated_memory - num - 1); + return static_cast(end_allocated_memory - num - 1); } // the new members will be initialised with the default constructor for T // but this should change in the future -template +template void -VectorWithOffset::reserve(const int new_capacity_min_index, const int new_capacity_max_index) +VectorWithOffset::reserve(const indexT new_capacity_min_index, const indexT new_capacity_max_index) { this->check_state(); - const int actual_capacity_min_index + const indexT actual_capacity_min_index = length == 0 ? new_capacity_min_index : std::min(this->get_capacity_min_index(), new_capacity_min_index); - const int actual_capacity_max_index + const indexT actual_capacity_max_index = length == 0 ? new_capacity_max_index : std::max(this->get_capacity_max_index(), new_capacity_max_index); if (actual_capacity_min_index > actual_capacity_max_index) return; - const unsigned int new_capacity = static_cast(actual_capacity_max_index - actual_capacity_min_index) + 1; + const size_type new_capacity = static_cast(actual_capacity_max_index - actual_capacity_min_index + 1); if (new_capacity <= this->capacity()) return; @@ -400,7 +413,9 @@ VectorWithOffset::reserve(const int new_capacity_min_index, const int new_cap assert(pointer_access == false); // TODO use allocator here instead of new shared_ptr new_allocated_memory_sptr(new T[new_capacity]); - const unsigned extra_at_the_left = length == 0 ? 0U : std::max(0, this->get_min_index() - actual_capacity_min_index); + + const size_type extra_at_the_left + = length == 0 ? 0U : std::max(indexT(0), indexT(this->get_min_index() - actual_capacity_min_index)); std::copy(this->begin(), this->end(), new_allocated_memory_sptr.get() + extra_at_the_left); this->_destruct_and_deallocate(); allocated_memory_sptr = std::move(new_allocated_memory_sptr); @@ -410,20 +425,20 @@ VectorWithOffset::reserve(const int new_capacity_min_index, const int new_cap this->check_state(); } -template +template void -VectorWithOffset::reserve(const unsigned int new_size) +VectorWithOffset::reserve(const size_type new_size) { // note: for 0 new_size, we avoid a wrap-around // otherwise we would be reserving quite a lot of memory! if (new_size != 0) - reserve(0, static_cast(new_size - 1)); + reserve(0, static_cast(new_size - 1)); } // the new members will be initialised with the default constructor for T -template +template void -VectorWithOffset::resize(const int min_index, const int max_index) +VectorWithOffset::resize(const indexT min_index, const indexT max_index) { this->check_state(); if (min_index > max_index) @@ -433,16 +448,16 @@ VectorWithOffset::resize(const int min_index, const int max_index) num = begin_allocated_memory; return; } - const unsigned old_length = length; + const size_type old_length = length; if (old_length > 0) { if (min_index == this->get_min_index() && max_index == this->get_max_index()) return; // determine overlapping range to avoid copying too much data when calling reserve() - const int overlap_min_index = std::max(this->get_min_index(), min_index); - const int overlap_max_index = std::min(this->get_max_index(), max_index); + const indexT overlap_min_index = std::max(this->get_min_index(), min_index); + const indexT overlap_max_index = std::min(this->get_max_index(), max_index); // TODO when using non-initialised memory, call delete here on elements that go out of range - length = overlap_max_index - overlap_min_index < 0 ? 0 : static_cast(overlap_max_index - overlap_min_index) + 1; + length = overlap_max_index - overlap_min_index < 0 ? 0 : static_cast(overlap_max_index - overlap_min_index) + 1; if (length == 0) { start = 0; @@ -454,11 +469,11 @@ VectorWithOffset::resize(const int min_index, const int max_index) start = overlap_min_index; } } // end if (length>0) - const unsigned overlapping_length = length; + const size_type overlapping_length = length; this->reserve(min_index, max_index); // TODO when using allocator, call default constructor for new elements here // (and delete the ones that go out of range!) - length = static_cast(max_index - min_index) + 1; + length = static_cast(max_index - min_index) + 1; start = min_index; if (overlapping_length > 0) { @@ -472,9 +487,9 @@ VectorWithOffset::resize(const int min_index, const int max_index) this->check_state(); } -template +template void -VectorWithOffset::resize(const unsigned new_size) +VectorWithOffset::resize(const size_type new_size) { if (new_size == 0) { @@ -483,26 +498,26 @@ VectorWithOffset::resize(const unsigned new_size) num = begin_allocated_memory; } else - this->resize(0, static_cast(new_size - 1)); + this->resize(0, static_cast(new_size - 1)); } -template +template void -VectorWithOffset::grow(const int min_index, const int max_index) +VectorWithOffset::grow(const indexT min_index, const indexT max_index) { this->resize(min_index, max_index); } -template +template void -VectorWithOffset::grow(const unsigned new_size) +VectorWithOffset::grow(const size_type new_size) { - this->grow(0, static_cast(new_size - 1)); + this->resize(new_size); } -template -VectorWithOffset& -VectorWithOffset::operator=(const VectorWithOffset& il) +template +VectorWithOffset& +VectorWithOffset::operator=(const VectorWithOffset& il) { this->check_state(); if (this == &il) @@ -525,33 +540,33 @@ VectorWithOffset::operator=(const VectorWithOffset& il) return *this; } -template -VectorWithOffset::VectorWithOffset(const VectorWithOffset& il) +template +VectorWithOffset::VectorWithOffset(const VectorWithOffset& il) : pointer_access(false) { this->init(); *this = il; // Uses assignment operator (above) } -template -int -VectorWithOffset::get_length() const +template +indexT +VectorWithOffset::get_length() const { this->check_state(); - return static_cast(length); + return static_cast(length); } -template +template size_t -VectorWithOffset::size() const +VectorWithOffset::size() const { this->check_state(); return size_t(length); } -template +template bool -VectorWithOffset::operator==(const VectorWithOffset& iv) const +VectorWithOffset::operator==(const VectorWithOffset& iv) const { this->check_state(); if (length != iv.length || start != iv.start) @@ -559,34 +574,34 @@ VectorWithOffset::operator==(const VectorWithOffset& iv) const return std::equal(this->begin(), this->end(), iv.begin()); } -template +template bool -VectorWithOffset::operator!=(const VectorWithOffset& iv) const +VectorWithOffset::operator!=(const VectorWithOffset& iv) const { return !(*this == iv); } -template +template void -VectorWithOffset::fill(const T& n) +VectorWithOffset::fill(const T& n) { this->check_state(); std::fill(this->begin(), this->end(), n); this->check_state(); } -template +template inline void -VectorWithOffset::apply_lower_threshold(const T& lower) +VectorWithOffset::apply_lower_threshold(const T& lower) { this->check_state(); threshold_lower(this->begin(), this->end(), lower); this->check_state(); } -template +template inline void -VectorWithOffset::apply_upper_threshold(const T& upper) +VectorWithOffset::apply_upper_threshold(const T& upper) { this->check_state(); threshold_upper(this->begin(), this->end(), upper); @@ -607,9 +622,9 @@ VectorWithOffset::apply_upper_threshold(const T& upper) get_const_data_ptr() and release_data_ptr(). (This is checked with assert() in DEBUG mode.) */ -template +template T* -VectorWithOffset::get_data_ptr() +VectorWithOffset::get_data_ptr() { assert(!pointer_access); @@ -629,9 +644,9 @@ VectorWithOffset::get_data_ptr() \see get_data_ptr() */ -template +template const T* -VectorWithOffset::get_const_data_ptr() const +VectorWithOffset::get_const_data_ptr() const { assert(!pointer_access); @@ -650,9 +665,9 @@ VectorWithOffset::get_const_data_ptr() const \see get_data_ptr() */ -template +template void -VectorWithOffset::release_data_ptr() +VectorWithOffset::release_data_ptr() { assert(pointer_access); @@ -667,9 +682,9 @@ VectorWithOffset::release_data_ptr() \see get_const_data_ptr() */ -template +template void -VectorWithOffset::release_const_data_ptr() const +VectorWithOffset::release_const_data_ptr() const { assert(pointer_access); @@ -677,9 +692,9 @@ VectorWithOffset::release_const_data_ptr() const } /********************** arithmetic operators ****************/ -template -inline VectorWithOffset& -VectorWithOffset::operator+=(const VectorWithOffset& v) +template +inline VectorWithOffset& +VectorWithOffset::operator+=(const VectorWithOffset& v) { this->check_state(); #if 1 @@ -693,15 +708,15 @@ VectorWithOffset::operator+=(const VectorWithOffset& v) } this->grow(std::min(get_min_index(), v.get_min_index()), std::max(get_max_index(), v.get_max_index())); #endif - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) num[i] += v.num[i]; this->check_state(); return *this; } -template -inline VectorWithOffset& -VectorWithOffset::operator-=(const VectorWithOffset& v) +template +inline VectorWithOffset& +VectorWithOffset::operator-=(const VectorWithOffset& v) { this->check_state(); #if 1 @@ -716,15 +731,15 @@ VectorWithOffset::operator-=(const VectorWithOffset& v) } grow(std::min(get_min_index(), v.get_min_index()), std::max(get_max_index(), v.get_max_index())); #endif - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) num[i] -= v.num[i]; this->check_state(); return *this; } -template -inline VectorWithOffset& -VectorWithOffset::operator*=(const VectorWithOffset& v) +template +inline VectorWithOffset& +VectorWithOffset::operator*=(const VectorWithOffset& v) { this->check_state(); #if 1 @@ -740,15 +755,15 @@ VectorWithOffset::operator*=(const VectorWithOffset& v) } grow(std::min(get_min_index(), v.get_min_index()), std::max(get_max_index(), v.get_max_index())); #endif - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) num[i] *= v.num[i]; this->check_state(); return *this; } -template -inline VectorWithOffset& -VectorWithOffset::operator/=(const VectorWithOffset& v) +template +inline VectorWithOffset& +VectorWithOffset::operator/=(const VectorWithOffset& v) { this->check_state(); #if 1 @@ -764,7 +779,7 @@ VectorWithOffset::operator/=(const VectorWithOffset& v) } grow(std::min(get_min_index(), v.get_min_index()), std::max(get_max_index(), v.get_max_index())); #endif - for (int i = v.get_min_index(); i <= v.get_max_index(); i++) + for (auto i = v.get_min_index(); i <= v.get_max_index(); i++) num[i] /= v.num[i]; this->check_state(); return *this; @@ -774,9 +789,9 @@ VectorWithOffset::operator/=(const VectorWithOffset& v) #if 0 // disabled for now // warning: not tested -template -inline VectorWithOffset& -VectorWithOffset::operator+= (const T &t) +template +inline VectorWithOffset& +VectorWithOffset::operator+= (const T &t) { typename iterator iter = this->begin(); const typename iterator end_iter = this->end(); @@ -784,27 +799,27 @@ VectorWithOffset::operator+= (const T &t) *iter++ += t; } -template -inline VectorWithOffset& -VectorWithOffset::operator-= (const T &t) +template +inline VectorWithOffset& +VectorWithOffset::operator-= (const T &t) { typename iterator iter = this->begin(); const typename iterator end_iter = this->end(); while (iter != end_iter) *iter++ -= t; } -template -inline VectorWithOffset& -VectorWithOffset::operator*= (const T &t) +template +inline VectorWithOffset& +VectorWithOffset::operator*= (const T &t) { typename iterator iter = this->begin(); const typename iterator end_iter = this->end(); while (iter != end_iter) *iter++ *= t; } -template -inline VectorWithOffset& -VectorWithOffset::operator/= (const T &t) +template +inline VectorWithOffset& +VectorWithOffset::operator/= (const T &t) { typename iterator iter = this->begin(); const typename iterator end_iter = this->end(); @@ -817,9 +832,9 @@ VectorWithOffset::operator/= (const T &t) /**** operator* etc ********/ // addition -template -inline VectorWithOffset -VectorWithOffset::operator+(const VectorWithOffset& v) const +template +inline VectorWithOffset +VectorWithOffset::operator+(const VectorWithOffset& v) const { this->check_state(); VectorWithOffset retval(*this); @@ -827,9 +842,9 @@ VectorWithOffset::operator+(const VectorWithOffset& v) const } // subtraction -template -inline VectorWithOffset -VectorWithOffset::operator-(const VectorWithOffset& v) const +template +inline VectorWithOffset +VectorWithOffset::operator-(const VectorWithOffset& v) const { this->check_state(); VectorWithOffset retval(*this); @@ -837,9 +852,9 @@ VectorWithOffset::operator-(const VectorWithOffset& v) const } // elem by elem multiplication -template -inline VectorWithOffset -VectorWithOffset::operator*(const VectorWithOffset& v) const +template +inline VectorWithOffset +VectorWithOffset::operator*(const VectorWithOffset& v) const { this->check_state(); VectorWithOffset retval(*this); @@ -847,9 +862,9 @@ VectorWithOffset::operator*(const VectorWithOffset& v) const } // elem by elem division -template -inline VectorWithOffset -VectorWithOffset::operator/(const VectorWithOffset& v) const +template +inline VectorWithOffset +VectorWithOffset::operator/(const VectorWithOffset& v) const { this->check_state(); VectorWithOffset retval(*this); diff --git a/src/include/stir/VectorWithOffsetFwd.h b/src/include/stir/VectorWithOffsetFwd.h new file mode 100644 index 0000000000..49174cedc2 --- /dev/null +++ b/src/include/stir/VectorWithOffsetFwd.h @@ -0,0 +1,26 @@ + +/* + Copyright (C) 2025, University College London + This file is part of STIR. + + SPDX-License-Identifier: Apache-2.0 + + See STIR/LICENSE.txt for details +*/ + +#ifndef ___stir_VectorWithOffsetFwd_H___ +#define ___stir_VectorWithOffsetFwd_H___ +/*! + \file + \ingroup Array + \brief forward declaration of stir::VectorWithOffset class +*/ +#include "stir/common.h" + +namespace stir +{ +template +class VectorWithOffset; +} // namespace stir + +#endif diff --git a/src/include/stir/array_index_functions.h b/src/include/stir/array_index_functions.h index 0ef8504331..1d5d684df3 100644 --- a/src/include/stir/array_index_functions.h +++ b/src/include/stir/array_index_functions.h @@ -34,21 +34,21 @@ START_NAMESPACE_STIR //! an alternative for array indexing using BasicCoordinate objects /*! Case where the index has lower dimension than the array*/ -template -inline const Array& get(const Array& a, - const BasicCoordinate& c); +template +inline const Array& get(const Array& a, + const BasicCoordinate& c); //! an alternative for array indexing using BasicCoordinate objects /*! Case where the index has the same dimension as the array*/ -template -inline const elemT& get(const Array& a, const BasicCoordinate& c); +template +inline const elemT& get(const Array& a, const BasicCoordinate& c); //! Get the first multi-dimensional index of the array /*! \todo If the array \arg a is empty, we return an object where all indices are 0. It would be better to throw an exception. */ -template -inline BasicCoordinate get_min_indices(const Array& a); +template +inline BasicCoordinate get_min_indices(const Array& a); //! Given an index into an array, increment it to the next one /*! @@ -57,16 +57,15 @@ inline BasicCoordinate get_min_indices(const Array array = ...; - BasicCoordinate indices = - get_min_indices(array); + auto indices = get_min_indices(array); do { something with indices } while (next(indices, array)); \endcode \warning The above loop will fail for empty arrays */ -template -inline bool next(BasicCoordinate& indices, const Array& a); +template +inline bool next(BasicCoordinate& indices, const Array& a); //@} diff --git a/src/include/stir/array_index_functions.inl b/src/include/stir/array_index_functions.inl index 835ee6bfa0..730ce2953c 100644 --- a/src/include/stir/array_index_functions.inl +++ b/src/include/stir/array_index_functions.inl @@ -35,9 +35,9 @@ START_NAMESPACE_STIR namespace detail { -template -inline BasicCoordinate -get_min_indices_help(is_not_1d, const Array& a) +template +inline BasicCoordinate +get_min_indices_help(is_not_1d, const Array& a) { if (a.get_min_index() <= a.get_max_index()) return join(a.get_min_index(), get_min_indices(*a.begin())); @@ -45,23 +45,23 @@ get_min_indices_help(is_not_1d, const Array& a) { // a is empty. Not clear what to return, so we just return 0 // It would be better to throw an exception. - BasicCoordinate tmp(0); + BasicCoordinate tmp(0); return tmp; } } -template -inline BasicCoordinate<1, int> -get_min_indices_help(is_1d, const Array<1, T>& a) +template +inline BasicCoordinate<1, indexT> +get_min_indices_help(is_1d, const Array<1, elemT, indexT>& a) { - BasicCoordinate<1, int> result; + BasicCoordinate<1, indexT> result; result[1] = a.get_min_index(); return result; } -template +template inline bool -next_help(is_1d, BasicCoordinate<1, int>& index, const Array& a) +next_help(is_1d, BasicCoordinate<1, indexT>& index, const Array& a) { if (a.get_min_index() > a.get_max_index()) return false; @@ -71,13 +71,13 @@ next_help(is_1d, BasicCoordinate<1, int>& index, const Array return index[1] <= a.get_max_index(); } -template +template inline bool -next_help(is_not_1d, BasicCoordinate& index, const Array& a) +next_help(is_not_1d, BasicCoordinate& index, const Array& a) { if (a.get_min_index() > a.get_max_index()) return false; - BasicCoordinate upper_index = cut_last_dimension(index); + auto upper_index = cut_last_dimension(index); assert(index[num_dimensions] >= get(a, upper_index).get_min_index()); assert(index[num_dimensions] <= get(a, upper_index).get_max_index()); index[num_dimensions]++; @@ -95,36 +95,36 @@ next_help(is_not_1d, BasicCoordinate& index, const Array -inline BasicCoordinate -get_min_indices(const Array& a) +template +inline BasicCoordinate +get_min_indices(const Array& a) { return detail::get_min_indices_help(detail::test_if_1d(), a); } -template +template inline bool -next(BasicCoordinate& index, const Array& a) +next(BasicCoordinate& index, const Array& a) { return detail::next_help(detail::test_if_1d(), index, a); } -template -inline const Array& -get(const Array& a, const BasicCoordinate& c) +template +inline const Array& +get(const Array& a, const BasicCoordinate& c) { return get(a[c[1]], cut_first_dimension(c)); } -template +template inline const elemT& -get(const Array& a, const BasicCoordinate& c) +get(const Array& a, const BasicCoordinate& c) { return a[c]; } -template -inline const Array& -get(const Array& a, const BasicCoordinate<1, int>& c) +template +inline const Array& +get(const Array& a, const BasicCoordinate<1, indexT>& c) { return a[c[1]]; } diff --git a/src/include/stir/centre_of_gravity.h b/src/include/stir/centre_of_gravity.h index 1cc103929c..8e2a1a02b8 100644 --- a/src/include/stir/centre_of_gravity.h +++ b/src/include/stir/centre_of_gravity.h @@ -22,10 +22,6 @@ START_NAMESPACE_STIR // predeclerations to avoid having to include the files and create unnecessary // dependencies -template -class BasicCoordinate; -template -class VectorWithOffset; template class CartesianCoordinate3D; template diff --git a/src/include/stir/convert_array.h b/src/include/stir/convert_array.h index 3e2f4e7336..5bce4c69e6 100644 --- a/src/include/stir/convert_array.h +++ b/src/include/stir/convert_array.h @@ -52,7 +52,7 @@ class NumericInfo; In that case, the same scale_factor is used as in the 0 case. \param data_in - some Array object, elements are of some numeric type \a T1 + some Array-like object, i.e. it needs to have \c begin_all(), \c end_all() methods \param info_for_out_type \a T2 is the desired output type @@ -61,9 +61,8 @@ class NumericInfo; \see convert_array */ -template -inline void -find_scale_factor(scaleT& scale_factor, const Array& data_in, const NumericInfo info_for_out_type); +template +inline void find_scale_factor(scaleT& scale_factor, const ArrayLike& data_in, const NumericInfo info_for_out_type); /*! \ingroup Array @@ -97,13 +96,13 @@ find_scale_factor(scaleT& scale_factor, const Array& data_in numbers are cut out) when T2 is an unsigned type. */ -template -inline Array -convert_array(scaleT& scale_factor, const Array& data_in, const NumericInfo info2); +template +inline Array +convert_array(scaleT& scale_factor, const Array& data_in, const NumericInfo info2); /*! \ingroup Array - \brief Converts the \c data_in Array to \c data_out (with elements of different types) such that \c data_in == \c data_out * \c - scale_factor + \brief Converts the \c data_in Array-like object to \c data_out (potentially with elements of different types) such that \c + data_in == \c data_out * \c scale_factor \par example \code @@ -112,9 +111,10 @@ convert_array(scaleT& scale_factor, const Array& data_in, co \see convert_array(scale_factor, data_in, info2) for more info */ - -template -inline void convert_array(Array& data_out, scaleT& scale_factor, const Array& data_in); +template +inline void convert_array(Array& data_out, + scaleT& scale_factor, + const Array& data_in); END_NAMESPACE_STIR diff --git a/src/include/stir/convert_array.inl b/src/include/stir/convert_array.inl index d424751a16..81f3fe3f3e 100644 --- a/src/include/stir/convert_array.inl +++ b/src/include/stir/convert_array.inl @@ -23,26 +23,28 @@ #include "stir/convert_range.h" START_NAMESPACE_STIR -template +template void -find_scale_factor(scaleT& scale_factor, const Array& data_in, const NumericInfo info_for_out_type) +find_scale_factor(scaleT& scale_factor, const ArrayLike& data_in, const NumericInfo info_for_out_type) { find_scale_factor(scale_factor, data_in.begin_all(), data_in.end_all(), info_for_out_type); } -template -Array -convert_array(scaleT& scale_factor, const Array& data_in, const NumericInfo info_for_out_type) +template +Array +convert_array(scaleT& scale_factor, const Array& data_in, const NumericInfo info_for_out_type) { - Array data_out(data_in.get_index_range()); + Array data_out(data_in.get_index_range()); convert_array(data_out, scale_factor, data_in); return data_out; } -template +template void -convert_array(Array& data_out, scaleT& scale_factor, const Array& data_in) +convert_array(Array& data_out, + scaleT& scale_factor, + const Array& data_in) { convert_range(data_out.begin_all(), scale_factor, data_in.begin_all(), data_in.end_all()); } diff --git a/src/include/stir/evaluation/compute_ROI_values.h b/src/include/stir/evaluation/compute_ROI_values.h index e4430f470c..58d2bf281e 100644 --- a/src/include/stir/evaluation/compute_ROI_values.h +++ b/src/include/stir/evaluation/compute_ROI_values.h @@ -20,6 +20,7 @@ #define __stir_evaluation_compute_ROI_values__H__ #include "stir/evaluation/ROIValues.h" +#include "stir/VectorWithOffsetFwd.h" #include "stir/VoxelsOnCartesianGrid.h" START_NAMESPACE_STIR @@ -28,8 +29,6 @@ template class CartesianCoordinate2D; template class CartesianCoordinate3D; -template -class VectorWithOffset; template class DiscretisedDensity; class Shape3D; diff --git a/src/include/stir/interpolate_projdata.h b/src/include/stir/interpolate_projdata.h index 0fe97060c7..80ee4bb9be 100644 --- a/src/include/stir/interpolate_projdata.h +++ b/src/include/stir/interpolate_projdata.h @@ -18,12 +18,11 @@ */ #include "stir/numerics/BSplines.h" +#include "stir/BasicCoordinateFwd.h" START_NAMESPACE_STIR class ProjData; -template -class BasicCoordinate; template class Sinogram; template diff --git a/src/include/stir/make_array.h b/src/include/stir/make_array.h index f5359a4811..8425191711 100644 --- a/src/include/stir/make_array.h +++ b/src/include/stir/make_array.h @@ -105,78 +105,80 @@ inline Array<1, T> make_1d_array(const T& a0, const T& a8, const T& a9); -template -inline Array make_array(const Array& a0); - -template -inline Array make_array(const Array& a0, const Array& a1); - -template -inline Array -make_array(const Array& a0, const Array& a1, const Array& a2); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7, - const Array& a8); - -template -inline Array make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7, - const Array& a8, - const Array& a9); +template +inline Array make_array(const Array& a0); + +template +inline Array make_array(const Array& a0, + const Array& a1); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7, + const Array& a8); + +template +inline Array make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7, + const Array& a8, + const Array& a9); END_NAMESPACE_STIR diff --git a/src/include/stir/make_array.inl b/src/include/stir/make_array.inl index dbb4568093..412a41e97d 100644 --- a/src/include/stir/make_array.inl +++ b/src/include/stir/make_array.inl @@ -2,6 +2,7 @@ // /* Copyright (C) 2005- 2005, Hammersmith Imanet Ltd + Copyright (C) 2025, University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 @@ -20,41 +21,41 @@ START_NAMESPACE_STIR -template -VectorWithOffset -make_vector(const T& a0) +template +VectorWithOffset +make_vector_with_index_type(const T& a0) { - VectorWithOffset a(1); + VectorWithOffset a(1); a[0] = a0; return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1) { - VectorWithOffset a(2); + VectorWithOffset a(2); a[0] = a0; a[1] = a1; return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1, const T& a2) { - VectorWithOffset a(3); + VectorWithOffset a(3); a[0] = a0; a[1] = a1; a[2] = a2; return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3) { - VectorWithOffset a(4); + VectorWithOffset a(4); a[0] = a0; a[1] = a1; a[2] = a2; @@ -62,11 +63,11 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3) return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) { - VectorWithOffset a(5); + VectorWithOffset a(5); a[0] = a0; a[1] = a1; a[2] = a2; @@ -75,11 +76,11 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5) { - VectorWithOffset a(6); + VectorWithOffset a(6); a[0] = a0; a[1] = a1; a[2] = a2; @@ -89,11 +90,11 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, con return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6) +template +VectorWithOffset +make_vector_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6) { - VectorWithOffset a(7); + VectorWithOffset a(7); a[0] = a0; a[1] = a1; a[2] = a2; @@ -104,11 +105,12 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, con return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7) +template +VectorWithOffset +make_vector_with_index_type( + const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7) { - VectorWithOffset a(8); + VectorWithOffset a(8); a[0] = a0; a[1] = a1; a[2] = a2; @@ -120,11 +122,38 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, con return a; } -template -VectorWithOffset -make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7, const T& a8) +template +VectorWithOffset +make_vector_with_index_type( + const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7, const T& a8) +{ + VectorWithOffset a(9); + a[0] = a0; + a[1] = a1; + a[2] = a2; + a[3] = a3; + a[4] = a4; + a[5] = a5; + a[6] = a6; + a[7] = a7; + a[8] = a8; + return a; +} + +template +VectorWithOffset +make_vector_with_index_type(const T& a0, + const T& a1, + const T& a2, + const T& a3, + const T& a4, + const T& a5, + const T& a6, + const T& a7, + const T& a8, + const T& a9) { - VectorWithOffset a(9); + VectorWithOffset a(10); a[0] = a0; a[1] = a1; a[2] = a2; @@ -134,9 +163,73 @@ make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, con a[6] = a6; a[7] = a7; a[8] = a8; + a[9] = a9; return a; } +template +VectorWithOffset +make_vector(const T& a0) +{ + return make_vector_with_index_type(a0); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1) +{ + return make_vector_with_index_type(a0, a1); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2) +{ + return make_vector_with_index_type(a0, a1, a2); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3) +{ + return make_vector_with_index_type(a0, a1, a2, a3); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) +{ + return make_vector_with_index_type(a0, a1, a2, a3, a4); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5) +{ + return make_vector_with_index_type(a0, a1, a2, a3, a4, a5); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6) +{ + return make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7) +{ + return make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7); +} + +template +VectorWithOffset +make_vector(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7, const T& a8) +{ + return make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8); +} + template VectorWithOffset make_vector(const T& a0, @@ -150,90 +243,166 @@ make_vector(const T& a0, const T& a8, const T& a9) { - VectorWithOffset a(10); - a[0] = a0; - a[1] = a1; - a[2] = a2; - a[3] = a3; - a[4] = a4; - a[5] = a5; - a[6] = a6; - a[7] = a7; - a[8] = a8; - a[9] = a9; + return make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9); +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0) +{ + const Array<1, T, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0)); return a; } -template +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1) +{ + const Array<1, T, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0, a1)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1, const T& a2) +{ + const Array<1, T, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3) +{ + const Array<1, T, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) +{ + const Array<1, T, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5) +{ + const Array<1, T, indexT> a + = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4, a5)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6) +{ + const Array<1, T, indexT> a + = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type( + const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7) +{ + const Array<1, T, indexT> a + = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type( + const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7, const T& a8) +{ + const Array<1, T, indexT> a + = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8)); + return a; +} + +template +Array<1, T, indexT> +make_1d_array_with_index_type(const T& a0, + const T& a1, + const T& a2, + const T& a3, + const T& a4, + const T& a5, + const T& a6, + const T& a7, + const T& a8, + const T& a9) +{ + const Array<1, T, indexT> a + = NumericVectorWithOffset(make_vector_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9)); + return a; +} + +template Array<1, T> make_1d_array(const T& a0) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0)); - return a; + return make_1d_array_with_index_type(a0); } template Array<1, T> make_1d_array(const T& a0, const T& a1) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1)); - return a; + return make_1d_array_with_index_type(a0, a1); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2)); - return a; + return make_1d_array_with_index_type(a0, a1, a2); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4, a5)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4, a5); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4, a5, a6)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4, a5, a6); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4, a5, a6, a7)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7); } template Array<1, T> make_1d_array(const T& a0, const T& a1, const T& a2, const T& a3, const T& a4, const T& a5, const T& a6, const T& a7, const T& a8) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4, a5, a6, a7, a8)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8); } template @@ -249,134 +418,139 @@ make_1d_array(const T& a0, const T& a8, const T& a9) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9)); - return a; + return make_1d_array_with_index_type(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9); } -template -Array -make_array(const Array& a0) +template +Array +make_array(const Array& a0) { - const Array<1, T> a = NumericVectorWithOffset(make_vector(a0)); + const Array<1, elemT, indexT> a = NumericVectorWithOffset(make_vector_with_index_type(a0)); return a; } -template -Array -make_array(const Array& a0, const Array& a1) +template +Array +make_array(const Array& a0, const Array& a1) { - const Array a = NumericVectorWithOffset, T>(make_vector(a0, a1)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1)); return a; } -template -Array -make_array(const Array& a0, const Array& a1, const Array& a2) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2) { - const Array a = NumericVectorWithOffset, T>(make_vector(a0, a1, a2)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3) { - const Array a = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4) { - const Array a = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5) { - const Array a - = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4, a5)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4, a5)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6) { - const Array a - = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4, a5, a6)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4, a5, a6)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7) +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7) { - const Array a - = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4, a5, a6, a7)); + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4, a5, a6, a7)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7, - const Array& a8) -{ - const Array a - = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4, a5, a6, a7, a8)); +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7, + const Array& a8) +{ + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4, a5, a6, a7, a8)); return a; } -template -Array -make_array(const Array& a0, - const Array& a1, - const Array& a2, - const Array& a3, - const Array& a4, - Array& a5, - const Array& a6, - const Array& a7, - const Array& a8, - const Array& a9) -{ - const Array a - = NumericVectorWithOffset, T>(make_vector(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9)); +template +Array +make_array(const Array& a0, + const Array& a1, + const Array& a2, + const Array& a3, + const Array& a4, + Array& a5, + const Array& a6, + const Array& a7, + const Array& a8, + const Array& a9) +{ + const Array a = NumericVectorWithOffset, elemT, indexT>( + make_vector_with_index_type, indexT>(a0, a1, a2, a3, a4, a5, a6, a7, a8, a9)); return a; } diff --git a/src/include/stir/numerics/norm.h b/src/include/stir/numerics/norm.h index 8c5ec51ed9..cfa44e3df2 100644 --- a/src/include/stir/numerics/norm.h +++ b/src/include/stir/numerics/norm.h @@ -95,7 +95,7 @@ inline double norm_squared(Iter begin, Iter end); /*! The l2-norm is defined as the sqrt of the sum of squares of the norms of each element in the sequence. - \see norm(const Array<1,elemT>&) for a convenience function for Array objects. + \see norm(const Array<1,elemT, indexT>&) for a convenience function for Array objects. */ template inline double norm(Iter begin, Iter end); @@ -105,16 +105,16 @@ inline double norm(Iter begin, Iter end); This returns the sqrt of the sum of the square of the absolute values of the elements of \a v1. */ -template -inline double norm(const Array<1, elemT>& v1); +template +inline double norm(const Array<1, elemT, indexT>& v1); //! square of the l2 norm of a 1D array /*! This returns the sum of the square of the absolute values of the elements of \a v1. */ -template -inline double norm_squared(const Array<1, elemT>& v1); +template +inline double norm_squared(const Array<1, elemT, indexT>& v1); //@} diff --git a/src/include/stir/numerics/norm.inl b/src/include/stir/numerics/norm.inl index 187daa589f..5760f47fae 100644 --- a/src/include/stir/numerics/norm.inl +++ b/src/include/stir/numerics/norm.inl @@ -53,16 +53,16 @@ norm(Iter begin, Iter end) return sqrt(norm_squared(begin, end)); } -template +template inline double -norm(const Array<1, elemT>& v1) +norm(const Array<1, elemT, indexT>& v1) { return norm(v1.begin(), v1.end()); } -template +template inline double -norm_squared(const Array<1, elemT>& v1) +norm_squared(const Array<1, elemT, indexT>& v1) { return norm_squared(v1.begin(), v1.end()); } diff --git a/src/include/stir/numerics/overlap_interpolate.h b/src/include/stir/numerics/overlap_interpolate.h index 8cbb34edf3..61f57718f7 100644 --- a/src/include/stir/numerics/overlap_interpolate.h +++ b/src/include/stir/numerics/overlap_interpolate.h @@ -23,13 +23,10 @@ #ifndef __stir_numerics_overlap_interpolate__H__ #define __stir_numerics_overlap_interpolate__H__ -#include "stir/common.h" +#include "stir/VectorWithOffsetFwd.h" START_NAMESPACE_STIR -template -class VectorWithOffset; - /*! \ingroup numerics @@ -40,9 +37,9 @@ class VectorWithOffset; const in_coord_iter_t in_coord_begin, const in_coord_iter_t in_coord_end, const bool only_add_to_output=false, const bool assign_rest_with_zeroes) */ -template -void overlap_interpolate(VectorWithOffset& out_data, - const VectorWithOffset& in_data, +template +void overlap_interpolate(VectorWithOffset& out_data, + const VectorWithOffset& in_data, const float zoom, const float offset, const bool assign_rest_with_zeroes = true); diff --git a/src/include/stir/recon_buildblock/ProjMatrixElemsForOneBinValue.h b/src/include/stir/recon_buildblock/ProjMatrixElemsForOneBinValue.h index a616551fee..5240987c71 100644 --- a/src/include/stir/recon_buildblock/ProjMatrixElemsForOneBinValue.h +++ b/src/include/stir/recon_buildblock/ProjMatrixElemsForOneBinValue.h @@ -24,13 +24,10 @@ #ifndef __ProjMatrixElemsForOneBinValue_H__ #define __ProjMatrixElemsForOneBinValue_H__ -#include "stir/common.h" +#include "stir/BasicCoordinateFwd.h" START_NAMESPACE_STIR -template -class BasicCoordinate; - /*! \ingroup projection \brief Stores voxel coordinates and the value of the matrix element. diff --git a/src/include/stir/recon_buildblock/SymmetryOperation.h b/src/include/stir/recon_buildblock/SymmetryOperation.h index 077874f5f0..238e91db2a 100644 --- a/src/include/stir/recon_buildblock/SymmetryOperation.h +++ b/src/include/stir/recon_buildblock/SymmetryOperation.h @@ -22,12 +22,10 @@ #ifndef __stir_recon_buildblock_SymmetryOperation_H__ #define __stir_recon_buildblock_SymmetryOperation_H__ -#include "stir/common.h" +#include "stir/BasicCoordinateFwd.h" START_NAMESPACE_STIR -template -class BasicCoordinate; class ViewSegmentNumbers; class ProjMatrixElemsForOneBin; class ProjMatrixElemsForOneDensel; diff --git a/src/include/stir/zoom.h b/src/include/stir/zoom.h index 1b8779d0ff..e1cb49254d 100644 --- a/src/include/stir/zoom.h +++ b/src/include/stir/zoom.h @@ -81,7 +81,7 @@ DiscretisedDensity\<3,float\>::get_physical_coordinates_for_indices. */ -#include "stir/common.h" +#include "stir/BasicCoordinateFwd.h" #include "stir/ZoomOptions.h" START_NAMESPACE_STIR @@ -96,8 +96,6 @@ template class PixelsOnCartesianGrid; template class CartesianCoordinate3D; -template -class BasicCoordinate; /*! \ingroup buildblock diff --git a/src/include/stir_experimental/Filter.h b/src/include/stir_experimental/Filter.h index 57fedf62c2..91016b0acc 100644 --- a/src/include/stir_experimental/Filter.h +++ b/src/include/stir_experimental/Filter.h @@ -25,9 +25,6 @@ START_NAMESPACE_STIR -template -class Array; - /*! \brief Preliminary class for 1D filtering using FFTs diff --git a/src/include/stir_experimental/fft.h b/src/include/stir_experimental/fft.h index 8f52831aaf..cea09efb04 100644 --- a/src/include/stir_experimental/fft.h +++ b/src/include/stir_experimental/fft.h @@ -25,13 +25,10 @@ #ifndef __FFT_H__ #define __FFT_H__ -#include "stir/common.h" +#include "stir/ArrayFwd.h" START_NAMESPACE_STIR -template -class Array; - //! 1-dimensional FFT /*! Replaces data by its discrete Fourier transform, if isign is input diff --git a/src/include/stir_experimental/multiply_plane_scale_factorsImageProcessor.h b/src/include/stir_experimental/multiply_plane_scale_factorsImageProcessor.h index e7df7abd1d..23447f7a47 100644 --- a/src/include/stir_experimental/multiply_plane_scale_factorsImageProcessor.h +++ b/src/include/stir_experimental/multiply_plane_scale_factorsImageProcessor.h @@ -21,13 +21,11 @@ #include "stir/RegisteredParsingObject.h" #include "stir/DataProcessor.h" #include "stir/DiscretisedDensity.h" +#include "stir/VectorWithOffsetFwd.h" #include START_NAMESPACE_STIR -template -class VectorWithOffset; - /*! \brief Simply multiplies each plane in an image with a scale factor. */ diff --git a/src/swig/stir_array.i b/src/swig/stir_array.i index 02b9fa1524..5f42d35d9d 100644 --- a/src/swig/stir_array.i +++ b/src/swig/stir_array.i @@ -115,19 +115,19 @@ namespace stir { - %template(IndexRange1D) IndexRange<1>; + %template(IndexRange1D) IndexRange<1, int>; // %template(IndexRange1DVectorWithOffset) VectorWithOffset >; - %template(IndexRange2D) IndexRange<2>; + %template(IndexRange2D) IndexRange<2, int>; //%template(IndexRange2DVectorWithOffset) VectorWithOffset >; - %template(IndexRange3D) IndexRange<3>; - %template(IndexRange4D) IndexRange<4>; + %template(IndexRange3D) IndexRange<3, int>; + %template(IndexRange4D) IndexRange<4, int>; %ADD_indexaccess(int,T,VectorWithOffset); - %template(FloatVectorWithOffset) VectorWithOffset; - %template(IntVectorWithOffset) VectorWithOffset; + %template(FloatVectorWithOffset) VectorWithOffset; + %template(IntVectorWithOffset) VectorWithOffset; // TODO need to instantiate with name? - %template (FloatNumericVectorWithOffset) NumericVectorWithOffset; + %template (FloatNumericVectorWithOffset) NumericVectorWithOffset; #ifdef SWIGPYTHON // TODO this extends ND-Arrays, but apparently not 1D Arrays (because specialised template?) @@ -264,7 +264,7 @@ namespace stir %ADD_indexaccess(%arg(const BasicCoordinate&),elemT, Array); - %template(FloatArray1D) Array<1,float>; + %template(FloatArray1D) Array<1,float, int>; // this doesn't work because of bug in swig (incorrectly parses num_dimensions) //%ADD_indexaccess(int,%arg(Array), Array); @@ -282,18 +282,18 @@ namespace stir // Todo need to instantiate with name? // TODO Swig doesn't see that Array<2,float> is derived from it anyway becuse of num_dimensions bug -%template (FloatNumericVectorWithOffset2D) stir::NumericVectorWithOffset, float>; + %template (FloatNumericVectorWithOffset2D) stir::NumericVectorWithOffset, float, int>; - %template(FloatArray2D) stir::Array<2,float>; - %apply const stir::Array<2,float>& { const stir::ArrayType<2,float>& }; - %apply stir::Array<2,float>& { stir::ArrayType<2,float>& }; - %apply stir::Array<2,float>* { stir::ArrayType<2,float>* }; + %template(FloatArray2D) stir::Array<2,float, int>; + %apply const stir::Array<2,float, int>& { const stir::ArrayType<2,float, int>& }; + %apply stir::Array<2,float, int>& { stir::ArrayType<2,float, int>& }; + %apply stir::Array<2,float, int>* { stir::ArrayType<2,float, int>* }; // TODO name - %template (FloatNumericVectorWithOffset3D) stir::NumericVectorWithOffset, float>; - %template(FloatArray3D) stir::Array<3,float>; + %template (FloatNumericVectorWithOffset3D) stir::NumericVectorWithOffset, float, int>; + %template(FloatArray3D) stir::Array<3,float, int>; #if 0 - %ADD_indexaccess(int,%arg(stir::Array<2,float>),%arg(stir::Array<3,float>)); + %ADD_indexaccess(int,%arg(stir::Array<2,float, int>),%arg(stir::Array<3,float, int>)); #endif - %template (FloatNumericVectorWithOffset4D) stir::NumericVectorWithOffset, float>; - %template(FloatArray4D) stir::Array<4,float>; + %template (FloatNumericVectorWithOffset4D) stir::NumericVectorWithOffset, float, int>; + %template(FloatArray4D) stir::Array<4,float, int>; diff --git a/src/test/CMakeLists.txt b/src/test/CMakeLists.txt index d1d69a91e8..6c045c4194 100644 --- a/src/test/CMakeLists.txt +++ b/src/test/CMakeLists.txt @@ -16,6 +16,10 @@ set(dir_INVOLVED_TEST_EXE_SOURCES ${dir}_INVOLVED_TEST_EXE_SOURCES) set(buildblock_simple_tests test_Array.cxx test_VectorWithOffset.cxx + test_NestedIterator.cxx + test_convert_array.cxx + test_IndexRange.cxx + test_coordinates.cxx ) if (NOT MINI_STIR) @@ -53,10 +57,6 @@ list (APPEND buildblock_simple_tests test_ArrayFilter.cxx test_SeparableMetzArrayFilter.cxx test_SeparableGaussianArrayFilter.cxx - test_NestedIterator.cxx - test_convert_array.cxx - test_IndexRange.cxx - test_coordinates.cxx test_filename_functions.cxx test_KeyParser.cxx test_VoxelsOnCartesianGrid.cxx diff --git a/src/test/test_Array.cxx b/src/test/test_Array.cxx index 00c0ff9014..c7819bbb3e 100644 --- a/src/test/test_Array.cxx +++ b/src/test/test_Array.cxx @@ -2,7 +2,7 @@ Copyright (C) 2000 PARAPET partners Copyright (C) 2000-2011, Hammersmith Imanet Ltd Copyright (C) 2013 Kris Thielemans - Copyright (C) 2013, 2020, 2023, 2024, 2025 University College London + Copyright (C) 2013, 2020, 2023-2026 University College London This file is part of STIR. @@ -59,6 +59,7 @@ #include #include +#include using std::ofstream; using std::ifstream; using std::plus; @@ -80,10 +81,13 @@ static_assert(!has_iterator_and_no_full_iterator_v>); namespace detail { -static Array<2, float> +template +static Array<2, float, indexT> test_make_array() { - return make_array(make_1d_array(1.F, 0.F, 0.F), make_1d_array(0.F, 1.F, 1.F), make_1d_array(0.F, -2.F, 2.F)); + return make_array(make_1d_array_with_index_type(1.F, 0.F, 0.F), + make_1d_array_with_index_type(0.F, 1.F, 1.F), + make_1d_array_with_index_type(0.F, -2.F, 2.F)); } } // namespace detail @@ -94,20 +98,21 @@ test_make_array() output.flt and output.other. Existing files with these names will be overwritten. */ +template class ArrayTests : public RunTests { private: // this function tests the next() function and compare it to using full_iterators // sadly needs to be declared in the class for VC 6.0 template - void run_tests_on_next(const Array& test) + void run_tests_on_next(const Array& test) { // exit if empty array (as do..while() loop would fail) - if (test.size() == 0) + if (test.empty()) return; - BasicCoordinate index = get_min_indices(test); - typename Array::const_full_iterator iter = test.begin_all(); + auto index = get_min_indices(test); + auto iter = test.begin_all(); do { check(*iter == test[index], "test on next(): element out of sequence?"); @@ -119,7 +124,7 @@ class ArrayTests : public RunTests // functions that runs IO tests for an array of arbitrary dimension // sadly needs to be declared in the class for VC 6.0 template - void run_IO_tests(const Array& t1) + void run_IO_tests(const Array& t1) { std::fstream os; std::fstream is; @@ -130,14 +135,14 @@ class ArrayTests : public RunTests run_IO_tests_with_file_args(ofptr, ifptr, t1); } template - void run_IO_tests_with_file_args(OFSTREAM& os, IFSTREAM& is, const Array& t1) + void run_IO_tests_with_file_args(OFSTREAM& os, IFSTREAM& is, const Array& t1) { { open_write_binary(os, "output.flt"); check(write_data(os, t1) == Succeeded::yes, "write_data could not write array"); close_file(os); } - Array t2(t1.get_index_range()); + Array t2(t1.get_index_range()); { open_read_binary(is, "output.flt"); check(read_data(is, t2) == Succeeded::yes, "read_data could not read from output.flt"); @@ -148,7 +153,7 @@ class ArrayTests : public RunTests { open_write_binary(os, "output.flt"); - const Array copy = t1; + const Array copy = t1; check(write_data(os, t1, ByteOrder::swapped) == Succeeded::yes, "write_data could not write array with swapped byte order"); check_if_equal(t1, copy, "test out with byte-swapping didn't change the array"); close_file(os); @@ -173,7 +178,7 @@ class ArrayTests : public RunTests for ostream to be able to write again, but that's not defined for FILE*. */ { - const Array copy = t1; + const Array copy = t1; cerr << "\n\tYou should now see a warning that writing failed. That's by intention.\n"; check(write_data(os, t1, ByteOrder::swapped) != Succeeded::yes, "write_data with swapped byte order should have failed"); check_if_equal(t1, copy, "test out with byte-swapping didn't change the array even with failed IO"); @@ -185,7 +190,7 @@ class ArrayTests : public RunTests template void run_IO_tests_mixed(OFSTREAM& os, IFSTREAM& is, - const Array& orig, + const Array& orig, NumericInfo output_type_info) { { @@ -209,7 +214,7 @@ class ArrayTests : public RunTests if (write_data_ok) { // only do reading test if data was written - Array data_read_back(orig.get_index_range()); + Array data_read_back(orig.get_index_range()); { open_read_binary(is, "output.other"); check(read_data(is, data_read_back) == Succeeded::yes, "read_data could not read from output.other"); @@ -220,21 +225,21 @@ class ArrayTests : public RunTests // compare with convert() { float newscale = static_cast(scale); - Array origconverted = convert_array(newscale, orig, NumericInfo()); + Array origconverted = convert_array(newscale, orig, NumericInfo()); check_if_equal(newscale, scale, "test read_data <-> convert : scale factor "); check_if_equal(origconverted, data_read_back, "test read_data <-> convert : data"); } // compare orig/scale with data_read_back { - const Array orig_scaled(orig / scale); + const Array orig_scaled(orig / scale); this->check_array_equality_with_rounding( orig_scaled, data_read_back, "test out/in: data written as other_type, read as other_type"); } // compare data written as original, but read as other_type { - Array data_read_back2(orig.get_index_range()); + Array data_read_back2(orig.get_index_range()); ifstream is; open_read_binary(is, "output.orig"); @@ -243,7 +248,7 @@ class ArrayTests : public RunTests check(read_data(is, data_read_back2, NumericInfo(), in_scale) == Succeeded::yes, "read_data could not read from output.orig"); // compare orig/in_scale with data_read_back2 - const Array orig_scaled(orig / in_scale); + const Array orig_scaled(orig / in_scale); this->check_array_equality_with_rounding( orig_scaled, data_read_back2, "test out/in: data written as original_type, read as other_type"); } @@ -255,14 +260,14 @@ class ArrayTests : public RunTests /*! we check up to .5 if output_type is integer, and up to tolerance otherwise */ template - bool check_array_equality_with_rounding(const Array& orig, - const Array& data_read_back, + bool check_array_equality_with_rounding(const Array& orig, + const Array& data_read_back, const char* const message) { NumericInfo output_type_info; bool test_failed = false; - typename Array::const_full_iterator diff_iter = orig.begin_all(); - typename Array::const_full_iterator data_read_back_iter = data_read_back.begin_all_const(); + auto diff_iter = orig.begin_all(); + auto data_read_back_iter = data_read_back.begin_all_const(); while (diff_iter != orig.end_all()) { if (output_type_info.integer_type()) @@ -303,18 +308,18 @@ vec_to_shared(std::vector& v) return sptr; } +template void -ArrayTests::run_tests() +ArrayTests::run_tests() { - cerr << "Testing Array classes\n"; { cerr << "Testing 1D stuff" << endl; { { - Array<1, int> testint(IndexRange<1>(5)); + Array<1, int, indexT> testint(IndexRange<1, indexT>(5)); testint[0] = 2; check_if_equal(testint.size(), size_t(5), "test size()"); check_if_equal(testint.size_all(), size_t(5), "test size_all()"); @@ -322,66 +327,69 @@ ArrayTests::run_tests() check_if_equal(testint[0], 4, "test assign"); } - Array<1, float> test(IndexRange<1>(10)); + Array<1, float, indexT> test(IndexRange<1, indexT>(10)); check_if_zero(test, "Array1D not initialised to 0"); test[1] = (float)10.5; - test.set_offset(-1); + const indexT offset = std::is_signed_v ? -1 : 10; + test.set_offset(offset); check_if_equal(test.size(), size_t(10), "test size() with non-zero offset"); check_if_equal(test.size_all(), size_t(10), "test size_all() with non-zero offset"); - check_if_equal(test[0], 10.5F, "test indexing of Array1D"); + check_if_equal(test[offset + 1], 10.5F, "test indexing of Array1D"); test += 1; - check_if_equal(test[0], 11.5F, "test operator+=(float)"); + check_if_equal(test[offset + 1], 11.5F, "test operator+=(float)"); check_if_equal(test.sum(), 20.5F, "test operator+=(float) and sum()"); check_if_zero(test - test, "test operator-(Array1D)"); - BasicCoordinate<1, int> c; - c[1] = 0; - check_if_equal(test[c], 11.5F, "test operator[](BasicCoordinate)"); + BasicCoordinate<1, indexT> c; + c[1] = offset + 1; + check_if_equal(test[c], 11.5F, "test operator[](BasicCoordinate) 1"); test[c] = 12.5; - check_if_equal(test[c], 12.5F, "test operator[](BasicCoordinate)"); + check_if_equal(test[c], 12.5F, "test operator[](BasicCoordinate) 2"); { - Array<1, float> ref(-1, 2); - ref[-1] = 1.F; - ref[0] = 3.F; - ref[1] = 3.14F; - Array<1, float> test = ref; + const indexT offset + = std::is_signed_v ? -1 : 10; // note: for unsigned, has to be larger than 1 for tests below to work + Array<1, float, indexT> ref(offset, offset + 3); + ref[offset] = 1.F; + ref[offset + 1] = 3.F; + ref[offset + 2] = 3.14F; + Array<1, float, indexT> test = ref; test += 1; - for (int i = ref.get_min_index(); i <= ref.get_max_index(); ++i) + for (auto i = ref.get_min_index(); i <= ref.get_max_index(); ++i) check_if_equal(test[i], ref[i] + 1, "test operator+=(float)"); test = ref; test -= 4; - for (int i = ref.get_min_index(); i <= ref.get_max_index(); ++i) + for (auto i = ref.get_min_index(); i <= ref.get_max_index(); ++i) check_if_equal(test[i], ref[i] - 4, "test operator-=(float)"); test = ref; test *= 3; - for (int i = ref.get_min_index(); i <= ref.get_max_index(); ++i) + for (auto i = ref.get_min_index(); i <= ref.get_max_index(); ++i) check_if_equal(test[i], ref[i] * 3, "test operator*=(float)"); test = ref; test /= 3; - for (int i = ref.get_min_index(); i <= ref.get_max_index(); ++i) + for (auto i = ref.get_min_index(); i <= ref.get_max_index(); ++i) check_if_equal(test[i], ref[i] / 3, "test operator/=(float)"); } { - Array<1, float> test2; + Array<1, float, indexT> test2; test2 = test * 2; - check_if_equal(2 * test[0], test2[0], "test operator*(float)"); + check_if_equal(2 * test[offset + 1], test2[offset + 1], "test operator*(float)"); } { - Array<1, float> test2 = test; - test.grow(-2, test.get_max_index()); - Array<1, float> test3 = test2 + test; - check_if_zero(test3[-2], "test growing during operator+"); + Array<1, float, indexT> test2 = test; + test.grow(offset - 1, test.get_max_index()); + Array<1, float, indexT> test3 = test2 + test; + check_if_zero(test3[offset - 1], "test growing during operator+"); } // using preallocated memory { std::vector mem(test.get_index_range().size_all()); std::copy(test.begin_all_const(), test.end_all_const(), mem.begin()); - Array<1, float> preallocated(test.get_index_range(), vec_to_shared(mem)); + Array<1, float, indexT> preallocated(test.get_index_range(), vec_to_shared(mem)); // shared_ptr mem_sptr(new float [test.get_index_range().size_all()]); // auto mem = mem_sptr.get(); // std::copy(test.begin_all_const(), test.end_all_const(), mem); @@ -419,7 +427,7 @@ ArrayTests::run_tests() { std::vector mem(test.get_index_range().size_all()); std::copy(test.begin_all_const(), test.end_all_const(), mem.begin()); - Array<1, float> test_from_mem(test.get_index_range(), reinterpret_cast(mem.data())); + Array<1, float, indexT> test_from_mem(test.get_index_range(), reinterpret_cast(mem.data())); check(test_from_mem.owns_memory_for_data(), "test preallocated with copy: should own memory"); check_if_equal(test, test_from_mem, "test construct from mem: equality"); std::copy(test.begin_all_const(), test.end_all_const(), test_from_mem.begin_all()); @@ -434,18 +442,19 @@ ArrayTests::run_tests() #if 1 { // tests on log/exp - Array<1, float> test(-3, 10); + const indexT offset = std::is_signed_v ? -3 : 20; + Array<1, float, indexT> test(offset, offset + 13); test.fill(1.F); in_place_log(test); { - Array<1, float> testeq(-3, 10); + Array<1, float, indexT> testeq(test.get_index_range()); check_if_equal(test, testeq, "test in_place_log of Array1D"); } { - for (int i = test.get_min_index(); i <= test.get_max_index(); i++) + for (auto i = test.get_min_index(); i <= test.get_max_index(); i++) test[i] = 3.5F * i + 100; } - Array<1, float> test_copy = test; + Array<1, float, indexT> test_copy = test; in_place_log(test); in_place_exp(test); @@ -457,10 +466,11 @@ ArrayTests::run_tests() { cerr << "Testing 2D stuff" << endl; { - const IndexRange<2> range(Coordinate2D(0, 0), Coordinate2D(9, 9)); - Array<2, float> test2(range); - check_if_equal(test2.size(), size_t(10), "test size()"); - check_if_equal(test2.size_all(), size_t(100), "test size_all()"); + const IndexRange<2, indexT> range(Coordinate2D(2, 2), + Coordinate2D(9, 9)); // Note: for unsigned, min needs to be >= 2 for test below + Array<2, float, indexT> test2(range); + check_if_equal(test2.size(), size_t(8), "test size()"); + check_if_equal(test2.size_all(), size_t(64), "test size_all()"); // KT 17/03/98 added check on initialisation check_if_zero(test2, "test Array<2,float> not initialised to 0"); @@ -481,9 +491,9 @@ ArrayTests::run_tests() } { - IndexRange<2> range(Coordinate2D(0, 0), Coordinate2D(3, 3)); - Array<2, float> testfp(range); - Array<2, float> t2fp(range); + IndexRange<2, indexT> range(Coordinate2D(1, 0), Coordinate2D(3, 3)); + Array<2, float, indexT> testfp(range); + Array<2, float, indexT> t2fp(range); #if 0 // KT 06/04/98 removed operator() testfp(3,2) = 3.3F; @@ -500,7 +510,7 @@ ArrayTests::run_tests() # pragma GCC diagnostic push # pragma GCC diagnostic ignored "-Waddress" // disable gcc warning that the following will always be satisfied #endif - check(dynamic_cast const*>(&t2) != 0, "test operator +(Array2D) return correct type"); + check(dynamic_cast const*>(&t2) != 0, "test operator +(Array2D) return correct type"); #if defined __GNUC__ # pragma GCC diagnostic pop #endif @@ -525,7 +535,7 @@ ArrayTests::run_tests() } { - BasicCoordinate<2, int> c; + BasicCoordinate<2, indexT> c; c[1] = 3; c[2] = 2; check_if_equal(t2[c], 5.5F, "test on operator[](BasicCoordinate)"); @@ -534,13 +544,13 @@ ArrayTests::run_tests() } // assert should break on next line (in Debug build) if uncommented - // t2[-4][3]=1.F; + // t2[t2.get_max_index() + 3][3]=1.F; // at() should throw error { bool exception_thrown = false; try { - t2.at(-4).at(3); + t2.at(t2.get_max_index() + 3).at(3); } catch (...) { @@ -549,10 +559,10 @@ ArrayTests::run_tests() check(exception_thrown, "out-of-range index should throw an exception"); } - // t2.grow_height(-5,5); - IndexRange<2> larger_range(Coordinate2D(-5, 0), Coordinate2D(5, 3)); + const indexT offset = std::is_signed_v ? -5 : 1; // Note: for unsigned, needs to be >= 1 for test below + IndexRange<2, indexT> larger_range(Coordinate2D(offset, offset), Coordinate2D(5, 3)); t2.grow(larger_range); - t2[-4][3] = 1.F; + t2[5][3] = 1.F; check_if_equal(t2[3][2], 6.F, "test on grow"); // test assignment @@ -561,11 +571,11 @@ ArrayTests::run_tests() { // copy constructor; - Array<2, float> t21(t2); + Array<2, float, indexT> t21(t2); check_if_equal(t21[3][2], 6.F, "test Array2D copy consructor element"); check_if_equal(t21, t2, "test Array2D copy constructor all elements"); // 'assignment constructor' (this simply calls copy constructor) - Array<2, float> t22 = t2; + Array<2, float, indexT> t22 = t2; check_if_equal(t22, t2, "test Array2D copy constructor"); } @@ -581,7 +591,7 @@ ArrayTests::run_tests() "check row-major order in 2D"); } // Array<2,float> preallocated(t2.get_index_range(), &mem[0], false); - Array<2, float> preallocated(t2.get_index_range(), vec_to_shared(mem)); + Array<2, float, indexT> preallocated(t2.get_index_range(), vec_to_shared(mem)); // check(!preallocated.owns_memory_for_data(), "test preallocated without copy: should not own memory"); check_if_equal(t2, preallocated, "test preallocated: equality"); std::copy(t2.begin_all_const(), t2.end_all_const(), preallocated.begin_all()); @@ -599,14 +609,14 @@ ArrayTests::run_tests() check_if_equal(*(preallocated.begin_all()), mem[0], "test preallocated: indirect buffer mod"); // test resize { - BasicCoordinate<2, int> min, max; + BasicCoordinate<2, indexT> min, max; preallocated.get_regular_range(min, max); // resizing to smaller range will keep pointing to the same memory - preallocated.resize(IndexRange<2>(min + 1, max - 1)); + preallocated.resize(IndexRange<2, indexT>(min + 1, max - 1)); std::fill(mem.begin(), mem.end(), 12345.F); check_if_equal(preallocated[min + 1], 12345.F, "test preallocated: resize smaller uses same memory"); check(!preallocated.is_contiguous(), "test preallocated: no longer contiguous after resize"); - preallocated.resize(IndexRange<2>(min - 1, max - 1)); + preallocated.resize(IndexRange<2, indexT>(min - 1, max - 1)); std::fill(mem.begin(), mem.end(), 123456.F); check_if_equal(preallocated[min + 1], 12345.F, "test preallocated: grow uses different memory"); } @@ -614,30 +624,30 @@ ArrayTests::run_tests() } // size_all with irregular range { - const IndexRange<2> range(Coordinate2D(-1, 1), Coordinate2D(1, 2)); - Array<2, float> test2(range); + const IndexRange<2, indexT> range(Coordinate2D(1, 1), Coordinate2D(3, 2)); + Array<2, float, indexT> test2(range); check(test2.is_regular(), "test is_regular() with regular"); check_if_equal(test2.size(), size_t(3), "test size() with non-zero offset"); check_if_equal(test2.size_all(), size_t(6), "test size_all() with non-zero offset"); - test2[0].resize(-1, 2); + test2[2].resize(1, 5); check(!test2.is_regular(), "test is_regular() with irregular"); check_if_equal(test2.size(), size_t(3), "test size() with irregular range"); - check_if_equal(test2.size_all(), size_t(6 + 2), "test size_all() with irregular range"); + check_if_equal(test2.size_all(), size_t(6 + 3), "test size_all() with irregular range"); } // full iterator { - IndexRange<2> range(Coordinate2D(0, 0), Coordinate2D(2, 2)); - Array<2, float> test2(range); + IndexRange<2, indexT> range(Coordinate2D(0, 0), Coordinate2D(2, 2)); + Array<2, float, indexT> test2(range); { float value = 1.2F; - for (Array<2, float>::full_iterator iter = test2.begin_all(); iter != test2.end_all();) + for (auto iter = test2.begin_all(); iter != test2.end_all();) *iter++ = value++; } { float value = 1.2F; - Array<2, float>::const_full_iterator iter = test2.begin_all_const(); - for (int i = test2.get_min_index(); i <= test2.get_max_index(); ++i) - for (int j = test2[i].get_min_index(); j <= test2[i].get_max_index(); ++j) + auto iter = test2.begin_all_const(); + for (auto i = test2.get_min_index(); i <= test2.get_max_index(); ++i) + for (auto j = test2[i].get_min_index(); j <= test2[i].get_max_index(); ++j) { check(iter != test2.end_all_const(), "test on 2D full iterator"); check_if_equal(*iter++, test2[i][j], "test on 2D full iterator vs. index"); @@ -645,17 +655,17 @@ ArrayTests::run_tests() } } - const Array<2, float> empty; + const Array<2, float, indexT> empty; check(empty.begin_all() == empty.end_all(), "test on 2D full iterator for empty range"); } // tests for next() { - const IndexRange<2> range(Coordinate2D(-1, 1), Coordinate2D(1, 2)); - Array<2, int> test(range); + const IndexRange<2, indexT> range(Coordinate2D(1, 1), Coordinate2D(3, 2)); + Array<2, int, indexT> test(range); // fill array with numbers in sequence { - Array<2, int>::full_iterator iter = test.begin_all(); - for (int i = 0; iter != test.end_all(); ++iter, ++i) + auto iter = test.begin_all(); + for (auto i = 0; iter != test.end_all(); ++iter, ++i) { *iter = i; } @@ -663,16 +673,16 @@ ArrayTests::run_tests() std::cerr << "\tTest on next() with regular array\n"; this->run_tests_on_next(test); // now do test with irregular array - test[0].resize(0, 2); - test[0][2] = 10; + test[2].resize(0, 2); + test[2][2] = 10; std::cerr << "\tTest on next() with irregular array, case 1\n"; this->run_tests_on_next(test); - test[1].resize(-2, 2); - test[1][-2] = 20; + test[1].resize(0, 2); + test[1][0] = 20; std::cerr << "\tTest on next() with irregular array, case 2\n"; this->run_tests_on_next(test); - test[-1].resize(-2, 0); - test[-1][-2] = 30; + test[1].resize(0, 0); + test[1][0] = 30; std::cerr << "\tTest on next() with irregular array, case 3\n"; this->run_tests_on_next(test); } @@ -681,48 +691,48 @@ ArrayTests::run_tests() { cerr << "Testing 3D stuff" << endl; - IndexRange<3> range(Coordinate3D(0, -1, 1), Coordinate3D(3, 3, 3)); - Array<3, float> test3(range); + IndexRange<3, indexT> range(Coordinate3D(0, 1, 1), Coordinate3D(3, 5, 3)); + Array<3, float, indexT> test3(range); check_if_equal(test3.size(), size_t(4), "test size()"); check_if_equal(test3.size_all(), size_t(60), "test size_all() with non-zero offset"); // KT 06/04/98 removed operator() #if 0 - test3(1,2,1) = (float)6.6; + test3(1,4,1) = (float)6.6; #else - test3[1][2][1] = (float)6.6; + test3[1][4][1] = (float)6.6; #endif - test3[1][0][2] = (float)7.3; - test3[1][0][1] = -1; + test3[1][2][2] = (float)7.3; + test3[1][2][1] = -1; check_if_equal(test3.sum(), 12.9F, "test on sum"); check_if_equal(test3.find_max(), 7.3F, "test on find_max"); check_if_equal(test3.find_min(), -1.F, "test on find_min"); { - Array<3, float> test3copy(test3); - BasicCoordinate<3, int> c; + Array<3, float, indexT> test3copy(test3); + BasicCoordinate<3, indexT> c; c[1] = 1; - c[2] = 0; + c[2] = 2; c[3] = 2; check_if_equal(test3[c], 7.3F, "test on operator[](BasicCoordinate)"); test3copy[c] = 8.; - check_if_equal(test3copy[1][0][2], 8.F, "test on operator[](BasicCoordinate)"); + check_if_equal(test3copy[1][2][2], 8.F, "test on operator[](BasicCoordinate)"); } - Array<3, float> test3bis(range); - test3bis[1][2][1] = (float)6.6; - test3bis[1][0][1] = (float)1.3; - Array<3, float> test3ter = test3bis; + Array<3, float, indexT> test3bis(range); + test3bis[1][4][1] = (float)6.6; + test3bis[1][2][1] = (float)1.3; + Array<3, float, indexT> test3ter = test3bis; test3ter += test3; - check_if_equal(test3ter[1][0][1], .3F, "test on operator+=(Array3D)"); + check_if_equal(test3ter[1][2][1], .3F, "test on operator+=(Array3D)"); - Array<3, float> test3quat = test3 + test3bis; + Array<3, float, indexT> test3quat = test3 + test3bis; check_if_equal(test3quat, test3ter, "test summing Array3D"); { - Array<3, float> tmp = test3 - 2; - Array<3, float> tmp2 = test3; + Array<3, float, indexT> tmp = test3 - 2; + Array<3, float, indexT> tmp2 = test3; tmp2.fill(1.F); check_if_zero(test3.sum() - 2 * tmp2.sum() - tmp.sum(), "test operator-(float)"); @@ -734,31 +744,31 @@ ArrayTests::run_tests() // size_all with irregular range { - const IndexRange<3> range(Coordinate3D(-1, 1, 4), Coordinate3D(1, 2, 6)); - Array<3, float> test(range); + const IndexRange<3, indexT> range(Coordinate3D(1, 1, 4), Coordinate3D(3, 2, 6)); + Array<3, float, indexT> test(range); check(test.is_regular(), "test is_regular() with regular"); check_if_equal(test.size(), size_t(3), "test size() with non-zero offset"); check_if_equal(test.size_all(), size_t(3 * 2 * 3), "test size_all() with non-zero offset"); - test[0][1].resize(-1, 2); + test[2][1].resize(1, 4); check(!test.is_regular(), "test is_regular() with irregular"); check_if_equal(test.size(), size_t(3), "test size() with irregular range"); check_if_equal(test.size_all(), size_t(3 * 2 * 3 + 4 - 3), "test size_all() with irregular range"); } // full iterator { - IndexRange<3> range(Coordinate3D(0, 0, 1), Coordinate3D(2, 2, 3)); - Array<3, float> test(range); + IndexRange<3, indexT> range(Coordinate3D(0, 0, 1), Coordinate3D(2, 2, 3)); + Array<3, float, indexT> test(range); { float value = 1.2F; - for (Array<3, float>::full_iterator iter = test.begin_all(); iter != test.end_all();) + for (auto iter = test.begin_all(); iter != test.end_all();) *iter++ = value++; } { float value = 1.2F; - Array<3, float>::const_full_iterator iter = test.begin_all_const(); - for (int i = test.get_min_index(); i <= test.get_max_index(); ++i) - for (int j = test[i].get_min_index(); j <= test[i].get_max_index(); ++j) - for (int k = test[i][j].get_min_index(); k <= test[i][j].get_max_index(); ++k) + auto iter = test.begin_all_const(); + for (auto i = test.get_min_index(); i <= test.get_max_index(); ++i) + for (auto j = test[i].get_min_index(); j <= test[i].get_max_index(); ++j) + for (auto k = test[i][j].get_min_index(); k <= test[i][j].get_max_index(); ++k) { check(iter != test.end_all_const(), "test on 3D full iterator"); check_if_equal(*iter++, test[i][j][k], "test on 3D full iterator vs. index"); @@ -767,13 +777,13 @@ ArrayTests::run_tests() } // test empty container { - const Array<3, float> empty; + const Array<3, float, indexT> empty; check(empty.begin_all() == empty.end_all(), "test on 3D full iterator for empty range"); } // test conversion from full_iterator to const_full_iterator { - Array<3, float>::full_iterator titer = test.begin_all(); - Array<3, float>::const_full_iterator ctiter = titer; // this should compile + auto titer = test.begin_all(); + typename Array<3, float, indexT>::const_full_iterator ctiter = titer; // this should compile check_if_equal(*ctiter, *titer, "test on full/const_full_iterator"); } } @@ -783,7 +793,7 @@ ArrayTests::run_tests() std::iota(test3.begin_all(), test3.end_all(), 1.5F); // regular { - Array<3, float> data_to_fill(test3.get_index_range()); + Array<3, float, indexT> data_to_fill(test3.get_index_range()); fill_from(data_to_fill, test3.begin_all(), test3.end_all()); check_if_equal(test3, data_to_fill, "test on 3D fill_from"); copy_to(test3, data_to_fill.begin_all()); @@ -791,8 +801,8 @@ ArrayTests::run_tests() } // irregular { - test3[0][1].resize(-1, 2); - Array<3, float> data_to_fill(test3.get_index_range()); + test3[0][1].resize(1, 2); + Array<3, float, indexT> data_to_fill(test3.get_index_range()); fill_from(data_to_fill, test3.begin_all(), test3.end_all()); check_if_equal(test3, data_to_fill, "test on 3D fill_from, irregular range"); copy_to(test3, data_to_fill.begin_all()); @@ -803,64 +813,64 @@ ArrayTests::run_tests() { cerr << "Testing 4D stuff" << endl; - const IndexRange<4> range(Coordinate4D(-3, 0, -1, 1), Coordinate4D(-2, 3, 3, 3)); - Array<4, float> test4(range); + const IndexRange<4, indexT> range(Coordinate4D(3, 0, 6, 1), Coordinate4D(4, 3, 10, 3)); + Array<4, float, indexT> test4(range); test4.fill(1.); - test4[-3][1][2][1] = (float)6.6; + test4[3][1][9][1] = (float)6.6; #if 0 - test4(-2,1,0,2) = (float)7.3; + test4(4,1,7,2) = (float)7.3; #else - test4[-2][1][0][2] = (float)7.3; + test4[4][1][7][2] = (float)7.3; #endif { float sum = test4.sum(); check_if_equal(sum, 131.9F, "test on sum()"); } - const IndexRange<4> larger_range(Coordinate4D(-3, 0, -1, 1), Coordinate4D(-1, 3, 3, 5)); + const IndexRange<4, indexT> larger_range(Coordinate4D(3, 0, 6, 1), Coordinate4D(5, 3, 10, 5)); test4.grow(larger_range); check_if_equal(test4.get_index_range(), larger_range, "test Array4D grow index range"); check_if_equal(test4.sum(), 131.9F, "test Array4D grow sum"); { - const Array<4, float> test41 = test4; + const Array<4, float, indexT> test41 = test4; check_if_equal(test4, test41, "test Array4D copy constructor"); - check_if_equal(test41[-3][1][2][1], 6.6F, "test on indexing after grow"); + check_if_equal(test41[3][1][9][1], 6.6F, "test on indexing after grow"); } { - Array<4, float> test41 = test4; - const IndexRange<4> mixed_range(Coordinate4D(-4, 1, 0, 1), Coordinate4D(-2, 3, 3, 6)); + Array<4, float, indexT> test41 = test4; + const IndexRange<4, indexT> mixed_range(Coordinate4D(2, 1, 7, 1), Coordinate4D(4, 3, 10, 6)); test41.resize(mixed_range); check_if_equal(test41.get_index_range(), mixed_range, "test Array4D resize index range"); - check_if_equal(test41[-3][1][2][1], 6.6F, "test on indexing after resize"); + check_if_equal(test41[3][1][9][1], 6.6F, "test on indexing after resize"); } { - BasicCoordinate<4, int> c; - c[1] = -2; + BasicCoordinate<4, indexT> c; + c[1] = 4; c[2] = 1; - c[3] = 0; + c[3] = 7; c[4] = 2; check_if_equal(test4[c], 7.3F, "test on operator[](BasicCoordinate)"); test4[c] = 1.; check_if_equal(test4[c], 1.F, "test on operator[](BasicCoordinate)"); } { - Array<4, float> test4bis(range); - test4bis[-2][1][2][1] = (float)6.6; - test4bis[-3][1][0][1] = (float)1.3; - Array<4, float> test4ter = test4bis; + Array<4, float, indexT> test4bis(range); + test4bis[4][1][9][1] = (float)6.6; + test4bis[3][1][7][1] = (float)1.3; + Array<4, float, indexT> test4ter = test4bis; test4ter += test4; - check_if_equal(test4ter[-3][1][0][1], 2.3F, "test on operator+=(Array4D)"); + check_if_equal(test4ter[3][1][7][1], 2.3F, "test on operator+=(Array4D)"); check(test4ter.get_index_range() == larger_range, "test range for operator+=(Array4D) with grow"); // Note that test4 is bigger in size than test4bis. - Array<4, float> test4quat = test4bis + test4; + Array<4, float, indexT> test4quat = test4bis + test4; check_if_equal(test4quat, test4ter, "test summing Array4D with grow"); check(test4quat.get_index_range() == larger_range, "test range for operator+=(Array4D)"); } // test on scalar multiplication, division { - Array<4, float> test4bis = test4; + Array<4, float, indexT> test4bis = test4; test4bis *= 6.F; check_if_equal(test4bis.sum(), test4.sum() * 6, "test operator *=(float)"); test4bis /= 5.F; @@ -869,30 +879,30 @@ ArrayTests::run_tests() // test on element-wise multiplication, division { - Array<4, float> test4bis(range); + Array<4, float, indexT> test4bis(range); { - for (int i = test4bis.get_min_index(); i <= test4bis.get_max_index(); i++) + for (auto i = test4bis.get_min_index(); i <= test4bis.get_max_index(); i++) test4bis[i].fill(i + 10.F); } // save for comparison later on - Array<4, float> test4ter = test4bis; + Array<4, float, indexT> test4ter = test4bis; // Note that test4 is bigger than test4bis, so it will grow with the *= // new elements in test4bis will remain 0 because we're using multiplication - test4[-1].fill(666); + test4[5].fill(666); test4bis *= test4; - check_if_zero(test4bis[-1], "test operator *=(Array4D) grows ok"); + check_if_zero(test4bis[5], "test operator *=(Array4D) grows ok"); check(test4.get_index_range() == test4bis.get_index_range(), "test operator *=(Array4D) grows ok: range"); // compute the new sum. { float sum_check = 0; - for (int i = test4.get_min_index(); i <= -2; i++) + for (auto i = test4.get_min_index(); i <= indexT(4); i++) // note: up to 4, as that was the original size sum_check += test4[i].sum() * (i + 10.F); check_if_equal(test4bis.sum(), sum_check, "test operator *=(Array4D)"); } // divide test4, but add a tiny number to avoid division by zero - const Array<4, float> test4quat = test4bis / (test4 + .00001F); + const Array<4, float, indexT> test4quat = test4bis / (test4 + .00001F); test4ter.grow(test4.get_index_range()); check_if_equal(test4ter, test4quat, "test operator /(Array4D)"); } @@ -900,8 +910,8 @@ ArrayTests::run_tests() // test operator+(float) { // KT 31/01/2000 new - Array<4, float> tmp = test4 + 2; - Array<4, float> tmp2 = test4; + Array<4, float, indexT> tmp = test4 + 2; + Array<4, float, indexT> tmp2 = test4; tmp2.fill(1.F); // KT 20/12/2001 made check_if_zero compare relative to 1 by dividing @@ -910,8 +920,8 @@ ArrayTests::run_tests() // test axpby { - Array<4, float> tmp(test4.get_index_range()); - Array<4, float> tmp2(test4 + 2); + Array<4, float, indexT> tmp(test4.get_index_range()); + Array<4, float, indexT> tmp2(test4 + 2); #if defined __GNUC__ # pragma GCC diagnostic push # pragma GCC diagnostic ignored "-Wdeprecated-declarations" @@ -920,16 +930,16 @@ ArrayTests::run_tests() #if defined __GNUC__ # pragma GCC diagnostic pop #endif - const Array<4, float> by_hand = test4 * 2.F + (test4 + 2) * 3.3F; + const Array<4, float, indexT> by_hand = test4 * 2.F + (test4 + 2) * 3.3F; check_if_equal(tmp, by_hand, "test axpby (Array4D)"); } // test xapyb, a and b scalar { - Array<4, float> tmp(test4.get_index_range()); + Array<4, float, indexT> tmp(test4.get_index_range()); tmp.xapyb(test4, 2.F, test4 + 2, 3.3F); - const Array<4, float> by_hand = test4 * 2.F + (test4 + 2) * 3.3F; + const Array<4, float, indexT> by_hand = test4 * 2.F + (test4 + 2) * 3.3F; check_if_equal(tmp, by_hand, "test xapyb scalar (Array4D)"); tmp = test4; @@ -939,10 +949,10 @@ ArrayTests::run_tests() // test xapyb, a and b vector { - Array<4, float> tmp(test4.get_index_range()); + Array<4, float, indexT> tmp(test4.get_index_range()); tmp.xapyb(test4, test4 + 4, test4 + 2, test4 + 6); - const Array<4, float> by_hand = test4 * (test4 + 4) + (test4 + 2) * (test4 + 6); + const Array<4, float, indexT> by_hand = test4 * (test4 + 4) + (test4 + 2) * (test4 + 6); check_if_equal(tmp, by_hand, "test xapyb vector (Array4D)"); tmp = test4; @@ -951,13 +961,13 @@ ArrayTests::run_tests() } { - typedef NumericVectorWithOffset, float> NVecArr; - typedef NVecArr::iterator NVecArrIter; - NVecArr tmp(-1, 2); + typedef NumericVectorWithOffset, float, indexT> NVecArr; + typedef typename NVecArr::iterator NVecArrIter; + NVecArr tmp(5, 8); - NVecArr x(-1, 2); - NVecArr y(-1, 2); - NVecArr by_hand(-1, 2); + NVecArr x(5, 8); + NVecArr y(5, 8); + NVecArr by_hand(5, 8); NVecArrIter iter_tmp = tmp.begin(); NVecArrIter iter_x = x.begin(); @@ -984,15 +994,15 @@ ArrayTests::run_tests() check_if_equal(x, by_hand, "test sapyb scalar (NumericVectorWithOffset)"); } { - typedef NumericVectorWithOffset, float> NVecArr; - typedef NVecArr::iterator NVecArrIter; - NVecArr tmp(-1, 2); + typedef NumericVectorWithOffset, float, indexT> NVecArr; + typedef typename NVecArr::iterator NVecArrIter; + NVecArr tmp(5, 8); - NVecArr x(-1, 2); - NVecArr y(-1, 2); - NVecArr a(-1, 2); - NVecArr b(-1, 2); - NVecArr by_hand(-1, 2); + NVecArr x(5, 8); + NVecArr y(5, 8); + NVecArr a(5, 8); + NVecArr b(5, 8); + NVecArr by_hand(5, 8); NVecArrIter iter_tmp = tmp.begin(); NVecArrIter iter_x = x.begin(); @@ -1029,17 +1039,17 @@ ArrayTests::run_tests() #if 1 { cerr << "Testing 1D float IO" << endl; - Array<1, float> t1(IndexRange<1>(-1, 10)); - for (int i = -1; i <= 10; i++) + Array<1, float, indexT> t1(IndexRange<1, indexT>(1, 10)); + for (auto i = t1.get_min_index(); i <= t1.get_max_index(); i++) t1[i] = static_cast(sin(i * _PI / 15.)); run_IO_tests(t1); } { cerr << "Testing 2D double IO" << endl; - IndexRange<2> range(Coordinate2D(-1, 11), Coordinate2D(10, 20)); - Array<2, double> t1(range); - for (int i = -1; i <= 10; i++) - for (int j = 11; j <= 20; j++) + IndexRange<2, indexT> range(Coordinate2D(1, 11), Coordinate2D(10, 20)); + Array<2, double, indexT> t1(range); + for (auto i = t1.get_min_index(); i <= t1.get_max_index(); i++) + for (auto j = t1[i].get_min_index(); j <= t1[i].get_max_index(); j++) t1[i][j] = static_cast(sin(i * j * _PI / 15.)); run_IO_tests(t1); } @@ -1048,11 +1058,11 @@ ArrayTests::run_tests() // construct test array which has rows of very different magnitudes, // numbers in last rows do not fit into short integers - IndexRange<3> range(Coordinate3D(-1, 11, 21), Coordinate3D(10, 20, 30)); - Array<3, float> t1(range); - for (int i = -1; i <= 10; i++) - for (int j = 11; j <= 20; j++) - for (int k = 21; k <= 30; k++) + IndexRange<3, indexT> range(Coordinate3D(1, 11, 21), Coordinate3D(10, 20, 30)); + Array<3, float, indexT> t1(range); + for (auto i = t1.get_min_index(); i <= t1.get_max_index(); i++) + for (auto j = t1[i].get_min_index(); j <= t1[i].get_max_index(); j++) + for (auto k = t1[i][j].get_min_index(); k <= t1[i][j].get_max_index(); k++) t1[i][j][k] = static_cast(20000. * k * sin(i * j * k * _PI / 3000.)); run_IO_tests(t1); } @@ -1060,28 +1070,31 @@ ArrayTests::run_tests() { cerr << "Testing make_array" << endl; - const Array<2, float> arr1 - = make_array(make_1d_array(1.F, 0.F, 0.F), make_1d_array(0.F, 1.F, 1.F), make_1d_array(0.F, -2.F, 2.F)); + const Array<2, float, indexT> arr1 = make_array(make_1d_array_with_index_type(1.F, 0.F, 0.F), + make_1d_array_with_index_type(0.F, 1.F, 1.F), + make_1d_array_with_index_type(0.F, -2.F, 2.F)); - const Array<2, float> arr2( - make_array(make_1d_array(1.F, 0.F, 0.F), make_1d_array(0.F, 1.F, 1.F), make_1d_array(0.F, -2.F, 2.F))); + const Array<2, float, indexT> arr2(make_array(make_1d_array_with_index_type(1.F, 0.F, 0.F), + make_1d_array_with_index_type(0.F, 1.F, 1.F), + make_1d_array_with_index_type(0.F, -2.F, 2.F))); - const Array<2, float> arr3 = detail::test_make_array(); - const Array<2, float> arr4(detail::test_make_array()); + const Array<2, float, indexT> arr3 = detail::test_make_array(); + const Array<2, float, indexT> arr4(detail::test_make_array()); check_if_equal(arr1[2][1], -2.F, "make_array element comparison"); check_if_equal(arr1, arr2, "make_array inline assignment vs constructor"); check_if_equal(arr1, arr3, "make_array inline vs function with assignment"); check_if_equal(arr1, arr4, "make_array inline constructor from function"); } +#if 0 std::cerr << "timings\n"; { HighResWallClockTimer t; - IndexRange<4> range(Coordinate4D(20, 100, 400, 600)); + IndexRange<4, indexT> range(Coordinate4D(20, 100, 400, 600)); t.start(); double create_duration; { - Array<4, int> a1(range); + Array<4, int, indexT> a1(range); t.stop(); std::cerr << "creation of non-contiguous 4D Array " << t.value() * 1000 << "ms\n"; create_duration = t.value(); @@ -1101,7 +1114,7 @@ ArrayTests::run_tests() std::cerr << "vector creation " << t.value() * 1000 << "ms\n"; t.start(); // Array<4,int> a1(range, v.data(), false); - Array<4, int> a1(range, vec_to_shared(v)); + Array<4, int, indexT> a1(range, vec_to_shared(v)); t.stop(); // check(!a1.owns_memory_for_data(), "test preallocated without copy: should not own memory"); create_duration = t.value(); @@ -1111,6 +1124,7 @@ ArrayTests::run_tests() t.stop(); std::cerr << "deletion " << (t.value() - create_duration) * 1000 << " ms\n"; } +#endif } END_NAMESPACE_STIR @@ -1120,7 +1134,34 @@ USING_NAMESPACE_STIR int main() { - ArrayTests tests; - tests.run_tests(); - return tests.main_return_value(); + { + cerr << "Testing Array classes with int as indices\n"; + cerr << "=========================================\n"; + ArrayTests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing Array classes with short int as indices\n"; + ArrayTests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing Array classes with long long as indices\n"; + ArrayTests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing Array classes with unsigned int as indices\n"; + ArrayTests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + return EXIT_SUCCESS; } diff --git a/src/test/test_convert_array.cxx b/src/test/test_convert_array.cxx index e444402ee2..abe4af6db7 100644 --- a/src/test/test_convert_array.cxx +++ b/src/test/test_convert_array.cxx @@ -3,6 +3,7 @@ /* Copyright (C) 2000 PARAPET partners Copyright (C) 2000- 2008, Hammersmith Imanet Ltd + Copyright (C) 2025 University College London This file is part of STIR. SPDX-License-Identifier: Apache-2.0 AND License-ref-PARAPET-license @@ -24,39 +25,42 @@ #include "stir/NumericInfo.h" #include "stir/IndexRange3D.h" #include "stir/RunTests.h" +#include "stir/stream.h" #include #include -#include +#include +#include +#include using std::cerr; using std::endl; START_NAMESPACE_STIR //! tests convert_array functionality +template class convert_array_Tests : public RunTests { public: void run_tests() override; }; +template void -convert_array_Tests::run_tests() +convert_array_Tests::run_tests() { - cerr << "Test program for 'convert_array'." << endl << "Everything is fine when there is no output below." << endl; - // 1D { - Array<1, float> tf1(1, 20); + Array<1, float, indexT> tf1(1, 20); tf1.fill(100.F); - Array<1, short> ti1(1, 20); + Array<1, short, indexT> ti1(1, 20); ti1.fill(100); { // float -> short with a preferred scale factor float scale_factor = float(1); - Array<1, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<1, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(scale_factor == float(1), "test convert_array float->short 1D"); check_if_equal(ti1, ti2, "test convert_array float->short 1D"); @@ -65,10 +69,10 @@ convert_array_Tests::run_tests() { // float -> short with automatic scale factor float scale_factor = 0; - Array<1, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<1, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(fabs(NumericInfo().max_value() / 1.01 / ti2[1] - 1) < 1E-4); - for (int i = 1; i <= 20; i++) + for (indexT i = 1; i <= 20; i++) ti2[i] = short(double(ti2[i]) * scale_factor); check(ti1 == ti2); } @@ -77,22 +81,22 @@ convert_array_Tests::run_tests() { // float -> short with a preferred scale factor that needs to be adjusted float scale_factor = 1; - Array<1, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<1, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(fabs(NumericInfo().max_value() / 1.01 / ti2[1] - 1) < 1E-4); - for (int i = 1; i <= 20; i++) + for (indexT i = 1; i <= 20; i++) check(fabs(double(ti2[i]) * scale_factor / tf1[i] - 1) < 1E-4); } { // short -> float with a scale factor = 1 float scale_factor = 1; - Array<1, float> tf2 = convert_array(scale_factor, ti1, NumericInfo()); - Array<1, short> ti2(1, 20); + Array<1, float, indexT> tf2 = convert_array(scale_factor, ti1, NumericInfo()); + Array<1, short, indexT> ti2(1, 20); check(scale_factor == float(1)); check(tf2[1] == 100.F); - for (int i = 1; i <= 20; i++) + for (indexT i = 1; i <= 20; i++) ti2[i] = short(double(tf2[i]) * scale_factor); check(ti1 == ti2); } @@ -100,12 +104,12 @@ convert_array_Tests::run_tests() { // short -> float with a preferred scale factor = .01 float scale_factor = .01F; - Array<1, float> tf2 = convert_array(scale_factor, ti1, NumericInfo()); - Array<1, short> ti2(1, 20); + Array<1, float, indexT> tf2 = convert_array(scale_factor, ti1, NumericInfo()); + Array<1, short, indexT> ti2(1, 20); check(scale_factor == float(.01)); // TODO double->short - for (int i = 1; i <= 20; i++) + for (indexT i = 1; i <= 20; i++) ti2[i] = short(double(tf2[i]) * scale_factor + 0.5); check(ti1 == ti2); } @@ -115,19 +119,19 @@ convert_array_Tests::run_tests() { // positive float -> unsigned short with a preferred scale factor float scale_factor = 1; - Array<1, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<1, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(scale_factor == float(1)); check(ti1 == ti2); } { - Array<1, unsigned short> ti3(1, 20); + Array<1, unsigned short, indexT> ti3(1, 20); ti3.fill(0); // negative float -> unsigned short with a preferred scale factor float scale_factor = 1; - Array<1, unsigned short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<1, unsigned short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(scale_factor == float(1)); check(ti3 == ti2); @@ -136,16 +140,17 @@ convert_array_Tests::run_tests() // 3D { - Array<3, float> tf1(IndexRange3D(1, 30, 1, 182, -2, 182)); + const auto min_indices = make_coordinate(1, 1, std::is_signed_v ? -2 : 0); + Array<3, float, indexT> tf1(IndexRange<3, indexT>(min_indices, make_coordinate(30, 182, 182))); tf1.fill(100.F); - Array<3, short> ti1(tf1.get_index_range()); + Array<3, short, indexT> ti1(tf1.get_index_range()); ti1.fill(100); { // float -> short with a preferred scale factor float scale_factor = float(1); - Array<3, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<3, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); check(scale_factor == float(1)); check(ti1 == ti2); @@ -154,11 +159,11 @@ convert_array_Tests::run_tests() { // float -> short with automatic scale factor float scale_factor = 0; - Array<3, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<3, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); #ifndef DO_TIMING_ONLY check(fabs(NumericInfo().max_value() / 1.01 / (*ti2.begin_all()) - 1) < 1E-4); - const Array<3, short>::full_iterator iter_end = ti2.end_all(); - for (Array<3, short>::full_iterator iter = ti2.begin_all(); iter != iter_end; ++iter) + const auto iter_end = ti2.end_all(); + for (auto iter = ti2.begin_all(); iter != iter_end; ++iter) *iter = short(double((*iter)) * scale_factor); check(ti1 == ti2); #endif @@ -168,13 +173,13 @@ convert_array_Tests::run_tests() { // float -> short with a preferred scale factor that needs to be adjusted float scale_factor = 1; - Array<3, short> ti2 = convert_array(scale_factor, tf1, NumericInfo()); + Array<3, short, indexT> ti2 = convert_array(scale_factor, tf1, NumericInfo()); #ifndef DO_TIMING_ONLY check(fabs(NumericInfo().max_value() / 1.01 / (*ti2.begin_all()) - 1) < 1E-4); - Array<3, short>::full_iterator iter_ti2 = ti2.begin_all(); - const Array<3, short>::full_iterator iter_ti2_end = ti2.end_all(); - Array<3, float>::full_iterator iter_tf1 = tf1.begin_all(); + auto iter_ti2 = ti2.begin_all(); + const auto iter_ti2_end = ti2.end_all(); + auto iter_tf1 = tf1.begin_all(); for (; iter_ti2 != iter_ti2_end; ++iter_ti2, ++iter_tf1) check(fabs(double(*iter_ti2) * scale_factor / *iter_tf1 - 1) < 1E-4); #endif @@ -216,7 +221,35 @@ USING_NAMESPACE_STIR int main() { - convert_array_Tests tests; - tests.run_tests(); - return tests.main_return_value(); + { + cerr << "Testing convert_array with int as indices\n"; + cerr << "=========================================\n"; + convert_array_Tests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing convert_array with short int as indices\n"; + cerr << "=========================================\n"; + convert_array_Tests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing convert_array with long long as indices\n"; + cerr << "=========================================\n"; + convert_array_Tests tests; + tests.run_tests(); + if (!tests.is_everything_ok()) + return tests.main_return_value(); + } + { + cerr << "Testing convert_array with unsigned int as indices\n"; + cerr << "=========================================\n"; + convert_array_Tests tests; + tests.run_tests(); + return tests.main_return_value(); + } }