A牛牛牛牛 发表于 2014-12-1 10:35

三分之一倍频程程序

01.clear all

02.close all

03.%读取数据文件

04.importfile('250.txt');

05.%****************************

06.%参数的初始化

07.%****************************

08.%定义三分之一倍频程的中心频率

09. center_frequency= ;

10. %中心频率总长度

11.Num_Channel = length(center_frequency);

12.%中心频率与上、下限频率的比值

13.oc6=2^(1/6);

14.%定义频段的上下限频率

15.lower_bound = center_frequency./oc6;

16.upper_bound = center_frequency.*oc6;

17.bound=;

18.data_length = 25000;

19.%采样时间间隔

20.time_interval = 0.0039062;

21.%采样频率

22.fs = 1./time_interval;

23.

24.

25.%进行FFT计算的数据长度

26.Num_dot = data_length;

27.

28.% Num_dot = min(find(t>55));

29.Num_FFT = Num_dot;

30.% %大于并最接近Num_dot的2的幂次方长度

31.% Num_FFT = 2^nextpow2(Num_dot);

32.data_fft = data(1:Num_dot);

33.% t_fft = t(1:Num_dot);

34.f=(0:Num_FFT - 1)*fs/Num_FFT;

35.%傅里叶变换

36.data_freq = fft(data_fft,Num_FFT);

37.%取幅度值

38.%--------------

39.%对于FFT以后的序列,其频率分辨率 = 采样频率/采样点数;

40.%直流分量处abs/N,其他点处abs/(N/2),就可得到幅值谱(频谱)

41.%--------------

42.% amplitude_freq=zeros(1,length(data_freq));

43.

44.amplitude_freq = abs(data_freq)*2./Num_dot;

45.amplitude_freq(1) = abs(data_freq(1))/Num_dot;

46.

47.t=(0:data_length-1)*time_interval;

48.for i = 1:Num_Channel

49.    %找到处在频段内的数据

50.    data_temp = amplitude_freq(find(f>bound(1,i) & f<bound(2,i)));

51.    %求均方值

52.    total_amplitude = sqrt(sum(data_temp.^2));

53.    center_freq_data(i) = total_amplitude;

54.end

55.figure

56.

57.plot(center_frequency,center_freq_data);

58.grid on;

59.xlabel('中心频率(Hz)');

60.ylabel('RMS(m/s2)')
复制代码

A牛牛牛牛 发表于 2014-12-1 10:53

自己顶一下

宇安乙昂 发表于 2017-2-23 10:56

楼主能否来一个三分之一倍频程谱图,我这边计算的谱图不知对不对,想和你的对比一下看看

TestGuru 发表于 2017-2-23 14:58

提供两个方法校验1/3倍频程序所需要的数据:

1. 粉红噪音 (1秒长44100KHZ采样率)

粉红噪音文本文件:
粉红噪音WAV格式文件:

在1/3倍频分析下,频谱应该是水平的。



2. 1/3倍频带的各中心频率的多音合成 (1秒长44100KHZ采样率)
1/3倍频带的各中心频率的多音合成文本文件:
1/3倍频带的各中心频率的多音合成WAV格式文件:

在1/3倍频分析下,频谱应该是水平的。



这个多音文件由32个位于各1/3倍频带中心的单音频率合成。
1:1,16Hz,1,0D
2:1,20Hz,1,0D
3:1,25Hz,1,0D
4:1,32Hz,1,0D
5:1,40Hz,1,0D
6:1,50Hz,1,0D
7:1,63Hz,1,0D
8:1,80Hz,1,0D
9:1,100Hz,1,0D
10:1,125Hz,1,0D
11:1,160Hz,1,0D
12:1,200Hz,1,0D
13:1,250Hz,1,0D
14:1,315Hz,1,0D
15:1,400Hz,1,0D
16:1,500Hz,1,0D
17:1,630Hz,1,0D
18:1,800Hz,1,0D
19:1,1000Hz,1,0D
20:1,1250Hz,1,0D
21:1,1600Hz,1,0D
22:1,2000Hz,1,0D
23:1,2500Hz,1,0D
24:1,3150Hz,1,0D
25:1,4000Hz,1,0D
26:1,5000Hz,1,0D
27:1,6300Hz,1,0D
28:1,8000Hz,1,0D
29:1,10000Hz,1,0D
30:1,12500Hz,1,0D
31:1,16000Hz,1,0D
32:1,20000Hz,1,0D


ghostzer 发表于 2017-11-1 15:34

TestGuru 发表于 2017-2-23 14:58
提供两个方法校验1/3倍频程序所需要的数据:

1. 粉红噪音 (1秒长44100KHZ采样率)


大神你好,我在用pink noise的时间域数据作三分之一倍频程谱分析后,图形如下请问对吗,
2.5H处那个峰很高,这个pink noise数据中包含2.5Hz的频率吗?


TestGuru 发表于 2018-1-9 19:26

ghostzer 发表于 2017-11-1 15:34
大神你好,我在用pink noise的时间域数据作三分之一倍频程谱分析后,图形如下请问对吗,
2.5H处那个峰很 ...

Pink Noise中含有2.5Hz,不过这个频率比较低,有几点需要注意:

1. 任何噪音都需要平均后才更容易体现出其统计意义上的特点,比如Pink Noise在倍频分析下频谱是平坦的,平均次数越多越平坦。
2. 频率越低,需要的时间序列越长,加上前述的需要做平均,因此所需要的时间序列更长。
3. 倍频分析通常都是先做FFT,然后按倍频带的上下频率范围计算倍频带内的幅度/功率/能量的。如果要分辨以2.5Hz为中心的1/3倍频带, 要求FFT的分辨率应比此倍频带的宽度2.81-2.23=0.58Hz更细得多。

下面是一个500秒长的Pink Noise, 采用1/3倍频分析,FFT频率分辨率为0.18Hz, 平均次数183,有图可见,频谱低到1.3Hz也基本是平坦的。





TestGuru 发表于 2018-1-9 19:29

对了,地板楼提供的Pink Noise 文件可能不够长,若要分析比较低的频率,需要直接去下载Multi-Instrument软件,利用其中的信号发生器来产生更长的WAV文件。
页: [1]
查看完整版本: 三分之一倍频程程序