index.html

信号处理系列导航

options.n = n;
[x,fs] = load_sound('bird', n);

%%
% You can actually play a sound.

sound(x,fs);

%% 
% We can display the sound.

clf;
plot(1:n,x);
axis('tight');
set_graphic_sizes([], 20);
title('Signal');

%%
% _Exercice 1:_ (the solution is <../private/audio_processing/exo1.m exo1.m>)
% Compute the local Fourier transform around a point |t0| of |x|, which is the FFT (use the
% function |fft|) of the windowed signal |x.*h| where |h| is smooth
% windowing function located around |t0|. For instance you can use for |h|
% a Gaussian bump centered at |t0|. To center the FFT for display, use
% |fftshift|.

exo1;

%%
% A good windowing function should balance both time localization and
% frequency localization.

t = linspace(-10,10,2048);
eta = 1e-5;
vmin = -2;

%%
% The block window has a sharp transition and thus a poor frequency
% localization.

h = abs(t)<1;
hf = fftshift(abs(fft(h)));
hf = log10(eta+hf); hf = hf/max(hf);
clf;
subplot(2,1,1);
title('Block window');
plot(t, h); axis([-2 2, -.1, 1.1]);
subplot(2,1,2);
plot(t, hf); axis([-2 2, vmin, 1.1]);
title('Fourier transform');

%%
% A Hamming window is smoother.

h = cos(t*pi/2) .* (abs(t)<1);
hf = fftshift(abs(fft(h)));
hf = log10(eta+hf); hf = hf/max(hf);
clf;
subplot(2,1,1);
title('Hamming window');
plot(t, h); axis([-2 2, -.1, 1.1]);
subplot(2,1,2);
plot(t, hf); axis([-2 2, vmin, 1.1]);
title('Fourier transform');


%%
% A Haning window has continuous derivatives.

h = (cos(t*pi)+1)/2 .* (abs(t)<1);
hf = fftshift(abs(fft(h)));
hf = log10(eta+hf); hf = hf/max(hf);
clf;
subplot(2,1,1);
title('Haning window');
plot(t, h); axis([-2 2, -.1, 1.1]);
subplot(2,1,2);
plot(t, hf); axis([-2 2, vmin, 1.1]);
title('Fourier transform');

%%
% A normalized Haning window has a sharper transition. It has the advantage
% of generating a tight frame STFT, and is used in the following.

h = sqrt(2)/2 * (1+cos(t*pi)) ./ sqrt( 1+cos(t*pi).^2 ) .* (abs(t)<1);
hf = fftshift(abs(fft(h)));
hf = log10(eta+hf); hf = hf/max(hf);
clf;
subplot(2,1,1);
title('Normalized Haning window');
plot(t, h); axis([-2 2, -.1, 1.1]);
subplot(2,1,2);
plot(t, hf); axis([-2 2, vmin, 1.1]);
title('Fourier transform');


%% Short time Fourier transform.
% Gathering a local Fourier transform at equispaced point create a local
% Fourier transform, also called *spectrogram*. By carefully chosing the
% window, this transform corresponds to the decomposition of the signal in
% a redundant tight frame. The redundancy corresponds to the overlap of the
% windows, and the tight frame corresponds to the fact that the
% pseudo-inverse is simply the transposed of the transform (it means that
% the same window can be used for synthesis with a simple summation of the
% reconstructed signal over each window).

%%
% We can compute a spectrogram of the sound to see its local Fourier
% content. The number of windowed used is |(n-noverlap)/(w-noverlap)|

% size of the window
w = 64*2;  
% overlap of the window
q = w/2;
S = perform_stft(x,w,q, options);

%%
% To see more clearly the evolution of the harmonics, we can display the
% spectrogram in log coordinates. The top of the spectrogram corresponds to
% low frequencies.

% display the spectrogram
clf; imageplot(abs(S)); axis('on');
% display log spectrogram
plot_spectrogram(S,x);

%%
% The STFT transform is decomposing the signal in a redundant tight frame.
% This can be checked by measuring the energy conservation.

% energy of the signal
e = norm(x,'fro').^2;
% energy of the coefficients
eS = norm(abs(S),'fro').^2;
disp(strcat(['Energy conservation (should be 1)=' num2str(e/eS)]));

%%
% One can also
% check that the inverse transform (which is just the transposed operator -
% it implements exactly the pseudo inverse) is working fine.

% one must give the signal size for the reconstruction
x1 = perform_stft(S,w,q, options);
disp(strcat(['Reconstruction error (should be 0)=' num2str( norm(x-x1, 'fro')./norm(x,'fro') ) ]));

%% Audio Denoising
% One can perform denosing by a non-linear thresholding over the
% transfomede Fourier domain.

%%
% First we create a noisy signal
sigma = .2;
xn = x + randn(size(x))*sigma;

%%
% Play the noisy sound.

sound(xn,fs);


%%
% Display the Sounds.

clf;
subplot(2,1,1);
plot(x); axis([1 n -1.2 1.2]);
set_graphic_sizes([], 20);
title('Original signal');
subplot(2,1,2);
plot(xn); axis([1 n -1.2 1.2]);
set_graphic_sizes([], 20);
title('Noisy signal');

%%
% One can threshold the spectrogram.

% perform thresholding
Sn = perform_stft(xn,w,q, options);
SnT = perform_thresholding(Sn, 2*sigma, 'hard');
% display the results
subplot(2,1,1);
plot_spectrogram(Sn);
subplot(2,1,2);
plot_spectrogram(SnT);

%%
% _Exercice 2:_ (the solution is <../private/audio_processing/exo2.m exo2.m>)
% A denoising is performed by hard or soft thresholding the STFT of the
% noisy signal. Compute the denosing SNR with both soft and hard
% thresholding, and compute the threshold that minimize the SNR. Remember that a soft thresholding
% should be approximately twice smaller than a hard thresholding. Check the
% result by listening. What can you conclude about the quality of the
% denoised signal ?

exo2;

%%
% _Exercice 3:_ (the solution is <../private/audio_processing/exo3.m exo3.m>)
% Display and hear the results. What do you notice ?

exo3;

%% Audio Block Thresholding
% It is possible to remove musical noise by thresholding blocks of STFT
% coefficients.

%%
% Denoising is performed by block soft thresholding.

% perform thresholding
Sn = perform_stft(xn,w,q, options);
SnT = perform_thresholding(Sn, sigma, 'block');
% display the results
subplot(2,1,1);
plot_spectrogram(Sn);
subplot(2,1,2);
plot_spectrogram(SnT);

%%
% _Exercice 4:_ (the solution is <../private/audio_processing/exo4.m exo4.m>)
% Trie for various block sizes and report the best results.

exo4;


##### SOURCE END #####
-->
   </body>
</html>