Forums

Help with THD+N calculations

Started by rcoup May 31, 2006
Hi All,

My DSP knowledge is far too rusty and I've been hitting my head against
the wall for a little while now. Archives of the group have got
tantalisingly close but I'm still confused. The fun situation is
calculating THD+N of an audio signal.

What I'm doing:
a. 1KHz sine wave is my fundamental.
b. Notch Filter to remove the 1KHz
c. FFT (8Kpoint)
d. Get magnitudes (mag[] = sqrt(re[]^2 + im[]^2))
e. Calculate the total energy over the range 20Hz to 20KHz
     Ex = sum(sqr(mag[k])) with k in 3-3143 (the FFT bins i'm
interested in)
f. thdn = sqrt(Ex) / (3143-3+1)

Now I need to compare that with the amplitude of the fundamental to get
my %age... whats an easy way to do that? Have I messed up anything so
far (i'm thinking i should be removing DC...)?

Thanks in advance,

Robert

"rcoup" <robert.coup@gmail.com> wrote in message 
news:1149047348.326166.195250@h76g2000cwa.googlegroups.com...
> Hi All, > > My DSP knowledge is far too rusty and I've been hitting my head against > the wall for a little while now. Archives of the group have got > tantalisingly close but I'm still confused. The fun situation is > calculating THD+N of an audio signal. > > What I'm doing: > a. 1KHz sine wave is my fundamental. > b. Notch Filter to remove the 1KHz > c. FFT (8Kpoint) > d. Get magnitudes (mag[] = sqrt(re[]^2 + im[]^2)) > e. Calculate the total energy over the range 20Hz to 20KHz > Ex = sum(sqr(mag[k])) with k in 3-3143 (the FFT bins i'm > interested in) > f. thdn = sqrt(Ex) / (3143-3+1) > > Now I need to compare that with the amplitude of the fundamental to get > my %age... whats an easy way to do that? Have I messed up anything so > far (i'm thinking i should be removing DC...)?
You didn't say what the sample rate was / is. With 8k points and 20kHz bandwidth (?) and some reasonable margin in sampling rate, the sample rate might be 50kHz and 4K samples would be around 25kHz??? and 3143 would be (3143/4000)*25kHz = 19643.75Hz so maybe the sample rate is 3143/4000 = 20,000/fs/2 so fs= 20,000 *(4000/3143)*2 = 50906.78Hz?? The frequency resolution then is 50906.78/8000 = 6.37Hz and 1kHz falls quite closely to being an integer multiple of 6.37 - which is good. I'm not sure of the best procedure for separating THD and N if that's what you want or how accurately you need to separate them if at all. Anyway, the THD is the ratio of rms values. I don't see the mean being calculated above - or at least not where I'd expect it to be. I also don't see the magnitude of the signal at 1kHz being used as a reference. An estimate of the THD would be the sqrt of the mean of the sum of the squared magnitudes at 2kHz, 3kHz, 4kHz, etc. divided by the magnitude at 1kHz. An estimate of the noise would be the sqrt of the mean of the sum of the squared magnitudes at all other frequencies - so, divided by the magnitude at 1kHz would be the SNR. I think that's correct.... Fred
Hi Fred,

Thanks for helping!

Fred Marshall wrote:
> You didn't say what the sample rate was / is.
48K
> With 8k points and 20kHz bandwidth (?) and some reasonable margin in > sampling rate, the sample rate might be 50kHz and 4K samples would be around > 25kHz??? and 3143 would be (3143/4000)*25kHz = 19643.75Hz so maybe the > sample rate is 3143/4000 = 20,000/fs/2 so fs= 20,000 *(4000/3143)*2 = > 50906.78Hz?? > > The frequency resolution then is 50906.78/8000 = 6.37Hz and 1kHz falls quite > closely to being an integer multiple of 6.37 - which is good.
Bin size at 48K/8192 is 48000/8192 = 5.86Hz. The closest bins are 1002.0Hz and 996.1. So should I be using 1002Hz as my fundamental to get it closely aligned (seems sensible)?
> I'm not sure of the best procedure for separating THD and N if that's what > you want or how accurately you need to separate them if at all. Anyway, the > THD is the ratio of rms values.
I don't need to separate them.
> I don't see the mean being calculated above - or at least not where I'd > expect it to be. > I also don't see the magnitude of the signal at 1kHz being used as a > reference. > > An estimate of the THD would be the sqrt of the mean of the sum of the > squared magnitudes at 2kHz, 3kHz, 4kHz, etc. divided by the magnitude at > 1kHz. > > An estimate of the noise would be the sqrt of the mean of the sum of the > squared magnitudes at all other frequencies - so, divided by the magnitude > at 1kHz would be the SNR.
Okay, so for THD+N I need to get the "rms summation of all signal components (excluding the fundamental) over some prescribed bandwidth" [Rane]. Which is the harmonics plus all the other frequencies. So in practise I should: Apply a notch filter to notch out 1002Hz and therabouts. Do the FFT. Sum the squared magnitude of all FFT bins in 20Hz-20KHz. Divide it by the number of bins 20Hz-20KHz to get the mean Take the square root of the mean. Taking that over the magnitude of the 1002Hz bin, and I get 1.11%. This is about 100 times too large (another piece of software is returning numbers around 0.07-0.08%) What else have I missed? Thanks again, Rob :)
"rcoup" <robert.coup@gmail.com> wrote in message 
news:1149136959.555994.299960@c74g2000cwc.googlegroups.com...
> Hi Fred, > > Thanks for helping! > > Fred Marshall wrote: >> You didn't say what the sample rate was / is. > > 48K > >> With 8k points and 20kHz bandwidth (?) and some reasonable margin in >> sampling rate, the sample rate might be 50kHz and 4K samples would be >> around >> 25kHz??? and 3143 would be (3143/4000)*25kHz = 19643.75Hz so maybe the >> sample rate is 3143/4000 = 20,000/fs/2 so fs= 20,000 *(4000/3143)*2 = >> 50906.78Hz?? >> >> The frequency resolution then is 50906.78/8000 = 6.37Hz and 1kHz falls >> quite >> closely to being an integer multiple of 6.37 - which is good. > > Bin size at 48K/8192 is 48000/8192 = 5.86Hz. The closest bins are > 1002.0Hz and 996.1. So should I be using 1002Hz as my fundamental to > get it closely aligned (seems sensible)? > >> I'm not sure of the best procedure for separating THD and N if that's >> what >> you want or how accurately you need to separate them if at all. Anyway, >> the >> THD is the ratio of rms values. > > I don't need to separate them. > >> I don't see the mean being calculated above - or at least not where I'd >> expect it to be. >> I also don't see the magnitude of the signal at 1kHz being used as a >> reference. >> >> An estimate of the THD would be the sqrt of the mean of the sum of the >> squared magnitudes at 2kHz, 3kHz, 4kHz, etc. divided by the magnitude at >> 1kHz. >> >> An estimate of the noise would be the sqrt of the mean of the sum of the >> squared magnitudes at all other frequencies - so, divided by the >> magnitude >> at 1kHz would be the SNR. > > Okay, so for THD+N I need to get the "rms summation of all signal > components (excluding the fundamental) over some prescribed bandwidth" > [Rane]. Which is the harmonics plus all the other frequencies. So in > practise I should: > > Apply a notch filter to notch out 1002Hz and therabouts. > Do the FFT. > Sum the squared magnitude of all FFT bins in 20Hz-20KHz. > Divide it by the number of bins 20Hz-20KHz to get the mean > Take the square root of the mean. > > Taking that over the magnitude of the 1002Hz bin, and I get 1.11%. This > is about 100 times too large (another piece of software is returning > numbers around 0.07-0.08%) > > What else have I missed?
Maybe exactly how you are getting the rms value at 1002Hz....? Note that the fft will scale the magnitudes so you can't just assume a magnitude from the time domain amplitude of the signal. One method might be to use the sample at 1002Hz for the reference and forget about a notch filter, just don't include that one sample or the rms value of that one sample and the two adjacent ones in the sum across frequencies. Fred
Fred Marshall wrote:

> One method might be to use the sample at 1002Hz for the reference and > forget about a notch filter, just don't include that one sample or the rms > value of that one sample and the two adjacent ones in the sum across > frequencies.
Or even better use the integral of the peak at 1002 Hz, meaning the sum of the bins at and close to 1002 Hz. Include all those bins to the integral which are significantly higher than the noise in the surrounding. The remaining bins are use for the THD+N summation. bye Andreas -- Andreas H&#2013266172;nnebeck | email: acmh@gmx.de ----- privat ---- | www : http://www.huennebeck-online.de Fax/Anrufbeantworter: 0721/151-284301 GPG-Key: http://www.huennebeck-online.de/public_keys/andreas.asc PGP-Key: http://www.huennebeck-online.de/public_keys/pgp_andreas.asc
rcoup wrote:

  ...

> Okay, so for THD+N I need to get the "rms summation of all signal > components (excluding the fundamental) over some prescribed bandwidth" > [Rane]. Which is the harmonics plus all the other frequencies. So in > practise I should: > > Apply a notch filter to notch out 1002Hz and therabouts. > Do the FFT. > Sum the squared magnitude of all FFT bins in 20Hz-20KHz. > Divide it by the number of bins 20Hz-20KHz to get the mean > Take the square root of the mean. > > Taking that over the magnitude of the 1002Hz bin, and I get 1.11%. This > is about 100 times too large (another piece of software is returning > numbers around 0.07-0.08%) > > What else have I missed?
You can't measure the amplitude of the carrier after the notch filter. Jerry -- Engineering is the art of making what you want from things you can get. &#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;
Andreas Huennebeck wrote:
> Fred Marshall wrote: > > > One method might be to use the sample at 1002Hz for the reference and > > forget about a notch filter, just don't include that one sample or the rms > > value of that one sample and the two adjacent ones in the sum across > > frequencies. > > Or even better use the integral of the peak at 1002 Hz, meaning the sum > of the bins at and close to 1002 Hz. Include all those bins to the integral > which are significantly higher than the noise in the surrounding. The remaining > bins are use for the THD+N summation.
Thanks again everyone for your help. To explain, previously I was doing another FFT (no-notch) to get the freq of the fundamental. not ideal but at least it was clear. Going back to one FFT, no notch. okay, so looking at my FFT the peak width is approximately +/- 3 bins from 1002Hz. - RMS of FFT mag from bins except those (noise) = 3.48E10 - RMS of FFT mag from peak bins = 7.37E10 - THD+N = NOISErms / PKrms = 47.2% well, thats not right. If I just sum the FFT magnitudes for the peak bins (reading back thats what Andreas actually suggested) then I get 9.4%.... 100 times too big :) Expanding the peak width to 21 bins instead of 7 changes it to 46%/7%... hmm. At a closer look it was clipping, which explains the 9%. Lowering the levels and we get a nice 0.094%! Cheers! Rob :)