Computing analytical signal using FFT in C++

Question

I'm currently lost trying to figure out how to implement an equivalent version of MATLAB's hilbert() function in C++. I'm very new to signal processing, but, ultimately, I would like to figure out a way to phase shift any given signal by 90 degrees. I was attempting to follow the method suggested in this question on MATLAB central, which appears to work based on tests using GNU Octave.

I have what I believe to be a working implementation of both FFT and the inverse FFT, and I have tried implementing the method described in the answer to this post in order to compute the analytical signal. I have tried doing this by applying the FFT, setting the upper half of the array to zero, and then applying the inverse FFT, but, based on graphs I made of output from a test, there must be a problem with the way I have implemented finding the analytical signal.

What would be a suitable way to implement the hilbert() function from MATLAB in C++ given a working implementation of FFT and inverse FFT? Is there a better way of achieving the 90 degree phase shift?

I'm not sure about the hilbert function in MATLAB but given complex coefficient in your FFT `a + i*b`, you can perform a 90 degree phase shift by replacing it with `-b + i*a` — jodag, Aug 19 '16 at 03:46
Also, if you have MATLAB you can type `edit hilbert` to see the implementation they use. It looks pretty straight forward once you have the FFT and inverse FFT implemented. — jodag, Aug 19 '16 at 04:03
@jodag Unfortunately, replacing `a + i*b` with `-b + i*a` doesn't seem to be having the desired effect. Maybe the implementation of FFT or inverse FFT I'm using is wrong? I'm actually getting back the same data from before applying the FFT. Also, I do not have MATLAB. — fakedad, Aug 19 '16 at 04:17
Yeah it looks like its a little different than that. I'll look through it and get back to you — jodag, Aug 19 '16 at 04:18
Did you really mean "FTT" all the time, or is it a consequent typo of "FFT"? — Matthias W., Aug 19 '16 at 07:52

jodag · Accepted Answer · 2016-08-19T04:48:16.633

Checking the MATLAB implementation the following should return the same result as the hilbert function. You'll obviously have to modify it to match your specific implementation. I'm assuming a signal class of some sort exists.

signal hilbert(const signal &x)
{
    int limit1, limit2;
    signal xfreq = fft(x);
    if (x.numel % 2 == 0) {
        limit1 = x.numel/2;
        limit2 = limit1 + 1;
    } else {
        limit1 = (x.numel + 1)/2;
        limit2 = limit1;
    }
    // multiply the first half by 2 (except the first element)
    for (int i = 1; i < limit1; ++i) {
        xfreq[i].real *= 2;
        xfreq[i].imag *= 2;
    }
    for (int i = limit2; i < x.numel; ++i) {
        xfreq[i].real = 0;
        xfreq[i].imag = 0;
    }
    return ifft(xfreq);
}

Edit: Forgot to set the second half to zeros. Edit2: Fixed a logical error. I coded the following up in MATLAB which matches hilbert.

function h = hil(x)
    n = numel(x);
    if (mod(n,2) == 0)
        limit1 = n/2;
        limit2 = limit1 + 2;
    else
        limit1 = (n+1)/2;
        limit2 = limit1+1;
    end

    xfreq = fft(x);

    for i = 2:limit1
        xfreq(i) = xfreq(i)*2;
    end
    for i = limit2:n
        xfreq(i) = 0;
    end

    h = ifft(xfreq);
end

This is how to do it: implement the algorithm in a C style in Matlab, check that you’re getting the same answer as the built-in function, then port it to actual C (or Rust or whatever). — Ahmed Fasih, Aug 19 '16 at 05:20
For some reason I'm having trouble getting this to work after switching to using FFTW. Any suggestions? — fakedad, Sep 10 '16 at 04:57
If you had it working using a previous FFT implementation then you should probably create a simple case and see how FFTW differs from the previous implementation. MATLAB uses FFTW on the backend so I'm sure you can make it work using FFTW. — jodag, Sep 10 '16 at 05:10

Computing analytical signal using FFT in C++

1 Answers1