9

我有一个采样率为 16e3 的信号,它的频率范围是 125 到 1000 Hz。因此,如果我绘制一个频谱图,由于所有未使用的频率,我会得到一个非常小的颜色范围。

我试图通过设置斧头限制来修复它,但这不起作用。

有没有办法切断未使用的频率或用 NaN 替换它们?

将数据重新采样到 2e3 将不起作用,因为在 125 Hz 以下仍有一些未使用的频率。

谢谢你的帮助。

4

4 回答 4

17

specgram() 正在为您完成所有工作。如果你在轴.py 中查看 specgram 函数,你可以看到它是如何工作的。原来的功能在Python27\Lib\site-packages\matplotlib\axes.py我的电脑上。

    <snip>  
    Pxx, freqs, bins = mlab.specgram(x, NFFT, Fs, detrend,
         window, noverlap, pad_to, sides, scale_by_freq)

    Z = 10. * np.log10(Pxx)
    Z = np.flipud(Z)

    if xextent is None: xextent = 0, np.amax(bins)
    xmin, xmax = xextent
    freqs += Fc
    extent = xmin, xmax, freqs[0], freqs[-1]
    im = self.imshow(Z, cmap, extent=extent, **kwargs)
    self.axis('auto')

    return Pxx, freqs, bins, im

您可能必须在此基础上创建自己的函数并剪辑 Pxx 数据以满足您的需要。

    Pxx, freqs, bins = mlab.specgram(x, NFFT, Fs, detrend,
         window, noverlap, pad_to, sides, scale_by_freq)
    # ****************
    # create a new limited Pxx and freqs
    # 
    # ****************
    Z = 10. * np.log10(Pxx)
    Z = np.flipud(Z)

Pxx 是一个二维数组,形状为 (len(freqs),len(bins)

>>> Pxx.shape
(129, 311)
>>> freqs.shape
(129,)
>>> bins.shape
(311,)
>>> 

这将限制 Pxx 和频率

Pxx = Pxx[(freqs >= 125) & (freqs <= 1000)]
freqs = freqs[(freqs >= 125) & (freqs <= 1000)]

这是一个完整的解决方案 - my_specgram() - 与来自画廊的 specgram_demo 一起使用

from pylab import *
from matplotlib import *


# 100, 400 and 200 Hz sine 'wave'
dt = 0.0005
t = arange(0.0, 20.0, dt)
s1 = sin(2*pi*100*t)
s2 = 2*sin(2*pi*400*t)
s3 = 2*sin(2*pi*200*t)

# create a transient "chirp"
mask = where(logical_and(t>10, t<12), 1.0, 0.0)
s2 = s2 * mask

# add some noise into the mix
nse = 0.01*randn(len(t))

x = s1 + s2 + +s3 + nse # the signal
NFFT = 1024       # the length of the windowing segments
Fs = int(1.0/dt)  # the sampling frequency

# modified specgram()
def my_specgram(x, NFFT=256, Fs=2, Fc=0, detrend=mlab.detrend_none,
             window=mlab.window_hanning, noverlap=128,
             cmap=None, xextent=None, pad_to=None, sides='default',
             scale_by_freq=None, minfreq = None, maxfreq = None, **kwargs):
    """
    call signature::

      specgram(x, NFFT=256, Fs=2, Fc=0, detrend=mlab.detrend_none,
               window=mlab.window_hanning, noverlap=128,
               cmap=None, xextent=None, pad_to=None, sides='default',
               scale_by_freq=None, minfreq = None, maxfreq = None, **kwargs)

    Compute a spectrogram of data in *x*.  Data are split into
    *NFFT* length segments and the PSD of each section is
    computed.  The windowing function *window* is applied to each
    segment, and the amount of overlap of each segment is
    specified with *noverlap*.

    %(PSD)s

      *Fc*: integer
        The center frequency of *x* (defaults to 0), which offsets
        the y extents of the plot to reflect the frequency range used
        when a signal is acquired and then filtered and downsampled to
        baseband.

      *cmap*:
        A :class:`matplotlib.cm.Colormap` instance; if *None* use
        default determined by rc

      *xextent*:
        The image extent along the x-axis. xextent = (xmin,xmax)
        The default is (0,max(bins)), where bins is the return
        value from :func:`mlab.specgram`

      *minfreq, maxfreq*
        Limits y-axis. Both required

      *kwargs*:

        Additional kwargs are passed on to imshow which makes the
        specgram image

      Return value is (*Pxx*, *freqs*, *bins*, *im*):

      - *bins* are the time points the spectrogram is calculated over
      - *freqs* is an array of frequencies
      - *Pxx* is a len(times) x len(freqs) array of power
      - *im* is a :class:`matplotlib.image.AxesImage` instance

    Note: If *x* is real (i.e. non-complex), only the positive
    spectrum is shown.  If *x* is complex, both positive and
    negative parts of the spectrum are shown.  This can be
    overridden using the *sides* keyword argument.

    **Example:**

    .. plot:: mpl_examples/pylab_examples/specgram_demo.py

    """

    #####################################
    # modified  axes.specgram() to limit
    # the frequencies plotted
    #####################################

    # this will fail if there isn't a current axis in the global scope
    ax = gca()
    Pxx, freqs, bins = mlab.specgram(x, NFFT, Fs, detrend,
         window, noverlap, pad_to, sides, scale_by_freq)

    # modified here
    #####################################
    if minfreq is not None and maxfreq is not None:
        Pxx = Pxx[(freqs >= minfreq) & (freqs <= maxfreq)]
        freqs = freqs[(freqs >= minfreq) & (freqs <= maxfreq)]
    #####################################

    Z = 10. * np.log10(Pxx)
    Z = np.flipud(Z)

    if xextent is None: xextent = 0, np.amax(bins)
    xmin, xmax = xextent
    freqs += Fc
    extent = xmin, xmax, freqs[0], freqs[-1]
    im = ax.imshow(Z, cmap, extent=extent, **kwargs)
    ax.axis('auto')

    return Pxx, freqs, bins, im

# plot
ax1 = subplot(211)
plot(t, x)
subplot(212, sharex=ax1)

# the minfreq and maxfreq args will limit the frequencies 
Pxx, freqs, bins, im = my_specgram(x, NFFT=NFFT, Fs=Fs, noverlap=900, 
                                cmap=cm.Accent, minfreq = 180, maxfreq = 220)
show()
close()
于 2013-10-19T20:25:34.323 回答
4

如今,有一种比提出问题时更简单的方法:您可以使用它matplotlib.pyplot.axis来设置yminymax调整您想要的频率。这很容易;这是我的程序的一个片段:

plt.specgram(xmit, NFFT=65536, Fs=Fs)
plt.axis(ymin=Fc-Fa*10, ymax=Fc+Fa*10)
plt.show()
于 2017-06-15T08:58:50.663 回答
0

specgram() 返回 (Pxx, freqs, bins, im),其中 im 是 AxesImage 实例 [1]。您可以使用它来设置情节的限制:

Pxx, freqs, bins, im = plt.specgram(signal, Fs=fs)
im.set_ylim((125,1000))

[1] http://matplotlib.org/api/pyplot_api.html#matplotlib.pyplot.specgram

于 2013-10-20T11:57:51.370 回答
0

这是一个改编版本:http ://matplotlib.org/examples/pylab_examples/specgram_demo.html ,它改变了绘制的频率范围。

#!/usr/bin/env python
#### from the example
#### 
from pylab import *

dt = 0.0005
t = arange(0.0, 20.0, dt)
s1 = sin(2*pi*100*t)
s2 = 2*sin(2*pi*400*t)
mask = where(logical_and(t>10, t<12), 1.0, 0.0)
s2 = s2 * mask
nse = 0.01*randn(len(t))
x = s1 + s2 + nse # the signal

NFFT = 1024       # the length of the windowing segments
Fs = int(1.0/dt)  # the sampling frequency

ax1 = subplot(211)
plot(t, x)
subplot(212, sharex=ax1)
Pxx, freqs, bins, im = specgram(x, NFFT=NFFT, Fs=Fs, noverlap=900,
                                cmap=cm.gist_heat)

#### edited from the example
#### 
# here we get access to the axes
x1,x2,y1,y2 = axis()
# leave x range the same, change y (frequency) range
axis((x1,x2,25,500))

show()
于 2015-05-01T16:18:09.193 回答