2

Using PyWavelets and Matplotbib.Specgram on a signal gives more detailed plots with pywt.dwt then pywt.cwt. How can I get a pywt.cwt specgram in a similar way?

With dwt:

import pywt
import pywt.data
import matplotlib.pyplot as plot
from scipy import signal
from scipy.io import wavfile

bA, bD = pywt.dwt(datamean, 'db2')
powerSpectrum, freqenciesFound, time, imageAxis = plot.specgram(bA, NFFT =  387, Fs=100)

plot.xlabel('Time')
plot.ylabel('Frequency') 
plot.show()

with this spectrogram plot:

https://i.stack.imgur.com/0YDnC.jpg

With cwt:

widths = np.arange(1,5)
coef, freqs = pywt.cwt(datamean, widths,'morl')
powerSpectrum, freqenciesFound, time, imageAxis = plot.specgram(coef, NFFT = 129, Fs=100)

plot.xlabel('Time')
plot.ylabel('Frequency') 
plot.show()

with this spectrogram plot:

https://i.stack.imgur.com/RIxhz.jpg

and for better results:

sig  = datamean
widths = np.arange(1, 31)
cwtmatr = signal.cwt(sig, signal.ricker, widths)

plt.imshow(cwtmatr, extent=[-1, 1, 1, 5], cmap='PRGn', aspect='auto',
       vmax=abs(cwtmatr).max(), vmin=-abs(cwtmatr).max())
plt.show()

with this spectrogram plot:

https://i.stack.imgur.com/NqRej.jpg

How can I get for cwt (spectrogram plot 2 and 3) a similar spectogram plot and style like in the first one? It seems like the 1st spectrogram plot compared to the 3rd has much more details.

Luc Pingo
  • 23
  • 2
  • 5

1 Answers1

3

This would be better as a comment, but since I lack the Karma to do that:

You don't want to make a spectrogram with wavelets, but a scalogram instead. What it looks like you're doing above is projecting your data in a scale subspace (that correlates to frequency), then taking those scales and finding the frequency content of them which is not what you probably want.

The detail and approximation coefficients are what you would want to use directly. Unfortunately, PyWavelets doesn't have a simple plotting function to do this for you, AFAIK. Matlab does, and their help page may be illuminating if I fail.

def scalogram(data):
    wave='db4'
    coeff=pywt.wavedec(data,wave)
    levels=len(coeff)
    lengths=[len(co) for co in coeff]
    col=np.max(lengths)
    im=np.ones([levels,col])
    col=col.astype(float)
    for level in range(levels):
        #print [lengths[level],col]
        y=coeff[level]
        if lengths[1+level]<col:
            x=col/(lengths[1+level]+1)*np.arange(1,len(y)+1)
            xi=np.linspace(0,int(col),int(col))
            yi=griddata(points=x,values=y,xi=xi,method='nearest')
        else:
            yi=y
        im[level,:]=yi
    im[im==0]=np.nan
    tiles=sum(lengths)-lengths[0]
    return im,tiles
Wxx,tiles=scalogram(data)
IM=plt.imshow(np.log10(abs(Wxx)),aspect='auto')
plt.show()

There are better ways of doing that, but it works. This produces a square matrix similar to spectrogram in "Wxx", and tiles is simply a counter of the number of time-frequency tilings to compare to the number used in a SFFT.

I've attached a picture of what these tilings look like enter image description here

NRG
  • 91
  • 5