Matlab任务我定稿版

clear;
close all;
fid=fopen('internal-sinewave.dat','r')
a=fread(fid,'int32')
a=a'
t=0:0.002:4095*0.002
 
subplot(3,1,1);
plot(t,a)
xlabel('时间/s');
ylabel('采样值');
title('原始信号(时域)');
grid on;
%a=a-mean(a);
Fs=500;
N=4096;
z=fft(a);
Pyy =z.*conj(z)/N;
f = 500*(0:256)/N;
subplot(3,1,2)
plot(f,Pyy(1:257),'r') 
xlabel('频率 (Hz)')
ylabel('能量')
title('频谱图(频域)')
grid on;
h=0:2999
subplot(3,1,3)
plot(h,Pyy(1:3000),'r') 
xlabel('频率 (Hz)')
ylabel('能量')
title('频谱图(频域)')
grid on;
he=sum(Pyy)
maxp=max(Pyy)
thd=(sum(Pyy)-maxp)/maxp
[b c]=max(Pyy(1:N));
fprintf('\n数据基本信息：\n')
fprintf('        采样点数 = %7.0f \n',length(a))      
fprintf('        采样时间 = %7.3f ms\n',length(a)*2)                   
fprintf('        采样频率 = %7.1f Hz\n',1000/2)
fprintf('        标准方差 = %7.3f \n',std(a))
fprintf('         协 方 差 = %7.3f \n',cov(a))
fprintf('     自相关系数 = %7.3f \n\n',corrcoef(a)) 
fprintf('\n傅立叶变换结果：\n')
fprintf('           基波频率= %1.3f Hz\n',f(c))            
fprintf('           基波周期 = %1.3f s\n',1/f(c)) 
fprintf('           谐波畸变 = %1.3f \n',thd)