Visualising the convolution theorem.

In terms of understanding the convolution theorem i've created the
following MatLab code:

[ BEGIN CODE]
function checkconvolutiontheorem

[data Fs] = wavread('pop.wav'); % a 5 second sample of pop-music

signal = data(1:453); % Signal consists of 453 non-zero elements.
kernel = data(1:60);  % Kernel consists of 60 non-zero elements.

timesignal = conv(signal,kernel);   % Convolution in spatial domain.
plot(1:256, timesignal(1:256),'b');
hold on;

h = zeros(1,512);
h(1:453) = signal;
signal = h;

k = zeros(1,512);
k(1:60) = kernel;
kernel = k;

fs = fft(signal,512); % Convert signal to frequencydomain.
fk = fft(kernel,512); % Convert kernel to frequencydomain.

% Use polar notation to calculate Re and Im of the pointwise product.
ReY = real(fs(1,:)).*real(fk(1,:)) - imag(fs(1,:)).*imag(fk(1,:));
ImY = imag(fs(1,:)).*real(fk(1,:)) + real(fs(1,:)).*imag(fk(1,:));


convertedrealsignal = ifft(ReY,512); % Convert back to spatial domain
convertedimagsignal = ifft(ImY,512); % Convert back to spatial domain.

restructeredsignal = convertedrealsignal+abs(convertedimagsignal);

plot(1:256, restructeredsignal(1:256), 'r');
hold off;
end
[END CODE]

The signal and kernel are chosen so that they convolution fits into a
512 sample period. My problem is that the two plots arent identical.
Can anyone correct my mistake?

Thanks.

On 5 Des, 23:04, Casper.Pt...@gmail.com wrote:
> In terms of understanding the convolution theorem i've created the
> following MatLab code:
>
> [ BEGIN CODE]
> function checkconvolutiontheorem
>
> [data Fs] = wavread('pop.wav'); % a 5 second sample of pop-music

First, use random data or whatever, so others who don't have
this wav file can check your code.

> signal = data(1:453); % Signal consists of 453 non-zero elements.
> kernel = data(1:60); �% Kernel consists of 60 non-zero elements.
>
> timesignal = conv(signal,kernel); � % Convolution in spatial domain.
> plot(1:256, timesignal(1:256),'b');
> hold on;
>
> h = zeros(1,512);
> h(1:453) = signal;
> signal = h;
>
> k = zeros(1,512);
> k(1:60) = kernel;
> kernel = k;
>
> fs = fft(signal,512); % Convert signal to frequencydomain.
> fk = fft(kernel,512); % Convert kernel to frequencydomain.

The lines below are wrong:

> % Use polar notation to calculate Re and Im of the pointwise product.
> ReY = real(fs(1,:)).*real(fk(1,:)) - imag(fs(1,:)).*imag(fk(1,:));
> ImY = imag(fs(1,:)).*real(fk(1,:)) + real(fs(1,:)).*imag(fk(1,:));
>
> convertedrealsignal = ifft(ReY,512); % Convert back to spatial domain
> convertedimagsignal = ifft(ImY,512); % Convert back to spatial domain.

Instead use

Y = fs.*fk;
restructuredsignal = real(ifft(Y));

> % restructeredsignal = convertedrealsignal+abs(convertedimagsignal);
>
> plot(1:256, restructeredsignal(1:256), 'r');
> hold off;
> end
> [END CODE]
>
> The signal and kernel are chosen so that they convolution fits into a
> 512 sample period. My problem is that the two plots arent identical.
> Can anyone correct my mistake?

Rune

Casper,

Rune correctly told you how to fix it and write it efficiently.
However, he has identified a couple of lines as wrong that are not
wrong, but how you used them later is  wrong.

Since you appear to understand (or be trying to understand) the
operations in a particular manner, I have made minimal changes to fix
your code to make it correct. Again, Rune's way is correct and more
efficient, but I think it helps to understand where you went wrong.
Note comments denoted by %DB.

Dirk

==========================

function checkconvolutiontheorem

%DB Used data I could get other than 'pop.wav'
%DB [data Fs] = wavread('pop.wav'); % a 5 second sample of pop-music
data =randn(1,512);

signal = data(1:453); % Signal consists of 453 non-zero elements.
kernel = data(1:60);  % Kernel consists of 60 non-zero elements.

timesignal = conv(signal,kernel);   % Convolution in spatial domain.

plot(timesignal,'b');  %DB changed to plot whole thing, optional
hold on;

%DB h = zeros(1,512);      %DB Removed, already set, optional
%DB h(1:453) = signal;
%DB signal = h;

%DB k = zeros(1,512);      %DB Removed, already set, optional
%DB k(1:60) = kernel;
%DB kernel = k;

fs = fft(signal,512); % Convert signal to frequencydomain.
fk = fft(kernel,512); % Convert kernel to frequencydomain.

% Use polar notation to calculate Re and Im of the pointwise product.
%DB "Polar"? Lines are okay.
ReY = real(fs(1,:)).*real(fk(1,:)) - imag(fs(1,:)).*imag(fk(1,:));
ImY = imag(fs(1,:)).*real(fk(1,:)) + real(fs(1,:)).*imag(fk(1,:));

%DB Removed 3 lines, incorrect
%DB "restructeredsignal" was probably meant to be
"reconstructedsignal",
%DB it's only a name, but it means something.

%DB convertedrealsignal = ifft(ReY,512); % Convert back to spatial
domain
%DB convertedimagsignal = ifft(ImY,512); % Convert back to spatial
domain.
%DB restructeredsignal = convertedrealsignal+abs
(convertedimagsignal);

%DB Added, correct, !THIS IS WHAT YOU MISSED!
restructeredsignal=ifft(ReY+j*ImY,512);

plot(restructeredsignal, 'r'); %DB changed to plot whole thing,
optional

%DB Added to verify results, should be =0, optional
error=timesignal-restructeredsignal;
plot(error,'k')

hold off;
end

Thank you both for your responses (and for pointing out efficiency
issues). I greatly appreciate it.