Speech Processing, Proj06

//Purdue Cal - ECE 595C.
//Spring 2010.
//Project 06.
//Copyright @ 2010 antonio081014 ;
//All codes are in Matlab.

Report Link:

Task:
1. Speech enhancement by using spectral subtraction;




%% Proj06;
clear all;
clc;
close all;

fd = fopen('chj.16', 'rb', 'b');
[spch, K] = fread(fd, 'int16');
fclose(fd);

% Remove the explosive in the background;
spch(28080:28100) = 0;
sspch = (spch - mean(spch));
y = sspch /max(sspch);

fs = 8000;
FL = 20 * fs / 1000; % 20ms;
NF = 1024;
frame = zeros(4, FL);
frame(1,:) = y(2241:2400);
frame(2,:) = y(6001:6160);
frame(3,:) = y(14912:15071);
frame(4,:) = y(33226:33385);

STFT = zeros(4, NF); % short-time spectrum;
STauto = zeros(4, 2 * FL -1); % short-time autocorrelation;
STPow = zeros(4, NF); % short-time
for i=1:4
    STFT(i, :) = fft(frame(i,:).*hamming(FL)', NF);
    STauto(i, :) = xcorr(frame(i,:));
    STPow(i, :) = abs(fft(frame(i,:).*hamming(FL)', NF)).^2; %Using hamming window to reduce the error.
%     STPow(i, :) = abs(fft(frame(i,:), NF)).^2;
end

AvergSTPow = sum(STPow, 1) / 4;

Newspch = zeros(1, length(y));
Fstart = 1;
Fend = Fstart + FL - 1;
overlap = FL / 2;
Fcnt = 0;
while(Fend <= length(y))
  
    data = y(Fstart:Fend).*hamming(FL);
    FFTdata = fft(data, NF);
  
    Re = max(abs((FFTdata)').^2 - AvergSTPow, 0);  %Get the Power of the spectral differences;
    phase = angle((FFTdata)'); % Keep the phase of the original speech;
    Newdata = sqrt(Re) .* (cos(phase) + i * sin(phase));
%    Newdata = sqrt(Re) .* exp(i * phase); % This way doesn't work. I don't know why.
    tmp = real(ifft(Newdata));
    Newspch(1, Fstart: Fend) = Newspch(1, Fstart: Fend) + tmp(1, 1:160);
%  don't use specified length when using ifft.
%     Newspch(1, Fstart: Fend) = Newspch(1, Fstart: Fend) + real(ifft(Newdata, FL)); % Reconstruct the speech with 50% overlap.
    Fstart = Fstart + overlap;
    Fend = Fstart + FL - 1;
    Fcnt = Fcnt + 1;
end

% soundsc(y, fs);
% soundsc(Newspch, fs);
figure; hold on;
plot(y, 'r'); plot(Newspch, 'g');