singal_system_analysis.m

随机信号分析MATLAB的重要的仿真实例

function varargout = singal_system_analysis(varargin)

%===============================================================
% Begin initialization code - DO NOT EDIT
gui_Singleton = 1;
gui_State = struct('gui_Name',       mfilename, ...
                   'gui_Singleton',  gui_Singleton, ...
                   'gui_OpeningFcn', @singal_system_analysis_OpeningFcn, ...
                   'gui_OutputFcn',  @singal_system_analysis_OutputFcn, ...
                   'gui_LayoutFcn',  [] , ...
                   'gui_Callback',   []);
if nargin & isstr(varargin{1})
    gui_State.gui_Callback = str2func(varargin{1});
end

if nargout
    [varargout{1:nargout}] = gui_mainfcn(gui_State, varargin{:});
else
    gui_mainfcn(gui_State, varargin{:});
end
% End initialization code - DO NOT EDIT
%=================================================================

%=================================================================
% Executes just before singal_system_analysis is made visible.
function singal_system_analysis_OpeningFcn(hObject, eventdata, handles, varargin)
% Choose default command line output for singal_system_analysis
handles.output = hObject;
% Update handles structure
guidata(hObject, handles);
set(handles.insingal_listbox,'string','均匀分布白噪声|正态分布白噪声|鸟叫声|火车鸣笛声|Chi 2分布随机噪声|窄带正态过程');
set(handles.insingal_edit,'string','幅度服从均匀分布的白噪声');
set(handles.outsingal_edit,'string','输入信号经过系统后的输出信号');
set(handles.insignal_popupmenu,'string','时域波形|频域波形|自相关估计|功率谱估计|概率分布图');
set(handles.outsingal_popupmenu,'string','时域波形|频域波形|自相关估计|功率谱估计|概率分布图');
set(handles.system_listbox,'string','全通系统|低通系统|高通系统|带通系统|平方律检波|全波线性检波|半波线性检波');
set(handles.system_edit,'string','即输出等于输入');

% --- Outputs from this function are returned to the command line.
function varargout = singal_system_analysis_OutputFcn(hObject, eventdata, handles)
% Get default command line output from handles structure
varargout{1} = handles.output;
%=================================================================

%==input===========================================================
% --- Executes during object creation, after setting all properties.
function insingal_listbox_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end

% --- Executes on selection change in insingal_listbox.
function insingal_listbox_Callback(hObject, eventdata, handles)
inlistindex=get(handles.insingal_listbox,'Value');
inliststr={'幅度服从均匀分布的白噪声','幅度服从正态（高斯）分布的白噪声',...
              '鸟叫声,能量主要集中在2600~4700Hz,可看做宽带随机过程',...
              '火车鸣笛声,可看做宽带随机过程',...
              'Chi 2分布的宽带随机过程'...
              '窄带正态随机过程'};
set(handles.insingal_edit,'string',inliststr(inlistindex));

% --- Executes during object creation, after setting all properties.
function insignal_popupmenu_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end
%=================================================================
 
%==system=========================================================
% --- Executes during object creation, after setting all properties.
function system_listbox_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end
% --- Executes on selection change in system_listbox.
function system_listbox_Callback(hObject, eventdata, handles)
systemlistindex=get(handles.system_listbox,'Value');
systemliststr={'即输出等于输入',...
                       '输入信号通过具有低通特性的滤波器(低通截至频率2000Hz)',...
                       '输入信号通过具有高通特性的滤波器(高通截至频率2000Hz)',...
                       '输入信号通过具有带通特性的滤波器(通带3000~3300Hz)',...
                       '输入信号通过平方律检波器件',...  
                       '输入信号通过全波线性检波器件',...
                       '输入信号通过半波线性检波器件'};
set(handles.system_edit,'string',systemliststr(systemlistindex));

Fs=10000;   
%输入信号经过选择的系统生成输出信号
sysindex=get(handles.system_listbox,'Value');
if sysindex==1
    M=[1 1];
    f=[0 Fs/2]/(Fs/2);
    [B,A]=yulewalk(11,f,M);
elseif sysindex==2
    M=[1 1 0 0];
    f=[0 2000 2100 Fs/2]/(Fs/2);
    [B,A]=yulewalk(30,f,M);
elseif sysindex==3
    M=[0 0 1 1];
    f=[0 2000 2100 Fs/2]/(Fs/2);
    [B,A]=yulewalk(25,f,M);
elseif sysindex==4
    M=[0 0 1 1 0 0 ];
    f=[0 3000 3050 3350 3400 Fs/2]/(Fs/2);
    [B,A]=yulewalk(50,f,M);
end

%画出系统的特性曲线
axes(handles.system_axes)
if sysindex==5  %y=x^2
    x=-1:0.1:1;
    y=x.^2;
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}平方律检波');
elseif sysindex==6
    x=-1:0.1:1;
    y=abs(x);
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}全波线性检波');
elseif sysindex==7
    x=-1:0.1:1;
    y=(x+abs(x))/2;
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}半波线性检波');
else
    [H,W]=freqz(B,A,128);
    plot(f*Fs/2,M,':r',W*Fs/(2*pi),abs(H),'-');grid on;
    legend('\fontsize{9}理想曲线','\fontsize{9}实际曲线');xlabel('f/Hz');ylabel('\fontsize{9}幅频特性');
    title('\fontsize{9}系统');
end
%================================================================


%==output=========================================================
% --- Executes during object creation, after setting all properties.
function outsingal_edit_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end

% --- Executes during object creation, after setting all properties.
function outsingal_listbox_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end

%===============================================================
% --- Executes during object creation, after setting all properties.
function outsingal_popupmenu_CreateFcn(hObject, eventdata, handles)
if ispc
    set(hObject,'BackgroundColor','white');
else
    set(hObject,'BackgroundColor',get(0,'defaultUicontrolBackgroundColor'));
end

% --- Executes on button press in outsingal_refreshbutton.
function outsingal_refreshbutton_Callback(hObject, eventdata, handles)
    inlistindex=get(handles.insingal_listbox,'Value');
    %产生适当的信号向量
    Fs=10000;Ts=1/Fs;N=8096;t=[0:N-1]*Ts;
    %均值和方差
    Mu=0;Sigma=1;
   
    %改变种子，可以使结果更“精确”
    rand('state',1)
    Nwhite=rand(1,N);
    Nwhite=Nwhite-0.5+Mu;
    temp=sqrt(12*(Mu^2+Sigma));
    Nwhite=temp*Nwhite;
    %改变种子，可以使结果更“精确”
    rand('state',2)
    Ngauss=randn(1,N);
    temp=sqrt(Mu^2+Sigma);
    Ngauss=temp*Ngauss;
    
%===============================================================
%根据选择产生输入信号
if inlistindex==1
    X=Nwhite;
elseif inlistindex==2  
    X=Ngauss;
elseif inlistindex==3
    load('chirp','Fs','y')
    X=y(1:N);Fs=10000;
elseif inlistindex==4
    load('train','Fs','y')
    X=y(1:N);Fs=10000;
elseif inlistindex==5
    X=chi2rnd(4,1,N);
elseif inlistindex==6
    M=[0 0 1 1 0 0 ];
    f=[0 3000 3050 3350 3400 Fs/2]/(Fs/2);
    [B,A]=yulewalk(40,f,M);
    Nt=filter(B,A,Ngauss);
    X=Nt;
end

%===============================================================
%画出输入信号的图形
inpopindex=get(handles.insignal_popupmenu,'Value');
axes(handles.insingal_axes)
if inpopindex==1
    plot(t,X);grid on;
    axis([0 200*Ts min(X) max(X)]);
    xlabel('t/s');title('\fontsize{9}时域波形');
elseif inpopindex==2
    Xf=fft(X,512);
    semilogy((0:511)*Fs/512,abs(Xf));grid on;
    axis([0 Fs/2 log(min(Xf)) inf]);
    xlabel('f/Hz');title('\fontsize{9}频域波形');
elseif inpopindex==3
    R=xcorr(X,'biased');
    plot((-N:length(R)-N-1)*Ts,R);grid on;
    axis([-N*Ts (length(R)-N-1)*Ts min(R) max(R)]);
    xlabel('t/s');title('\fontsize{9}自相关函数');
elseif inpopindex==4
    [Pxx,f] = pwelch(X,[],[],1024,Fs);
    semilogy(f,Pxx);grid on;
    axis([min(f) max(f) log(min(Pxx)) inf]);
    xlabel('f/Hz');title('\fontsize{9}功率谱密度');
elseif inpopindex==5
    histfit(X,50);
    title('\fontsize{9}概率分布直方图');grid on;
    h = findobj(gca,'Type','line');
    set(h,'Color','r');
end

%输入信号经过选择的系统生成输出信号
sysindex=get(handles.system_listbox,'Value');
if sysindex==1
    M=[1 1];
    f=[0 Fs/2]/(Fs/2);
    [B,A]=yulewalk(11,f,M);
    Y=filter(B,A,X);
elseif sysindex==2
    M=[1 1 0 0];
    f=[0 2000 2100 Fs/2]/(Fs/2);
    [B,A]=yulewalk(30,f,M);
    Y=filter(B,A,X);
elseif sysindex==3
    M=[0 0 1 1];
    f=[0 2000 2100 Fs/2]/(Fs/2);
    [B,A]=yulewalk(25,f,M);
    Y=filter(B,A,X);
elseif sysindex==4
    M=[0 0 1 1 0 0 ];
    f=[0 3000 3050 3350 3400 Fs/2]/(Fs/2);
    [B,A]=yulewalk(50,f,M);
    Y=filter(B,A,X);
elseif sysindex==5
    Y=X.^2;
elseif sysindex==6
    Y=abs(X);
elseif sysindex==7
    Y=(X+abs(X))/2;
end

%画出系统的特性曲线
axes(handles.system_axes)
if sysindex==5  %y=x^2
    x=-1:0.1:1;
    y=x.^2;
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}平方律检波');
elseif sysindex==6
    x=-1:0.1:1;
    y=abs(x);
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}全波线性检波');
elseif sysindex==7
    x=-1:0.1:1;
    y=(x+abs(x))/2;
    plot(x,y);grid on;axis([-2 2 0 inf]);xlabel('X');ylabel('Y');title('\fontsize{9}半波线性检波');
else
    [H,W]=freqz(B,A,128);
    plot(f*Fs/2,M,':r',W*Fs/(2*pi),abs(H),'-');grid on;
    legend('\fontsize{9}理想曲线','\fontsize{9}实际曲线');xlabel('f/Hz');ylabel('\fontsize{9}幅频特性');
    title('\fontsize{9}系统');
end


%===============================================================
%画出输出信号图形
outpopindex=get(handles.outsingal_popupmenu,'Value');
axes(handles.outsingal_axes)
if outpopindex==1
    plot(t,Y);
    axis([0 200*Ts min(Y) max(Y)]);
    xlabel('t/s');title('\fontsize{9}时域波形');grid on;
elseif outpopindex==2
    Yf=fft(Y,512);
    semilogy((0:511)*Fs/512,abs(Yf));grid on;
    axis([0 Fs/2 log(min(Yf)) inf]);
    xlabel('f/Hz'); title('\fontsize{9}频域波形');
elseif outpopindex==3
    R=xcorr(Y,'biased');
    plot((-N:length(R)-N-1)*Ts,R);
    axis([-N*Ts (length(R)-N-1)*Ts min(R) max(R)]);
    xlabel('t/s');title('\fontsize{9}自相关函数');grid on;
elseif outpopindex==4
    [Pyy,f] = pwelch(Y,[],[],1024,Fs);
    semilogy(f,Pyy);
    axis([min(f) max(f) log(min(Pyy)) inf]);
    xlabel('f/Hz');title('\fontsize{9}功率谱密度');grid on;
elseif outpopindex==5
    histfit(Y,55);
    title('\fontsize{9}概率分布直方图');grid on;
    h = findobj(gca,'Type','patch');
    set(h,'FaceColor','r','EdgeColor','w');
    h = findobj(gca,'Type','line');
    set(h,'Color','b');
end
%=================================================================