title: Using FFTW in C# to compute HilbertTransform
body: <p>I want to implement the Hilbert Transform in C#. From <a href="http://www.codeproject.com/Articles/1095473/Comparison-of-FFT-implementations-for-NET" rel="nofollow noreferrer">this</a> article I saw that the fastest FFT open source implementation seems to be the FFTW, so I downloaded that example and used it to learn how to use the fftw wrapper for C#.</p>
  <p>I have a current signal of 200.000 points which I'm using for testing. Getting the Hilbert transform through the fft is relatively simple:</p>
  <ol>
  <li>Compute the fft.</li>
  <li>Multiply by 2 all positive frequencies except for the DC and Nyquist components (0 and n/2 + 1, if the sample size is even).</li>
  <li>Multiply by 0 all the negative frequencies ([n/2 + 1, n]).</li>
  <li>Compute the inverse fft.</li>
  </ol>
  <p>This far, I've done all of it. The only problem is the inverse fft. I'm not able to get the same results with fftw than with the ifft from Matlab.</p>
  <h2>My code</h2>
  <pre><code>    RealArray _input;
      ComplexArray _fft;
  
      void ComputeFFT()
      {
          _fft = new ComplexArray(_length / 2 + 1);
          _input.Set(Data);
          _plan = Plan.Create1(_length, _input, _fft, Options.Estimate);
          _plan.Execute();
      }
  </code></pre>
  <p>This far, I've a fft with only the positive frequencies. So I don't need to multiply by zero the negative frequencies: they don't even exist. With the following code, I can get my original signal back:</p>
  <pre><code>    double[] ComputeIFFT(ComplexArray input)
      {
          double[] temp = new double[_length];
          RealArray output = new RealArray(_length);
  
          _plan = Plan.Create1(_length, input, output, Options.Estimate);
          _plan.Execute();
          temp = output.ToArray();
  
          for (int i = 0; i &lt; _length; ++i)
          {
              temp[i] /= _length;
          }
  
          return temp;
      }
  </code></pre>
  <p>The problem comes when I try to get a complex inverse from the signal.</p>
  <pre><code>    void ComputeHilbert()
      {
          double[] fft = FFT.ToArray();
          double[] h = new double[_length / 2 + 1];
          double[] temp = new double[_length * 2];
  
          bool fftLengthIsOdd = (_length | 1) == 1;
  
          h[0] = 1;
  
          for (int i = 1; i &lt; _length / 2; i++) h[i] = 2;
  
          if (!fftLengthIsOdd) h[(_length / 2)] = 1;
  
          for (int i = 0; i &lt;= _length / 2; i++)
          {
              temp[2 * i] = fft[2*i] * h[i];
              temp[2 * i + 1] = fft[2*i + 1] * h[i];
          }
  
          ComplexArray _tempHilbert = new ComplexArray(_length);
  
          _tempHilbert.Set(temp);
  
          _hilbert = ComputeIFFT(_tempHilbert);
          _hilbertComputed = true;
      }
  </code></pre>
  <p>It's important to note that, when I do apply the <code>ToArray()</code> method on a <code>ComplexArray</code> object, I get as result a <code>double[]</code> with twice as length as the original array, having the real and imaginary parts consecutive. That's it, for a <code>ComplexArray</code> object containing &quot;3 + 1i&quot;, I would get a double vector with [3, 1].</p>
  <p>So, at this moment, what I have is something like:</p>
  <p><code>[DC Frequency, 2*positive frequencies, Nyquist Frequency, zeros]</code>
  If I export this data to Matlab and compute the IFFT, I get the same result as its hilbert(signal).</p>
  <p>However, if I try to apply the IFFT provided by fftw, I get weird values from Nyquist Frequency to the end (that is to say, the zeros mess with fftw).</p>
  <p>This is the ifft I'm using to do this:</p>
  <pre><code>double[] ComputeIFFT(ComplexArray input)
  {
      double[] temp;
      ComplexArray output = new ComplexArray(_length);
  
      _plan = Plan.Create1(_length, input, output, Direction.Backward, Options.Estimate);
      _plan.Execute();
      temp = output.ToArray();
  
      for (int i = 0; i &lt; _length; ++i)
      {
          temp[i] /= _length;
      }
  
      return temp;
  }
  </code></pre>
  <p>So, just to sum it up, my problem is the way I'm using to calculate the ifft. It doesn't seems to work well with zeros. Or maybe Matlab is capable to understand that it has to apply some different approach and I should do it manually, but I don't know how.</p>
  <p>Thank you very much for your help in advance, much appreciated!</p>