精华内容
下载资源
问答
  • Untitledcohere:两信号互相关函数互功率谱untitledxcov:自相关函数Untitledwelch:welch谱估计方法
  • 自功率谱于互功率谱

    千次阅读 2019-06-11 15:04:09
    周期图法是直接将信号的采样数据x(n)进行Fourier变换求取功率谱密度估计的方法。假定有限长随机信号序列为x(n)。它的Fourier变换和功率谱密度估计存在下面的关系: 式中,N为随机信号序列x(n)的长度。在离散的频率点...

    1. 基本方法
    周期图法是直接将信号的采样数据x(n)进行Fourier变换求取功率谱密度估计的方法。假定有限长随机信号序列为x(n)。它的Fourier变换和功率谱密度估计存在下面的关系:

    式中,N为随机信号序列x(n)的长度。在离散的频率点f=kΔf,有:

    其中,FFT[x(n)]为对序列x(n)的Fourier变换,由于FFT[x(n)]的周期为N,求得的功率谱估计以N为周期,因此这种方法称为周期图法。下面用例子说明如何采用这种方法进行功率谱
    用有限长样本序列的Fourier变换来表示随机序列的功率谱,只是一种估计或近似,不可避免存在误差。为了减少误差,使功率谱估计更加平滑,可采用分段平均周期图法(Bartlett法)、加窗平均周期图法(Welch法)等方法加以改进。
    2. 分段平均周期图法(Bartlett法)
    将信号序列x(n),n=0,1,…,N-1,分成互不重叠的P个小段,每小段由m个采样值,则P*m=N。对每个小段信号序列进行功率谱估计,然后再取平均作为整个序列x(n)的功率谱估计。
    平均周期图法还可以对信号x(n)进行重叠分段,如按2:1重叠分段,即前一段信号和后一段信号有一半是重叠的。对每一小段信号序列进行功率谱估计,然后再取平均值作为整个序列x(n)的功率谱估计。这两种方法都称为平均周期图法,一般后者比前者好。程序运行结果为图9-5,上图采用不重叠分段法的功率谱估计,下图为2:1重叠分段的功率谱估计,可见后者估计曲线较为平滑。与上例比较,平均周期图法功率谱估计具有明显效果(涨落曲线靠近0dB)。
    3.加窗平均周期图法
    加窗平均周期图法是对分段平均周期图法的改进。在信号序列x(n)分段后,用非矩形窗口对每一小段信号序列进行预处理,再采用前述分段平均周期图法进行整个信号序列x(n)的功率谱估计。由窗函数的基本知识(第7章)可知,采用合适的非矩形窗口对信号进行处理可减小“频谱泄露”,同时可增加频峰的宽度,从而提高频谱分辨率。
    其中上图采用无重叠数据分段的加窗平均周期图法进行功率谱估计,而下图采用重叠数据分段的加窗平均周期图法进行功率谱估计,显然后者是更佳的,信号谱峰加宽,而噪声谱均在0dB附近,更为平坦(注意采用无重叠数据分段噪声的最大的下降分贝数大于5dB,而重叠数据分段周期图法噪声的最大下降分贝数小于5dB)。
    4. Welch法估计及其MATLAB函数
    Welch功率谱密度就是用改进的平均周期图法来求取随机信号的功率谱密度估计的。Welch法采用信号重叠分段、加窗函数和FFT算法等计算一个信号序列的自功率谱估计(PSD如上例中的下半部分的求法)和两个信号序列的互功率谱估计(CSD)。
    MATLAB信号处理工具箱函数提供了专门的函数PSD和CSD自动实现Welch法估计,而不需要自己编程。
    (1) 函数psd利用Welch法估计一个信号自功率谱密度,函数调用格式为:
    [Pxx[,f]]=psd(x[,Nfft,Fs,window,Noverlap,’dflag’])
    式中,x为信号序列;Nfft为采用的FFT长度。这一值决定了功率谱估计速度,当Nfft采用2的幂时,程序采用快速算法;Fs为采样频率;Window定义窗函数和x分段序列的长度。窗函数长度必须小于或等于Nfft,否则会给出错误信息;Noverlap为分段序列重叠的采样点数(长度),它应小于Nfft;dflag为去除信号趋势分量的选择项:’linear’,去除线性趋势分量,’mean’去除均值分量,’none’不做去除趋势处理。Pxx为信号x的自功率谱密度估计。f为返回的频率向量,它和Pxx对应,并且有相同长度。
    在psd函数调用格式中,缺省值为:Nfft=min(256,length(x)),Fs=2Hz, window=hanning(Nfft),noverlap=0. 若x是实序列,函数psd仅计算频率为正的功率
    注意程序前半部分中频率向量f的创建方法。它与函数psd的输出Pxx长度的关系如下:若x为实序列,当Nfft为奇数时,f=(0:(Nfft+1)/2-1)/Nfft;当Nfft为偶数时,f=(0:Nfft/2)/Nfft。
    函数还有一种缺省返回值的调用格式,用于直接绘制信号序列x的功率谱估计曲线。
    函数还可以计算带有置信区间的功率谱估计,调用格式为:
    [Pxx,Pxxc,f]=psd(x,Nfft,Fs,window,Noverlap,p)
    式中,p为置信区间,0<=p<=1。

    由此可知,滤波器输入白噪声序列的输出信号的功率谱或自相关可以确定滤波器的频率特性。
    (2)函数csd利用welch法估计两个信号的互功率谱密度,函数调用格式为:
    [Pxy[,f]]=csd(x,y,Nfft,Fs,window,Noverlap,’dflag’)
    [Pxy,Pxyc[,f]]=csd(x,y,Nfft,Fs,window,Noverlap,p)
    这里,x,y为两个信号序列;Pxy为x,y的互功率谱估计;其他参数的意义同自功率谱函数psd。
    可以看到,两个白噪声信号的互功率谱(上图)杂乱无章,看不出周期成分,大部分功率谱在-5dB以下。然而白噪声与带有噪声的周期信号的功率谱在其周期(频率为1000Hz)处有一峰值,清楚地表明了周期信号的周期或频率。因此,利用未知信号与白噪声信号的互功率谱也可以检测未知信号中所含有的频率成分。

    5 多 窗 口 法
    多窗口法(Multitaper method,简称MTM法)利用多个正交窗口(Tapers)获得各自独立的近似功率谱估计,然后综合这些估计得到一个序列的功率谱估计。相对于普通的周期图法,这种功率谱估计具有更大的自由度,并在估计精度和估计波动方面均有较好的效果。普通的功率谱估计只利用单一窗口,因此在序列始端和末端均会丢失相关信息,而且无法找回。而MTM法估计增加窗口使得丢失的信息尽量减少。
    MTM法简单地采用一个参数:时间带宽积(Time-bandwidth product)NW,这个参数用以定义计算功率谱所用窗的数目,为2*NW-1。NW越大,功率谱计算次数越多,时间域分辨率越高,而频率域分辨率降低,使得功率谱估计的波动减小。随着NW增大,每次估计中谱泄漏增多,总功率谱估计的偏差增大。对于每一个数据组,通常有一个最优的NW使得在估计偏差和估计波动两方面求得折中,这需要在程序中反复调试来获得。
    MATLAB信号处理工具箱中函数PMTM就是采用MTM法估计功率谱密度。函数调用格式为:
    [Pxx[,f]]=pmtm(x[,nw,Nfft,Fs])
    式中,x为信号序列;nw为时间带宽积,缺省值为4。通常可取2,5/2,3,7/2;Nfft为FFT长度;Fs为采样频率
    上面的函数还可以通过无返回值而绘出置信区间,如pmtm(x,nw,Nfft,Fs,’option’,p)绘制带置信区间的功率谱密度估计曲线,0<=p<=1。
    6 最 大 熵 法(Maxmum entropy method, MEM法)
    如上所述,周期图法功率谱估计需要对信号序列“截断”或加窗处理,其结果是使估计的功率谱密度为信号序列真实谱和窗谱的卷积,导致误差的产生。
    最大熵功率谱估计的目的是最大限度地保留截断后丢失的“窗口”以外信号的信息,使估计谱的熵最大。主要方法是以已知的自相关序列rxx(0),rxx(1),…,rxx(p)为基础,外推自相关序列rxx(p+1),rxx(p+2),…,保证信息熵最大。
    最大熵功率谱估计法假定随机过程是平稳高斯过程,可以证明,随机信号的最大熵谱与AR自回归(全极点滤波器)模型谱是等价的。
    MATLAB信号处理工具箱提供最大熵功率谱估计函数pmem,其调用格式为:
    [Pxx,f,a]=pmem(x,p,Nfft,Fs,’xcorr’)
    式中,x为输入信号序列或输入相关矩阵;p为全极点滤波器阶次;a为全极点滤波器模型系数向量;’xcorr’是把x认为是相关矩阵。

    比较最大熵功率谱估计(MEM)和改进的平均周期图功率谱估计,可见,MEM法估计的功率谱曲线较光滑。在这一方法中,MEM法选定全极点滤波器的阶数取得越大,能够获得的窗口外的信息越多,但计算量也越大,需要根据情况折中考虑。

    7 多信号分类法
    MATLAB信号处理工具箱还提供另一种功率谱估计函数pmusic。该函数执行多信号分类法(multiple signal classification, Music法)。将数据自相关矩阵看成由信号自相关矩阵和噪声自相关矩阵两部分组成,即数据自相关矩阵R包含有两个子空间信息:信号子空间和噪声子空间。这样,矩阵特征值向量(Eigen vector)也可分为两个子空间:信号子空间和噪声子空间。为了求得功率谱估计,函数pmusic计算信号子空间和噪声子空间的特征值向量函数,使得在周期信号频率处函数值最大,功率谱估计出现峰值,而在其他频率处函数值最小。其调用格式为:
    [Pxx,f,a]=pmusic(x,p[[,thresh],Nfft,Fs,window,Noverlap])
    式中,x为输入信号的向量或矩阵;p为信号子空间维数;thresh为阈值,其他参数的意义与函数psd相同。

    功率谱密度相关方法的MATLAB实现
    %%
    %分段平均周期图法(Bartlett法)
    %运用信号不重叠分段估计功率谱
    Nsec=256;n=0:sigLength-1;t=n/Fs; %数据点数,分段间隔,时间序列
    pxx1=abs(fft(y(1:256),Nsec).^2)/Nsec; %第一段功率谱
    pxx2=abs(fft(y(257:512),Nsec).^2)/Nsec; %第二段功率谱
    pxx3=abs(fft(y(515:768),Nsec).^2)/Nsec; %第三段功率谱
    pxx4=abs(fft(y(769:1024),Nsec).^2)/Nsec; %第四段功率谱
    Pxx=10*log10((pxx1+pxx2+pxx3+pxx4)/4); %平均得到整个序列功率谱
    f=(0:length(Pxx)-1)*Fs/length(Pxx); %给出功率谱对应的频率
    %%plot(f(1:Nsec/2),Pxx(1:Nsec/2)); %绘制功率谱曲线
    figure,plot(f(1:Nsec),Pxx(1:Nsec)); %绘制功率谱曲线
    xlabel('频率/Hz');ylabel('功率谱 /dB');
    title('平均周期图(无重叠) N=4*256');
    grid on
    %%
    %运用信号重叠分段估计功率谱
    pxx1=abs(fft(y(1:256),Nsec).^2)/Nsec; %第一段功率谱
    pxx2=abs(fft(y(129:384),Nsec).^2)/Nsec; %第二段功率谱
    pxx3=abs(fft(y(257:512),Nsec).^2)/Nsec; %第三段功率谱
    pxx4=abs(fft(y(385:640),Nsec).^2)/Nsec; %第四段功率谱
    pxx5=abs(fft(y(513:768),Nsec).^2)/Nsec; %第五段功率谱
    pxx6=abs(fft(y(641:896),Nsec).^2)/Nsec; %第六段功率谱
    pxx7=abs(fft(y(769:1024),Nsec).^2)/Nsec; %第七段功率谱
    Pxx=10*log10((pxx1+pxx2+pxx3+pxx4+pxx5+pxx6+pxx7)/7); %功率谱平均并转化为dB
    f=(0:length(Pxx)-1)*Fs/length(Pxx); %频率序列
    figure,plot(f(1:Nsec/2),Pxx(1:Nsec/2)); %绘制功率谱曲线
    xlabel('频率/Hz'); ylabel('功率谱/dB');
    title('平均周期图(重叠一半) N=1024');
    grid on

    %%
    Nsec=256;n=0:sigLength-1;t=n/Fs; %数据长度,分段数据长度、时间序列
    w=hanning(256); %采用的窗口数据
    %采用不重叠加窗方法的功率谱估计
    pxx1=abs(fft(w.*y(1:256),Nsec).^2)/norm(w)^2; %第一段加窗振幅谱平方
    pxx2=abs(fft(w.*y(257:512),Nsec).^2)/norm(w)^2; %第二段加窗振幅谱平方
    pxx3=abs(fft(w.*y(513:768),Nsec).^2)/norm(w)^2; %第三段加窗振幅谱平方
    pxx4=abs(fft(w.*y(769:1024),Nsec).^2)/norm(w)^2; %第四段加窗振幅谱平方
    Pxx=10*log10((pxx1+pxx2+pxx3+pxx4)/4); %求得平均功率谱,转换为dB
    f=(0:length(Pxx)-1)*Fs/length(Pxx); %求得频率序列
    figure
    subplot(2,1,1),plot(f(1:Nsec/2),Pxx(1:Nsec/2)); %绘制功率谱曲线
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('加窗平均周期图(无重叠) N=4*256');
    grid on
    %采用重叠加窗方法的功率谱估计
    pxx1=abs(fft(w.*y(1:256),Nsec).^2)/norm(w)^2; %第一段加窗振幅谱平方
    pxx2=abs(fft(w.*y(129:384),Nsec).^2)/norm(w)^2; %第二段加窗振幅谱平方
    pxx3=abs(fft(w.*y(257:512),Nsec).^2)/norm(w)^2; %第三段加窗振幅谱平方
    pxx4=abs(fft(w.*y(385:640),Nsec).^2)/norm(w)^2; %第四段加窗振幅谱平方
    pxx5=abs(fft(w.*y(513:768),Nsec).^2)/norm(w)^2; %第五段加窗振幅谱平方
    pxx6=abs(fft(w.*y(641:896),Nsec).^2)/norm(w)^2; %第六段加窗振幅谱平方
    pxx7=abs(fft(w.*y(769:1024),Nsec).^2)/norm(w)^2; %第七段加窗振幅谱平方
    Pxx=10*log10((pxx1+pxx2+pxx3+pxx4+pxx5+pxx6+pxx7)/7);%平均功率谱转换为dB
    f=(0:length(Pxx)-1)*Fs/length(Pxx); %频率序列
    subplot(2,1,2),plot(f(1:Nsec/2),Pxx(1:Nsec/2)); %绘制功率谱曲线
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('加窗平均周期图(重叠一半)N=1024');
    grid on
    %%
    %4分段平均周期图法(hanning窗)
    Nsec=256;
    n=0:sigLength-1;
    t=n/Fs;
    w=hanning(256);
    Pxx1=abs(fft(w.*y(1:256),Nsec).^2)/Nsec;
    Pxx2=abs(fft(w.*y(257:512),Nsec).^2)/Nsec;
    Pxx3=abs(fft(w.*y(513:768),Nsec).^2)/Nsec;
    Pxx4=abs(fft(w.*y(769:1024),Nsec).^2)/Nsec;
    Pxx=10*log10((Pxx1+Pxx2+Pxx3+Pxx4)/4);
    f=(0:length(Pxx)-1)*Fs/length(Pxx);
    figure
    subplot(2,1,1)
    plot(f,Pxx);
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('Averaged Modified Periodogram (none overlap) N=4*256');
    grid
    %4分段(2:1重叠)平均周期图法(hanning窗)
    Nsec=256;
    n=0:sigLength-1;
    t=n/Fs;
    w=hanning(256);
    Pxx1=abs(fft(w.*y(1:256),Nsec).^2)/Nsec;
    Pxx2=abs(fft(w.*y(129:384),Nsec).^2)/Nsec;
    Pxx3=abs(fft(w.*y(257:512),Nsec).^2)/Nsec;
    Pxx4=abs(fft(w.*y(385:640),Nsec).^2)/Nsec;
    Pxx5=abs(fft(w.*y(513:768),Nsec).^2)/Nsec;
    Pxx6=abs(fft(w.*y(641:896),Nsec).^2)/Nsec;
    Pxx7=abs(fft(w.*y(769:1024),Nsec).^2)/Nsec;
    Pxx=10*log10((Pxx1+Pxx2+Pxx3+Pxx4+Pxx5+Pxx6+Pxx7)/7);
    f=(0:length(Pxx)-1)*Fs/length(Pxx);
    subplot(2,1,2);plot(f,Pxx);
    xlabel('频率/Hz');ylabel('Power Spectrum (dB)');
    title('Averaged Modified Periodogram (half overlap) N=1024');
    grid
    %%
    %PSD_WELCH方法
    %采样频率
    Nfft=256;n=0:sigLength-1;t=n/Fs; %数据长度、时间序列
    window=hanning(256); %选用的窗口
    noverlap=128; %分段序列重叠的采样点数(长度)
    dflag='none'; %不做趋势处理
    [Pxx,Pxxc,f]=psd(y,Nfft,Fs,window,noverlap,0.95); %功率谱估计,并以0.95的置信度给出置信区间,无返回值是绘制出置信区间
    figure
    plot(f,10*log10(Pxx)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('PSD—Welch方法'); grid on

    %%
    %最大熵法(MEM法)
    Nfft=256;n=0:sigLength-1;t=n/Fs; %数据长度、分段长度和时间序列
    window=hanning(256); %采用窗口
    [Pxx1,f]=pmem(x,20,Nfft,Fs); %采用最大熵法,采用滤波器阶数14,估计功率谱
    figure,subplot(2,1,1),plot(f,10*log10(Pxx1)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('最大熵法 Order=20原始信号功率谱');grid on
    [Pxx1,f]=pmem(y0,20,Nfft,Fs); %采用最大熵法,采用滤波器阶数14,估计功率谱
    subplot(2,1,2),plot(f,10*log10(Pxx1)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -20 0])
    title('最大熵法 Order=20滤波后的信号功率谱');grid on
    %%
    %%功率谱密度
    %PSD_WELCH方法
    Nfft=512;n=0:sigLength-1;t=n/Fs; %数据长度、时间序列
    window=hanning(256); %选用的窗口
    noverlap=128; %分段序列重叠的采样点数(长度)
    dflag='none'; %不做趋势处理
    [Pxx,Pxxc,f]=psd(x,Nfft,Fs,window,noverlap,0.95); %功率谱估计,并以0.95的置信度给出置信区间,无返回值是绘制出置信区间
    figure;subplot(211);
    plot(f,10*log10(Pxx)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');
    grid on;title('PSD—Welch方法的原始信号功率谱')
    subplot(212)
    [Pxx,Pxxc,f]=psd(y0,Nfft,Fs,window,noverlap,0.95); %功率谱估计,并以0.95的置信度给出置信区间,无返回值是绘制出置信区间
    plot(f,10*log10(Pxx)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -30 0])
    grid on;title('PSD—Welch方法的滤波后的信号功率谱')
    %%
    %用多窗口法(MTM)
    n=0:sigLength-1;t=n/Fs; %数据长度、分段数据长度,时间序列
    [Pxx1,f]=pmtm(x,2,Nfft,Fs); %用多窗口法(NW=4)估计功率谱
    figure;subplot(2,1,1),plot(f,10*log10(Pxx1)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');
    title('多窗口法(MTM) nw=2原始信号功率谱');grid on
    [Pxx,f]=pmtm(y0,2,Nfft,Fs); %用多窗口法(NW=2)估计功率谱
    subplot(2,1,2),plot(f,10*log10(Pxx)); %绘制功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB');axis([0 4000 -80 -20])
    title('多窗口法(MTM) nw=2滤波后的信号功率谱');grid on

    %%
    %采用Welch方法估计功率谱
    noverlap=128; %重叠数据
    dflag='none'; Nfft=1024;
    figure;subplot(2,1,1)
    psd(x,Nfft,Fs,window,noverlap,dflag); %采用Welch方法估计功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB')
    title('Welch方法原始信号功率谱');grid on
    subplot(2,1,2)
    psd(y0,Nfft,Fs,window,noverlap,dflag); %采用Welch方法估计功率谱
    xlabel('频率/Hz');ylabel('功率谱/dB'),axis([0 4000 -30 0])
    title('Welch方法滤波后的信号功率谱');grid on
    pmusic(y0,[7,1.1],Nfft,Fs,32,16); %采用多信号分类法估计功率谱
    xlabel('频率/Hz'); ylabel('功率谱/dB')
    title('通过MUSIC法估计的伪谱')

    展开全文
  • 利用相关函数求信号功率谱,信号的自相关函数及不同信号的相关函数-Use of correlation function for the signal power spectrum, the signal of the autocorrelation function and different signal cross-...
  • csd函数估计互功率谱

    2021-08-09 20:48:46
    csd函数估计互功率谱
  • matlab cpsd互功率谱

    千次阅读 2020-06-24 08:27:46
    clear; clc; fs = 1024; len = 1024; t = 0:1/fs:(len-1)/fs; x = 3.7 * sin(2pi16t); y = 1.9 * sin(2pi32t); w = blackman(1024, ‘periodic’); [p, f] = cpsd(x, y, w, 0, 1024, fs...cnx = z(len:end) + [0 z(1:le

    clear;
    clc;

    fs = 1024;
    len = 1024;
    t = 0:1/fs:(len-1)/fs;
    x = 3.7 * sin(2pi16t);
    y = 1.9 * sin(2
    pi32t);

    w = blackman(1024, ‘periodic’);

    [p, f] = cpsd(x, y, w, 0, 1024, fs);
    plot(f, abs§)
    title(‘cpsd’)
    max(abs§)

    z = xcorr(x, y);
    cnx = z(len:end) + [0 z(1:len-1)];
    w = w’;
    x1 = cnx .* w;
    cxk1 = fft(x1, 1024);
    p1=abs(cxk1) * 2 / (1024 * fs);

    w1 = w .* w;
    u = sum(w1) / 1024;
    %p = p1/u;
    p = p1;

    f = 0:fs/1024:fs/2;
    figure
    plot(f, p(1:1024/2+1))
    title(‘xcorr’)
    max§

    展开全文
  • https://blog.csdn.net/weixin_42612518/article/details/93491303

    https://blog.csdn.net/weixin_42612518/article/details/93491303

    展开全文
  • 基于互功率谱(相位相关)的全局运动检测方法,可以对存在平移、旋转、缩放情况下图像运动检测,实现对图像的快速配准。基于互功率谱(相位相关)的理论基础是傅里叶变换,目前在傅里叶变换领域有了快速算法fft,...

    基于互功率谱(相位相关)的全局运动检测方法,可以对存在平移、旋转、缩放情况下图像运动检测,实现对图像的快速配准。基于互功率谱(相位相关)的理论基础是傅里叶变换,目前在傅里叶变换领域有了快速算法fft,因此速度较快,在图像配准、模式识别特征匹配等有着广泛应用。
      1)图像间有平移变换。
            图像f2(x,y)是图像f1(x,y)经平移(x0,y0)后得到的图像,即:
                f2(x,y)=f1(x-x0,y-y0)
         由傅里叶时移性质对应傅里叶变换F1和F2的关系如下:
               F2(u,v)=exp(-j*2*pi(u*x0+v*y0))*F1(u,v)
         计算互功率谱可得:
      exp(j2pi(u*x0+v*y0))=F1(u,v)*conj(F2)/F2(u,v)*conj(F2)
        当对上式进行傅里叶反变换后,将会得到一个冲击函数,该函数在两幅图像的相对位移(x0,y0)即“匹配点”处取到极大值,其他地方几乎为零。
      2)针对图像间有平移旋转缩放变换关系:
         若图像f2(x,y)是图像f1(x,y)经平移(x0,y0)、旋转a角度、缩放尺度σ后得到的图像,用下面公式表示为:
    f2(x,y)=f1(σ(x*cosa+y*sina)-x0,σ(-x*sina+y*cosa))-y0))
         由傅里叶旋转平移特性,fft变换后两图像间的关系如下:
    F2(u,v)=exp(-j2pi(u*x0+v*y0))/(σ*σ)*F1((u*cosa+v*sina)/σ,(-u*sina+v*cosa)/σ)
    用M1、M2分别表示F1、F2的能量,则:
    M2(u,v)=1/(σ*σ)*M1(u*cosa+v*sina,-u*sina+v*cosa);
    由上式看出F1、F2能量是相同的。把直角坐标转到极坐标可表示如下:
           M1(θ,logρ)=1/(σ*σ)*M2(θ-a,logρ-σ)
          再由互所述方法,在极坐标系下用相位相关可求出旋转角度a和尺度缩放σ,最后对图像以角度a做旋转,旋转得到图像与原图再次相位相关就可求出图像间的平移参数,完成图像的配准。
    实验验证:
    2、下面给出一组实验样本,以检验算法有效性,样本f1和f2是发生平移、旋转、缩放倍情况下的图像配准。
    1)样本图

    [转载]基于互功率谱或相位相关的图像配准、拼接  [转载]基于互功率谱或相位相关的图像配准、拼接

    两幅发生了平移、旋转、缩放情况的图像配准

    2)幅度谱
    [转载]基于互功率谱或相位相关的图像配准、拼接  [转载]基于互功率谱或相位相关的图像配准、拼接

     3)对极转换,极平面谱

    [转载]基于互功率谱或相位相关的图像配准、拼接  [转载]基于互功率谱或相位相关的图像配准、拼接

    4)第一次相位相关,计算旋转参数和尺度缩放参数
    [转载]基于互功率谱或相位相关的图像配准、拼接  [转载]基于互功率谱或相位相关的图像配准、拼接
       相关平面的平面图和立体图

    5)图像进行解旋和尺度变换

    [转载]基于互功率谱或相位相关的图像配准、拼接
    解旋后图片

     

    5)第二次计算互功率谱相关平面,求取平移参数
    [转载]基于互功率谱或相位相关的图像配准、拼接   [转载]基于互功率谱或相位相关的图像配准、拼接
            相关平面的平面图和立体图
     
    6)根据参数完成图像配准,拼接
    [转载]基于互功率谱或相位相关的图像配准、拼接

    展开全文
  • 随机信号及其自相关函数功率谱密度的MATLAB实现
  • 3种经典功率谱估计方法的MATLABA代码-功率谱代码.doc 3种MATLAB的经典谱估计方法 希望对大家有用~ 附件所有代码: 直接法: 直接法又称周期图法,它是把随机序列x的N个观测数据视为一能量有限的序列,直接计算x...
  • 基于MATLAB软件仿真分析输出信号的自相关函数功率谱密度,并画出图形。
  • 用matlab中求自()相关的xcorr函数,参考帮助文档 [r,lags]=xcorr(z,'biased'); 其中z是上述高斯白噪声,r是自相关函数,lags是时间偏移量(索引),尤其注意’biased’参数,这是调试了半天才发现的问题。 帮助...
  • 信号处理之功率谱原理与python实现

    千次阅读 2019-12-05 19:45:44
    目录功率谱简介功率谱、能量谱、幅值谱之间的关系功率谱python...功率谱功率谱密度函数的简称,它定义为单位频带内的信号功率。它表示了信号功率随着频率的变化情况,即信号功率在频域的分布状况。 功率谱表示了...
  • 自相关、相关函数学习笔记

    千次阅读 2019-07-17 15:56:50
    为什么相关性有效 参考文献 【1】自相关的物理意义 【2】自相关和相关函数计算方法总结及心得体会
  • 基于matlab,计算信号功率谱,调用了高阶谱分析工具箱的函数 bispecd 和 bispeci计算信号的双谱
  • 自相关函数相关函数的matlab 计算和作图 1. 首先说说自相关和相关的概念 这个是信号分析里的概念他们分别表示的是两个时间序列之间和同一个时 间序列在任意两个不同时刻的取值之间的相关程度即相关函数是...
  • 随机信号也称为随机过程、随机函数或随机序列。 确定性信号:如果序列{s(t)}在每个时刻的取值不是随机的,而是服从某种固定函数的关系,则称之为确定性信号。如阶跃信号、符号信号或矩形脉冲等。 随机信号:若序列{...
  • 我利用膨胀函数对原图进行了【50,140】的平移,然后利用互功率谱计算,再逆变换,得到峰值的坐标,但是下一步需要怎么换算才能得到我的偏移量呢,请求大神的解答。Thanks♪(・ω・)ノ
  • 对话功率谱与自相关函数

    万次阅读 多人点赞 2018-06-18 13:55:31
    引入:能量有限而平均功率为0的信号称之为能量信号;能量无限而平均功率有限的信号称之为功率信号;...计算能量谱密度以及功率谱密度:能量谱密度:功率谱密度:信号的自相关函数:同卷积的比较:由上图可以看出...
  • 功率谱密度相关方法的MATLAB实现

    万次阅读 2016-07-11 19:30:32
    1. 基本方法 ...周期图法是直接将信号的采样数据x(n)进行Fourier变换求取功率谱密度估计的方法。假定有限长随机信号序列为x(n)。它的Fourier变换和功率谱密度估计存在下面的关系:   式中,N
  • 功率谱学习及matlab代码

    千次阅读 2021-01-27 21:41:08
    只接触很少信号处理的问题,该篇是查阅资料总结的,先对概念等内容进行介绍,最后附matlab的功率谱代码。 看了很多资料,没有说明白为啥可以有这么多种方法计算,也不清楚具体这种方法计算出来的是否正确,就写了一...
  • 采用自功率谱、互功率谱和相干函数分析,识别了高桩码头物理模型前4阶固有频率和阻尼比。同时运用有限元软件ANSYS分析得到相应的固有频率,试验数据与ANSYS分析结果的对比表明,2种方法所得结果比较接近。
  • 此提交提供了使用 Welch 方法计算功率谱密度 (PSD) 的可能性。 该文件基于使用信号处理工具箱的 Matlab 实现。 我排除了计算机密间隔的可能性。 如果需要,请发表评论,我会更新必要的依赖项。 ...
  • 功率信号:能量无限,不能用能量表示,所以用平均功率表示; 能量信号:能量有限,平均功率为0; 二、功率信号的分析 频谱(离散):C(nf0)=1T0∫−T/2T/2s(t)e−j2πnf0tdtC\left( n{{f}_{0}} \r
  • 使用python 实现时间序列信号的频谱、倒频谱以及功率谱 ,资源中以振动信号为例,并封装了相应的函数,可以见我的博客,马上就能理解
  •   实际上,如果引入冲击函数 δ(·),功率对频率微分也可得到周期信号的功率谱密度,功率谱密度在基频整数倍为脉冲形式,即 G(ω) = ΣPnδ(ω-nΩ₀),同样满足功率谱密度的积分就等于平均功率 P。   以三角...
  • Matlab实现经典功率谱分析和估计

    千次阅读 2020-05-20 15:33:11
    功率谱功率谱密度函数的简称,它定义为单位频带内的信号功率。它表示了信号功率随着频率的变化情况,即信号功率在频域的分布状况。功率谱表示了信号功率随着频率的变化关系 。 常用于功率信号(区别于能量信号)的...
  • 深入理解傅里叶变换的性质:实函数、卷积、相关、功率谱、频响函数 1实函数傅里叶变换的性质 1.1实函数傅里叶变换的性质 所以,实函数x(t)的傅里叶变换X(w)的共轭 X*(w)=X(-w) 1.2实偶函数傅里叶变换的性质 1.2实...
  • [Matlab科学计算] 功率谱一点介绍

    千次阅读 2019-09-06 20:36:59
    利用给定的N个样本数据估计一个平稳随机信号的功率谱密度叫做谱估计,功率谱密度一般简称功率谱。谱估计方法分为参数化方法和非参数化方法。非参数化方法又叫经典谱估计,如周期图法、自相关法等,其主要缺点是描述...
  • 功率谱密度

    万次阅读 多人点赞 2017-05-13 20:59:12
    功率谱密度 缩写:PSD 定义:单位频率间隔的光功率或者噪声功率 在光学中,功率谱密度(有时称为功率密度)会以下面两种形式出现: 光功率谱密度,定义为单位频率(或者波长)间隔的光功率,例如,单位为 ...

空空如也

空空如也

1 2 3 4 5 ... 20
收藏数 775
精华内容 310
关键字:

互功率谱函数