clear all
clc

[filename,pathname] = uigetfile( ...
{'*.wav', 'Zvukové soubory (*.wav)';
'*.wav', 'Wave (*.wav)';...
'*.*', 'Vsechny soubory (*.*)'},...
'Vyber soubor');
[pathstr,name,ext,versn] = fileparts(filename);   
[x,fsx]=wavread(fullfile(pathname,filename));
disp(' ');
disp( strcat('Původní vzorkovací kmitočet: ' ,num2str( fsx), ' Hz'));
disp('---------------------------------------------');

fsy = input('Zadejte nový vzorkovací kmitočet [Hz]: ');

Nx = length(x);
nx = 0:Nx-1;
tx = nx(:)/fsx;
Xf = fft( x)/Nx;
% oboustranne spektrum
Xf = [Xf( ceil(Nx/2)+1:end); Xf(1:ceil( Nx/2))];
fx = ((0:Nx-1)-floor( Nx/2))/Nx * fsx;



K = gcd(fsx,fsy);
L = fsy/K;
M = fsx/K;

%Nadvzorkovani s cinitelem L
Ny = Nx*L;
ny = 0:Ny-1;
ty = ny/(fsx*L);
% provedeni nadvzorkovani
y = zeros(Ny, 1);
y (nx*L+1) = x;


%Antialiasingovy filtr
if 1/L < 1/M
    [N,wn] = cheb1ord( 1/L, 1.2/L, 2, 40 );
else
    [N,wn] = cheb1ord( 0.8/M, 1/M, 2, 40 );
end
[b, a] = cheby1( N, 2, wn);
y = L*(filter(b, a, y));


%Podvzorkovani s cinitelem M
Nz = floor(Ny/M);   
nz = 0:Nz-1;
tz = nz/(fsy/M);
% +1 z duvodu indexovani od 1
z = y(nz*M+1);
z=z./max(z);
% vypocet spektra podvzorkovaneho signalu
Zf = fft( z)/Nz;
Zf = [Zf( ceil( Nz/2)+1:end); Zf( 1:ceil( Nz/2))];
fz = (( 0:Nz-1)-floor( Nz/2))/Nz * fsy;

%Spektrum
figure( 1);
% spektrum vstupniho signalu
plot( fx, 20*log10( abs( Xf)));
% spektrum vystupniho signalu
hold on; 
plot( fz, 20*log10( abs( Zf)),'r'); 
hold off;
grid on;
xlabel( '\rightarrow {\it f} [Hz]');
ylabel( '\uparrow |{\it S}({\it f})| [dB]');
legend( {'vstupní signál', ...
  'podvzorkovaný signál'});

[save_file,save_path] = uiputfile('*.wav','Uložit jako',[name,'_resampl']);
savefile = [ save_path,save_file];
wavwrite(z,fsy,savefile); 


