4种方法合1的matlab数据谱分析作图及自动筛选谱的峰值函数包

研究过程中实现的，分享出来。 是一种多种分析图绘制合一的完整函数并对生成谱进行自动的peaks筛选输出。已经形成自定义函数，放到特定路径下即可直接调用，配有较详细的使用说明 style=1,plot amplitude value spectrumï¼›%style=2,plot power spectrum;%style=3,plot loglog PSD %style=4,pwelch with segements number 此绘图函数可以输入多种调用参数，如'FrequencyBan'，'Findpeaks'，'MinPeakHeight','MinPeakProminence'，'MinPeakDistance','Threshold','WindowNum','Overlap','Nfft', 'ustar'等参数控制

function [pkslocspk_wpk_p] = hua_fft（Xfsstylevarargin）
% fs Sampling fequency of signal
% hua_fft（Xfsstylevarargin）
%style=1plot amplitude value spectrum；
%style=2plot power spectrum;
%style=3plot loglog
%style=4pwelch with segements number equals fsplot loglog
%style=other plot both
% ===========varargin======================
%@@‘FrequencyBan‘[f1f2]
%When style=1 varargin input parameters can be used to control the frequency bands
%f1 is start frequency
%f2 is end frequency
%@@‘Findpeaks‘ is used to control the peaks finding
%@@‘MinPeakHeight‘value
% please refer help ‘findpeaks‘ with the same parameter
%@@‘MinPeakProminence‘value
% please refer help ‘findpeaks‘ with the same parameter
%@@‘MinPeakDistance‘value
% please refer help ‘findpeaks‘ with the same parameter
%@@‘Threshold‘value
% please refer help ‘findpeaks‘ with the same parameter
% value = 1e-4 can exclude the flat peak
%@@‘WindowNum‘value
%windows number equal fs as default if no input
%@@‘Overlap‘value
%input overlap from 33%~50% here take 50% overlap as default
%@@‘Nfft‘ value
%input number of DFT points no less than 256
%@@‘ustar‘ value
%input ustar and normalize the PSD to PSD*fs/ustar^2

%http://blog.sina.com.cn/s/blog_574d08530100qu18.html

nfft= 2^nextpow2（length（X））;%set step length automatically
%nfft=1024;%set step length manually
%y=X-mean（X）;
dtr_X=detrend（X）;
y_ft_0=fft（dtr_Xnfft）;%DFTget amplitude vaule distribution
y_ft = y_ft_0（1:nfft/2）; % throwing away half of X : Nyquist criteria

Amp = 2*abs（y_ft）./nfft;       %Amplitude
mx = （abs（y_ft））.^2;    %one-sided PSD
y_E=（y_ft.*conj（y_ft））./nfft;%compute the energy spectrum
y_p=y_ft.*conj（y_ft）./nfft;      %Power
fn=fs*（0:nfft/2-1）./nfft;  %frequency serie

ArgStr = varargin;

switch style
    case 1
        if nargin==3
            Xvar = fn;
            Yvar = Amp;
            plot（XvarYvar）;%use matlab help method plot FFT
            ylabel（‘Amplitude‘）;xlabel（‘Frequency‘）;%title（‘Spectrum of amplitude value‘）;
            %plot（y_fabs（y_ft（1:nfft/2）））;%another method to plot FFT
        else
            index= find（strcmp（‘FrequencyBan‘ArgStr））;
            if ~isempty（index）
                v_f1=ArgStr{index+1};
                v_fn=ArgStr{index+2};
                ni=round（v_f1 * nfft/fs+1）;
                na=round（v_fn * nfft/fs+1）;
                Xvar = fn（ni:na）;
                Yvar = abs（y_ft_0（ni:na）*2/nfft）;
                loglog（XvarYvar）;
            end
        end
        if nargin > 3
            [arg_strcom4]= SetArgStr（ArgStr）;
            eval_str= [‘[pkslocspk_wpk_p]=findpeaks（YvarXvar‘arg_strcom4‘‘‘Annotate‘‘‘‘extents‘‘‘‘WidthReference‘‘‘‘halfheight‘‘）;‘];
            eval（eval_str