pythonaudiohistogramrmsamplitude

Python histogram of RMS amplitude of audio file


I would like to make a histogram that bins RMS amplitude for an audio file. The goal is to show, for the full duration of the file, how much of the signal is at each amplitude level. Is there a Python package with a function for this? If not, how can it be coded?

I would also like to set the frequency range over which the analysis will be computed, e.g. between 1 and 6 kHz.

I have the following as a crude start, though I don't yet understand what it represents, and it certainly does not use RMS:

import numpy as np
import matplotlib.pyplot as plt
   
Fs, data = wavfile.read('file') 
print('data =',data)
print('number of samples in data =',len(data))

subset = data[0:44100] 
subset = abs(subset)
print('number of samples in subset =',len(subset))

plt.hist(subset, bins='auto')  
plt.show()

Solution

  • As far as I know, there is no special function in numpy for RMS, but you can do it like this

    RMS = np.sqrt(np.mean(x**2))
    

    And the question is, for which data (for which x) do you want to calculate the RMS. For example, you can apply RMS for each sample, than assuming you only have one channel in your wav-file:

    length = data.shape[0] / Fs
    print(f"length = {length}s")
    
    RMS = lambda x: np.sqrt(np.mean(x**2))
    
    sample = np.arange(int(length))
    RMS_of_sample = np.zeros(sample.shape)
    for ns in sample:
        # here you can apply the frequency window for the sample 
        RMS_of_sample[ns] = RMS(data[ns*Fs:(ns+1)*Fs])
    
    plt.hist(RMS_of_sample, label="Left channel")
    plt.show()
    

    here you can also apply some signal windows. This code gives you something like this

    enter image description here

    for incoming signal:

    enter image description here


    ADDITION to the question in the comment regarding full/partial frequency range

    If you want to analyze a complete signal in a certain frequency domain, you can apply, for example, simple filter (rectangular frequency window) for the frequency range [filter_freq_min, filter_freq_max] like this:

    from scipy.fft import fft, ifft, fftfreq
    
    filter_freq_min = 1000 # Hz
    filter_freq_max = 2000 # Hz
    
    freq = fftfreq(len(data), 1 / Fs)
    data_fft = fft(data)
    
    condition = np.logical_or(abs(freq) <= filter_freq_min, abs(freq) >= filter_freq_max)
    filtered_data_fft = np.copy(data_fft)
    filtered_data_fft [condition] = 0
    filtered_data = np.real(ifft(filtered_data_fft ))
    
    # show fft for incoming signal (blue) and filtered signal (orange) 
    plt.plot(freq, np.abs(data_fft),'.')
    plt.plot(freq, np.abs(filtered_data_fft ),'.')
    plt.xlim( [10, Fs/2] )
    plt.xlabel( 'Frequency (Hz)' )
    plt.show()
    
    # check RMS for filtered and unfiltered signal
    print(RMS(filtered_data),RMS(data))
    

    enter image description here

    In this way, you can cycle through the required frequency ranges.

    To play sound directly in Python, you can use

    import sounddevice as sd # For playing/recording audio
    sd.play(data, Fs)
    sd.play(filtered_data, Fs)