I'm transferring some code I have from Python/numpy to C++/Eigen. One feature that's pretty handy in numpy is the argsort method, which returns a vector of the indices that would sort a sequence of numbers. The following code works, but since I'm not experienced with C++ or Eigen I'm sure it could be more efficient. The code below can take a row or column of a matrix and return a vector such that the asPermutation() method will provide a permutation that can be used to sort the matrix by the row or column. If anyone has a better way of doing this, I'd love to see it.

A few other questions:

Is there any interest in some kind of wiki space where Eigen "recipes" could be shared?

Coefficient-wise functions all appear to be methods in Eigen. Out of curiosity, why not just use functions that accept arrays as arguments (sin(a), cos(a), etc)? Would there be interest in creating a namespace of coefficient-wise functions that, in usage, would be more akin to Matlab/numpy/Fortran?

Don't take my questions as complaints. Eigen is an amazing library, and once I get a bit more practice I'd love to contribute where I can.

Code: Select all

#include <iostream>
#include <vector>

#include <Eigen/Dense>

using namespace Eigen;

typedef std::pair<int, double> argsort_pair;

bool argsort_comp(const argsort_pair& left, const argsort_pair& right) {
    return left.second < right.second;
}

template<typename Derived>
VectorXi argsort(const MatrixBase<Derived> &x) {
    VectorXi indices(x.size());
    std::vector<argsort_pair> data(x.size());
    for(int i=0;i<x.size();i++) {
        data[i].first = i;
        data[i].second = x(i);
    }
    std::sort(data.begin(), data.end(), argsort_comp);
    for(int i=0;i<data.size();i++) {
        indices(data[i].first) = i;
    }   
    return indices;
}



int main() {
    MatrixXd A = MatrixXd::Random(4, 4);
    VectorXi sort_order;

    std::cout << "BEFORE SORTING" << std::endl << A << std::endl;

    sort_order = argsort(A.col(2));

    std::cout << "AFTER SORTING BY 3rd COLUMN" << std::endl << sort_order.asPermutation()*A << std::endl;
}

/Edited because asPermutation() expects ( i -> x[i] ), whereas numpy expects (x[i] -> i)