I have array data (128 channels) that I'd like to filter using an FIR filter. For any set of samples, if my data is arranged like this:

Code: Select all

   x1(t-3)   x1(t-2)   x1(t-1)   x1(t)
   x2(t-3)   x2(t-2)   x2(t-1)   x2(t)
D=                ...
                  ...
   x128(t-3) x128(t-2) x128(t-1) x128)t)

and I have an FIR filter

Code: Select all

f = [f1 f2 f3 f4]

, I can compute the filtered value for time y for all 128 channels at once as

Code: Select all

y = D * f'

Then, x(t) becomes x(t-1), and a new value is copied for x(t).

And since eigen wraps all my nice BLAS/linpack matrix vector stuff in C++ code, I can do this using the SSE stuff without having to delve into that myself.

Here's the problem: I have to use a circular buffer for the data - at some point, x(t) has to be copied into element one of the matrix, and the addressing becomes wanky. I feel like the aforementioned primitives incorporate modulo addressing - does Eigen have any facility for accessing that?

thanks!

If I understod correctly, you would like to be able to view, e.g.,[x(t-3) x(t) x(t-1) x(t-2)] as [x(t) x(t-1) x(t-2) x(t-3)], right? If so, that's not really possible at the moment. BLAS and Lapack cannot do that neither and that would significantly reduce the performance of the matrix-vector product. I see two solutions:
- you shift (copy) the entire matrix such that you always have [x(t) x(t-1) x(t-2) x(t-3)]
- you only shift the coefficient of the vector f and use an indirect table lookup to access the columns of X.

In a near future you will be able to do something like X(All,Vector4i(1, 2, 3, 0)) to facilitate this.

In the past with FIR code I have found the following to be useful:
http://www.mikrocontroller.net/wikifiles/1/1e/FirAlgs.c
http://www.dspguru.com/dsp/faqs/fir/implementation

Some of the methods avoid using a circular buffer.

I have not yet implemented in Eigen.

To be a bit more clear, in response to ggael, what I was hoping for was a circular shift as part of a vector multiply so the ability to do the product

Code: Select all

z = x * y'

where in memory, x is stored as:

Code: Select all

x ~ [x4 x0 x1 x2 x3]

but why is stored in order as:

Code: Select all

y ~ [y0 y1 y2 y3 y4]

So the addressing of x would be modulo 5 in this example. I believe some of the primitives allow for modulo addressing (which is easy when vector lengths are a power of two).


Thanks kp0987 - I hadn't thought about duplicating the coefficients. I think that will work with minimal effort. What's awesome about Eigen is that since I can get a pointer to the memory storage, I can make my matrix slide along the circular data buffer as new samples come in. When the data wraps, all I need to do is to slide the coefficient vector in the same way...

thanks!

I have a similar problem where I would like to keep a moving cache of values in phi_. I am not convinced that my method is particularly fast. Is it possible to be recommended code that could be faster. Thanks!

Code: Select all

chi_= MatrixXd::Constant(100, phi_.cols(),0);
tChi_= MatrixXd::Constant100, phi_.cols(),0);

updateChi()
{
      tChi_.row(0)= phi_;
      tChi_.bottomRows(tChi_.rows()-1)= chi_.topRows(chi_.rows()-1);
      chi_= tChi_;               
}