在python中计算thd

Question

我正在尝试计算提供的交流电压的总谐波失真值。我正在使用 Arduino 以超过 8 KHz 的速率对电压数据进行采样，并将这些数据存储到文本文件中。然后我尝试使用以下用 python 编写的代码 sn-p 计算 thd：

    import numpy as np
    import scipy.fftpack
    from scipy.fftpack import fft
    from numpy import genfromtxt

    sampled_data = genfromtxt('/../file.txt',delimiter=',')
    abs_yf=np.abs(fft(sampled_data))

    #As far as I know, THD=sqrt(sum of square magnitude of
    #harmonics+noise)/Fundamental value (Is it correct?)So I'm
    #just summing up square of all frequency data obtained from FFT,
    #sqrt() them and dividing them with fundamental frequecy value.

    def thd(abs_data):
    sq_sum=0.0
    for r in range(len(abs_data)):
       sq_sum=sq_sum+(abs_data[r])**2
    sq_harmonics=sq_sum-(max(abs_data))**2.0
    thd=100*sq_harmonics**0.5/max(abs_data)
    return thd

    print "Total Harmonic Distortion(in percent):"
    print thd(abs_yf)

问题是，在我的情况下，获得的 Thd 值在 5% 到 25% 之间变化。 （实际上它实际上不超过 5%）。我究竟做错了什么？还有其他方法可以找出 thd 吗？

Answer 1

虽然这很安静，但对于像我这样遇到这篇文章的人来说：OP 方法存在一些问题。

1) FFT 返回的幅度包括 0 频率 bin 的幅度，因此如果信号中存在任何 DC 偏置，则假设 max(abs_data) 是对应于基频的幅度是不正确的。这是行中的问题

thd = 100*sq_harmonics**0.5 / max(abs_data)

与 0 频率相关的幅度可以忽略，作为一种快速解决方案。

2) abs_data 的后半部分应该被扔掉，它是前半部分的“镜像”反映。这是由于傅里叶变换的性质。

这两个问题都可以通过更改函数的输入来解决，即替换

打印 thd(abs_yf)

与

打印(thd(abs_yf[1:int(len(abs_yf)/2)]))

我们已将输入更改为不包含第一个或最后 N/2 个元素。

结果仍然不理想，因为正如上面提到的先前答案，窗口需要恰好是整数个周期。使用带偏移量和调整窗口的纯正弦进行测试表明该方法运行良好，但如果出现明显的窗口错误，则会严重失败。

t0=0
tf = 0.02  # integer number of cycles
dt = 1e-4
offset = 0.5
N = int((tf-t0)/dt)
time = np.linspace(0.0,tf,N )    #;

commandSigFreq = 100
Amplitude = 2

waveOfSin = Amplitude*np.sin(2.0*pi*commandSigFreq*time) + offset

abs_yf = np.abs(fft(waveOfSin))
#print("freq is" + str(scipy.fftpack.fftfreq(sampled_data, dt )  ))
#As far as I know, THD=sqrt(sum of square magnitude of
#harmonics+noise)/Fundamental value (Is it correct?)So I'm
#just summing up square of all frequency data obtained from FFT,
#sqrt() them and dividing them with fundamental frequency value.

def thd(abs_data):
    sq_sum=0.0
    for r in range( len(abs_data)):
       sq_sum = sq_sum + (abs_data[r])**2

    sq_harmonics = sq_sum -(max(abs_data))**2.0
    thd = 100*sq_harmonics**0.5 / max(abs_data)

    return thd

print("Total Harmonic Distortion(in percent):")
print(thd(abs_yf[1:int(len(abs_yf)/2) ]))

Answer 2

您很可能会通过测量过程本身添加额外的失真。

如果您将 Arduino ADC 与高级测量设备进行比较，Arduino 的值很可能会更差。至少你需要一个非常稳定且无抖动的时钟。

此外，数据的输出（我猜是通过 UART）可能会干扰 ADC 测量的时序。