Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Fastest method for calculating convolution

Anybody know about the fastest method for calculating convolution? Unfortunately the matrix which I deal with is very large (500x500x200) and if I use convn in MATLAB it takes a long time (I have to iterate this calculation in a nested loop). So, I used convolution with FFT and it is faster now. But, I am still looking for a faster method. Any idea?

like image 368
Nicole Avatar asked Dec 12 '13 21:12

Nicole


1 Answers

i have 2 way to calc fastconv

and 2 betther than 1

1- armadillo you can use armadillo library for calcing conv with this code

cx_vec signal(1024,fill::randn);
cx_vec code(300,fill::randn);
cx_vec ans = conv(signal,code);

2-use fftw ans sigpack and armadillo library for calcing fast conv in this way you must init fft of your code in constructor

FastConvolution::FastConvolution(cx_vec inpCode)
{
    filterCode = inpCode;
    fft_w = NULL;
}


cx_vec FastConvolution::filter(cx_vec inpData)
{
int length = inpData.size()+filterCode.size();
    if((length & (length - 1)) == 0)
    {

    }
    else
    {
        length = pow(2 , (int)log2(length) + 1);
    }
    if(length != fftCode.size())
        initCode(length);

    static cx_vec zeroPadedData;
    if(length!= zeroPadedData.size())
    {
        zeroPadedData.resize(length);
    }
    zeroPadedData.fill(0);
    zeroPadedData.subvec(0,inpData.size()-1) = inpData;


    cx_vec fftSignal = fft_w->fft_cx(zeroPadedData);
    cx_vec mullAns = fftSignal % fftCode;
    cx_vec ans = fft_w->ifft_cx(mullAns);
    return ans.subvec(filterCode.size(),inpData.size()+filterCode.size()-1);
}

void FastConvolution::initCode(int length)
{
    if(fft_w != NULL)
    {
        delete fft_w;
    }
    fft_w = new sp::FFTW(length,FFTW_ESTIMATE);
    cx_vec conjCode(length,fill::zeros);
    fftCode.resize(length);
    for(int i = 0; i < filterCode.size();i++)
    {
        conjCode.at(i) = filterCode.at(filterCode.size() - i - 1);
    }
    conjCode = conj(conjCode);
    fftCode = fft_w->fft_cx(conjCode);
}
like image 55
javad Avatar answered Oct 06 '22 04:10

javad