I am trying to calculate total harmonic distortion plus noise in python. When I input the same waveform to my code and another closed source total harmonic distortion plus noise calculator, I get different answers. What’s wrong with my code?
# Imports
from __future__ import division
import soundfile as sf
import numpy as np
# Load the audio files
wave_file = sf.SoundFile("Recording.wav")
signal = wave_file.read()
sample_rate = wave_file.samplerate
# Compute FFT of signal
fft = np.fft.fft(signal, n=sample_rate)
# Compute magnitude of FFT
mag = np.abs(fft)
print(len(mag))
# Compute fundamental frequency
fundamental = mag[1:].argmax() + 1
print(fundamental)
# Compute harmonic frequencies
harmonics = [fundamental * (i+2) for i in range(int((sample_rate/fundamental)-1))]
print(harmonics)
# Compute THD and THD+N
thd = np.sqrt(np.sum(mag[harmonics] ** 2)) / mag[fundamental]
print(thd)
thdn = np.sqrt(np.sum(mag[harmonics] ** 2)) / np.sqrt(np.sum(mag ** 2))
print(thdn)