Convolution

 

Follow these steps to gain understanding of the CONV (convolution) function.

 

  1. Set up an empty Worksheet with 9 windows.

  2. Move to W1 and generate a square wave with the following command:

     W1: gsqrwave(64, 1/64, 2)


    The first argument to GSQRWAVE is the number of points to be generated. The number 64 was chosen because a series length equal to a power of 2 will require less calculation time for a Fourier Transform. The second argument is the delta-x interval, which is 1/64 to make sure there are complete periods in the series. The third argument, frequency, is 2 so that the series contains two periods.

  3. Now move to W2 and convolve window 1 with itself. The command is:

     W2: conv(W1, W1)

  4. CONV computes convolution directly in the time domain. Let's try to perform the same function using the discrete Fourier transform. Take the FFT of the series, multiply it by itself, and then isolate the real part of the inverse transform.

  5. Move to window 3. Place a copy of the original square wave here:

     W3: W1

  6. Now move to window 4. Take the discrete Fourier transform of W4 and square it:

     W4: fft(w3)*fft(w3)

  7. The horizontal units of the resultant series is frequency. This series is actually the frequency spectrum of the original series.

  8. Now, move to window 5 and take the real part of the inverse transform.

     W5: real(ifft(w4))

 

But wait! This does not look like the original convolution. What went wrong?

 

We did not perform the calculation correctly. Not only do the two convolutions look different, their lengths are different. Use LENGTH or press the I button or [F2] key to find the number of points in a series.

The correct way to compare the FFT to the convolution is to extend the original square wave to the proper length before taking the FFT. To do this, use EXTRACT to extend the square wave to 128 points by adding zeros to the end of the series.

 

  1. Put a zero padded copy of the W1 series in window 6:

     W6: extract(w1, 1, 2*length(w1))

  2. This creates a 128-point series; since the original series had only 64 points, the remainder is padded with zeros. Now repeat the operation we performed above in windows 7 and 8.

  3. In window 7:

     W7: fft(w6)*fft(w6)

  4. In window 8:

     W8: real(ifft(w7))


     Now the result looks like it should.

  5. To check the difference between the two calculations, move to window 9:

     W9: W8 - W2

 

 

The error is quite small.

 

The FFT function accepts a length argument to automatically zero pad the input. The convolution procedure could be accomplished with:

 

real(ifft(fft(w1, 2*length(w1)) * fft(w1, 2*length(w1))))

 

The FCONV function performs all of the above steps to perform frequency domain convolution:

 

fconv(w1, w1)

 

Moving averages (and the MOVAVG function) are based on the concept of convolution. With convolution and proper scaling, any variation on moving averages can be achieved.