三分之一倍频程程序

A牛牛牛牛

01.clear all

02.close all

03.%读取数据文件

04.importfile('250.txt');

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

06.%参数的初始化

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

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

09. center_frequency  = [0.02 0.025 0.0315 0.04 0.05 0.063 0.08 0.1 0.125 0.160 0.20 0.25 0.315 0.40 0.50 0.63 0.8 1.0 1.25 1.6 2.0 2.50 3.15 4.0];

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=[lower_bound;upper_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牛牛牛牛

自己顶一下

宇安乙昂

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

TestGuru

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

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

粉红噪音文本文件： PinkNoise.txt (1.32 MB, 下载次数: 14)

2017-2-23 14:24 上传

点击文件名下载附件

粉红噪音WAV格式文件： PinkNoise.rar (128.66 KB, 下载次数: 8)

2017-2-23 14:42 上传

点击文件名下载附件

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



2. 1/3倍频带的各中心频率的多音合成 (1秒长44100KHZ采样率)
1/3倍频带的各中心频率的多音合成文本文件： Octave3CentreFrequencies.txt (1.32 MB, 下载次数: 5)

2017-2-23 14:43 上传

点击文件名下载附件

1/3倍频带的各中心频率的多音合成WAV格式文件： Octave3CentreFrequencies.rar (129.44 KB, 下载次数: 6)

2017-2-23 14:46 上传

点击文件名下载附件

在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

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

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

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

TestGuru

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

对了，地板楼提供的Pink Noise 文件可能不够长，若要分析比较低的频率，需要直接去下载Multi-Instrument软件，利用其中的信号发生器来产生更长的WAV文件。