基于谐波次数生成直方图的GNU Octave技术问询

我明白你现在的需求：把FFT得到的频率-幅值数据，转换成以谐波次数为横轴的直方图，就像你给出的示例图那样。咱们一步步来调整你的代码，实现这个目标：

步骤1：提取正频率部分的FFT数据 #

你的代码里用了fftshift，这会把频谱分成负频率、0频率、正频率三部分。我们只需要正频率部分，而且要对应基波（60Hz）和它的倍数谐波（120Hz、180Hz...直到600Hz，对应1到10次谐波）。

步骤2：匹配谐波频率对应的幅值 #

对于每个谐波次数n（1到num_harmonics），对应的频率是n*freq，我们需要在FFT的频率数组里找到最接近这个值的点，提取它的幅值。

修改后的完整代码 #

我把你的代码做了调整，加入了谐波直方图的绘制部分：

clear(); 
vrms = 120; 
freq = 60; 
nCycles = 2; 
level = 25; 
vpeak = sqrt(2) * vrms; 
sampleinterval = 0.00001; 
num_harmonics = 10; 

disp("Start"); 

% 生成波形
x = 0 : sampleinterval : nCycles * 1 / freq; % 时间序列
dimmed_wave = []; 
undimmed_wave = []; 
for i = 1 : columns(x) 
    rad_value = x(i) * 2 * pi * freq; 
    off_time = mod(rad_value, pi); 
    on_time = pi*(100-level)/100; 
    if (off_time < on_time) 
        dimmed_wave = [dimmed_wave, 0]; % 调光阶段输出0
    else 
        dimmed_wave = [dimmed_wave, sin(rad_value)]; % 非调光阶段输出正弦波
    endif 
    undimmed_wave = [undimmed_wave, sin(rad_value)]; 
endfor 

y = dimmed_wave * vpeak; % 实际电压值
undimmed = undimmed_wave * vpeak; 

% 绘制原始波形
subplot(3,1,1) 
plot(x*1000, y, '-', x*1000, undimmed, '--'); 
xlabel ("Time (ms)"); 
ylabel ("Voltage"); 
legend("Dimmed Wave", "Undimmed Wave");

% 傅里叶变换
N = length(dimmed_wave); 
fft_vals = abs(fftshift(fft(dimmed_wave))); 
frequency = [ -(ceil((N-1)/2):-1:1) ,0 ,(1:floor((N-1)/2)) ] * 1 / (N * sampleinterval);

% 绘制全频谱图
subplot(3,1,2) 
plot(frequency, fft_vals); 
axis([0, num_harmonics*freq]); 
xlabel ("Frequency (Hz)"); 
ylabel ("Amplitude"); 
title("Full Frequency Spectrum");

% 提取谐波幅值并绘制直方图
harmonic_nums = 1:num_harmonics;
harmonic_amps = zeros(1, num_harmonics);

for n = harmonic_nums
    target_freq = n * freq;
    % 找到最接近目标频率的索引
    [~, idx] = min(abs(frequency - target_freq));
    harmonic_amps(n) = fft_vals(idx);
endfor

% 绘制谐波直方图
subplot(3,1,3)
bar(harmonic_nums, harmonic_amps);
xlabel("Harmonic Number");
ylabel("Amplitude");
title("Harmonic Amplitude Histogram");
set(gca, 'XTick', harmonic_nums); % 确保横轴显示所有谐波次数

关键细节解释 #

• 我们用min(abs(frequency - target_freq))来精准定位每个谐波频率对应的FFT数据点，避免因为采样间隔导致的频率偏移问题。
• 用bar函数代替plot来绘制直方图，这样更符合你想要的柱状效果。
• 调整了子布局为3行1列，同时展示原始波形、全频谱和谐波直方图，方便对比。

这样运行代码后，你就能得到和示例图类似的谐波次数-幅值直方图啦！

内容的提问来源于stack exchange，提问作者docsteer