Search code examples
matlabeigeneigen3

Eigen: how can I substitute matrix positive values with 1 and 0 otherwise?


I want to write the following matlab code in Eigen (where K is pxp and W is pxb):

H = (K*W)>0;

However the only thing that I came up so far is:

H = ((K*W.array() > 0).select(1,0)); 

This code doesn't work as explained here, but replacing 0 with VectorXd::Constant(p,0) (as suggested in the link question) generates a runtime error:

Eigen::internal::variable_if_dynamic<T, Value>::variable_if_dynamic(T) [with T = long int; int Value = 1]: Assertion `v == T(Value)' failed.

How can I solve this?


Solution

  • You don't need .select(). You just need to cast an array of bool to an array of H's component type.

    H = ((K * W).array() > 0.0).cast<double>();
    

    Your original attempt failed because the size of your constant 1/0 array is not match with the size of H. Using VectorXd::Constant is not a good choice when H is MatrixXd. You also have a problem with parentheses. I think you want * rather than .* in matlab notation.

    #include <iostream>
    #include <Eigen/Eigen>
    using namespace Eigen;
    
    int main() {
      const int p = 5;
      const int b = 10;
      MatrixXd H(p, b), K(p, p), W(p, b);
      K.setRandom();
      W.setRandom();
      H = ((K * W).array() > 0.0).cast<double>();
      std::cout << H << std::endl << std::endl;
    
      H = ((K * W).array() > 0).select(MatrixXd::Constant(p, b, 1),
                                       MatrixXd::Constant(p, b, 0));
      std::cout << H << std::endl;
      return 0;
    }
    

    When calling a template member function in a template, you need to use the template keyword.

    #include <iostream>
    #include <Eigen/Eigen>
    using namespace Eigen;
    
    template<typename Mat, typename Vec>
    void createHashTable(const Mat &K, Eigen::MatrixXi &H, Mat &W, int b) {
      Mat CK = K;
      H = ((CK * W).array() > 0.0).template cast<int>();
    }
    
    int main() {
      const int p = 5;
      const int b = 10;
      Eigen::MatrixXi H(p, b);
      Eigen::MatrixXf W(p, b), K(p, p);
      K.setRandom();
      W.setRandom();
      createHashTable<Eigen::MatrixXf, Eigen::VectorXf>(K, H, W, b);
      std::cout << H << std::endl;
      return 0;
    }
    

    See this for some explanation.

    Issue casting C++ Eigen::Matrix types via templates