/ / numpy [closed]による単純なたたみ込みによる予期しない結果 - python、numpy、scipy、信号処理

numpy [閉じた] - 単純な畳み込みによる予期しない結果 - python、numpy、scipy、シグナル処理

私はしばらくの間DSPをやっていないが、私は手渡さない私は基本の私の把握がこれまでに滑ることを期待した。 複雑な指数関数を使ってトーンを畳み込むスクリプトがあります。その結果、私はシフトしたトーンになると期待しています。私の結果は非常に予想外です - 私は3トーンを得ています、そして、どれも私が期待する頻度にありません。誰かが私がなぜこれらの結果を得ているのか説明できますか?

ここにスクリプトがあります。

import sys
import numpy
import math
import scipy
from pylab import *

def gen_tone(f, fs, length):
t = linspace(0, length, length * fs)
return cos(2.0 * pi * f * t)

def gen_exp(f, fs, length):
t = linspace(0, length, length * fs)
return numpy.exp(1.0j * 2 * pi * f * t)

def plot_fft(f, fs):
FFT = abs(scipy.fft(f, 1024)) / f.size
figure()
plot(FFT)

f100 = gen_tone(8000, 44100, 1)
f200j = gen_exp(1000, 44100, 1)

res = scipy.signal.fftconvolve(f100, f200j, "full")

plot_fft(f100, 44100)
plot_fft(f200j, 44100)
plot_fft(res, 44100)

show()

回答:

回答№1は1

あなたは 周波数シフト特性。 (例えば、 http://ocw.usu.edu/Electrical_and_Computer_Engineering/Signals_and_Systems/5_6node6.html;ラベルの付いたセクションまで少し下にスクロールします。 周波数シフト特性つまり、フーリエ変換を f(t)F(w)その後、のフーリエ変換 f(t)*exp(j*w0*t)F(w - w0)。表現 f(t)*exp(j*w0*t) の点ごとの乗算 f(t) そして exp(j*w0*t), ない たたみ込み。

予想した結果を見るには、これを置き換えます。

res = scipy.signal.fftconvolve(f100, f200j, "full")

〜と

res = f100 * f200j

プロット関数を次のように変更すると、結果が見やすくなります。

def plot_fft(f, fs):
FFT = abs(fft(f, 1024)) / f.size
freq = fftfreq(1024, 1.0/fs)
ndx = freq.argsort()
figure()
plot(freq[ndx], FFT[ndx])
grid(True)

追加する

from scipy.fftpack import fft, fftfreq

あなたのスクリプトの一番上に。

のFFTプロットの-8000と8000のピークがわかります。 f100 のFFTプロットでは、-7000と9000にシフトします。 res.


回答№2の場合は0

必要なのは周波数領域での乗算(シフトトーン)なので、それを実行する必要があります。

res = f100 * f200j

の代わりに

res = scipy.signal.fftconvolve(f100, f200j, "full")

下の図は、FFTではなくスペクトルの写真です。 input (黒)そして res (青い)信号。シフトされた信号は9kHzでディラックを示します(1000の8000アップシフト)。 ここに画像の説明を入力

そしてこれが私のplotSpectrum関数です。

def plotSpectrum(y,Fs):
"""
Plots a Single-Sided Amplitude Spectrum of y(t)
"""
n = len(y) # length of the signal
k = arange(n)
T = n/Fs
frq = k/T # two sides frequency range
frq = frq[arange(int(n/2))] # one side frequency range

Y = fft(y)/n # fft computing and normalization
Y = Y[arange(int(n/2))]

plot(frq,abs(Y),"r") # plotting the spectrum
xlabel("Freq (Hz)")
ylabel("|Y(freq)|")