0
$\begingroup$

There are two related questions with good answers; for calculating the forward FFT in N/2 values, and calculating inverse FFT from N real values, now I'm in search for the missing link of packing the N/2+1 element complex signal X back to the N/2 element complex signal that it was derived in the first place for memory and arithmetic efficient inverse FFT.

Having calculated N/2+1 valued complex FFT from N real values (code from the first link):

x = [ 1 2 3 4 5 6 7 8 ]; n = 8; n2 = 4;
Z = fft(x(1:2:end) + i * x(2:2:end);
Ze = .5*( Z + conj([Z(1),Z(n2:-1:2)]));
Zo = -.5*j*( Z - conj([Z(1),Z(n2:-1:2)]));
X = [Ze,Ze(1)] + exp(-j*2*pi/n*(0:n2)).*[Zo,Zo(1)];

I have actually packed the real valued element X(n2+1) to imaginary part of X(1), saving some memory. Then I've point wise multiplied X with G representing the FFT of a same sized real signal g.

In order to finish the (circular) convolution of conv(x, g), I'd need to calculate the packed N/2 sized complex signal Z from (X.*G) for IFFT.

The first coefficient can be reverted as Z(1) = real(X(1)+X(n2+1))/2 + i*(imag(X(1)-X(n2+1)))/2;

As the other values X(2:n2-1) are originally produced from a bufferfly operation of x(k) +- exp(j*k/N) * x(n2-n), I would hope that all the elements would be recoverable. The question is how.

EDIT

I need to recheck the formulas, but since for a pair A=Z(k+1), B=Z(n2+1-k):

                                         [ 1-s 1+s   c     c ] [Real(A)]
 [Real(A') Imag(A') Real(B') Imag(B')] = [ 1+s 1-s  -c    -c ]*[Imag(A)]
                                         [ -c   c    1-s -1-s] [Real(B)]
                                         [ -c   c   -1-s  1-s] [Imag(B)]

with invertible matrix, there really is a solution.

Improve this question
$\endgroup$
5
  • $\begingroup$ Your code would be lot easier you to read if you add a few comments. $\endgroup$ – Hilmar Apr 8 '20 at 13:56
  • $\begingroup$ Let me fully understand your question: Are you asking how to compute the circular convolution of two real sequences efficiently using the relationship $xcorr = ifft(fft(a)*fft(b))$ where a and b are real sequences and specifically fft(a) and fft(b) are computed using N/2 complex values as given by the first link you gave? Also in your case is N = 64 (It's confusing how the index 32 and 33 show up as it's no where else in the code you gave or your description)? $\endgroup$ – Dan Boschen Apr 8 '20 at 14:55
  • $\begingroup$ @DanBoschen: Yes, that's what I'm after. And additionally I want to use the N/2 IFFT, which is just one element shorter than the X in the Matlab script. (My N might indeed be 64, but it's now edited, thanks for pointing out...) $\endgroup$ – Aki Suihkonen Apr 8 '20 at 16:57
  • $\begingroup$ please learn to use $\LaTeX$ in your mathematical expression. $\endgroup$ – robert bristow-johnson Apr 8 '20 at 19:02
  • $\begingroup$ just to clarify in my comment above it should be $xcorr = ifft(fft(a)*conj(fft(b))$ I was thinking they were real so showing the conjugate was not needed but the fft's are indeed complex! $\endgroup$ – Dan Boschen Apr 9 '20 at 3:08
1
$\begingroup$

This is the general solution to recover N real samples using an N/2 Length IFFT. N/2+1 samples of the FFT are required for an even sequence N.

First generating a sample sequence of the N/2+1 samples of an FFT of a real signal:

x = [1 2 3 4 5 6 7 8];
N = length(x);
N2 = N/2;
fx = fft(x);
X = fx(1:N2+1);             % starting result: positive spectrum of real sequence

So above is the same result as can be alternatively computed with an N/2 FFT as already demonstrated.

Given such an result, below shows how to get the original sequence x using an N/2 IFFT. The code is self-explantory:

W = @(k,N) exp(-j*2*pi/N*[0:k-1]);                   % FFT Twiddle factors
Xe =  (X(1:N2) + conj(X(N2+1:-1:2)))/2;              % FFT of even samples
Xo = (X(1:N2) - conj(X(N2+1:-1:2)))./(2*W(N2,N));    % FFT of odd samples

newX = Xe + j*Xo;                      

xr = ifft(newX);                      % N/2 IFFT 
newxo =  real(xr)                     % recovered odd x samples
newxe =  imag(xr)                     % recovered even x samples

Since each FFT in your formula for cross correlation is the FFT of a real signal, the product will also be the FFT of a real signal, hence the technique outlined above can be used.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.