<p>In order to demonstrate the relation between noise spectral density, noise bandwidth and rms estimation I offer a slightly modified example code</p>
<pre><code>library(tuneR)
library(soundgen)
library(seewave)
#
noise_power=1e-3                  # spectral noise power V^2/Hz
samp_frq=24000                    # sampling frequency
noise_BW = samp_frq/2             # noise bandwidth
noise_var=noise_power*noise_BW    # noise variance V^2/Hz * Hz 
 
#
# noise time series
noise = rnorm(5*samp_frq,sd=sqrt(noise_var))
#
# convert to wav object
w = Wave(left = noise, samp.rate = samp_frq, bit=32)

# custom spectral averaging
f &lt;- function(x){  return (mean(abs(x)^2)) }
#
# processing
window_length_in_ms=100
window_length_in_samples=(w@samp.rate * window_length_in_ms) / 1000
#
rms = getRMS(w, windowLength = window_length_in_ms, 
             normalize = F,  overlap = 50)
power_spec = meanspec(w, wn=&quot;rectangle&quot;, wl = window_length_in_samples, 
                      FUN=f, 
                      norm=F, ovlp = 50, plot = F)
#
# print results
print(10*log10(noise_power))
print(20*log10(sd(noise))               -10*log10(noise_BW))
print(20*log10(rms<span class="math-container">$summary$</span>ampl_mean)   -10*log10(noise_BW))
print(10*log10(sum(power_spec[,&quot;y&quot;])/2) -10*log10(noise_BW))
</code></pre>
<p>the results of which are</p>
<pre><code>[1] -30
[1] -30.01798
[1] -30.01745
[1] -30.01693
</code></pre>
<p>The users are invited to play with sampling frequency and window length to note that the estimated (rms-based) noise power levels remain constant at -30dB or 1e-3</p>

CC BY-SA 4.0
2023-12-12T16:53:37.457Z
2023-12-14T06:21:39.423Z
2023-12-14T06:21:39.423Z
WMXZ
<p>The differences are due to different scalings.</p>
<p>the following modifications give all the same answer (or similar due to windowing)</p>
<pre><code># now everything relative to amplitude unit 1
f &lt;- function(x)
{
  return(mean(abs(x)^2))
}

std = sqrt(mean(w@left^2))

rms = getRMS(w, windowLength = window_length_in_ms, normalize = F,  overlap = 50)
power_spec = meanspec(w, wn=&quot;rectangle&quot;, wl = window_length_in_samples, norm=F, FUN=f, ovlp = 50, plot = T)

print(20*log10(1000))
print(20*log10(std))
print(20*log10(rms<span class="math-container">$summary$</span>ampl_mean))
print(10*log10(sum(power_spec[,&quot;y&quot;])/2))
</code></pre>
<p>generating as output</p>
<pre><code>[1] 60
[1] 59.98897
[1] 59.97903
[1] 59.99018
</code></pre>
<p>The obtained level is 60 dB (or similar due to random samples) relative to 1 amplitude unit as it should be.</p>
<p>The are multiple differences:</p>
<ol>
<li>the meanSpec routine seems not to correct for the windowing shape (using different windows gives different results) So, rectangle window should be used to compare</li>
<li>the SeeWave meanSpec function uses as default average the mean of the modulus of the spectrum. This is physically not correct. The power (or modulus square) should be used for averaging. Therefore I added a custom function.</li>
<li>To convert the power spectrum to dB, 10*log10(.) must be used</li>
<li>to convert power spectrum to broadband rms values, the sum over all
frequencies should be used and not the mean. This can easily be verified by changing the window length.The factor 2 is still unclear and most likely due to the stdft function.</li>
</ol>
<p>I would strongly recommend to use scaling to max value (in data or spectrum) ONLY for visualization purpose, and otherwise scale to 1 (or acoustic reference 1uPa or 20uPa if sensitivity is known).</p>

CC BY-SA 4.0
2023-12-12T12:03:02.467Z
2023-12-12T12:03:02.467Z
WMXZ
<p>I want to calculate the sound pressure level (in dB) of a wave and then display the mean power spectrum.  The issue I face is that the SPL value I get from getRMS() is very different from the mean amplitude values from meanspec().</p>
<pre><code>library(tuneR)
library(soundgen)
library(seewave)

#make random noise
noise = rnorm(24000*5, sd = 1000)
w = Wave(left = noise, samp.rate = 24000)

#set window lengths (ms for getRMS, samples for meanspec)
window_length_in_ms = 50
window_length_in_samples = w@samp.rate * (window_length_in_ms / 1000)

#calculate SPL
rms = getRMS(w, windowLength = window_length_in_ms, overlap = 50)
SPLdB = ((20 * log10(rms<span class="math-container">$summary$</span>ampl_mean)))
SPLdB  #this is the same RMS value I get from Audacity

#graph power spectrum
power_spec = meanspec(w, wl = window_length_in_samples, dB = &quot;max0&quot;, ovlp = 50, plot = T)
mean(power_spec[,&quot;y&quot;]). #this value is VERY different from my SPL
</code></pre>
<p>Obviously, I am missing something fundamental.  I appreciate any feedback.</p>

<r><rms><seewave><soundgen>
Can't reconcile amplitude from getRMS (soundgen) with output of meanspec (seewave) in R
# Can't reconcile amplitude from getRMS (soundgen) with output of meanspec (seewave) in R

(Asked by: username_0 on 2023-12-11)

I want to calculate the sound pressure level (in dB) of a wave and then display the mean power spectrum. The issue I face is that the SPL value I get from getRMS() is very different from the mean amplitude values from meanspec().

```
library(tuneR)
library(soundgen)
library(seewave)

#make random noise
noise = rnorm(24000*5, sd = 1000)
w = Wave(left = noise, samp.rate = 24000)

#set window lengths (ms for getRMS, samples for meanspec)
window_length_in_ms = 50
window_length_in_samples = w@samp.rate * (window_length_in_ms / 1000)

#calculate SPL
rms = getRMS(w, windowLength = window_length_in_ms, overlap = 50)
SPLdB = ((20 * log10(rms$summary$ampl_mean)))
SPLdB  #this is the same RMS value I get from Audacity

#graph power spectrum
power_spec = meanspec(w, wl = window_length_in_samples, dB = "max0", ovlp = 50, plot = T)
mean(power_spec[,"y"]). #this value is VERY different from my SPL

```

Obviously, I am missing something fundamental. I appreciate any feedback.

**Tags:** r, rms, seewave, soundgen


---

## 2 answers

### Answer by username_1 on 2023-12-12. Score: 0 | Accepted answer

In order to demonstrate the relation between noise spectral density, noise bandwidth and rms estimation I offer a slightly modified example code

```
library(tuneR)
library(soundgen)
library(seewave)
#
noise_power=1e-3                  # spectral noise power V^2/Hz
samp_frq=24000                    # sampling frequency
noise_BW = samp_frq/2             # noise bandwidth
noise_var=noise_power*noise_BW    # noise variance V^2/Hz * Hz 
 
#
# noise time series
noise = rnorm(5*samp_frq,sd=sqrt(noise_var))
#
# convert to wav object
w = Wave(left = noise, samp.rate = samp_frq, bit=32)

# custom spectral averaging
f <- function(x){  return (mean(abs(x)^2)) }
#
# processing
window_length_in_ms=100
window_length_in_samples=(w@samp.rate * window_length_in_ms) / 1000
#
rms = getRMS(w, windowLength = window_length_in_ms, 
             normalize = F,  overlap = 50)
power_spec = meanspec(w, wn="rectangle", wl = window_length_in_samples, 
                      FUN=f, 
                      norm=F, ovlp = 50, plot = F)
#
# print results
print(10*log10(noise_power))
print(20*log10(sd(noise))               -10*log10(noise_BW))
print(20*log10(rms$summary$ampl_mean)   -10*log10(noise_BW))
print(10*log10(sum(power_spec[,"y"])/2) -10*log10(noise_BW))

```

the results of which are

```
[1] -30
[1] -30.01798
[1] -30.01745
[1] -30.01693

```

The users are invited to play with sampling frequency and window length to note that the estimated (rms-based) noise power levels remain constant at -30dB or 1e-3


---

### Answer by username_1 on 2023-12-12. Score: 0

The differences are due to different scalings.

the following modifications give all the same answer (or similar due to windowing)

```
# now everything relative to amplitude unit 1
f <- function(x)
{
  return(mean(abs(x)^2))
}

std = sqrt(mean(w@left^2))

rms = getRMS(w, windowLength = window_length_in_ms, normalize = F,  overlap = 50)
power_spec = meanspec(w, wn="rectangle", wl = window_length_in_samples, norm=F, FUN=f, ovlp = 50, plot = T)

print(20*log10(1000))
print(20*log10(std))
print(20*log10(rms$summary$ampl_mean))
print(10*log10(sum(power_spec[,"y"])/2))

```

generating as output

```
[1] 60
[1] 59.98897
[1] 59.97903
[1] 59.99018

```

The obtained level is 60 dB (or similar due to random samples) relative to 1 amplitude unit as it should be.

The are multiple differences:

1. the meanSpec routine seems not to correct for the windowing shape (using different windows gives different results) So, rectangle window should be used to compare
2. the SeeWave meanSpec function uses as default average the mean of the modulus of the spectrum. This is physically not correct. The power (or modulus square) should be used for averaging. Therefore I added a custom function.
3. To convert the power spectrum to dB, 10\*log10(.) must be used
4. to convert power spectrum to broadband rms values, the sum over all
   frequencies should be used and not the mean. This can easily be verified by changing the window length.The factor 2 is still unclear and most likely due to the stdft function.

I would strongly recommend to use scaling to max value (in data or spectrum) ONLY for visualization purpose, and otherwise scale to 1 (or acoustic reference 1uPa or 20uPa if sensitivity is known).


---

