% Section of the single-sided amplitude spectrum of the "noise", which
% is brought into the audio file by watermarking and the new quantization
% of the samples as a consequence to it, i.e. the spectrum of the
% difference signal
% R. Brigola, April 2010
clear all;
start_time=cputime;
% File names
audio_input1='your_file.wav'; % original audio file
audio_input2='your_file_marked.wav'; % watermarked version
% Choice of the frequency section, which shall be considered
% Since we know from watermark embedding the frequency section, where
% we can possibly find a trace of the watermark in the spectra
Hz1=0;
Hz2=200;
% Choice of a blocksize and number of blocks for the fft
blocksize=8820;
no_of_blocks=100;
no_of_samples=blocksize*no_of_blocks;
% Choice of a start block and an audio channel for the fft
start_block=31;
channel=1;
[y1,Fs,nbits]=wavread(audio_input1);
data1=y1(start_block:start_block+no_of_samples-1,channel);
y2=wavread(audio_input2);
data2=y2(start_block:start_block+no_of_samples-1,channel);
data=data2-data1;
% We use the Hann window to reduce spectral leakage in the fft
for k=1:no_of_samples
data(k)=data(k)*(0.5-0.5*cos(2*pi*(k-1)/no_of_samples));
end
duration=no_of_samples/Fs;
peak1=floor(Hz1*duration);
peak2=floor(Hz2*duration);
% see the index and frequency correspondance in Matlab's fft routine
f=1/duration*linspace(peak1,peak2,peak2-peak1+1);
% Now the FFT
Data=fft(data);
% Plot single-sided amplitude spectrum of the difference signal
figure(1);plot(f,2*abs(Data(peak1+1:peak2+1))/no_of_samples);
title('Single-Sided Amplitude Spectrum of the Difference Signal')
xlabel('Frequency in Hz')
ylabel('Amplitude')
computation_time=cputime - start_time;