Search code examples
c++eigen

Static reshape of Eigen matrix

I am experimenting with doing bicubic interpolation of some gridded data using Eigen, and I can't figure out how to reshape the 16x1 column vector of coefficients into a 4x4 matrix. Ideally I would like to do something along the lines of https://bitbucket.org/eigen/eigen/pull-request/41/reshape/diff without any copying, but I can't make heads or tails of the docs. Alternatively, a map would be fine as well, but I can't figure out how to use a map on an already existing matrix.

More here: http://en.wikipedia.org/wiki/Bicubic_interpolation

/// The inverse of the A matrix for the bicubic interpolation 
/// (http://en.wikipedia.org/wiki/Bicubic_interpolation)
static const double Ainv_data[16*16] = {
     1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
     0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
    -3,  3,  0,  0, -2, -1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
     2, -2,  0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,
     0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,
     0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0,
     0,  0,  0,  0,  0,  0,  0,  0, -3,  3,  0,  0, -2, -1,  0,  0,
     0,  0,  0,  0,  0,  0,  0,  0,  2, -2,  0,  0,  1,  1,  0,  0,
    -3,  0,  3,  0,  0,  0,  0,  0, -2,  0, -1,  0,  0,  0,  0,  0,
     0,  0,  0,  0, -3,  0,  3,  0,  0,  0,  0,  0, -2,  0, -1,  0,
     9, -9, -9,  9,  6,  3, -6, -3,  6, -6,  3, -3,  4,  2,  2,  1,
    -6,  6,  6, -6, -3, -3,  3,  3, -4,  4, -2,  2, -2, -2, -1, -1,
     2,  0, -2,  0,  0,  0,  0,  0,  1,  0,  1,  0,  0,  0,  0,  0,
     0,  0,  0,  0,  2,  0, -2,  0,  0,  0,  0,  0,  1,  0,  1,  0,
    -6,  6,  6, -6, -4, -2,  4,  2, -3,  3, -3,  3, -2, -1, -2, -1,
     4, -4, -4,  4,  2,  2, -2, -2,  2, -2,  2, -2,  1,  1,  1,  1};

Eigen::Matrix<double, 16, 16> Ainv(Ainv_data);

Eigen::Matrix<double, 16, 1> f;
f.setRandom();
Eigen::Matrix<double, 16, 1> alpha = Ainv*f;
// This next line works, but it is making a copy, right?
Eigen::Matrix<double, 4, 4> a(alpha.data());
Solution

  • The last line is indeed doing a copy, so you can use a Map as follow:

    Map<Matrix4d,Eigen::Aligned> a(alpha.data());

    a behaves like a Matrix4d and it is read-write. The Eigen::Aligned flag tells Eigen that the pointer you pass to Map is properly aligned for vectorization. The only difference with a pure Matrix4d is that the C++ type is not the same.