• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1@gmail.com>
5 // Copyright (C) 2009 Ricard Marxer <email@ricardmarxer.com>
6 // Copyright (C) 2009-2010 Gael Guennebaud <gael.guennebaud@inria.fr>
7 //
8 // This Source Code Form is subject to the terms of the Mozilla
9 // Public License v. 2.0. If a copy of the MPL was not distributed
10 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
11 
12 #ifndef EIGEN_REVERSE_H
13 #define EIGEN_REVERSE_H
14 
15 namespace Eigen {
16 
17 namespace internal {
18 
19 template<typename MatrixType, int Direction>
20 struct traits<Reverse<MatrixType, Direction> >
21  : traits<MatrixType>
22 {
23   typedef typename MatrixType::Scalar Scalar;
24   typedef typename traits<MatrixType>::StorageKind StorageKind;
25   typedef typename traits<MatrixType>::XprKind XprKind;
26   typedef typename ref_selector<MatrixType>::type MatrixTypeNested;
27   typedef typename remove_reference<MatrixTypeNested>::type _MatrixTypeNested;
28   enum {
29     RowsAtCompileTime = MatrixType::RowsAtCompileTime,
30     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
31     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
32     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
33     Flags = _MatrixTypeNested::Flags & (RowMajorBit | LvalueBit)
34   };
35 };
36 
37 template<typename PacketType, bool ReversePacket> struct reverse_packet_cond
38 {
39   static inline PacketType run(const PacketType& x) { return preverse(x); }
40 };
41 
42 template<typename PacketType> struct reverse_packet_cond<PacketType,false>
43 {
44   static inline PacketType run(const PacketType& x) { return x; }
45 };
46 
47 } // end namespace internal
48 
49 /** \class Reverse
50   * \ingroup Core_Module
51   *
52   * \brief Expression of the reverse of a vector or matrix
53   *
54   * \tparam MatrixType the type of the object of which we are taking the reverse
55   * \tparam Direction defines the direction of the reverse operation, can be Vertical, Horizontal, or BothDirections
56   *
57   * This class represents an expression of the reverse of a vector.
58   * It is the return type of MatrixBase::reverse() and VectorwiseOp::reverse()
59   * and most of the time this is the only way it is used.
60   *
61   * \sa MatrixBase::reverse(), VectorwiseOp::reverse()
62   */
63 template<typename MatrixType, int Direction> class Reverse
64   : public internal::dense_xpr_base< Reverse<MatrixType, Direction> >::type
65 {
66   public:
67 
68     typedef typename internal::dense_xpr_base<Reverse>::type Base;
69     EIGEN_DENSE_PUBLIC_INTERFACE(Reverse)
70     typedef typename internal::remove_all<MatrixType>::type NestedExpression;
71     using Base::IsRowMajor;
72 
73   protected:
74     enum {
75       PacketSize = internal::packet_traits<Scalar>::size,
76       IsColMajor = !IsRowMajor,
77       ReverseRow = (Direction == Vertical)   || (Direction == BothDirections),
78       ReverseCol = (Direction == Horizontal) || (Direction == BothDirections),
79       OffsetRow  = ReverseRow && IsColMajor ? PacketSize : 1,
80       OffsetCol  = ReverseCol && IsRowMajor ? PacketSize : 1,
81       ReversePacket = (Direction == BothDirections)
82                     || ((Direction == Vertical)   && IsColMajor)
83                     || ((Direction == Horizontal) && IsRowMajor)
84     };
85     typedef internal::reverse_packet_cond<PacketScalar,ReversePacket> reverse_packet;
86   public:
87 
88     EIGEN_DEVICE_FUNC explicit inline Reverse(const MatrixType& matrix) : m_matrix(matrix) { }
89 
90     EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Reverse)
91 
92     EIGEN_DEVICE_FUNC inline Index rows() const { return m_matrix.rows(); }
93     EIGEN_DEVICE_FUNC inline Index cols() const { return m_matrix.cols(); }
94 
95     EIGEN_DEVICE_FUNC inline Index innerStride() const
96     {
97       return -m_matrix.innerStride();
98     }
99 
100     EIGEN_DEVICE_FUNC const typename internal::remove_all<typename MatrixType::Nested>::type&
101     nestedExpression() const
102     {
103       return m_matrix;
104     }
105 
106   protected:
107     typename MatrixType::Nested m_matrix;
108 };
109 
110 /** \returns an expression of the reverse of *this.
111   *
112   * Example: \include MatrixBase_reverse.cpp
113   * Output: \verbinclude MatrixBase_reverse.out
114   *
115   */
116 template<typename Derived>
117 inline typename DenseBase<Derived>::ReverseReturnType
118 DenseBase<Derived>::reverse()
119 {
120   return ReverseReturnType(derived());
121 }
122 
123 
124 //reverse const overload moved DenseBase.h due to a CUDA compiler bug
125 
126 /** This is the "in place" version of reverse: it reverses \c *this.
127   *
128   * In most cases it is probably better to simply use the reversed expression
129   * of a matrix. However, when reversing the matrix data itself is really needed,
130   * then this "in-place" version is probably the right choice because it provides
131   * the following additional benefits:
132   *  - less error prone: doing the same operation with .reverse() requires special care:
133   *    \code m = m.reverse().eval(); \endcode
134   *  - this API enables reverse operations without the need for a temporary
135   *  - it allows future optimizations (cache friendliness, etc.)
136   *
137   * \sa VectorwiseOp::reverseInPlace(), reverse() */
138 template<typename Derived>
139 inline void DenseBase<Derived>::reverseInPlace()
140 {
141   if(cols()>rows())
142   {
143     Index half = cols()/2;
144     leftCols(half).swap(rightCols(half).reverse());
145     if((cols()%2)==1)
146     {
147       Index half2 = rows()/2;
148       col(half).head(half2).swap(col(half).tail(half2).reverse());
149     }
150   }
151   else
152   {
153     Index half = rows()/2;
154     topRows(half).swap(bottomRows(half).reverse());
155     if((rows()%2)==1)
156     {
157       Index half2 = cols()/2;
158       row(half).head(half2).swap(row(half).tail(half2).reverse());
159     }
160   }
161 }
162 
163 namespace internal {
164 
165 template<int Direction>
166 struct vectorwise_reverse_inplace_impl;
167 
168 template<>
169 struct vectorwise_reverse_inplace_impl<Vertical>
170 {
171   template<typename ExpressionType>
172   static void run(ExpressionType &xpr)
173   {
174     Index half = xpr.rows()/2;
175     xpr.topRows(half).swap(xpr.bottomRows(half).colwise().reverse());
176   }
177 };
178 
179 template<>
180 struct vectorwise_reverse_inplace_impl<Horizontal>
181 {
182   template<typename ExpressionType>
183   static void run(ExpressionType &xpr)
184   {
185     Index half = xpr.cols()/2;
186     xpr.leftCols(half).swap(xpr.rightCols(half).rowwise().reverse());
187   }
188 };
189 
190 } // end namespace internal
191 
192 /** This is the "in place" version of VectorwiseOp::reverse: it reverses each column or row of \c *this.
193   *
194   * In most cases it is probably better to simply use the reversed expression
195   * of a matrix. However, when reversing the matrix data itself is really needed,
196   * then this "in-place" version is probably the right choice because it provides
197   * the following additional benefits:
198   *  - less error prone: doing the same operation with .reverse() requires special care:
199   *    \code m = m.reverse().eval(); \endcode
200   *  - this API enables reverse operations without the need for a temporary
201   *
202   * \sa DenseBase::reverseInPlace(), reverse() */
203 template<typename ExpressionType, int Direction>
204 void VectorwiseOp<ExpressionType,Direction>::reverseInPlace()
205 {
206   internal::vectorwise_reverse_inplace_impl<Direction>::run(_expression().const_cast_derived());
207 }
208 
209 } // end namespace Eigen
210 
211 #endif // EIGEN_REVERSE_H
212