1
// This file is part of Eigen, a lightweight C++ template library
4
// Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
6
// Eigen is free software; you can redistribute it and/or
7
// modify it under the terms of the GNU Lesser General Public
8
// License as published by the Free Software Foundation; either
9
// version 3 of the License, or (at your option) any later version.
11
// Alternatively, you can redistribute it and/or
12
// modify it under the terms of the GNU General Public License as
13
// published by the Free Software Foundation; either version 2 of
14
// the License, or (at your option) any later version.
16
// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
17
// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18
// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
19
// GNU General Public License for more details.
21
// You should have received a copy of the GNU Lesser General Public
22
// License and a copy of the GNU General Public License along with
23
// Eigen. If not, see <http://www.gnu.org/licenses/>.
25
#ifndef EIGEN_HOMOGENEOUS_H
26
#define EIGEN_HOMOGENEOUS_H
28
/** \geometry_module \ingroup Geometry_Module
32
* \brief Expression of one (or a set of) homogeneous vector(s)
34
* \param MatrixType the type of the object in which we are making homogeneous
36
* This class represents an expression of one (or a set of) homogeneous vector(s).
37
* It is the return type of MatrixBase::homogeneous() and most of the time
38
* this is the only way it is used.
40
* \sa MatrixBase::homogeneous()
45
template<typename MatrixType,int Direction>
46
struct traits<Homogeneous<MatrixType,Direction> >
49
typedef typename traits<MatrixType>::StorageKind StorageKind;
50
typedef typename nested<MatrixType>::type MatrixTypeNested;
51
typedef typename remove_reference<MatrixTypeNested>::type _MatrixTypeNested;
53
RowsPlusOne = (MatrixType::RowsAtCompileTime != Dynamic) ?
54
int(MatrixType::RowsAtCompileTime) + 1 : Dynamic,
55
ColsPlusOne = (MatrixType::ColsAtCompileTime != Dynamic) ?
56
int(MatrixType::ColsAtCompileTime) + 1 : Dynamic,
57
RowsAtCompileTime = Direction==Vertical ? RowsPlusOne : MatrixType::RowsAtCompileTime,
58
ColsAtCompileTime = Direction==Horizontal ? ColsPlusOne : MatrixType::ColsAtCompileTime,
59
MaxRowsAtCompileTime = RowsAtCompileTime,
60
MaxColsAtCompileTime = ColsAtCompileTime,
61
TmpFlags = _MatrixTypeNested::Flags & HereditaryBits,
62
Flags = ColsAtCompileTime==1 ? (TmpFlags & ~RowMajorBit)
63
: RowsAtCompileTime==1 ? (TmpFlags | RowMajorBit)
65
CoeffReadCost = _MatrixTypeNested::CoeffReadCost
69
template<typename MatrixType,typename Lhs> struct homogeneous_left_product_impl;
70
template<typename MatrixType,typename Rhs> struct homogeneous_right_product_impl;
72
} // end namespace internal
74
template<typename MatrixType,int _Direction> class Homogeneous
75
: public MatrixBase<Homogeneous<MatrixType,_Direction> >
79
enum { Direction = _Direction };
81
typedef MatrixBase<Homogeneous> Base;
82
EIGEN_DENSE_PUBLIC_INTERFACE(Homogeneous)
84
inline Homogeneous(const MatrixType& matrix)
88
inline Index rows() const { return m_matrix.rows() + (int(Direction)==Vertical ? 1 : 0); }
89
inline Index cols() const { return m_matrix.cols() + (int(Direction)==Horizontal ? 1 : 0); }
91
inline Scalar coeff(Index row, Index col) const
93
if( (int(Direction)==Vertical && row==m_matrix.rows())
94
|| (int(Direction)==Horizontal && col==m_matrix.cols()))
96
return m_matrix.coeff(row, col);
99
template<typename Rhs>
100
inline const internal::homogeneous_right_product_impl<Homogeneous,Rhs>
101
operator* (const MatrixBase<Rhs>& rhs) const
103
eigen_assert(int(Direction)==Horizontal);
104
return internal::homogeneous_right_product_impl<Homogeneous,Rhs>(m_matrix,rhs.derived());
107
template<typename Lhs> friend
108
inline const internal::homogeneous_left_product_impl<Homogeneous,Lhs>
109
operator* (const MatrixBase<Lhs>& lhs, const Homogeneous& rhs)
111
eigen_assert(int(Direction)==Vertical);
112
return internal::homogeneous_left_product_impl<Homogeneous,Lhs>(lhs.derived(),rhs.m_matrix);
115
template<typename Scalar, int Dim, int Mode, int Options> friend
116
inline const internal::homogeneous_left_product_impl<Homogeneous,Transform<Scalar,Dim,Mode,Options> >
117
operator* (const Transform<Scalar,Dim,Mode,Options>& lhs, const Homogeneous& rhs)
119
eigen_assert(int(Direction)==Vertical);
120
return internal::homogeneous_left_product_impl<Homogeneous,Transform<Scalar,Dim,Mode,Options> >(lhs,rhs.m_matrix);
124
const typename MatrixType::Nested m_matrix;
129
* \return an expression of the equivalent homogeneous vector
133
* Example: \include MatrixBase_homogeneous.cpp
134
* Output: \verbinclude MatrixBase_homogeneous.out
136
* \sa class Homogeneous
138
template<typename Derived>
139
inline typename MatrixBase<Derived>::HomogeneousReturnType
140
MatrixBase<Derived>::homogeneous() const
142
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived);
148
* \returns a matrix expression of homogeneous column (or row) vectors
150
* Example: \include VectorwiseOp_homogeneous.cpp
151
* Output: \verbinclude VectorwiseOp_homogeneous.out
153
* \sa MatrixBase::homogeneous() */
154
template<typename ExpressionType, int Direction>
155
inline Homogeneous<ExpressionType,Direction>
156
VectorwiseOp<ExpressionType,Direction>::homogeneous() const
158
return _expression();
163
* \returns an expression of the homogeneous normalized vector of \c *this
165
* Example: \include MatrixBase_hnormalized.cpp
166
* Output: \verbinclude MatrixBase_hnormalized.out
168
* \sa VectorwiseOp::hnormalized() */
169
template<typename Derived>
170
inline const typename MatrixBase<Derived>::HNormalizedReturnType
171
MatrixBase<Derived>::hnormalized() const
173
EIGEN_STATIC_ASSERT_VECTOR_ONLY(Derived);
174
return ConstStartMinusOne(derived(),0,0,
175
ColsAtCompileTime==1?size()-1:1,
176
ColsAtCompileTime==1?1:size()-1) / coeff(size()-1);
181
* \returns an expression of the homogeneous normalized vector of \c *this
183
* Example: \include DirectionWise_hnormalized.cpp
184
* Output: \verbinclude DirectionWise_hnormalized.out
186
* \sa MatrixBase::hnormalized() */
187
template<typename ExpressionType, int Direction>
188
inline const typename VectorwiseOp<ExpressionType,Direction>::HNormalizedReturnType
189
VectorwiseOp<ExpressionType,Direction>::hnormalized() const
191
return HNormalized_Block(_expression(),0,0,
192
Direction==Vertical ? _expression().rows()-1 : _expression().rows(),
193
Direction==Horizontal ? _expression().cols()-1 : _expression().cols()).cwiseQuotient(
194
Replicate<HNormalized_Factors,
195
Direction==Vertical ? HNormalized_SizeMinusOne : 1,
196
Direction==Horizontal ? HNormalized_SizeMinusOne : 1>
197
(HNormalized_Factors(_expression(),
198
Direction==Vertical ? _expression().rows()-1:0,
199
Direction==Horizontal ? _expression().cols()-1:0,
200
Direction==Vertical ? 1 : _expression().rows(),
201
Direction==Horizontal ? 1 : _expression().cols()),
202
Direction==Vertical ? _expression().rows()-1 : 1,
203
Direction==Horizontal ? _expression().cols()-1 : 1));
208
template<typename MatrixOrTransformType>
209
struct take_matrix_for_product
211
typedef MatrixOrTransformType type;
212
static const type& run(const type &x) { return x; }
215
template<typename Scalar, int Dim, int Mode,int Options>
216
struct take_matrix_for_product<Transform<Scalar, Dim, Mode, Options> >
218
typedef Transform<Scalar, Dim, Mode, Options> TransformType;
219
typedef typename TransformType::ConstAffinePart type;
220
static const type run (const TransformType& x) { return x.affine(); }
223
template<typename Scalar, int Dim, int Options>
224
struct take_matrix_for_product<Transform<Scalar, Dim, Projective, Options> >
226
typedef Transform<Scalar, Dim, Projective, Options> TransformType;
227
typedef typename TransformType::MatrixType type;
228
static const type& run (const TransformType& x) { return x.matrix(); }
231
template<typename MatrixType,typename Lhs>
232
struct traits<homogeneous_left_product_impl<Homogeneous<MatrixType,Vertical>,Lhs> >
234
typedef typename take_matrix_for_product<Lhs>::type LhsMatrixType;
235
typedef typename remove_all<MatrixType>::type MatrixTypeCleaned;
236
typedef typename remove_all<LhsMatrixType>::type LhsMatrixTypeCleaned;
237
typedef typename make_proper_matrix_type<
238
typename traits<MatrixTypeCleaned>::Scalar,
239
LhsMatrixTypeCleaned::RowsAtCompileTime,
240
MatrixTypeCleaned::ColsAtCompileTime,
241
MatrixTypeCleaned::PlainObject::Options,
242
LhsMatrixTypeCleaned::MaxRowsAtCompileTime,
243
MatrixTypeCleaned::MaxColsAtCompileTime>::type ReturnType;
246
template<typename MatrixType,typename Lhs>
247
struct homogeneous_left_product_impl<Homogeneous<MatrixType,Vertical>,Lhs>
248
: public ReturnByValue<homogeneous_left_product_impl<Homogeneous<MatrixType,Vertical>,Lhs> >
250
typedef typename traits<homogeneous_left_product_impl>::LhsMatrixType LhsMatrixType;
251
typedef typename remove_all<LhsMatrixType>::type LhsMatrixTypeCleaned;
252
typedef typename remove_all<typename LhsMatrixTypeCleaned::Nested>::type LhsMatrixTypeNested;
253
typedef typename MatrixType::Index Index;
254
homogeneous_left_product_impl(const Lhs& lhs, const MatrixType& rhs)
255
: m_lhs(take_matrix_for_product<Lhs>::run(lhs)),
259
inline Index rows() const { return m_lhs.rows(); }
260
inline Index cols() const { return m_rhs.cols(); }
262
template<typename Dest> void evalTo(Dest& dst) const
264
// FIXME investigate how to allow lazy evaluation of this product when possible
265
dst = Block<const LhsMatrixTypeNested,
266
LhsMatrixTypeNested::RowsAtCompileTime,
267
LhsMatrixTypeNested::ColsAtCompileTime==Dynamic?Dynamic:LhsMatrixTypeNested::ColsAtCompileTime-1>
268
(m_lhs,0,0,m_lhs.rows(),m_lhs.cols()-1) * m_rhs;
269
dst += m_lhs.col(m_lhs.cols()-1).rowwise()
270
.template replicate<MatrixType::ColsAtCompileTime>(m_rhs.cols());
273
const typename LhsMatrixTypeCleaned::Nested m_lhs;
274
const typename MatrixType::Nested m_rhs;
277
template<typename MatrixType,typename Rhs>
278
struct traits<homogeneous_right_product_impl<Homogeneous<MatrixType,Horizontal>,Rhs> >
280
typedef typename make_proper_matrix_type<typename traits<MatrixType>::Scalar,
281
MatrixType::RowsAtCompileTime,
282
Rhs::ColsAtCompileTime,
283
MatrixType::PlainObject::Options,
284
MatrixType::MaxRowsAtCompileTime,
285
Rhs::MaxColsAtCompileTime>::type ReturnType;
288
template<typename MatrixType,typename Rhs>
289
struct homogeneous_right_product_impl<Homogeneous<MatrixType,Horizontal>,Rhs>
290
: public ReturnByValue<homogeneous_right_product_impl<Homogeneous<MatrixType,Horizontal>,Rhs> >
292
typedef typename remove_all<typename Rhs::Nested>::type RhsNested;
293
typedef typename MatrixType::Index Index;
294
homogeneous_right_product_impl(const MatrixType& lhs, const Rhs& rhs)
295
: m_lhs(lhs), m_rhs(rhs)
298
inline Index rows() const { return m_lhs.rows(); }
299
inline Index cols() const { return m_rhs.cols(); }
301
template<typename Dest> void evalTo(Dest& dst) const
303
// FIXME investigate how to allow lazy evaluation of this product when possible
304
dst = m_lhs * Block<const RhsNested,
305
RhsNested::RowsAtCompileTime==Dynamic?Dynamic:RhsNested::RowsAtCompileTime-1,
306
RhsNested::ColsAtCompileTime>
307
(m_rhs,0,0,m_rhs.rows()-1,m_rhs.cols());
308
dst += m_rhs.row(m_rhs.rows()-1).colwise()
309
.template replicate<MatrixType::RowsAtCompileTime>(m_lhs.rows());
312
const typename MatrixType::Nested m_lhs;
313
const typename Rhs::Nested m_rhs;
316
} // end namespace internal
318
#endif // EIGEN_HOMOGENEOUS_H