www.gusucode.com > 声音的处理有：LPC，FFT，共振峰，频谱源码程序 > siganlandsystemusingMatlab/SSUM/filters/iir/iirexpofn.m
function iirexpofn(action,datastruct)
if nargin < 1
action='init';
end
% Fix freq labels
% Apply filter shouldn't know about Undo
% Add decibel constraints for elliptical and chebyshev
name = mfilename;
figname = [name(1:end-2) '_fig'];
f=findobj('Tag',figname);
handles = get(f,'UserData');
switch action
case 'help'
display_help(figname);
case 'init'
setdefaults;
reset(handles.specplot);
reset(handles.timeplot);
set(handles.rp,'Visible','off');
set(handles.rs,'Visible','off');
handles.exponum = 0;
case 'loadsound'
handles.audiodata = load_audiodata;
handles.original = handles.audiodata;
if ~isfield(handles.audiodata, 'filenamepath')
return;
end
if (size(handles.audiodata.data,2) > 1)
handles.audiodata.data = to_mono(handles.audiodata.data);
end
set(handles.freqzoom,'Value',1.0);
contents = get(handles.fftsize,'String');
fftsize = str2double(contents{get(handles.fftsize,'Value')});
contents = get(handles.Window,'String');
shape = contents{get(handles.Window,'Value')};
[handles.spectrum, handles.bin, handles.st] = ...
spectrogram(handles.audiodata, fftsize, shape);
handles = makeTimePlot(handles);
handles = makeSpecPlot(handles);
linkedzoom([handles.timeplot, handles.specplot],'onx');
place_header(f,handles);
case 'readinput'
handles.audiodata = datastruct;
clear datastruct;
if (size(handles.audiodata.data,2) > 1)
handles.audiodata.data = to_mono(handles.audiodata.data);
end
set(handles.freqzoom,'Value',1.0);
handles.original = handles.audiodata;
contents = get(handles.fftsize,'String');
fftsize = str2double(contents{get(handles.fftsize,'Value')});
contents = get(handles.Window,'String');
shape = contents{get(handles.Window,'Value')};
[handles.spectrum, handles.bin, handles.st] = ...
spectrogram(handles.audiodata, fftsize, shape);
handles = makeTimePlot(handles);
handles = makeSpecPlot(handles);
linkedzoom([handles.timeplot, handles.specplot],'onx');
place_header(f,handles);
case {'playsound', 'playoriginal'}
if isfield(handles,'audiodata')
times = get(handles.specplot,'XLim');
switch action
case 'playsound'
audiodata = handles.audiodata;
button = handles.play;
case 'playoriginal'
audiodata = handles.original;
button = handles.playoriginal;
end
samples = floor(times*audiodata.Fs);
if (samples(1) <= 0)
samples(1) = 1;
end
if (samples(2) > length(audiodata.data))
samples(2) = length(audiodata.data);
end
audiodata.data = audiodata.data(samples(1):samples(2));
if (max(abs(audiodata.data)) > 1.0)
audiodata.data = normalize(audiodata.data);
end
play_audiodata(audiodata, button);
end
case {'db','inverse','interpolate','colormap','fftdisplay'}
if isfield(handles,'audiodata')
updateSpecPlot(handles);
set(handles.specplot,'XLim',get(handles.timeplot,'XLim'));
iirexpofn 'freqzoom';
end
case {'fftsize','window'}
if isfield(handles,'audiodata')
contents = get(handles.fftsize,'String');
fftsize = str2double(contents{get(handles.fftsize,'Value')});
contents = get(handles.Window,'String');
shape = contents{get(handles.Window,'Value')};
[handles.spectrum, handles.bin, handles.st] = ...
spectrogram(handles.audiodata, fftsize, shape);
updateSpecPlot(handles);
set(handles.specplot,'XLim',get(handles.timeplot,'XLim'));
iirexpofn 'freqzoom';
end
case 'normalize'
if isfield(handles,'audiodata')
handles.audiodata.data = normalize(handles.audiodata.data);
xlim = get(handles.timeplot,'XLim');
handles = makeTimePlot(handles);
contents = get(handles.fftsize,'String');
fftsize = str2double(contents{get(handles.fftsize,'Value')});
contents = get(handles.Window,'String');
shape = contents{get(handles.Window,'Value')};
[handles.spectrum, handles.bin, handles.st] = ...
spectrogram(handles.audiodata, fftsize, shape);
updateSpecPlot(handles);
set(handles.timeplot,'XLim',xlim);
set(handles.specplot,'XLim',xlim);
iirexpofn 'freqzoom';
end
case 'explore_data'
explore_data(handles);
case 'plotfreqz'
plotfreqz(handles);
case 'apply_filter'
if isfield(handles, 'audiodata')
handles = apply_filter(handles);
contents = get(handles.fftsize,'String');
fftsize = str2double(contents{get(handles.fftsize,'Value')});
contents = get(handles.Window,'String');
shape = contents{get(handles.Window,'Value')};
[handles.spectrum, handles.bin, handles.st] = ...
spectrogram(handles.audiodata, fftsize, shape);
updateSpecPlot(handles);
Xlim = get(handles.timeplot,'XLim');
handles = makeTimePlot(handles);
set(handles.specplot,'XLim',Xlim);
set(handles.timeplot,'XLim',Xlim);
iirexpofn 'freqzoom';
end
case 'freqzoom'
if isfield(handles,'audiodata')
val = get(handles.freqzoom,'Value');
if val == 0,
val = 0.001;
end
Fs = handles.audiodata.Fs;
set(handles.specplot,'YLim',[0 val*Fs/2])
end
case 'zoomreset'
if isfield(handles,'audiodata')
handles = makeTimePlot(handles);
set(handles.specplot,'XLim',get(handles.timeplot,'XLim'));
linkedzoom([handles.timeplot, handles.specplot],'onx');
end
case 'save'
if isfield(handles, 'audiodata')
save_audiodata(handles.audiodata);
end
case 'undo'
if isfield(handles,'audiodata')
handles = undo(handles);
updateSpecPlot(handles);
Xlim = get(handles.timeplot,'XLim');
handles = makeTimePlot(handles);
set(handles.specplot,'XLim',Xlim);
set(handles.timeplot,'XLim',Xlim);
iirexpofn 'freqzoom';
end
case {'pzexpo','pzfilterexpo'}
[a_k,b_k] = designfilter(handles);
[data.zeropos,data.polepos] = tf2zpk(b_k,a_k);
switch action
case 'pzexpo'
pzexpogui(data);
case 'pzfilterexpo'
pzfilterexpogui(data);
end
case 'convexpo'
[a_k,b_k] = designfilter(handles);
data.impulse = filter(b_k, a_k, [1, zeros(1,99)]);
convexpogui(data);
case {'fourier'}
if isfield(handles,'audiodata')
audiodata = handles.audiodata;
fourierexpogui(audiodata);
end
case 'print'
print_figure(f);
case 'close'
close_figure(f,figname(1:end-4));
return;
end
set(f,'UserData',handles);
% --------------------------------------------------------------------
function handles = makeTimePlot(handles);
handles.t = [0:1/handles.audiodata.Fs:(length(handles.audiodata.data)-1)/...
handles.audiodata.Fs];
axes(handles.timeplot)
plot(handles.t,handles.audiodata.data)
maxtime = length(handles.t)/handles.audiodata.Fs;
set(handles.timeplot,'XLim',[0 maxtime]);
set(handles.timeplot,'YLim',[-1.0 1.0]);
grid;
xlabel('time (s)');
function handles = makeSpecPlot(handles);
axes(handles.specplot);
handles.pos = get(gca,'Position'); % Save axes position
if (get(handles.interpolate,'Value'))
if (get(handles.dBcheckbox,'Value'))
pcolor(handles.st,handles.bin,20*log10(abs(handles.spectrum)));
else
pcolor(handles.st,handles.bin,abs(handles.spectrum));
end
shading interp;
else
if (get(handles.dBcheckbox,'Value'))
imagesc(handles.t,handles.bin,20*log10(abs(handles.spectrum)));
else
imagesc(handles.t,handles.bin,abs(handles.spectrum));
end
end
axis xy;
ylabel('Frequency (Hz)');
colormap('Jet');
set(handles.specplot,'XTickLabel',['']);
h = colorbar('vert');
if (get(handles.dBcheckbox,'Value'))
set(get(h, 'YLabel'), 'String', 'dB');
else
set(get(h, 'YLabel'), 'String', 'Amplitude');
end
set(handles.specplot,'Units','Characters');
sonopos = get(handles.specplot,'Position')
set(h,'units','Characters', ...
'Position',[sonopos(1)+sonopos(3)+18 sonopos(2) 7 sonopos(4)])
switch computer
case {'GLNX86','MAC'}
%set(gca,'Position',handles.pos + [0 0.07 0 0])
set(gca,'Position',handles.pos)
case 'PCWIN'
set(gca,'Position',handles.pos)
end
handles.xlim_specplot_orig = get(handles.specplot,'XLim');
handles.ylim_specplot_orig = get(handles.specplot,'YLim');
handles.xlim_timeplot_orig = get(handles.timeplot,'XLim');
handles.ylim_timeplot_orig = get(handles.timeplot,'YLim');
handles.xlim_specplot = handles.xlim_specplot_orig;
handles.ylim_specplot = handles.ylim_specplot_orig;
handles.xlim_timeplot = handles.xlim_timeplot_orig;
handles.ylim_timeplot = handles.ylim_timeplot_orig;
handles.clim = get(handles.specplot,'CLim');
% --------------------------------------------------------------------
function updateTimePlot(handles)
axes(handles.timeplot)
plot(handles.t,handles.audiodata)
set(handles.timeplot,'XLim',handles.xlim_timeplot);
set(handles.timeplot,'YLim',handles.ylim_timeplot);
xlabel('time (s)');
grid;
function updateSpecPlot(handles)
axes(handles.specplot);
plotype = get(handles.fftdisplay,'String');
switch lower(plotype{get(handles.fftdisplay,'Value')})
case 'sonogram'
% Check Interpolate and dB
if (get(handles.interpolate,'Value'))
if (get(handles.dBcheckbox,'Value'))
pcolor(handles.st,handles.bin,20*log10(abs(handles.spectrum)));
else
pcolor(handles.st,handles.bin,abs(handles.spectrum));
end
shading interp;
else
if (get(handles.dBcheckbox,'Value'))
imagesc(handles.t,handles.bin,20*log10(abs(handles.spectrum)));
else
imagesc(handles.t,handles.bin,abs(handles.spectrum));
end
end
set(handles.specplot,'XTickLabel',['']);
axis xy;
ylabel('Frequency (Hz)');
% Get Colormap
contents = get(handles.colormap,'String');
cmap = colormap(lower(contents{get(handles.colormap,'Value')}));
% Handle Inverse
if (get(handles.inverse,'Value'))
colormap(flipud(cmap));
else
colormap(cmap);
end
set(handles.specplot,'XLim',handles.xlim_specplot);
set(handles.specplot,'YLim',handles.ylim_specplot);
set(handles.specplot,'CLim',handles.clim);
h = colorbar('vert');
if (get(handles.dBcheckbox,'Value'))
set(get(h, 'YLabel'), 'String', 'dB');
else
set(get(h, 'YLabel'), 'String', 'Amplitude');
end
set(handles.specplot,'Units','Characters');
sonopos = get(handles.specplot,'Position')
set(h,'units','Characters', ...
'Position',[sonopos(1)+sonopos(3)+18 sonopos(2) 7 sonopos(4)])
switch computer
case {'GLNX86','MAC'}
%set(gca,'Position',handles.pos + [0 0.07 0 0])
set(gca,'Position',handles.pos)
case 'PCWIN'
set(gca,'Position',handles.pos)
end
case 'waterfall'
% Waterfall does not work on Linux, Mac, or Windows!
%waterfall(handles.st,handles.bin,abs(handles.spectrum));
%if (get(handles.dBcheckbox,'Value'))
%me = surf(handles.st,handles.bin,20*log10(abs(handles.spectrum)));
%else
% perhaps downsample
%me = surf(handles.st,handles.bin,abs(handles.spectrum));
% This would be nice but unecessary
%explore3d(gca);
%end
%zoom off;
%xlabel('Time')
%ylabel('Frequency')
%if (get(handles.interpolate,'Value'))
% shading interp;
%end
end
% --------------------------------------------------------------------
function plotfreqz(handles)
% Find a way to keep original xlims
[a_k,b_k] = designfilter(handles);
%[H, W] = freqz(b_k, a_k, 512, 'whole');
[H, W] = freqz(b_k, a_k);
if isfield(handles,'audiodata')
Fs = handles.audiodata.Fs;
else
Fs = 44100;
end
f=findobj('Tag','freqzplot');
if isempty(f)
f = figure('units','normal','paperunits','normal',...
'position',[.1,.1,.6,.6], 'Tag', 'freqzplot','MenuBar','none');
uicontrol('String','dB','Style','checkbox','Value', 0, ...
'Position',[5 400 40 20],'Tag','dB',...
'Callback','iirexpofn plotfreqz');
uicontrol('String','unwrap','Style','checkbox','Value', 1, ...
'Position',[5 200 80 20],'Tag','unwrap',...
'Callback','iirexpofn plotfreqz');
end
dB = findobj(f,'Tag','dB');
unwrapbox = findobj(f,'Tag','unwrap');
% Plot frequency response of H(w)
freqplot = subplot(3,1,1);
W = W.*Fs/(2*pi);
if get(dB,'Value')
plot(W,20*log10(abs(H)) - max(20*log10(abs(H))));
axis([0.0 Fs/2 -200 10.1]);
ylabel('Magnitude (dB)');
else
plot(W,abs(H)./max(abs(H)));
axis([0.0 Fs/2 -0.1 1.1]);
ylabel('Magnitude');
end
set(freqplot,'XTickLabel','');
t_h = title('Frequency, Phase and Impulse Response');
set(t_h, 'fontsize',15, 'fontweight','bold');
grid;
% Plot phase response of H(w)
phasplot = subplot(3,1,2);
if get(unwrapbox,'Value')
plot(W,unwrap(angle(H)));
axis([0.0 Fs/2 1.1*min(unwrap(angle(H))) 1.1*max(unwrap(angle(H)))]);
else
plot(W,angle(H));
axis([0.0 Fs/2 -pi +pi]);
end
xlabel('Frequency (Hz)');
ylabel('Phase (radians)');
grid;
linkedzoom([freqplot,phasplot],'onx');
% Plot impulse response
implot = subplot(3,1,3);
impulse = filter(b_k, a_k, [1, zeros(1,99)]);
plot(impulse);
axis([0 length(impulse) -max(abs(impulse))*1.1 max(abs(impulse))*1.1]);
xlabel('Sample');
ylabel('Amplitude');
grid;
function [a,b] = designfilter(handles)
order = str2double(get(handles.filtorder,'String'));
contents = get(handles.filtermenu,'String');
typecontents = get(handles.type,'String');
cut1 = str2double(get(handles.cut1,'String'));
cut2 = str2double(get(handles.cut2,'String'));
Rp = str2num(get(handles.rp,'String'));
Rs = str2num(get(handles.rs,'String'));
if isfield(handles,'audiodata')
Fs = handles.audiodata.Fs;
else
Fs = 44100;
end
cut1 = cut1/(Fs/2);
cut2 = cut2/(Fs/2);
switch lower(typecontents{get(handles.type,'Value')})
case 'butterworth'
switch (lower(contents{get(handles.filtermenu,'Value')}))
case 'allpass'
b = 1;
a = 1;
case 'lowpass'
[b,a] = butter(order, cut1);
case 'highpass'
[b,a] = butter(order, cut1,'high');
case 'bandpass'
[b,a] = butter(order, [cut1 cut2]);
case 'notch'
[b,a] = butter(order, [cut1 cut2], 'stop');
end
case 'elliptical'
switch (lower(contents{get(handles.filtermenu,'Value')}))
case 'allpass'
b = 1;
a = 1;
case 'lowpass'
[b,a] = ellip(order,Rp,Rs, cut1);
case 'highpass'
[b,a] = ellip(order,Rp,Rs, cut1,'high');
case 'bandpass'
[b,a] = ellip(order,Rp,Rs, [cut1 cut2]);
case 'notch'
[b,a] = ellip(order,Rp,Rs, [cut1 cut2], 'stop');
end
case 'chebyshev'
switch (lower(contents{get(handles.filtermenu,'Value')}))
case 'allpass'
b = 1;
a = 1;
case 'lowpass'
[b,a] = cheby1(order,Rp, cut1);
case 'highpass'
[b,a] = cheby1(order,Rp, cut1,'high');
case 'bandpass'
[b,a] = cheby1(order,Rp, [cut1 cut2]);
case 'notch'
[b,a] = cheby1(order,Rp,[cut1 cut2], 'stop');
end
end
function handles = apply_filter(handles)
[a_k,b_k] = designfilter(handles);
handles.audiodata2 = handles.audiodata;
handles.audiodata.data = filter(b_k, a_k, handles.audiodata.data);
% Store old data for undo
handles.spectrum2 = handles.spectrum;
handles.bin2 = handles.bin;
handles.st2 = handles.st;
set(handles.undo,'String','Undo');
function handles = undo(handles)
hObject = handles.undo;
if (strcmp(get(hObject,'String'),'Undo'))
% Store old data for redo
handles.spectrum3 = handles.spectrum;
handles.bin3 = handles.bin;
handles.st3 = handles.st;
handles.audiodata3 = handles.audiodata;
% Reload old data
handles.spectrum = handles.spectrum2;
handles.bin = handles.bin2;
handles.st = handles.st2;
handles.audiodata = handles.audiodata2;
set(hObject,'String','Redo');
else
handles.spectrum = handles.spectrum3;
handles.bin = handles.bin3;
handles.st = handles.st3;
handles.audiodata = handles.audiodata3;
set(hObject,'String','Undo');
end