Bonjour afin de compléter le sujet voici le code source , tout d'abord la fonction principale puis wiener.m qui est utilisée ds le fichier principal
clear all;
close all;
[ y, fs, nbits ] = wavread( 'Blue20s.wav' );
hist(y(:,2),1000);
title('gaussienne distribution');
Ts = 1/fs
t=0:Ts:0.5;
left = y(:,1)';
right = y(:,2)';
figure;
subplot(2,1,1);
plot( t, left(1:length(t)) );
title( 'left ear' );
subplot(2,1,2);
plot( t, right(1:length(t)) );
title( 'right ear' );
%La s´eparation de sources s’effectue dans le domaine frequentiel. On represente les
%signaux par leur TFCT dont le calcul peut se decomposer en 3 etapes:
%decoupage avec un recouvrement de 1024 points soit 50% des points de la FFT
fs = 1/Ts;
N=2048; % nombre de points pris pour les fonctions à étuder
t=[1:N]/fs;
%parametre fct
NFFT = 2048; %nombre de points pris pour la FFT
fe = 8000; %frequence d'echantillonnage du signal
NBR = 1024; % nombre de points de recouvrement
fenetre = hamming(NFFT); % fenetre de recouvrement
s =[left' right']; %signal a etudie
t=1:length(s);
figure;
subplot(211);
plot(t, s);
title('concatenation de left et right');
subplot(212);
plot(fenetre);
figure;
i = 1;
IndInf = 1;
IndSup = NFFT;
while IndSup<length(s)
ss = s(IndInf : IndSup);
ss = ss.*fenetre'; %multiplication des T trames par une fenetre de ponderation
X(i,:) = abs( fft(ss) ); %calcul de la Transformee de Fourier Discrete (TFD) sur N points de chacune
%d’elles
IndInf = IndInf + NFFT - NBR;
IndSup = IndSup + NFFT - NBR;
i=i+1;
end
%surf(X(:, 1: NFFT/2));
[L, c] = size(X);
mesh( [1:L], [1:NFFT/2]*fe/NFFT, X(:, 1:NFFT/2)');
title('representation de la TFCT ');
wiener(y,4,'carre',0.01,0.5);
while IndSup<length(s)
ss = y(IndInf : IndSup);
ss = ss.*fenetre'; %multiplication des T trames par une fenetre de ponderation
X(i,:) = abs( ifft(ss) ); %calcul de la Transformee de Fourier Discrete (TFD) sur N points de chacune
%d’elles
IndInf = IndInf + NFFT - NBR;
IndSup = IndSup + NFFT - NBR;
i=i+1;
end
wavwrite(X,'Blue20s.wav' );
////////////// wienr.m
function wiener(y,sbox,type,sigma,alpha);
[L,c]=size(y); % Taille de l'image
% Convolution
box = zeros (L,c);
switch type
case 'carre'
box(1:sbox,1:sbox) = 1;
end
box = box / sum(box(:));
Yy = fft2(y);
Hh = fft2(box,L,c);
Cc = Yy.*Hh;
co = real(ifft2(Cc));
% Calcul du bruit
sigma = sqrt(sigma);
noisy = (co + sigma*randn(size(co)) + 0);
noisy = max(0,min(noisy,1));
Gg = fft2(noisy);
% Filtre inverse
%Ffinv=Gg./([abs(Hh)<1e-1].*1e-1+Hh);
%finv=im2uint8(abs(ifft2(Ffinv)));
% Filtre de Wiener
H2=abs(Hh).^2;
Ffwin=H2.*Gg./([(H2+sigma).*Hh<1e-14].*1e-14+((H2+sigma).*Hh));
fwin=abs(ifft(Ffwin));
figure;
plot(fwin);
title('filtre de wiener');
MERCI
