声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 2476|回复: 10

[FFT] matlab如何实现功率谱幅值曲线

[复制链接]
发表于 2008-12-15 17:13 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
我有一列响应数据,4096行。
怎么用matlab实现功率谱曲线。谢谢!

如果有几列数据,互功率谱曲线怎么画。

多谢!!!急!
回复
分享到:

使用道具 举报

 楼主| 发表于 2008-12-15 18:10 | 显示全部楼层
期待有人帮忙
发表于 2008-12-15 22:57 | 显示全部楼层
可以用pwelch命令
发表于 2008-12-16 09:35 | 显示全部楼层
本帖最后由 wdhd 于 2016-9-13 14:04 编辑

PWELCH Power Spectral Density estimate via Welch's method.
    Pxx = PWELCH(X) returns the Power Spectral Density (PSD) estimate, Pxx, of a discrete-time signal vector X using Welch's averaged,  modified periodogram method.  By default, X is divided into eight sections with 50% overlap, each section is windowed with a Hamming window and eight modified periodograms are computed and averaged.

    If the length of X is such that it cannot be divided exactly into eight sections with 50% overlap, X will be truncated accordingly.

    Pxx is the distribution of power per unit frequency. For real signals, PWELCH returns the one-sided PSD by default; for complex signals, it returns the two-sided PSD.  Note that a one-sided PSD contains the total power of the input signal.
发表于 2008-12-16 09:36 | 显示全部楼层
本帖最后由 wdhd 于 2016-9-13 14:05 编辑

See also periodogram, pcov, pmcov, pburg, pyulear, peig, pmtm, pmusic, spectrum/welch, dspdata/psd, dspdata/msspectrum.
 楼主| 发表于 2008-12-19 09:51 | 显示全部楼层
多谢回复!
有人能给出详细的中文回答吗?英文不好。
发表于 2008-12-19 10:45 | 显示全部楼层
nfft=4096;
p=abs(fft(x,nfft));%x为你的数列名称
plot((0:nfft/2-1)/nfft*fs,p(1:nfft/2));% 即可画出功率谱曲线
发表于 2008-12-19 20:18 | 显示全部楼层
本帖最后由 wdhd 于 2016-9-13 14:05 编辑
原帖由 camel2018 于 2008-12-19 10:45 发表
nfft=4096;
p=abs(fft(x,nfft));%x为你的数列名称
plot((0:nfft/2-1)/nfft*fs,p(1:nfft/2));% 即可画出功率谱曲线

好像不对吧,这样画出来的是FFT变换后的幅值谱,不是功率谱
发表于 2008-12-19 22:38 | 显示全部楼层
这个论坛中好像有啊
多搜索下啊
我把我以前找到的贴一下
你参考下
把load data换成你的数据就行了

clear;

load x1.mat;
load y.mat;

y = y - mean(y);

plot(x,y);

%% 经典功率谱估计
Fs=416.6667; %采样频率
window=boxcar(length(y)); %矩形窗
nfft=512;
[Pxx,f]=periodogram(y,window,nfft,Fs); %直接法

figure
plot(f,Pxx);

%% 间接法
% Fs=100; %采样频率
% nfft=416.6667;
% cxn=xcorr(y,'unbiased'); %计算序列的自相关函数
% CXk=fft(cxn,nfft);
% Pxx=abs(CXk);
% index=0:round(nfft/2-1);
% k=index*Fs/nfft;
% plot_Pxx=10*log10(Pxx(index+1));
% plot(k,plot_Pxx);

%% Bartlett法(改进方法一)
% Fs=100;
% nfft=416.6667;
% window=boxcar(256); %矩形窗
% noverlap=0; %数据无重叠
% p=0.9; %置信概率
% [Pxx,Pxxc]=psd(y,nfft,Fs,window,noverlap,p);
% index=0:round(nfft/2-1);
% k=index*Fs/nfft;
% plot_Pxx=(Pxx(index+1));
% plot_Pxxc=(Pxxc(index+1));
% figure(1)
% plot(k,plot_Pxx);
% pause;
% figure(2)
% plot(k,[plot_Pxx plot_Pxx-plot_Pxxc plot_Pxx+plot_Pxxc]);

%% 韦尔奇法

% Fs=100;
% nfft=416.6667;
% window=boxcar(100); %矩形窗
% window1=hamming(100); %海明窗
% window2=blackman(100); %blackman窗
% noverlap=20; %数据无重叠
% range='half'; %频率间隔为[0 Fs/2],只计算一半的频率
% [Pxx,f]=pwelch(y,window,noverlap,nfft,Fs,range);
% [Pxx1,f]=pwelch(y,window1,noverlap,nfft,Fs,range);
% [Pxx2,f]=pwelch(y,window2,noverlap,nfft,Fs,range);
% plot_Pxx=(Pxx);
% plot_Pxx1=(Pxx1);
% plot_Pxx2=(Pxx2);
%
% figure(1)
% plot(f,plot_Pxx);
%
% pause;
%
% figure(2)
% plot(f,plot_Pxx1);
%
% pause;
%
% figure(3)
% plot(f,plot_Pxx2);

评分

1

查看全部评分

 楼主| 发表于 2008-12-20 13:53 | 显示全部楼层
多谢!问题已解决!
 楼主| 发表于 2008-12-20 13:56 | 显示全部楼层
另外一个问题,我有一个3*4096的矩阵,想用它来构造hankel矩阵,比如说构造为30*30的hankel矩阵,改怎么做?
直接hankel不行。初学者,多谢!
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-11-11 18:21 , Processed in 0.062350 second(s), 19 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表