www.gusucode.com > 声音的处理有:LPC,FFT,共振峰,频谱源码程序 > siganlandsystemusingMatlab/SSUM/filters/fir/firexpofn.m
function firexpofn(action,datastruct) if nargin < 1 action='init'; end % Fix freq labels % Apply filter shouldn't know about Undo 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); 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')); firexpofn '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')); firexpofn '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); firexpofn '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); firexpofn '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); firexpofn '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 = b_k; 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); a_k = 1; [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', 1, ... 'Position',[5 430 40 20],'Tag','dB',... 'Callback','firexpofn plotfreqz'); uicontrol('String','unwrap','Style','checkbox','Value', 0, ... 'Position',[5 150 80 20],'Tag','unwrap',... 'Callback','firexpofn 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') h = plot(W,20*log10(abs(H)) - max(20*log10(abs(H)))); axis([0.0 Fs/2 -100 1.1]); ylabel('Magnitude (dB)'); else h = plot(W,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 = b_k; plot(b_k); axis([0 length(b_k) -max(abs(b_k))*1.1 max(abs(b_k))*1.1]); xlabel('Sample'); ylabel('Amplitude'); grid; function [a,b] = designfilter(handles) order = str2double(get(handles.filtorder,'String')); contents = get(handles.filtermenu,'String'); cut1 = str2double(get(handles.cut1,'String')); cut2 = str2double(get(handles.cut2,'String')); if isfield(handles,'audiodata') Fs = handles.audiodata.Fs; else Fs = 44100; end cut1 = cut1/(Fs/2); cut2 = cut2/(Fs/2); switch (lower(contents{get(handles.filtermenu,'Value')})) case 'allpass' b = 1; case 'lowpass' b = fir1(order, cut1, 'low'); case 'highpass' b = fir1(order, cut1, 'high'); case 'bandpass' b = fir1(order, [cut1 cut2], 'bandpass'); case 'notch' b = fir1(order, [cut1 cut2], 'stop'); case 'averager' b = ones(1,order)./order; end a = 1; 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