7

My question is the following:

I have all the values that I need for a spectrogram (scipy.fftpack.fft). I would like to create a 3D spectrogram in python.

In MATLAB this is a very simple task, while in python it seems much more complicated. I tried mayavi, 3D plotting matplotlib but I have not managed to do this.

Thanks

My code:

import numpy as np
import pandas as pd
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
from matplotlib.collections import PolyCollection

fs = 11240.
t = 10
time = np.arange(fs*t) / fs
frequency = 1000.
mysignal = np.sin(2.0 * np.pi * frequency * time)

nperseg = 2**14
noverlap = 2**13
f, t, Sxx = signal.spectrogram(mysignal, fs, nperseg=nperseg,noverlap=noverlap)

myfilter = (f>800) & (f<1200)

fig,ax = plt.subplots()

plt.pcolormesh(t, f[myfilter], 10*np.log10(Sxx[myfilter, :]), cmap='jet')
plt.show()

fig = plt.figure()
ax = fig.gca(projection='3d')
x = []
y = []

for counter,i in enumerate(f):
    x.append(np.array([i for k in t]))
    y.append(t)

ax.plot_surface(np.array(x), np.array(y), 10.0*np.log10(Sxx), cmap=cm.coolwarm)
plt.show()

Similar unanswered question: How to convert a spectrogram to 3d plot. Python

Desired plot in python like Matlab's figure (last plot here: https://www.mathworks.com/help/signal/ref/spectrogram.html)

enter image description here

seralouk
  • 30,938
  • 9
  • 118
  • 133
  • 3
    please provide a [Minimal, Reproducible Example](https://stackoverflow.com/help/minimal-reproducible-example) of what you have done so far – DrBwts Jun 27 '19 at 11:30
  • Do you have the corresponding MATLAB code? – Dan Jun 27 '19 at 12:01
  • Nice question by the OP. The matlab function is just a build-in function (https://www.mathworks.com/help/signal/ref/spectrogram.html). I want to do the same in python like the OP asked. – seralouk Jun 27 '19 at 12:02
  • 1
    There is also this: https://stackoverflow.com/questions/48598829/how-to-convert-a-spectrogram-to-3d-plot-python but it is not really clear –  Jun 27 '19 at 12:05
  • What are you trying to make `x` and `y` look like? In particular, what are you trying to do with `[i for k in t]`? Have you had a look at `np.meshdrid` in order to generate `x` and `y`? – Dan Jun 27 '19 at 12:05
  • Have you tried [`specgram`](https://matplotlib.org/3.1.0/gallery/images_contours_and_fields/specgram_demo.html#sphx-glr-gallery-images-contours-and-fields-specgram-demo-py)? – Dan Jun 27 '19 at 12:07

2 Answers2

7

You just need to get your arrays in the right shape:

fs = 11240.
t = 10
time = np.arange(fs*t) / fs
frequency = 1000.
mysignal = np.sin(2.0 * np.pi * frequency * time)

nperseg = 2**14
noverlap = 2**13
f, t, Sxx = signal.spectrogram(mysignal, fs, nperseg=nperseg,noverlap=noverlap)

myfilter = (f>800) & (f<1200)

f = f[myfilter]
Sxx = Sxx[myfilter, ...]

fig = plt.figure()
ax = fig.gca(projection='3d')

ax.plot_surface(f[:, None], t[None, :], 10.0*np.log10(Sxx), cmap=cm.coolwarm)
plt.show()
Nils Werner
  • 34,832
  • 7
  • 76
  • 98
  • What if t is a list of timestamps? I get an error: TypeError: list indices must be integers, not tuple. Example list: t = [Timestamp('2016-09-14 01:43:13'), Timestamp('2016-09-14 01:43:14')...] –  Jun 27 '19 at 13:32
4

Here is an adapted version of @Nils Werner's answer with different variable names and a complete set of imports.

import numpy as np
import matplotlib.pyplot as plt 
from scipy import signal # spectrogram function
from matplotlib import cm # colour map

# basic config
sample_rate = 11240.  # 
sig_len_secs = 10
frequency = 2000.

# generate the signal
timestamps_secs = np.arange(sample_rate*sig_len_secs) / sample_rate
mysignal = np.sin(2.0 * np.pi * frequency * timestamps_secs) 

# extract the spectrum
freq_bins, timestamps, spec = signal.spectrogram(mysignal, sample_rate)

# 3d plot
fig = plt.figure()
ax = plt.axes(projection='3d')
ax.plot_surface(freq_bins[:, None], timestamps[None, :], 10.0*np.log10(spec), cmap=cm.coolwarm)
plt.show()

3D plot of a spectrum

yeeking
  • 958
  • 8
  • 11