关于谱能量，有这样一种解释，你可以试着去算一算
信号可以分成能量信号与
功率信号，非周期能量信号具有能量谱密度，是傅立叶变换的平方，功率信号具有功率谱密度，其与自相关
函数是一对傅立叶变换对，等于傅立叶变换的平方/区间长度。不能混淆。能量信号是没有功率谱的。胡广书老师的书上找到这么一段话，“随机信号在时间上是无限的，在样本上也是无穷多，因此随机信号的能量是无限的，它应是功率信号。功率信号不满足付里叶变换的绝对可积的条件，因此其付里叶变换是不存在的。如确定性的正弦函数的付里叶变换是不存在，只有引入了冲激函数才求得其付里叶变换。因此，对随机信号的频谱分析，不再简单的是频谱，而是功率谱。”对于确定性信号而言，里面存在能量信号，是没有功率谱密度的，也存在功率信号，是有功率谱密度的。所以信号的频谱与是否是确定性信号没有必然联系。以下论点来源于研学论坛：频谱是信号的傅立叶变换。它描述了信号在各个频率上的分布大小。频谱的平方(当能量有限，平均功率为0时称为能量谱)描述了信号能量在各个频率上的分布大小。计算过程中，都是通过样本数据的快速傅立叶变换来计算。但不同的是，信号的频谱是复数，包含幅频响应和相频响应，重复计算时的结果基本相同。而随机信号的功率谱也可以对数据进行FFT，但必须计算模值的平方，因为功率谱是实数。而且换一组样本后，计算的结果略有不同，因为随机信号的样本取值不同。要得到真实的功率谱必须进行多次平均，次数越多越好。根据parseval定理，信号傅氏变换模平方被定义为能量谱，即单位频率范围内包含的信号能量。自然，能量跟功率有一个时间平均的关系，所以，能量谱密度在时间上平均就得到了功率谱。matlab实现经典功率谱估计fft做出来是频谱，psd做出来是功率谱；功率谱丢失了频谱的相位信息；频谱不同的信号其功率谱是可能相同的；功率谱是幅度取模后平方，结果是个实数matlab中自功率谱密度直接用psd函数就可以求，按照matlab的说法，psd能实现Welch法估计，即相当于用改进的平均周期图法来求取随机信号的功率谱密度估计。psd求出的结果应该更光滑吧。1、直接法：直接法又称周期图法，它是把随机序列x(n)的N个观测数据视为一能量有限的序列，直接计算x(n)的离散傅立叶变换，得X(k)，然后再取其幅值的平方，并除以N，作为序列x(n)真实功率谱的估计。Matlab代码示例：clearFs=1000%采样频率n=0:1/Fs:1%产生含有噪声的序列xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n))window=boxcar(length(xn))%矩形窗nfft=1024[Pxx,f]=periodogram(xn,window,nfft,Fs)%直接法plot(f,10*log10(Pxx))2、间接法：间接法先由序列x(n)估计出自相关函数R(n)，然后对R(n)进行傅立叶变换，便得到x(n)的功率谱估计。Matlab代码示例：clearFs=1000%采样频率n=0:1/Fs:1%产生含有噪声的序列xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n))nfft=1024cxn=xcorr(xn,'unbiased')%计算序列的自相关函数CXk=fft(cxn,nfft)Pxx=abs(CXk)index=0:round(nfft/2-1)k=index*Fs/nfftplot_Pxx=10*log10(Pxx(index+1))plot(k,plot_Pxx)3、改进的直接法：对于直接法的功率谱估计，当数据长度N太大时，谱曲线起伏加剧，若N太小，谱的分辨率又不好，因此需要改进。3.1、Bartlett法Bartlett平均周期图的方法是将N点的有限长序列x(n)分段求周期图再平均。Matlab代码示例：clear；Fs=1000n=0:1/Fs:1xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n))nfft=1024window=boxcar(length(n))%矩形窗noverlap=0%数据无重叠p=0.9%置信概率[Pxx,Pxxc]=psd(xn,nfft,Fs,window,noverlap,p)index=0:round(nfft/2-1)k=index*Fs/nfftplot_Pxx=10*log10(Pxx(index+1))plot_Pxxc=10*log10(Pxxc(index+1))figure(1)plot(k,plot_Pxx)pausefigure(2)plot(k,[plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc])3.2、Welch法Welch法对Bartlett法进行了两方面的修正，一是选择适当的窗函数w(n)，并再周期图计算前直接加进去，加窗的优点是无论什么样的窗函数均可使谱估计非负。二是在分段时，可使各段之间有重叠，这样会使方差减小。Matlab代码示例：clearFs=1000n=0:1/Fs:1xn=cos(2*pi*40*n)+3*cos(2*pi*100*n)+randn(size(n))nfft=1024window=boxcar(100)%矩形窗window1=hamming(100)%海明窗window2=blackman(100)%blackman窗noverlap=20%数据无重叠range='half'%频率间隔为[0 Fs/2]，只计算一半的频率[Pxx,f]=pwelch(xn,window,noverlap,nfft,Fs,range)[Pxx1,f]=pwelch(xn,window1,noverlap,nfft,Fs,range)[Pxx2,f]=pwelch(xn,window2,noverlap,nfft,Fs,range)plot_Pxx=10*log10(Pxx)plot_Pxx1=10*log10(Pxx1)plot_Pxx2=10*log10(Pxx2)figure(1)plot(f,plot_Pxx)pausefigure(2)plot(f,plot_Pxx1)pausefigure(3)plot(f,plot_Pxx2)自相关函数和互相关函数的matlab计算和作图
1. 首先说说自相关和互相关的概念。
这个是信号分析里的概念，他们分别表示的是两个时间序列之间和同一个时间序列在任意两个不同时刻的取值之间的相关程度，即互相关函数是描述随机信号x(t),y(t)在任意两个不同时刻t1，t2的取值之间的相关程度，自相关函数是描述随机信号x(t)在任意两个不同时刻t1，t2的取值之间的相关程度。互相关函数给出了在频域内两个信号是否相关的一个判断指标，把两测点之间信号的互谱与各自的自谱联系了起来。它能用来确定输出信号有多大程度来自输入信号，对修正测量中接入噪声源而产生的误差非常有效.
事实上，在图象处理中，自相关和互相关函数的定义如下：设原函数是f(t)，则自相关函数定义为R(u)=f(t)*f(-t)，其中*表示卷积；设两个函数分别是f(t)和g(t)，则互相关函数定义为R(u)=f(t)*g(-t)，它反映的是两个函数在不同的相对位置上互相匹配的程度。那么，如何在matlab中实现这两个相关并用图像显示出来呢？
dt=.1
t=[0:dt:100]
x=cos(t)
[a,b ]=xcorr(x,'unbiased')
plot(b*dt,a)
上面代码是求自相关函数并作图，对于互相关函数，稍微修改一下就可以了，即把[a,b]=xcorr(x,'unbiased')改为[a,b]=xcorr(x,y,'unbiased')便可。
2. 实现过程：
在Matalb中，求解xcorr的过程事实上是利用Fourier变换中的卷积定理进行的，即R(u)=ifft(fft(f)×fft(g))，其中×表示乘法，注：此公式仅表示形式计算，并非实际计算所用的公式。当然也可以直接采用卷积进行计算，但是结果会与xcorr的不同。事实上，两者既然有定理保证，那么结果一定是相同的，只是没有用对公式而已。下面是检验两者结果相同的代码：
dt=.1
t=[0:dt:100]
x=3*sin(t)
y=cos(3*t)
subplot(3,1,1)
plot(t,x)
subplot(3,1,2)
plot(t,y)
[a,b]=xcorr(x,y)
subplot(3,1,3)
plot(b*dt,a)
yy=cos(3*fliplr(t))% or use: yy=fliplr(y)
z=conv(x,yy)
pause
subplot(3,1,3)
plot(b*dt,z,'r')
即在xcorr中不使用scaling。
3. 其他相关问题：
1) 相关程度与相关函数的取值有什么联系？
相关系数只是一个比率，不是等单位量度，无什么单位名称，也不是相关的百分数，一般取小数点后两位来表示。相关系数的正负号只表示相关的方向，绝对值表示相关的程度。因为不是等单位的度量，因而不能说相关系数0.7是0.35两倍，只能说相关系数为0.7的二列变量相关程度比相关系数为0.35的二列变量相关程度更为密切和更高。也不能说相关系数从0.70到0.80与相关系数从0.30到0.40增加的程度一样大。对于相关系数的大小所表示的意义目前在统计学界尚不一致，但通常按下是这样认为的：
相关系数 相关程度
0.00-±0.30 微相关
±0.30-±0.50 实相关
±0.50-±0.80 显著相关
±0.80-±1.00 高度相关
自相关函数是描述随机信号X(t)在任意两个不同时刻t1，t2的取值之间的相关程度.设原函数是f(t)，则自相关函数定义为R(u)=f(t)*f(-t)，其中*表示卷积.
给个例子:
dt=.1
t=[0:dt:100]
x=cos(t)
[a,b]=xcorr(x,'unbiased')
plot(b*dt,a)
上面代码是求自相关函数并作图，
matlab中查看帮助时,
help xcorr 解释其意思是：
C(m) = E[A(n+m)*conj(B(n))] = E[A(n)*conj(B(n-m))]；
但是,在调用xcorr函数求自相关时,有 scaleopt参数
r=xcorr(s,SCALEOPT)
SCALEOPT有
'biased' - scales the raw cross-correlation by 1/M.
'unbiased' - scales the raw correlation by 1/(M-abs(lags)).
'coeff'- normalizes the sequence so that the auto-correlations
at zero lag are identically 1.0.
'none' - no scaling (this is the default).
注意观察下面的测试:
s = [1 2 3]
r = xcorr(s)
r =
3.00008.0000 14.00008.00003.0000
当用r=xcorr(s,'unbiased')时就能得到
r =3.00004.00004.66674.00003.0000
