信息处理 - 使用 Scipy 检查高斯信号的 Parseval 定理 - 吾爱随笔录

使用 Scipy 检查高斯信号的 Parseval 定理

信息处理 matlab fft 傅里叶变换 scipy 解析

2022-01-12 09:59:34

我正在尝试检查高斯信号的 Parseval 定理。众所周知，的傅立叶变换是。所以我通过使用 quad 和 simps 来实现它。我认为一个是连续积分，另一个是使用样本的积分。 $\exp(-t^2)$ $\sqrt{\pi}\exp(-\pi^2 k^2)$

代码

from scipy import integrate
import numpy as np

def f(x):
    return (np.exp(-x**2))**2
def F(k):
    return ((np.pi**0.5)*np.exp((-np.pi**2)*(k**2)))**2

a=integrate.quad(f,-np.inf,np.inf)
print(a)
#>>>(1.2533141373155017, 4.4674503165883495e-09)
b=integrate.quad(F,-np.inf,np.inf)
print(b)
#>>>(1.2533141373155066, 1.592743224014847e-08)

a=b 所以上面的代码遵循 Parseval 定理。问题出在以下代码上。

import scipy.fftpack as fft

N=50000
t  = np.linspace(-1000,1000, N)
h=np.exp(-t**2)
H=2*np.abs(fft.fftshift(fft.fft(h)/N))
freq=fft.fftshift(fft.fftfreq(H.shape[0],t[1]-t[0]))

S_h=integrate.simps(h**2,t)
print(S_h)
>>>1.25331413732

S_H=integrate.simps(H**2,freq)
print(S_H)
>>>1.25326400525e-06

S_h不等于S_H。问题是什么？我该如何解决？

2个回答

Parseval 定理背后有一个简单的想法。
酉变换保留范数。当正确的因子到位时，傅里叶变换（所有它的变体）确实是一元的。 ${L}_{2}$

问题在于对元素求和并分解它们的方式。
您需要通过考虑的是该系列的绝对值平方和。 $\frac{1}{N}$

查看以下 MATLAB 代码：

numSamples = 50000;
vT = linspace(-1000, 1000, numSamples);
vH = exp(-(vT .* vT));

vK = fft(vH);

% Summation by Parseval Theorem

sum(vH .* vH)
sum(vK .* conj(vK)) / numSamples

在 MATLAB 中，它产生

ans =

   31.3322

ans =

   31.3322

这是为什么？注意fft()MATLAB 和 Python 应用的不是 Unitary 形式。
对于单一形式，它应该由分解。 $\frac{1}{\sqrt{N}}$

然而，在 Parseval 中，我们看到：

\sum_{n = 0}^{N - 1} x^{2} [n] = \frac{1}{N} \sum_{k = 0}^{N - 1} {| X [k] |}^{2}

$\sum_{n = 0}^{N - 1} {x}^{2} \left[ n \right] = \frac{1}{N} \sum_{k = 0}^{N - 1} {\left| X \left[ k \right] \right|}^{2}$

因子正好是。 $\frac{1}{N}$ ${\left( \frac{1}{\sqrt{N}} \right)}^{2}$

我不知道 SciPy。但是这两个值在我看来非常接近。当一个人在时域或频域中进行高斯运算时，请记住，一旦您进行采样，这些单个高斯函数就会在另一个域中复制，并且您有一个高斯脉冲序列。

傅里叶积分与具有有限限制并使用黎曼求和的近似值并不完全相同，那么您将得到 DFT（或fft）。

对于任何离散时间、离散频域，有限极限积分的黎曼求和可能与积分不同。

其它你可能感兴趣的问题

上一篇实时从真实信号中获取 I/Q 分量 - 希尔伯特变换还是数字下变频？下一篇具有非线性测量的线性系统的扩展卡尔曼滤波器