声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 2353|回复: 18

[FFT] 请高手将记录的一段波形进行傅立叶分析并画图

[复制链接]
发表于 2009-9-24 21:28 | 显示全部楼层 |阅读模式

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

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

x
小弟我记录到一段波形,并保存为*.mat格式,请问如何在MATLAB将其转为可进行分析的函数,并进行傅立叶分析?
请高手给一个程序学习下.

附件为记录到的波形

medal2.0.mat

189.17 KB, 下载次数: 36

记录到的波形

回复
分享到:

使用道具 举报

发表于 2009-9-24 23:20 | 显示全部楼层
help load, fft
 楼主| 发表于 2009-9-25 11:21 | 显示全部楼层
load之后怎么把这个文件变成可处理的数据呢?
 楼主| 发表于 2009-9-25 11:47 | 显示全部楼层
我按书上例子做了以下程序,可以画出原波形,想问下,X轴怎么弄成跟原始图一样?

但傅立叶分析却有问题.
想问下,傅立叶分析中频率轴跟采样率有什么关系,我的采样率很高,有三万.FFT坐标轴该如何设?
那个f=1000*(0:2048)/4096; 这为什么要这样算频率?

d=load('D:\MATLAB\R2008a\work\medal2.0.mat');
y=d.medal_Ch1.values;
t=0:0.00003:25.6;
figure(1);
plot(y)
Y=fft(y,4096);
p=Y.*conj(Y)/4096;
f=1000*(0:2048)/4096;    %原来采样率是30000,
figure(2);
plot(f,p(1:2049))

[ 本帖最后由 laoxuzi 于 2009-9-25 11:49 编辑 ]

原波形图

原波形图
fft.jpg
发表于 2009-9-25 14:34 | 显示全部楼层
load medal2.0.mat;
TimeData=medal_Ch1.values; SampleRate=1/medal_Ch1.interval;
Time=[1:length(TimeData)]'/SampleRate;
subplot(2,1,1); plot(Time,TimeData);
nfft=8192; noverlap=nfft/2; window=ones(nfft,1);
[P,f]=spectrum(TimeData,nfft,noverlap,window,SampleRate);
subplot(2,1,2); plot(f,P(:,1))
发表于 2009-9-25 14:44 | 显示全部楼层
两点注意下
1.应该是f=30000*(1:2048)/4096;
2.LZ未去除趋势线, help detrend
LZ参考下楼上的东西
 楼主| 发表于 2009-9-25 17:06 | 显示全部楼层
感谢楼上热心帮助.不太看得懂你的程序,看了趋势线的帮助,没理解.能把我的程序修改下吗?
另外,运行你的程序出现的结果好象不太对.

用自带的记录软件处理后,出现5HZ 10hz的二个波峰.

[ 本帖最后由 laoxuzi 于 2009-9-25 17:09 编辑 ]
ad.jpg
发表于 2009-9-25 18:17 | 显示全部楼层
好奇?
有兴趣的频率(5/10Hz), 不高为何取样率取那麽大?
 楼主| 发表于 2009-9-25 20:43 | 显示全部楼层
这个波形是通过电生理放大器记录的, 这个设备本是记录神经元放电的,所以采样频率大.
我现在记录的是一个金属棒放在水里的波形,用这个配套的软件的功率普功能处理得到10/5HZ的峰,但坐标轴也是到1500,不拉坐标轴也看到这个二个峰值的.
不知道能不能在MATLAB里处理,只显示0-100r 的坐标范围,这样是不是可以看出来?

谢谢帮助~

[ 本帖最后由 laoxuzi 于 2009-9-25 20:51 编辑 ]
100.bmp
发表于 2009-9-25 22:18 | 显示全部楼层
我将频域的解析度提高, 即改为nfft=length(TimeData)
亦找不到5/10HZ, 个人水平专业有限, 待高人路过!:@L
还有先前用的spectrum是v5.6版的函数, 新版好像改了
axis([0,100,-inf,inf]) 可以改变坐标范围,
 楼主| 发表于 2009-9-25 23:27 | 显示全部楼层
难道不同的处理方法显示的频率会有差异?我用SPIKE 2处理的.

不过,几千至上万赫兹的信号,似乎也不太可能.
我们通常记录的交流电干扰也就是50/60HZ这样子.
我按我上面的程序处理时,坐标轴左(零点附近)有一个细线,但拉不开.

[ 本帖最后由 laoxuzi 于 2009-9-25 23:36 编辑 ]
12.jpg
 楼主| 发表于 2009-9-25 23:37 | 显示全部楼层
load medal_2.0.mat;
y=Data1__Stopped__Ch1.values;
Fs=1/Data1__Stopped__Ch1.interval;
t=[1:length(y)];
subplot(2,1,1);
plot(t,y)
Y=fft(y,512);
p=Y.*conj(Y)/512;
f=Fs*(0:256)/512;
subplot(2,1,2);
plot(f,p(1:257));
axis([0,100,-inf,inf])

[ 本帖最后由 laoxuzi 于 2009-9-25 23:38 编辑 ]
121.jpg

想做成这样的图

想做成这样的图
发表于 2009-9-26 01:00 | 显示全部楼层

回复 12楼 laoxuzi 的帖子

看看楼主贴的图与程序, 有些问题
1. f=Fs*(0:256)/512; 表示f(2)=Fs/512= 58.12Hz, 但图121.jpg第二点约为2Hz??
2. untitled.jpg中在50ms中约有50个点, 所以取样率约为1000Hz??
楼主something is wrong!!
 楼主| 发表于 2009-9-26 09:22 | 显示全部楼层
121图是采样率改成1000记录的.坐标拉开了.
我想做成图像untitled.jpg,untitled是书中的例子,也是1000采样.我也不知道哪里出错?我是一个初学者,还得请教.
发表于 2009-9-26 10:39 | 显示全部楼层
LZ的信号中有很大的直流分量,在FFT分析前要去除直流分量。又按LZ要求画出500HZ以下的图
load medal2.0.mat
fs=1/medal_Ch1.interval;
x=medal_Ch1.values;
N=length(x);
t=(0:N-1)/fs;
subplot 211; plot(t,x); title('Signal');
xlabel('Time');
x=x-mean(x);
S=fft(x);
N2=N/2;
n2=1:N2;
df=fs/N;
F=(n2-1)*df;
subplot 212; plot(F,abs(S(n2))); grid;
xlim([0 500]); title('Signal Spectrum');
xlabel('Frequency (Hz)');
得到的图如下
lz11a.jpg

评分

1

查看全部评分

您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-5-16 10:17 , Processed in 0.054621 second(s), 22 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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