From 93216c6c10a85b7cb00bb258ce22e5c34c7689f1 Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Fri, 14 Nov 2025 08:29:28 +0000 Subject: [PATCH 1/5] add indexT template argument to Array, VectorWithOffset etc Adding an indexT (defaulting to int) to various classes related to arrays, including BasicCoordinate, IndexRange, VectorWithOffset and Array. This will allow to use larger (or smaller) ranges for the indices. Most code using these classes has not been changed, therefore still using "int". This was a bit more work than anticipated, as I had to moved the forward declarations to separate files, such that - they are consistent - there is only one place where the default is defined (as required by C++) --- src/buildblock/CMakeLists.txt | 1 - src/buildblock/IndexRange.cxx | 88 --- src/buildblock/overlap_interpolate.cxx | 36 +- src/include/stir/Array.h | 106 ++-- src/include/stir/Array.inl | 504 +++++++++--------- .../stir/ArrayFilter1DUsingConvolution.h | 3 - ...yFilter1DUsingConvolutionSymmetricKernel.h | 3 - .../stir/ArrayFilter2DUsingConvolution.h | 4 - .../stir/ArrayFilter3DUsingConvolution.h | 3 - src/include/stir/ArrayFunction.h | 84 +-- src/include/stir/ArrayFunction.inl | 189 ++++--- src/include/stir/ArrayFunctionObject.h | 17 +- src/include/stir/ArrayFwd.h | 15 +- src/include/stir/BasicCoordinate.h | 2 +- src/include/stir/BasicCoordinateFwd.h | 26 + src/include/stir/IO/interfile.h | 6 - src/include/stir/IO/read_data.h | 13 +- src/include/stir/IO/read_data.inl | 27 +- src/include/stir/IO/read_data_1d.h | 10 +- src/include/stir/IO/read_data_1d.inl | 8 +- src/include/stir/IO/write_data.h | 16 +- src/include/stir/IO/write_data.inl | 27 +- src/include/stir/IO/write_data_1d.h | 14 +- src/include/stir/IO/write_data_1d.inl | 26 +- src/include/stir/IndexRange.h | 54 +- src/include/stir/IndexRange.inl | 165 ++++-- src/include/stir/IndexRangeFwd.h | 26 + src/include/stir/NumericVectorWithOffset.h | 8 +- src/include/stir/NumericVectorWithOffset.inl | 162 +++--- src/include/stir/RunTests.h | 16 +- src/include/stir/VectorWithOffset.h | 64 +-- src/include/stir/VectorWithOffset.inl | 413 +++++++------- src/include/stir/VectorWithOffsetFwd.h | 26 + src/include/stir/array_index_functions.h | 21 +- src/include/stir/array_index_functions.inl | 52 +- src/include/stir/centre_of_gravity.h | 4 - src/include/stir/convert_array.h | 24 +- src/include/stir/convert_array.inl | 18 +- .../stir/evaluation/compute_ROI_values.h | 3 +- src/include/stir/interpolate_projdata.h | 3 +- src/include/stir/make_array.h | 146 ++--- src/include/stir/make_array.inl | 488 +++++++++++------ src/include/stir/numerics/norm.h | 10 +- src/include/stir/numerics/norm.inl | 8 +- .../stir/numerics/overlap_interpolate.h | 11 +- .../ProjMatrixElemsForOneBinValue.h | 5 +- .../stir/recon_buildblock/SymmetryOperation.h | 4 +- src/include/stir/zoom.h | 4 +- src/include/stir_experimental/Filter.h | 3 - src/include/stir_experimental/fft.h | 5 +- ...ltiply_plane_scale_factorsImageProcessor.h | 4 +- src/test/CMakeLists.txt | 8 +- src/test/test_Array.cxx | 320 ++++++----- src/test/test_convert_array.cxx | 99 ++-- 54 files changed, 1843 insertions(+), 1559 deletions(-) delete mode 100644 src/buildblock/IndexRange.cxx create mode 100644 src/include/stir/BasicCoordinateFwd.h create mode 100644 src/include/stir/IndexRangeFwd.h create mode 100644 src/include/stir/VectorWithOffsetFwd.h 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/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..531dfc250f 100644 --- a/src/include/stir/Array.h +++ b/src/include/stir/Array.h @@ -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,7 @@ 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; //! return the total number of elements in this array inline size_t size_all() const; @@ -215,10 +215,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 +254,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 +413,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 +424,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 +438,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 +474,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 +537,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 +603,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 +632,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..0efc9944b9 100644 --- a/src/include/stir/Array.inl +++ b/src/include/stir/Array.inl @@ -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,20 +57,20 @@ 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) { base_type::resize(range.get_min_index(), range.get_max_index()); auto iter = this->begin(); @@ -83,9 +83,9 @@ 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) { base_type::resize(range.get_min_index(), range.get_max_index()); auto iter = this->begin(); @@ -98,21 +98,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 +125,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 +144,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 +153,8 @@ Array::Array(const base_type& t) } #endif -template -Array::~Array() +template +Array::~Array() { if (this->_allocated_full_data_ptr) { @@ -163,54 +163,54 @@ 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()) { @@ -221,9 +221,9 @@ Array::begin_all() 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()) { @@ -234,32 +234,32 @@ Array::begin_all_const() const 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,7 +268,7 @@ 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; } @@ -287,9 +287,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 +306,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 +324,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 +340,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 +360,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 +377,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 +395,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 +408,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 +421,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 +434,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 +558,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,37 +579,37 @@ 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); @@ -622,165 +622,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++) + for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) assign(this->num[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 +790,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 +807,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 +815,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 +832,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 +849,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 +873,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..d2d9654a7a 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 @@ -79,16 +80,16 @@ class IndexRange : public VectorWithOffset> inline IndexRange(); //! Make an IndexRange from the base type - inline IndexRange(const VectorWithOffset>& range); + inline IndexRange(const base_type& range); //! 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 +97,22 @@ 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; //! 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 +128,38 @@ 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; //! return the total number of elements in this range inline size_t size_all() 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; }; END_NAMESPACE_STIR diff --git a/src/include/stir/IndexRange.inl b/src/include/stir/IndexRange.inl index 8e02a65a26..ff63c6efd6 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,130 @@ 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 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) 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::is_regular() const +IndexRange::get_regular_range(BasicCoordinate& min, + BasicCoordinate& max) const +{ + // check if empty range + if (base_type::begin() == base_type::end()) + { + constexpr bool signed_type = std::is_signed_v; + std::fill(min.begin(), min.end(), signed_type ? 0 : 1); + std::fill(max.begin(), max.end(), signed_type ? -1 : 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 { switch (is_regular_range) { @@ -106,8 +163,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 +178,84 @@ IndexRange::is_regular() const 1D version ***************************************/ -IndexRange<1>::IndexRange() +template +IndexRange<1, indexT>::IndexRange() : min(0), max(0) {} -IndexRange<1>::IndexRange(const int min_v, const int max_v) +template +IndexRange<1, indexT>::IndexRange(const indexT min_v, const indexT max_v) : min(min_v), max(max_v) {} -IndexRange<1>::IndexRange(const BasicCoordinate<1, int>& min_v, const BasicCoordinate<1, int>& max_v) +template +IndexRange<1, indexT>::IndexRange(const BasicCoordinate<1, indexT>& min_v, const BasicCoordinate<1, indexT>& max_v) : min(min_v[1]), max(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 - min + 1); } +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,8 +263,9 @@ 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) { 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/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..80da41416f 100644 --- a/src/include/stir/VectorWithOffset.h +++ b/src/include/stir/VectorWithOffset.h @@ -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,20 +244,20 @@ 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; @@ -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..b91c1fcdb6 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/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..93f51b1719 100644 --- a/src/test/test_Array.cxx +++ b/src/test/test_Array.cxx @@ -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) 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,7 +327,7 @@ 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; @@ -335,45 +340,45 @@ ArrayTests::run_tests() check_if_equal(test.sum(), 20.5F, "test operator+=(float) and sum()"); check_if_zero(test - test, "test operator-(Array1D)"); - BasicCoordinate<1, int> c; + BasicCoordinate<1, indexT> c; c[1] = 0; check_if_equal(test[c], 11.5F, "test operator[](BasicCoordinate)"); test[c] = 12.5; check_if_equal(test[c], 12.5F, "test operator[](BasicCoordinate)"); { - Array<1, float> ref(-1, 2); + Array<1, float, indexT> ref(-1, 2); ref[-1] = 1.F; ref[0] = 3.F; ref[1] = 3.14F; - Array<1, float> test = ref; + 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)"); } { - Array<1, float> test2 = test; + Array<1, float, indexT> test2 = test; test.grow(-2, test.get_max_index()); - Array<1, float> test3 = test2 + test; + Array<1, float, indexT> test3 = test2 + test; check_if_zero(test3[-2], "test growing during operator+"); } @@ -381,7 +386,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> 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 +424,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 +439,18 @@ ArrayTests::run_tests() #if 1 { // tests on log/exp - Array<1, float> test(-3, 10); + Array<1, float, indexT> test(-3, 10); test.fill(1.F); in_place_log(test); { - Array<1, float> testeq(-3, 10); + Array<1, float, indexT> testeq(-3, 10); 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,8 +462,8 @@ ArrayTests::run_tests() { cerr << "Testing 2D stuff" << endl; { - const IndexRange<2> range(Coordinate2D(0, 0), Coordinate2D(9, 9)); - Array<2, float> test2(range); + const IndexRange<2, indexT> range(Coordinate2D(0, 0), Coordinate2D(9, 9)); + Array<2, float, indexT> test2(range); check_if_equal(test2.size(), size_t(10), "test size()"); check_if_equal(test2.size_all(), size_t(100), "test size_all()"); // KT 17/03/98 added check on initialisation @@ -481,9 +486,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(0, 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 +505,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 +530,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)"); @@ -550,7 +555,7 @@ ArrayTests::run_tests() } // t2.grow_height(-5,5); - IndexRange<2> larger_range(Coordinate2D(-5, 0), Coordinate2D(5, 3)); + IndexRange<2, indexT> larger_range(Coordinate2D(-5, 0), Coordinate2D(5, 3)); t2.grow(larger_range); t2[-4][3] = 1.F; check_if_equal(t2[3][2], 6.F, "test on grow"); @@ -561,11 +566,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 +586,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 +604,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,8 +619,8 @@ 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(1, 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"); @@ -626,18 +631,18 @@ ArrayTests::run_tests() } // 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 +650,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(1, 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; } @@ -681,8 +686,8 @@ 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, 3, 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() @@ -699,8 +704,8 @@ ArrayTests::run_tests() 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[3] = 2; @@ -709,20 +714,20 @@ ArrayTests::run_tests() check_if_equal(test3copy[1][0][2], 8.F, "test on operator[](BasicCoordinate)"); } - Array<3, float> test3bis(range); + Array<3, float, indexT> 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> test3ter = test3bis; test3ter += test3; check_if_equal(test3ter[1][0][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,8 +739,8 @@ 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(1, 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"); @@ -746,19 +751,19 @@ ArrayTests::run_tests() } // 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 +772,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 +788,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()); @@ -792,7 +797,7 @@ ArrayTests::run_tests() // irregular { test3[0][1].resize(-1, 2); - 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, irregular range"); copy_to(test3, data_to_fill.begin_all()); @@ -803,8 +808,8 @@ 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, -1, 1), Coordinate4D(-2, 3, 3, 3)); + Array<4, float, indexT> test4(range); test4.fill(1.); test4[-3][1][2][1] = (float)6.6; #if 0 @@ -816,24 +821,24 @@ ArrayTests::run_tests() 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, -1, 1), Coordinate4D(-1, 3, 3, 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"); } { - 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(-4, 1, 0, 1), Coordinate4D(-2, 3, 3, 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"); } { - BasicCoordinate<4, int> c; + BasicCoordinate<4, indexT> c; c[1] = -2; c[2] = 1; c[3] = 0; @@ -843,24 +848,24 @@ ArrayTests::run_tests() check_if_equal(test4[c], 1.F, "test on operator[](BasicCoordinate)"); } { - Array<4, float> test4bis(range); + Array<4, float, indexT> 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> test4ter = test4bis; test4ter += test4; check_if_equal(test4ter[-3][1][0][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,13 +874,13 @@ 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 @@ -887,12 +892,12 @@ ArrayTests::run_tests() // 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(-2); i++) // note: up to -2, 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 +905,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 +915,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 +925,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 +944,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,8 +956,8 @@ ArrayTests::run_tests() } { - typedef NumericVectorWithOffset, float> NVecArr; - typedef NVecArr::iterator NVecArrIter; + typedef NumericVectorWithOffset, float, indexT> NVecArr; + typedef typename NVecArr::iterator NVecArrIter; NVecArr tmp(-1, 2); NVecArr x(-1, 2); @@ -984,8 +989,8 @@ ArrayTests::run_tests() check_if_equal(x, by_hand, "test sapyb scalar (NumericVectorWithOffset)"); } { - typedef NumericVectorWithOffset, float> NVecArr; - typedef NVecArr::iterator NVecArrIter; + typedef NumericVectorWithOffset, float, indexT> NVecArr; + typedef typename NVecArr::iterator NVecArrIter; NVecArr tmp(-1, 2); NVecArr x(-1, 2); @@ -1029,17 +1034,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 +1053,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 +1065,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 +1109,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 +1119,7 @@ ArrayTests::run_tests() t.stop(); std::cerr << "deletion " << (t.value() - create_duration) * 1000 << " ms\n"; } +#endif } END_NAMESPACE_STIR @@ -1120,7 +1129,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(); + } } From bf6f69d85cfaf9c9f497bb6cf5c59cccb938a247 Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Fri, 14 Nov 2025 08:34:20 +0000 Subject: [PATCH 2/5] use Array<1,float,long long> for ProjDataInMemory buffer This removes the limitation on the number of elements in the proj-data. Fixes https://github.com/UCL/STIR/issues/1505 --- src/buildblock/ProjDataInMemory.cxx | 4 ++-- src/include/stir/ProjDataInMemory.h | 10 +++++----- 2 files changed, 7 insertions(+), 7 deletions(-) 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/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); From 8eda22811b1493478cf898d472aabbd8e6391730 Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Mon, 2 Feb 2026 22:52:38 +0000 Subject: [PATCH 3/5] [SWIG] fixes related to new indexT - Looks like SWIG doesn't understand default templates arguments in %template unfortunately. - cope with num_dimenssions swig bug --- src/include/stir/IndexRange.h | 3 +++ src/swig/stir_array.i | 22 +++++++++++----------- 2 files changed, 14 insertions(+), 11 deletions(-) diff --git a/src/include/stir/IndexRange.h b/src/include/stir/IndexRange.h index d2d9654a7a..8f6d574234 100644 --- a/src/include/stir/IndexRange.h +++ b/src/include/stir/IndexRange.h @@ -79,8 +79,11 @@ class IndexRange : public VectorWithOffset& range); diff --git a/src/swig/stir_array.i b/src/swig/stir_array.i index 02b9fa1524..2acf8a6bf7 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); @@ -284,16 +284,16 @@ namespace stir // 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(FloatArray2D) stir::Array<2,float>; + %template(FloatArray2D) stir::Array<2,float, int>; %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>* }; // TODO name %template (FloatNumericVectorWithOffset3D) stir::NumericVectorWithOffset, float>; - %template(FloatArray3D) stir::Array<3,float>; + %template(FloatArray3D) stir::Array<3,float, int>; #if 0 %ADD_indexaccess(int,%arg(stir::Array<2,float>),%arg(stir::Array<3,float>)); #endif %template (FloatNumericVectorWithOffset4D) stir::NumericVectorWithOffset, float>; - %template(FloatArray4D) stir::Array<4,float>; + %template(FloatArray4D) stir::Array<4,float, int>; From 757ff71dcaf42702ef5c9accc0a1091296d053e2 Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Sun, 8 Mar 2026 01:42:00 +0000 Subject: [PATCH 4/5] Fixes for VectorWithOffset/Array indexing - handle max_index < min_index in more places - add/override empty() for a faster check if an Array is empty - earlier exit for size_all() if empty - fix test_Array for unsigned indexT (we had various arrays with negative indices, which made tests fail) --- src/include/stir/Array.h | 4 +- src/include/stir/Array.inl | 46 +++++--- src/include/stir/IndexRange.h | 4 + src/include/stir/IndexRange.inl | 66 ++++++++--- src/include/stir/VectorWithOffset.h | 4 +- src/include/stir/VectorWithOffset.inl | 2 +- src/test/test_Array.cxx | 161 +++++++++++++------------- 7 files changed, 178 insertions(+), 109 deletions(-) diff --git a/src/include/stir/Array.h b/src/include/stir/Array.h index 531dfc250f..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 @@ -202,6 +202,8 @@ class Array : public NumericVectorWithOffset get_index_range() const; + inline bool empty() const override; + //! return the total number of elements in this array inline size_t size_all() const; diff --git a/src/include/stir/Array.inl b/src/include/stir/Array.inl index 0efc9944b9..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 @@ -72,6 +72,11 @@ template void 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(); @@ -87,6 +92,11 @@ template void 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(); @@ -212,11 +222,8 @@ 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()); } @@ -225,11 +232,8 @@ 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()); } @@ -273,6 +277,20 @@ Array::size_all() const 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. @@ -614,7 +632,7 @@ Array<1, elemT, indexT>::resize(const indexT min_index, const indexT max_index, base_type::resize(min_index, max_index); - if (!initialise_with_0) + if (!initialise_with_0 || this->size() == 0) { this->check_state(); return; @@ -622,8 +640,8 @@ Array<1, elemT, indexT>::resize(const indexT min_index, const indexT max_index, if (oldlength == 0) { - for (auto i = this->get_min_index(); i <= this->get_max_index(); i++) - assign(this->num[i], 0); + for (auto& i : *this) + assign(i, 0); } else { diff --git a/src/include/stir/IndexRange.h b/src/include/stir/IndexRange.h index 8f6d574234..9c2213666f 100644 --- a/src/include/stir/IndexRange.h +++ b/src/include/stir/IndexRange.h @@ -113,6 +113,7 @@ class IndexRange : public VectorWithOffset& min, BasicCoordinate& max) const; @@ -147,8 +148,10 @@ class IndexRange<1, indexT> 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, indexT>& range2) const; @@ -163,6 +166,7 @@ class IndexRange<1, indexT> private: 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 ff63c6efd6..0595d56c1c 100644 --- a/src/include/stir/IndexRange.inl +++ b/src/include/stir/IndexRange.inl @@ -68,12 +68,29 @@ IndexRange::IndexRange(const BasicCoordinatefill(lower_dims); } +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 { 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; @@ -102,12 +119,11 @@ bool IndexRange::get_regular_range(BasicCoordinate& min, BasicCoordinate& max) const { - // check if empty range - if (base_type::begin() == base_type::end()) + if (base_type::empty()) { - constexpr bool signed_type = std::is_signed_v; - std::fill(min.begin(), min.end(), signed_type ? 0 : 1); - std::fill(max.begin(), max.end(), signed_type ? -1 : 0); + // use empty range (suitable for unsigned indexT) + std::fill(min.begin(), min.end(), 1); + std::fill(max.begin(), max.end(), 0); return true; } @@ -178,22 +194,41 @@ IndexRange::is_regular() const 1D version ***************************************/ +template +bool +IndexRange<1, indexT>::empty() const +{ + return max < min; +} + +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() - : min(0), - max(0) -{} +{ + 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(); +} template IndexRange<1, indexT>::IndexRange(const BasicCoordinate<1, indexT>& min_v, const BasicCoordinate<1, indexT>& max_v) - : min(min_v[1]), - max(max_v[1]) + : IndexRange<1, indexT>(min_v[1], max_v[1]) {} template @@ -225,7 +260,7 @@ template typename IndexRange<1, indexT>::size_type IndexRange<1, indexT>::get_length() const { - return static_cast(max - min + 1); + return static_cast((max + 1) - min); // note: this order of calculation avoids problems with unsigned indexT } template @@ -267,6 +302,11 @@ template void 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/VectorWithOffset.h b/src/include/stir/VectorWithOffset.h index 80da41416f..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 @@ -260,7 +260,7 @@ class VectorWithOffset 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 //@{ diff --git a/src/include/stir/VectorWithOffset.inl b/src/include/stir/VectorWithOffset.inl index b91c1fcdb6..19a2d7ed6a 100644 --- a/src/include/stir/VectorWithOffset.inl +++ b/src/include/stir/VectorWithOffset.inl @@ -132,7 +132,7 @@ template indexT VectorWithOffset::get_max_index() const { - assert(std::is_signed_v || length > 0); + assert(std::is_signed_v || (length > 0)); return start + length - 1; } diff --git a/src/test/test_Array.cxx b/src/test/test_Array.cxx index 93f51b1719..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. @@ -108,7 +108,7 @@ class ArrayTests : public RunTests 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; auto index = get_min_indices(test); @@ -331,26 +331,29 @@ ArrayTests::run_tests() 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, indexT> c; - c[1] = 0; - check_if_equal(test[c], 11.5F, "test operator[](BasicCoordinate)"); + 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, indexT> ref(-1, 2); - ref[-1] = 1.F; - ref[0] = 3.F; - ref[1] = 3.14F; + 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; @@ -372,14 +375,14 @@ ArrayTests::run_tests() { 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, indexT> test2 = test; - test.grow(-2, test.get_max_index()); + test.grow(offset - 1, test.get_max_index()); Array<1, float, indexT> test3 = test2 + test; - check_if_zero(test3[-2], "test growing during operator+"); + check_if_zero(test3[offset - 1], "test growing during operator+"); } // using preallocated memory @@ -439,11 +442,12 @@ ArrayTests::run_tests() #if 1 { // tests on log/exp - Array<1, float, indexT> 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, indexT> testeq(-3, 10); + Array<1, float, indexT> testeq(test.get_index_range()); check_if_equal(test, testeq, "test in_place_log of Array1D"); } { @@ -462,10 +466,11 @@ ArrayTests::run_tests() { cerr << "Testing 2D stuff" << endl; { - const IndexRange<2, indexT> range(Coordinate2D(0, 0), Coordinate2D(9, 9)); + 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(10), "test size()"); - check_if_equal(test2.size_all(), size_t(100), "test size_all()"); + 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"); @@ -486,7 +491,7 @@ ArrayTests::run_tests() } { - IndexRange<2, indexT> range(Coordinate2D(0, 0), Coordinate2D(3, 3)); + IndexRange<2, indexT> range(Coordinate2D(1, 0), Coordinate2D(3, 3)); Array<2, float, indexT> testfp(range); Array<2, float, indexT> t2fp(range); #if 0 @@ -539,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 (...) { @@ -554,10 +559,10 @@ ArrayTests::run_tests() check(exception_thrown, "out-of-range index should throw an exception"); } - // t2.grow_height(-5,5); - IndexRange<2, indexT> 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 @@ -619,15 +624,15 @@ ArrayTests::run_tests() } // size_all with irregular range { - const IndexRange<2, indexT> range(Coordinate2D(-1, 1), Coordinate2D(1, 2)); + 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 { @@ -655,7 +660,7 @@ ArrayTests::run_tests() } // tests for next() { - const IndexRange<2, indexT> range(Coordinate2D(-1, 1), Coordinate2D(1, 2)); + const IndexRange<2, indexT> range(Coordinate2D(1, 1), Coordinate2D(3, 2)); Array<2, int, indexT> test(range); // fill array with numbers in sequence { @@ -668,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); } @@ -686,18 +691,18 @@ ArrayTests::run_tests() { cerr << "Testing 3D stuff" << endl; - IndexRange<3, indexT> range(Coordinate3D(0, -1, 1), Coordinate3D(3, 3, 3)); + 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"); @@ -707,20 +712,20 @@ ArrayTests::run_tests() 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, indexT> test3bis(range); - test3bis[1][2][1] = (float)6.6; - test3bis[1][0][1] = (float)1.3; + 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, indexT> test3quat = test3 + test3bis; check_if_equal(test3quat, test3ter, "test summing Array3D"); @@ -739,12 +744,12 @@ ArrayTests::run_tests() // size_all with irregular range { - const IndexRange<3, indexT> range(Coordinate3D(-1, 1, 4), Coordinate3D(1, 2, 6)); + 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"); @@ -796,7 +801,7 @@ ArrayTests::run_tests() } // irregular { - test3[0][1].resize(-1, 2); + 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"); @@ -808,40 +813,40 @@ ArrayTests::run_tests() { cerr << "Testing 4D stuff" << endl; - const IndexRange<4, indexT> range(Coordinate4D(-3, 0, -1, 1), Coordinate4D(-2, 3, 3, 3)); + 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, indexT> 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, 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, indexT> test41 = test4; - const IndexRange<4, indexT> mixed_range(Coordinate4D(-4, 1, 0, 1), Coordinate4D(-2, 3, 3, 6)); + 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, indexT> c; - c[1] = -2; + 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.; @@ -849,12 +854,12 @@ ArrayTests::run_tests() } { Array<4, float, indexT> test4bis(range); - test4bis[-2][1][2][1] = (float)6.6; - test4bis[-3][1][0][1] = (float)1.3; + 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. @@ -884,15 +889,15 @@ ArrayTests::run_tests() // 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 (auto i = test4.get_min_index(); i <= indexT(-2); i++) // note: up to -2, as that was the original size + 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)"); } @@ -958,11 +963,11 @@ ArrayTests::run_tests() { typedef NumericVectorWithOffset, float, indexT> NVecArr; typedef typename NVecArr::iterator NVecArrIter; - NVecArr tmp(-1, 2); + 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(); @@ -991,13 +996,13 @@ ArrayTests::run_tests() { typedef NumericVectorWithOffset, float, indexT> NVecArr; typedef typename NVecArr::iterator NVecArrIter; - NVecArr tmp(-1, 2); + 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(); @@ -1034,14 +1039,14 @@ ArrayTests::run_tests() #if 1 { cerr << "Testing 1D float IO" << endl; - Array<1, float, indexT> t1(IndexRange<1, indexT>(-1, 10)); + 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, indexT> range(Coordinate2D(-1, 11), Coordinate2D(10, 20)); + 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++) @@ -1053,7 +1058,7 @@ 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, indexT> range(Coordinate3D(-1, 11, 21), Coordinate3D(10, 20, 30)); + 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++) From ab693201fdaea0f13b751f9f92c633094134c63f Mon Sep 17 00:00:00 2001 From: Kris Thielemans Date: Mon, 9 Mar 2026 02:22:17 +0000 Subject: [PATCH 5/5] [SWIG] fix Array for index type update --- src/swig/stir_array.i | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/src/swig/stir_array.i b/src/swig/stir_array.i index 2acf8a6bf7..5f42d35d9d 100644 --- a/src/swig/stir_array.i +++ b/src/swig/stir_array.i @@ -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, int>; - %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>* }; + %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 (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 (FloatNumericVectorWithOffset4D) stir::NumericVectorWithOffset, float, int>; %template(FloatArray4D) stir::Array<4,float, int>;