compresie.m

audio compression matlab code

%*******************************************************
% Compresia surselor audio prin DPCM si ADPCM;
% DPCM = Differential Pulse Code Modulation;
% ADPCM = Adaptive Differential Pulse Code Modulation;
%*******************************************************
clc; clear;clf;

% definirea tipului de intrare ..
ver_set = {'vorbire', 'audio'};
ver = 'audio';    

% definirea metodei de compresie..
compression_set ={'fix', 'song', 'jayant'};
compression = 'fix';

switch ver,
    case 'audio',   file = 'sound.wav';
    case 'vorbire', file = 'speech.wav';
end;

% prelucrari primare:
[ys,fs,nbits]=wavread(file);
[m,m1] = size(ys);
if (m1 == 2),           % stereo
    y1=ys(:,1); y2=ys(:,2); ym=(y1+y2)/2;
else 
    ym = ys;
end;
N=length(ym); ts=1/fs; t=0:ts:(N-1)*ts;

% defineste vectorii scurti:
n = 10000;
x = ym(1:n);     % intrarea;

% initializari:
b=[-1];
q0 = 0.01;      % valoarea minima a pasului de cuantizare;
y=zeros(n,1);   % vectorul de iesire, aproximat/refacut local; 
q = q0;	    % valoarea pasului de cuantizare;
b_new = -1;     % valoarea bitului nou;
b_old = -1;     % valoarea bitului vechi;
p=1.5;          % parametru algoritm Jayant;

for k = 2:n,
    if  x(k) > y(k-1),  b_new = 1; else b_new = -1; end;
    b = [b b_new];
    switch compression,   % selectie metoda..
        case 'fix',         % pas fix
            q = q0;
        case 'song',
              if b_new == b_old, q = q + q0; else q = q - q0; end;
          case 'jayant',              
              if b_new == b_old, q = q * p; else q = q / p; end;
    end; % switch           
    y(k) = y(k-1) + q * b_new;       % reconstructie:
    b_old = b_new;             
end;
    
% analiza spectru:
[f,df,spec_in] =  f_fft(fs,x(1:n));
[f,df,spec_out] = f_fft(fs,y(1:n));
    
%vizualizare:
t = t(1:n);
[dt,dy]=stairs(t,y);
[dt,dx]=stairs(t,x);
subplot(211), plot(dt,dx,'b',dt,dy,'r-');
subplot(223), stairs(b(100:150)); axis([0 50 -1.1 1.1]);
nf = 1000;          % numarul de puncte depe axa frecventei;
f = (0:nf-1).*df;
subplot(224), plot(f,spec_in(1:nf), f, spec_out(1:nf),'r');  xlabel('Hz');

% calcul raport compresie:
file_input_size = n * nbits / 8 ; % mono channal;
file_output_size = n / 8 ; % mono channal;
compression_ratio = file_input_size / file_output_size;

% salvarea semnalului refacut dupa compresie:
file = strcat(compression,'_compr_',ver, '.wav');
wavwrite(y,fs,nbits, file);

%salvare fisier in format binar...
% conversie valoare -1 in valoare 0:
for i=1:length(b), if b(i) < 0, b(i) = 0; end; end;
fid = fopen('compressed.bin','wb')
fwrite(fid, b,'bit1');
fclose(fid);