【问题标题】:What is the source of discrepancy in 2D interpolated spectrogram with matplotlib?matplotlib 的 2D 插值频谱图差异的根源是什么?
【发布时间】:2015-12-10 11:40:57
【问题描述】:

我正在尝试使用 scipy 的 inetrp2d 函数对从 matplotlib 获得的频谱图进行插值,但不知何故未能获得相同的频谱图。数据可用here

实际的频谱图是:

插值频谱图为:

代码看起来不错,但即使这样也有问题。使用的代码是:

from __future__ import division
from matplotlib import ticker as mtick
from matplotlib.backends.backend_pdf import PdfPages
import matplotlib.pyplot as plt
import numpy as np
from bisect import bisect
from scipy import interpolate
from matplotlib.ticker import MaxNLocator
data = np.genfromtxt('spectrogram.dat', skiprows = 2, delimiter = ',')
pressure = data[:, 1] * 0.065
time = data[:, 0]
cax = plt.specgram(pressure * 100000, NFFT = 256, Fs = 50000, noverlap=4, cmap=plt.cm.gist_heat, zorder = 1)

f = interpolate.interp2d(cax[2], cax[1], cax[0], kind='cubic')
xnew = np.linspace(cax[2][0], cax[2][-1], 100)
ynew = np.linspace(cax[1][0], cax[1][-1], 100)
znew = 10 * np.log10(f(xnew, ynew))

fig = plt.figure(figsize=(6, 3.2))
ax = fig.add_subplot(111)
ax.set_title('colorMap')
plt.pcolormesh(xnew, ynew, znew, cmap=plt.cm.gist_heat)
# plt.colorbar()
plt.title('Interpolated spectrogram')
plt.colorbar(orientation='vertical')
plt.savefig('interp_spectrogram.pdf')

如何使用 Python 正确插入频谱图?

【问题讨论】:

    标签: python matplotlib scipy interpolation


    【解决方案1】:

    解决方案的关键在于此警告,您可能见过也可能没见过:

    RuntimeWarning: invalid value encountered in log10
        znew = 10 * np.log10(f(xnew, ynew))
    

    如果您的数据实际上是一种幂,您希望将其日志明确地视为分贝幂,请先获取对数,然后再拟合到样条曲线:

    spectrum, freqs, t, im = cax
    dB = 10*np.log10(spectrum)
    #f = interpolate.interp2d(t, freqs, dB, kind='cubic') # docs for this recommend next line
    f = interpolate.RectBivariateSpline(t, freqs,  dB.T) # but this uses xy not ij, hence the .T
    
    xnew = np.linspace(t[0], t[-1], 10*len(t))
    ynew = np.linspace(freqs[0], freqs[-1], 10*len(freqs)) # was it wider spaced than freqs on purpose?
    znew = f(xnew, ynew).T
    

    然后按你的方式绘制:


    上一个答案:

    如果您只想plot on logscale,请使用matplotlib.colors.LogNorm

    znew = f(xnew, ynew) # Don't take the log here
    
    plt.figure(figsize=(6, 3.2))
    plt.pcolormesh(xnew, ynew, znew, cmap=plt.cm.gist_heat, norm=colors.LogNorm())
    

    看起来像这样:

    当然,在对数刻度上绘制时,其值为负数时仍然存在差距。当值为负时,您的数据对您意味着什么应该决定您如何填写。一个简单的解决方案是将这些值设置为最小的正值,然后将它们填充为黑色:

    【讨论】:

    • 如果我错了请纠正我 cax[0] 是权力吗?必须通过 np.log10(cax[0]) 转换为 dB。
    • 好的,在这种情况下,我建议在执行样条之前记录输出日志。我会更新我的答案。
    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 2015-02-15
    • 2021-01-16
    • 1970-01-01
    • 1970-01-01
    • 2022-06-10
    • 1970-01-01
    • 2020-10-27
    相关资源
    最近更新 更多