1 year ago

#247353

test-img

JayJona

plot power spectrum of a signal with wavelet transform

Using this dataset https://philharmonia.co.uk/resources/sound-samples/ I'm trying to plot the power spectrum of note played by a specific instrument.

I'm using librosa to load the audio file and get some information with this code

import librosa
import numpy as np
import matplotlib.pyplot as plt
import pywt

y, sr = librosa.load(file_path)
duration = librosa.get_duration(y=y, sr=sr)
delta_t = duration / len(y)
t0=0
time = np.arange(0, len(y)) * delta_t + t0

I'm also following this https://ataspinar.com/2018/12/21/a-guide-for-using-the-wavelet-transform-in-machine-learning/ guide to plot the power spectrum and I'm using pywavelet library.

The problem that I have with this code is RuntimeWarning: divide by zero encountered in log2 and the plot is not shown.

def plot_wavelet(time, signal, scales, 
                 waveletname = 'cmor', 
                 cmap = plt.cm.seismic, 
                 title = 'Wavelet Transform (Power Spectrum) of signal', 
                 ylabel = 'Period (years)', 
                 xlabel = 'Time'):
    
    dt = time[1] - time[0]
    print("dt ", dt)
    [coefficients, frequencies] = pywt.cwt(signal, scales, waveletname, dt)
    power = (abs(coefficients)) ** 2
    period = 1. / frequencies
    levels = [0.0625, 0.125, 0.25, 0.5, 1, 2, 4, 8]
    contourlevels = np.log2(levels)
    
    fig, ax = plt.subplots(figsize=(15, 10))
    im = ax.contourf(time, np.log2(period), np.log2(power), contourlevels, extend='both',cmap=cmap)
    

    ax.set_title(title, fontsize=20)
    ax.set_ylabel(ylabel, fontsize=18)
    ax.set_xlabel(xlabel, fontsize=18)
    
    yticks = 2**np.arange(np.ceil(np.log2(period.min())), np.ceil(np.log2(period.max())))
    ax.set_yticks(np.log2(yticks))
    ax.set_yticklabels(yticks)
    ax.invert_yaxis()
    ylim = ax.get_ylim()
    ax.set_ylim(ylim[0], -1)
    
    cbar_ax = fig.add_axes([0.95, 0.5, 0.03, 0.25])
    fig.colorbar(im, cax=cbar_ax, orientation="vertical")
    plt.show()

scales = np.arange(1, 128)
plot_wavelet(time=time, signal=y, scales=scales, waveletname='gaus5')

Noting that some values in the power array are at -inf. How can I solve?

python

librosa

wavelet

wavelet-transform

pywavelets

0 Answers

Your Answer

Accepted video resources