标签: dft

可能重复:
如何从FFT结果中获取频率

我在某种程度上与Matlab中的fft(DFT)命令的x轴混淆.当我们对在n点采样的信号执行fft命令时,我们得到x轴为0到n-1的图.这是指以Hz为单位的频率吗？例如,当我在n = 2时出现尖峰时,是否意味着信号频率为2Hz？

这可能是一个非常天真的问题,但在这里.

我想计算函数f(x)的傅里叶变换.所以我定义了一个numpy数组X并通过向量化函数f.现在,如果我计算这个数组f(X)的FFT,它就不会像f(x)那样在一张纸上进行傅里叶变换.例如,如果我计算高斯的FFT,我应该得到高斯或数组,其实部非常接近地类似于高斯.

这是代码.请让我知道我需要改变什么来获得通常的傅立叶变换.

import matplotlib.pyplot as plt
import numpy as np

N = 128
x = np.linspace(-5, 5, N)
y = np.exp(-x**2)

y_fft = np.fft.fftshift(np.fft.fft(y).real)
plt.plot(x, y_fft)

plt.show()

让我重申一下.我想计算任何函数的傅里叶变换(例如高斯).FFT是计算数字数组的傅立叶变换的方法,但这与连续傅立叶变换公式的简单离散化不同.

我试图了解fftand ifft函数如何在python中工作。我给出了一个虚数奇函数的简单示例，以计算傅里叶逆变换，以期获得一个实数奇函数（应如此）。下面是我的代码：

v = np.array([-1,-2,0,2,1]) * 1j
t = [-2,-1,0,1,2]
V = ifft(fftshift(v))

显然，被采样的函数v是一个奇数虚函数，因此当我计算傅立叶逆变换并进行平移后，我应该得到一个实奇函数。但这种情况并非如此。我对傅立叶变换有什么误解？谢谢！

我想使用离散傅立叶变换来识别销售动态，然后聚类相似的模式。但是，我是使用 R 的新手，在搜索解决方案后，我找到了一个 prodecure fft()，但不确定我是否得到与 DFT 相同的结果。我想在图上呈现波浪，然后使用算法来聚类类似的销售动态。更重要的是，我想知道我是否可以使用过程fft来转换所有时间序列，而不是一个一个（所以建议R：26周后转换新时间序列-查看数据库）

http://imageshack.com/a/img854/1958/zlco.jpg我的数据库；三栏：Product - 展示产品组 Week - 自推出产品以来的时间（周），前 26 周 Sales_gain - 产品的销售额如何按周变化

http://imageshack.com/a/img703/6726/sru7.jpg这就是我的时间序列的样子

我相信我可以使用 fft() 来最终实现这个目标，但是从 fft() 的输出到我的目标的飞跃有点不清楚。

请注意，我对时间序列分析比较陌生（这就是为什么我不能把我的代码放在这里）所以你可以提供任何清晰度，把 fft() 的输出放在上下文中，或者你可以推荐的任何可以有效地完成这个任务的包将不胜感激

我是数字信号处理的新手.我有以下传感器样本数据

Time(milliseconds)            data
------------------    -------------------
0                     0.30865225195884705   
60                    0.14355185627937317   
100                  -0.16846869885921478   
156                  -0.2458019256591797    
198                  -0.19664153456687927
258                   0.27148059010505676   
305                  -0.16949564218521118   
350                  -0.227480947971344 
397                   0.23532353341579437   
458                   0.20740140974521637

这意味着在时间上0我有价值,0.30865225195884705并且在时间上60我有价值0.14355185627937317等等.

从每个传感器获取数据10 milliseconds.所以,我假设采样率应该设置为100 Hz.

我想计算时域信号的总能量.

我读到它可以使用Parseval定理计算如下:

其中X[k]是DFT的x[n],二者的长度N.

有任何建议,如何使用MATLAB计算总能量？

有人告诉我，将平均池化应用于矩阵 M 相当于丢弃 M 的傅里叶表示的高频分量。对于平均池化，我的意思是 2 x 2 平均池化，如下图所示：

我想验证这一点并看看它是如何使用 numpy 工作的。因此，我编写了平均池的简单实现，并复制了一个函数来从这里整齐地显示矩阵：

def prettyPrintMatrix(m):
    s = [['{:.3f}'.format(e) for e in row] for row in m]
    lens = [max(map(len, col)) for col in zip(*s)]
    fmt = '\t'.join('{{:{}}}'.format(x) for x in lens)
    table = [fmt.format(*row) for row in s]
    print '\n'.join(table)

def averagePool(im):
    imNew = np.empty((im.shape[0] /2, im.shape[1]/2))
    for i in range(imNew.shape[0]):
        for j in range(imNew.shape[1]):
            imNew[i,j] = np.average(im[(2*i):(2*i+2), (2*j):(2*j+2)])
    return imNew

现在为了测试傅里叶系数的变化，我运行了以下代码：

M = np.random.random((8,8))
Mpooled = averagePool(M)

# …

我正在努力解决功率谱密度（及其倒数）的正确归一化问题。

我遇到了一个实际问题，假设加速度计的读数以功率谱密度（psd）的形式（以幅度^2/Hz为单位）。我想将其转换回随机时间序列。然而，首先我想了解 PSD 的“前进”方向，即时间序列。

根据[1]，时间序列x(t)的PSD可以通过以下公式计算：

PSD(w) = 1/T * abs(F(w))^2 = df * abs(F(w))^2

其中T是x(t)的采样时间，F(w)是x(t)的傅里叶变换，df=1/T是傅里叶空间中的频率分辨率。然而，我得到的结果并不等于我使用 scipy Welch 方法得到的结果，请参见下面的代码。

第一个代码块取自 scipy.welch 纪录片：

from scipy import signal
import matplotlib.pyplot as plt

fs = 10e3
N = 1e5
amp = 2*np.sqrt(2)
freq = 1234.0
noise_power = 0.001 * fs / 2
time = np.arange(N) / fs
x = amp*np.sin(2*np.pi*freq*time)
x += np.random.normal(scale=np.sqrt(noise_power), size=time.shape)

f, Pxx_den = signal.welch(x, fs, nperseg=1024)
plt.semilogy(f, Pxx_den)
plt.ylim(\[0.5e-3, 1\])
plt.xlabel('frequency \[Hz\]')
plt.ylabel('PSD \[V**2/Hz\]')
plt.show()

我注意到的第一件事是绘制的 psd 随变量 fs 变化，这对我来说似乎很奇怪。（也许我需要相应地调整 …

我对这一切都很陌生，我想从图像中获得幅度谱，然后从修改后的幅度谱重建图像..但现在我得到了一个非常暗的重建。

import numpy as np
import cv2
from matplotlib import pyplot as plt

img = cv2.imread('IMG.jpg',0)

dft = cv2.dft(np.float32(img),flags = cv2.DFT_COMPLEX_OUTPUT)
dft_shift = np.fft.fftshift(dft)

m, a = np.log(cv2.cartToPolar(dft_shift[:,:,0],dft_shift[:,:,1]))

# do somthing with m

x, y = cv2.polarToCart(np.exp(m), a)


back = cv2.merge([x, y])


f_ishift = np.fft.ifftshift(back)
img_back = cv2.idft(f_ishift)
img_back = cv2.magnitude(img_back[:,:,0],img_back[:,:,1])

plt.subplot(131),plt.imshow(img, cmap = 'gray')
plt.title('Input Image'), plt.xticks([]), plt.yticks([])
plt.subplot(132),plt.imshow(m, cmap = 'gray')
plt.title('Magnitude Spectrum'), plt.xticks([]), plt.yticks([])
plt.subplot(133),plt.imshow(img_back, cmap = 'gray')
plt.title('result'), plt.xticks([]), plt.yticks([])
plt.show()

结果

你们能帮我弄清楚为什么这么黑吗？

预先感谢 ：）

编辑 …

在我的程序中，我有一个执行快速傅立叶变换的函数。我知道有非常好的免费实现，但这是一个学习的东西，所以我不想使用它们。我最终通过以下实现找到了此评论（它源自意大利语的 FFT 条目）：

void transform(complex<double>* f, int N) //
{
  ordina(f, N);    //first: reverse order
  complex<double> *W;
  W = (complex<double> *)malloc(N / 2 * sizeof(complex<double>));
  W[1] = polar(1., -2. * M_PI / N);
  W[0] = 1;
  for(int i = 2; i < N / 2; i++)
    W[i] = pow(W[1], i);
  int n = 1;
  int a = N / 2;
  for(int j = 0; j < log2(N); j++) {
    for(int k = 0; k < N; k++) { …

我对在 NumPy 中创建 2D hanning、hamming、Blackman 等窗口感兴趣。我知道 NumPy 中存在一维版本的现成函数，例如np.blackman(51)、np.hamming(51)、np.kaiser(51)、np.hanning(51)等。

如何创建它们的 2D 版本？我不确定以下解决方案是否是正确的方法。

window1d = np.blackman(51)
window2d = np.sqrt(np.outer(window1d,window1d)) 

- -编辑

令人担忧的是，np.sqrt只期望正值，而np.outer(window1d,window1d)肯定会有一些负值。一种解决方案是放弃np.sqrt

有什么建议如何将这些 1d 函数扩展到 2d 吗？