I am new to Eigen (using Eigen3) and would like to ask for some advice in designing my application.

I need to work with a new data type consisting of 2 complex numbers (say A and B) which can be stored in fairly large arrays (say 1024-32k rows). I need to perform multiple operations on these, including FFTs. I was hoping to use the unsupported FFT module (with FFTW backend).

My need for FFTs initially led me to use 2 Arrays for A and B so that the FFT can be efficiently performed on each (there is no need for 2D FFT or any couping from A to B). However for other operations this leads to either passing A&B separately into functions for efficient working, or combining by wrapping in a class, which means all the nice Eigen efficiency is lost. (Incidentally is there a good reference for best practices to wrap Eigen Arrays/Matrices in a class/struct?)

I notice though in FFTW docs that there is support for strided storage which suggests to me that I should instead develop a new Scalar class (much along the lines of "complex" but with A&B instead of Real&Imag). By using striding I should still be able to have an efficient FFT and also be able to use Arrays of the new Scalar in function arguments.

Does this seem like a good approach?

Some comments/questions on the unsupported FFT module (with FFTW backend). I notice that it does not seem to support Array (only Matrix) and in particular it looks like I need to define a temporary matrix for the inverse FFT and then convert back into an array. Presumably this is because it is based on MatrixBase rather than DenseBase? Also, it seems to only support striding via a temporary matrix, which I guess is a requirement for the kissfft since it is not for fftw.

I have read (among others):
http://eigen.tuxfamily.org/dox-devel/Cu ... Eigen.html
http://eigen.tuxfamily.org/dox-devel/To ... Types.html

Thanks for your help.

If you need a pair of complex numbers as scalar type, have you tried Array2cf / Array2cd ? It should be readily usable as scalar type.

kp0987 wrote:Some comments/questions on the unsupported FFT module (with FFTW backend). I notice that it does not seem to support Array (only Matrix) and in particular it looks like I need to define a temporary matrix for the inverse FFT and then convert back into an array. Presumably this is because it is based on MatrixBase rather than DenseBase?

Exactly. We should fix this.

Also, it seems to only support striding via a temporary matrix, which I guess is a requirement for the kissfft since it is not for fftw.

Interesting.

You should mail the FFT module author since I don't think he reads the forum. Find his email address in the source code.

Thanks for the rapid response. I had not seen that it is possible for another matrix to act like a scalar. What are the limits on this?

I tried the following:

Code: Select all

  Array<Array2cd, 2, 1> tmp;
  cout << tmp << endl;

The first line works fine, but the second errors with error: 'epsilon' is not a member of 'Eigen::NumTraits<Eigen::Array<double, 2, 1, 0, 2, 1> >'. I assume this is just due to the printing via cout needing ceil/epsilon etc which are not defined, but it is not a real problem. Or should I be worried?

I will try both options, but right now am erring to the side of a custom scalar. If I use Array2cd then I'll probably want to allow easy access to A and B (rather than tmp(0) and tmp(1)). I would use MatrixBaseAddons.h or DenseBaseAddons.h I presume to do this? Do you think the work involved in implementing either option will be much different in the end?

I will contact the FFT author.

Thanks.

Only Array can act as scalar type, Matrix can't.

The error you're getting while printing is a bug on our side. There should be no such restriction on using arrays as scalar type.

Indeed you can add accessors using EIGEN_MATRIXBASE_PLUGIN or EIGEN_DENSEBASE_PLUGIN. Though note that you already have x() and y() accessors for the 2 first coeffs. That's the solution with the fewest lines of code.

Thanks. The x(), y() accessors are useful here.

I have a question about how to efficiently access parts of the scalar. Say I have:

Code: Select all

Array<Array2cd, Dynamic, 1> tmp;

where I use Array2cd to store 2 complex numbers A&B, and then form an array of these.

How can I efficiently do a math operation on all of A or B (eg .abs())? I originally thought something like tmp.x().abs() would work, but of course that is just the same as tmp(0).abs() and not tmp(0...n).x().abs() which is what I need.

In pseudocode I am really wanting to do A.abs().square() + B.abs().square(). Maybe I am missing something obvious.

It seems that the use-Array-as-scalar-type feature is not quite polished enough for your needs. So here's another idea instead. A Vector of Array2cd is like a 2xN matrix of complex numbers. Since you want each pair of complex numbers to be stored contiguously, and default Eigen storage order is column-major, let's make each pair be its own column. So use this type:

Matrix<complex<double>, 2, Dynamic>

as your "vector of pairs" type.

Then "vec[index].A" is just matrix(0, index) and "vec[index].B" is matrix(1, index).

You can then address the expression of all A's as matrix.row(0), and the expression of all B's as matrix.row(1).

If you want |A|^2 + |B|^2 you get it as matrix.col(i).squaredNorm().

If you want the sum of all of these, it's just matrix.squaredNorm()

Of course you can use Array instead of Matrix. It's just that Array doesn't have a squaredNorm() method, but you can convert at no cost with the .matrix() method:

array.matrix().squaredNorm()

Thanks. I am trying a modification of your suggestion by just using a Nx2 Array of A and B. I don't really need each pair of complex numbers to be stored contiguously (in fact it causes a complication on the FFT which would need to be solved by adding stride support), but originally liked the idea of using a scalar. I added read/write accessors for A&B in a DenseBaseAddons file using a copy of the col function in BlockMethods. This all seems to be working OK.