blob: 95d3dd31a153c8f2cbc3ad4ab85feb7a8b5dd97c (
plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
|
#include <Eigen/Core>
#include <iostream>
using namespace Eigen;
// [circulant_func]
template<class ArgType>
class circulant_functor {
const ArgType &m_vec;
public:
circulant_functor(const ArgType& arg) : m_vec(arg) {}
const typename ArgType::Scalar& operator() (Index row, Index col) const {
Index index = row - col;
if (index < 0) index += m_vec.size();
return m_vec(index);
}
};
// [circulant_func]
// [square]
template<class ArgType>
struct circulant_helper {
typedef Matrix<typename ArgType::Scalar,
ArgType::SizeAtCompileTime,
ArgType::SizeAtCompileTime,
ColMajor,
ArgType::MaxSizeAtCompileTime,
ArgType::MaxSizeAtCompileTime> MatrixType;
};
// [square]
// [makeCirculant]
template <class ArgType>
CwiseNullaryOp<circulant_functor<ArgType>, typename circulant_helper<ArgType>::MatrixType>
makeCirculant(const Eigen::MatrixBase<ArgType>& arg)
{
typedef typename circulant_helper<ArgType>::MatrixType MatrixType;
return MatrixType::NullaryExpr(arg.size(), arg.size(), circulant_functor<ArgType>(arg.derived()));
}
// [makeCirculant]
// [main]
int main()
{
Eigen::VectorXd vec(4);
vec << 1, 2, 4, 8;
Eigen::MatrixXd mat;
mat = makeCirculant(vec);
std::cout << mat << std::endl;
}
// [main]
|