ITK  6.0.0
Insight Toolkit
itkSymmetricSecondRankTensor.h
Go to the documentation of this file.
1/*=========================================================================
2 *
3 * Copyright NumFOCUS
4 *
5 * Licensed under the Apache License, Version 2.0 (the "License");
6 * you may not use this file except in compliance with the License.
7 * You may obtain a copy of the License at
8 *
9 * https://www.apache.org/licenses/LICENSE-2.0.txt
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 *
17 *=========================================================================*/
18#ifndef itkSymmetricSecondRankTensor_h
19#define itkSymmetricSecondRankTensor_h
20
21// Undefine an eventual SymmetricSecondRankTensor macro
22#ifdef SymmetricSecondRankTensor
23# undef SymmetricSecondRankTensor
24#endif
25
26#include "itkIndent.h"
27#include "itkFixedArray.h"
28#include "itkMatrix.h"
30
31namespace itk
32{
74template <typename TComponent, unsigned int VDimension = 3>
75class ITK_TEMPLATE_EXPORT SymmetricSecondRankTensor : public FixedArray<TComponent, VDimension *(VDimension + 1) / 2>
76{
77public:
80 using Superclass = FixedArray<TComponent, VDimension *(VDimension + 1) / 2>;
81
83 static constexpr unsigned int Dimension = VDimension;
84 static constexpr unsigned int InternalDimension = VDimension * (VDimension + 1) / 2;
85
88
91
95
97 using ComponentType = TComponent;
98 using typename Superclass::ValueType;
101
104
107#ifdef ITK_FUTURE_LEGACY_REMOVE
109#else
110 SymmetricSecondRankTensor() { this->Fill(0); }
111#endif
114 SymmetricSecondRankTensor(const ComponentType & r) { this->Fill(r); }
115
117 template <typename TCoordRepB>
119 : BaseArray(pa)
120 {}
121
122 using ComponentArrayType = ComponentType[Self::InternalDimension];
123
126 : BaseArray(r)
127 {}
128
130 template <typename TCoordRepB>
131 Self &
133 {
134 BaseArray::operator=(pa);
135 return *this;
136 }
140 Self &
142
143 Self &
145
148 Self
149 operator+(const Self & r) const;
150
151 Self
152 operator-(const Self & r) const;
153
154 const Self &
155 operator+=(const Self & r);
156
157 const Self &
158 operator-=(const Self & r);
159
161 Self operator*(const RealValueType & r) const;
162
163 Self
164 operator/(const RealValueType & r) const;
165
166 const Self &
168
169 const Self &
171
173 static unsigned int
175 {
176 return Self::InternalDimension;
177 }
178
180 ComponentType
181 GetNthComponent(int c) const
182 {
183 return this->operator[](c);
184 }
185
187 void
189 {
190 this->operator[](c) = v;
191 }
192
194 ValueType &
195 operator()(unsigned int row, unsigned int col);
196
197 const ValueType &
198 operator()(unsigned int row, unsigned int col) const;
199
202 void
204
207 GetTrace() const;
208
210 void
212
215 void
217
221 template <typename TMatrixValueType>
222 Self
224 template <typename TMatrixValueType>
225 Self
226 Rotate(const vnl_matrix_fixed<TMatrixValueType, VDimension, VDimension> & m) const
227 {
228 return this->Rotate(static_cast<Matrix<TMatrixValueType, VDimension, VDimension>>(m));
229 }
230 template <typename TMatrixValueType>
231 Self
232 Rotate(const vnl_matrix<TMatrixValueType> & m) const
233 {
234 return this->Rotate(static_cast<Matrix<TMatrixValueType>>(m));
235 }
239 MatrixType
240 PreMultiply(const MatrixType & m) const;
241
244 PostMultiply(const MatrixType & m) const;
245
246private:
247};
248
251using OutputStreamType = std::ostream;
252using InputStreamType = std::istream;
253
254template <typename TComponent, unsigned int VDimension>
257
258template <typename TComponent, unsigned int VDimension>
261
262template <typename T>
263inline void
265{
266 a.swap(b);
267}
268} // end namespace itk
269
271
272#ifndef ITK_MANUAL_INSTANTIATION
273# include "itkSymmetricSecondRankTensor.hxx"
274#endif
275
276#endif
Pixel-wise addition of two images.
Simulate a standard C array with copy semantics.
Definition: itkFixedArray.h:54
void swap(FixedArray &other)
A templated class holding a M x N size Matrix.
Definition: itkMatrix.h:53
Represent a symmetric tensor of second rank.
void ComputeEigenAnalysis(EigenValuesArrayType &eigenValues, EigenVectorsMatrixType &eigenVectors) const
const Self & operator*=(const RealValueType &r)
Self operator/(const RealValueType &r) const
MatrixType PostMultiply(const MatrixType &m) const
Self operator*(const RealValueType &r) const
typename NumericTraits< ValueType >::RealType RealValueType
Self Rotate(const vnl_matrix< TMatrixValueType > &m) const
SymmetricSecondRankTensor(const ComponentType &r)
Self & operator=(const ComponentArrayType r)
Self operator+(const Self &r) const
SymmetricSecondRankTensor(const ComponentArrayType r)
Self operator-(const Self &r) const
ComponentType[Self::InternalDimension] ComponentArrayType
Self & operator=(const SymmetricSecondRankTensor< TCoordRepB, VDimension > &pa)
AccumulateValueType GetTrace() const
Self Rotate(const Matrix< TMatrixValueType, VDimension, VDimension > &m) const
Self & operator=(const ComponentType &r)
MatrixType PreMultiply(const MatrixType &m) const
ComponentType GetNthComponent(int c) const
const Self & operator+=(const Self &r)
void ComputeEigenValues(EigenValuesArrayType &eigenValues) const
const Self & operator-=(const Self &r)
Self Rotate(const vnl_matrix_fixed< TMatrixValueType, VDimension, VDimension > &m) const
const ValueType & operator()(unsigned int row, unsigned int col) const
typename NumericTraits< ValueType >::RealType AccumulateValueType
const Self & operator/=(const RealValueType &r)
void SetNthComponent(int c, const ComponentType &v)
ValueType & operator()(unsigned int row, unsigned int col)
SymmetricSecondRankTensor(const SymmetricSecondRankTensor< TCoordRepB, VDimension > &pa)
The "itk" namespace contains all Insight Segmentation and Registration Toolkit (ITK) classes....
void swap(Array< T > &a, Array< T > &b) noexcept
Definition: itkArray.h:242
std::istream & operator>>(std::istream &is, Point< T, VPointDimension > &vct)
ITKCommon_EXPORT std::ostream & operator<<(std::ostream &out, typename AnatomicalOrientation::CoordinateEnum value)
std::istream InputStreamType
std::ostream OutputStreamType