Reverse.h
Go to the documentation of this file.
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 // Eigen is free software; you can redistribute it and/or
9 // modify it under the terms of the GNU Lesser General Public
10 // License as published by the Free Software Foundation; either
11 // version 3 of the License, or (at your option) any later version.
12 //
13 // Alternatively, you can redistribute it and/or
14 // modify it under the terms of the GNU General Public License as
15 // published by the Free Software Foundation; either version 2 of
16 // the License, or (at your option) any later version.
17 //
18 // Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
19 // WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
21 // GNU General Public License for more details.
22 //
23 // You should have received a copy of the GNU Lesser General Public
24 // License and a copy of the GNU General Public License along with
25 // Eigen. If not, see <http://www.gnu.org/licenses/>.
26 
27 #ifndef EIGEN_REVERSE_H
28 #define EIGEN_REVERSE_H
29 
30 namespace Eigen {
31 
46 namespace internal {
47 
48 template<typename MatrixType, int Direction>
49 struct traits<Reverse<MatrixType, Direction> >
50  : traits<MatrixType>
51 {
52  typedef typename MatrixType::Scalar Scalar;
53  typedef typename traits<MatrixType>::StorageKind StorageKind;
54  typedef typename traits<MatrixType>::XprKind XprKind;
55  typedef typename nested<MatrixType>::type MatrixTypeNested;
56  typedef typename remove_reference<MatrixTypeNested>::type _MatrixTypeNested;
57  enum {
58  RowsAtCompileTime = MatrixType::RowsAtCompileTime,
59  ColsAtCompileTime = MatrixType::ColsAtCompileTime,
60  MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
61  MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
62 
63  // let's enable LinearAccess only with vectorization because of the product overhead
64  LinearAccess = ( (Direction==BothDirections) && (int(_MatrixTypeNested::Flags)&PacketAccessBit) )
65  ? LinearAccessBit : 0,
66 
67  Flags = int(_MatrixTypeNested::Flags) & (HereditaryBits | LvalueBit | PacketAccessBit | LinearAccess),
68 
69  CoeffReadCost = _MatrixTypeNested::CoeffReadCost
70  };
71 };
72 
73 template<typename PacketScalar, bool ReversePacket> struct reverse_packet_cond
74 {
75  static inline PacketScalar run(const PacketScalar& x) { return preverse(x); }
76 };
77 
78 template<typename PacketScalar> struct reverse_packet_cond<PacketScalar,false>
79 {
80  static inline PacketScalar run(const PacketScalar& x) { return x; }
81 };
82 
83 } // end namespace internal
84 
85 template<typename MatrixType, int Direction> class Reverse
86  : public internal::dense_xpr_base< Reverse<MatrixType, Direction> >::type
87 {
88  public:
89 
90  typedef typename internal::dense_xpr_base<Reverse>::type Base;
92  using Base::IsRowMajor;
93 
94  // next line is necessary because otherwise const version of operator()
95  // is hidden by non-const version defined in this file
96  using Base::operator();
97 
98  protected:
99  enum {
100  PacketSize = internal::packet_traits<Scalar>::size,
101  IsColMajor = !IsRowMajor,
102  ReverseRow = (Direction == Vertical) || (Direction == BothDirections),
103  ReverseCol = (Direction == Horizontal) || (Direction == BothDirections),
104  OffsetRow = ReverseRow && IsColMajor ? PacketSize : 1,
105  OffsetCol = ReverseCol && IsRowMajor ? PacketSize : 1,
106  ReversePacket = (Direction == BothDirections)
107  || ((Direction == Vertical) && IsColMajor)
108  || ((Direction == Horizontal) && IsRowMajor)
109  };
110  typedef internal::reverse_packet_cond<PacketScalar,ReversePacket> reverse_packet;
111  public:
112 
113  inline Reverse(const MatrixType& matrix) : m_matrix(matrix) { }
114 
116 
117  inline Index rows() const { return m_matrix.rows(); }
118  inline Index cols() const { return m_matrix.cols(); }
119 
120  inline Index innerStride() const
121  {
122  return -m_matrix.innerStride();
123  }
124 
125  inline Scalar& operator()(Index row, Index col)
126  {
127  eigen_assert(row >= 0 && row < rows() && col >= 0 && col < cols());
128  return coeffRef(row, col);
129  }
130 
131  inline Scalar& coeffRef(Index row, Index col)
132  {
133  return m_matrix.const_cast_derived().coeffRef(ReverseRow ? m_matrix.rows() - row - 1 : row,
134  ReverseCol ? m_matrix.cols() - col - 1 : col);
135  }
136 
137  inline CoeffReturnType coeff(Index row, Index col) const
138  {
139  return m_matrix.coeff(ReverseRow ? m_matrix.rows() - row - 1 : row,
140  ReverseCol ? m_matrix.cols() - col - 1 : col);
141  }
142 
143  inline CoeffReturnType coeff(Index index) const
144  {
145  return m_matrix.coeff(m_matrix.size() - index - 1);
146  }
147 
148  inline Scalar& coeffRef(Index index)
149  {
150  return m_matrix.const_cast_derived().coeffRef(m_matrix.size() - index - 1);
151  }
152 
153  inline Scalar& operator()(Index index)
154  {
155  eigen_assert(index >= 0 && index < m_matrix.size());
156  return coeffRef(index);
157  }
158 
159  template<int LoadMode>
160  inline const PacketScalar packet(Index row, Index col) const
161  {
162  return reverse_packet::run(m_matrix.template packet<LoadMode>(
163  ReverseRow ? m_matrix.rows() - row - OffsetRow : row,
164  ReverseCol ? m_matrix.cols() - col - OffsetCol : col));
165  }
166 
167  template<int LoadMode>
168  inline void writePacket(Index row, Index col, const PacketScalar& x)
169  {
170  m_matrix.const_cast_derived().template writePacket<LoadMode>(
171  ReverseRow ? m_matrix.rows() - row - OffsetRow : row,
172  ReverseCol ? m_matrix.cols() - col - OffsetCol : col,
173  reverse_packet::run(x));
174  }
175 
176  template<int LoadMode>
177  inline const PacketScalar packet(Index index) const
178  {
179  return internal::preverse(m_matrix.template packet<LoadMode>( m_matrix.size() - index - PacketSize ));
180  }
181 
182  template<int LoadMode>
183  inline void writePacket(Index index, const PacketScalar& x)
184  {
185  m_matrix.const_cast_derived().template writePacket<LoadMode>(m_matrix.size() - index - PacketSize, internal::preverse(x));
186  }
187 
188  const typename internal::remove_all<typename MatrixType::Nested>::type&
189  nestedExpression() const
190  {
191  return m_matrix;
192  }
193 
194  protected:
195  typename MatrixType::Nested m_matrix;
196 };
197 
204 template<typename Derived>
207 {
208  return derived();
209 }
210 
212 template<typename Derived>
213 inline const typename DenseBase<Derived>::ConstReverseReturnType
215 {
216  return derived();
217 }
218 
231 template<typename Derived>
233 {
234  derived() = derived().reverse().eval();
235 }
236 
237 } // end namespace Eigen
238 
239 #endif // EIGEN_REVERSE_H