MATLAB-轴心轨迹图、幅值谱、相位谱、功率谱(转)

MATLAB编程，题目如下x=6*sin(60*pi*t)+8*sin(120*pi*t+pi/6);y=8*sin(60*pi*t+pi/2)+3*sin(120*pi*t+pi/4);程序一：两个波形图和轴心图clear;closeall;sf=800;%采样频率;N=512;%采样点数t=0:1/sf:1;x=6*sin(60*pi*t)+8*sin(120*pi*t+pi/6);y=8*sin(60*pi*t+pi/2)+3*sin(120*pi*t+pi/4);subplot(2,1,1);%生成m*n个图，当前激活第一个图plot(t,x);xlabel('时间t单位s');%x轴注解ylabel('幅值x单位cm');%y轴注解title('A的波形图');%图形标题gridon;%显示格线subplot(2,1,2);plot(t,y,'g');xlabel('时间t单位s');%x轴注解ylabel('幅值y单位cm');%y轴注解title('B的波形图');%图形标题gridon;%显示格线%轴心轨迹图figure(2)plot(x,y,'b')gridon;程序二：计算X、Y的波形系数和波峰系数%xy的波形系数和波峰系数clc;closeall;sf=800;t=0:1/sf:1;x=6*sin(60*pi*t)+8*sin(120*pi*t+pi/6);xav=mean(abs(x));%平均绝对值xp=max(x);%幅值也就是最大值z=x.^2;d=mean(z);xrms=sqrtm(d);%有效值Ft=xrms/xav;Fc=xp/xrms;disp(['x的波形系数是']);disp(num2str(Ft));disp(['x的波峰系数是']);disp(num2str(Fc));%y的波形系数和波峰系数y=8*sin(60*pi*t+pi/2)+3*sin(120*pi*t+pi/4);x=y;xav=mean(abs(x));%平均绝对值xp=max(x);%幅值也就是最大值z=x.^2;d=mean(z);xrms=sqrtm(d);%有效值Ft=xrms/xav;Fc=xp/xrms;disp(['y的波形系数是']);disp(num2str(Ft));disp(['y的波峰系数是']);disp(num2str(Fc));程序三：求解幅频谱、相频谱、功率谱%幅频谱、相频谱、功率谱clc;fs=800;%设定采样频率N=512;%采样点数n=0:N-1;t=n/fs;y=8*sin(60*pi*t+pi/2)+3*sin(120*pi*t+pi/4);%进行FFT变换并做幅频谱yy=fft(y,N);%进行fft变换py1=yy/512;py2=py1*2;%单边谱乘以2mag=abs(py2);%求幅值f=(0:256)*800/512;subplot(2,1,1);plot(f,mag(1:257));%做幅频谱图xlabel('频率(Hz)');ylabel('幅值');title('B的幅频谱');grid;%功率谱power=mag(1:257).^2;subplot(2,1,2);plot(f,power);gridon;title('B的功率谱');xlabel('频率(Hz)');ylabel('功率');%相频图pyy=angle(yy);Pyy=180/pi*pyy;%弧度变角度figure(2);plot(f,Pyy(1:257));gridon;title('B的相频谱');xlabel('频率(Hz)');ylabel('相位（度）');程序四：%相关系数和相干性clc;closeall;dt=1/800;sf=800;t=0:1/800:1;x=6*sin(60*pi*t)+8*sin(120*pi*t+pi/6);y=8*sin(60*pi*t+pi/2)+3*sin(120*pi*t+pi/4);%相关系数N=128;a=xcorr(x,y,128);a=a(N+1:2*N+1);plot(t(1:N),a(1:N));gridon;%求A、B的相干函数figure(2);N=512;[cxy,w]=cohere(x,y,N/4,sf,hamming(N/4),0,'mean');plot(w,cxy);grid;title('相干分析');xlabel('频率(Hz)');以下是matlabhelp的例子t=0:0.001:0.6;x=sin(2*pi*50*t)+sin(2*pi*120*t);y=x+2*randn(size(t));plot(1000*t(1:50),y(1:50))title('SignalCorruptedwithZero-MeanRandomNoise')xlabel('time(milliseconds)'Y=fft(y,512);Pyy=Y.*conj(Y)/512;f=1000*(0:256)/512;plot(f,Pyy(1:257))title('Frequencycontentofy')xlabel('frequency(Hz)')