Ćwiczenia 2 2

Zadania

Współczynniki obliczane przez FFT

Najprostsza sytuacja: Badamy współczynniki zwracane przez fft dla sinusoid o różnych częstościach. Wykorzystujemy funkcje do generacji sygnałów napisane na poprzednich zajęciach.

• Proszę wygenerować sygnał delta położony w sekundzie 0,5 na odcinku czasu o długości 1s próbkowany 128Hz. Dla takiego sygnału proszę policzyć transformatę Fouriera i wykreślić zarówno sygnały jak i wartość bezwzględne otrzymanych współczynników.
  • Jak wygląda transformata funkcji delta? Jakie częstości w sobie zawiera?

*

Rekonstrukcja sygnału ze współczynników

Reprezentacja sygnałów w dziedzinie częstości jest dualna do reprezentacji w dziedzinie czasu. To znaczy, że jedną reprezentację można przekształcić w drugą. Do przejścia z dziedziny czasu do częstości używaliśmy transformaty Fouriera (zaimplemantowanej w fft). Przejścia z dziedziny częstości do czasu dokonujemy przy pomocy odwrotnej transformaty Fouriera (zaimplementowanej jako ifft. Mając (zespolone) współczynniki w dziedzinie częstości dla pewnego sygnału,możemy odzyskać jego przebieg czasowy.

• Proszę wygenerować sygnał [math]s(t) = \sin(2\pi t \cdot 1)+\sin\left(2 \pi t \cdot 3+\frac{\pi}{5}\right) [/math] o długości 2,5 s próbkowany 10 Hz, obliczyć jego transformatę Fouriera, a następnie zrekonstruować przebieg czasowy. Sygnał oryginalny i zrekonstruowany wykreślić na jednym rysunku. Korzystając z naszych funkcji generację sygnału można zapisać tak:

(s1,t) = sin(f=1, T=2.5, Fs=10, phi=0)
(s2,t) = sin(f=3, T=2.5, Fs=10, phi=np.pi/5)
s = s1 + s2

• Dla porównania proszę zrekonstruować sygnał korzystając jawnie z postaci odwrotnej transformaty Fouriera danej wzorem:

[math]x(n)=\frac{1}{N} \sum_{k=0}^{N-1} S[k] \cdot \exp\left(2j \pi n\frac{ k}{N}\right)[/math]

Zwróćmy uwagę, że n w powyższym wzorze to indeksy punktów czasu, zatem wiążą się one z czasem zwracanym przez nasze funkcje następująco:

n = np.arange(0,len(t),1).

zatem kod rekonstruujący może być np. taki:

n = np.arange(0,len(t),1)
s_rekonstrukcja = np.zeros(len(t))
for k in range(0,N):
	s_rekonstrukcja += 1.0/N * S[k]*np.exp((2j*np.pi*n*k)/N)

*

Efekt nieciągłości funkcji

• Wygenerować sinusoidę o następujących własnościach: f=10 Hz, T=1, Fs=100 Hz, i fazie = 1;
• Przy pomocy subplotów proszę sporządzić rysunek zgodnie z ponższym opisem:
  • subplot(2,2,1): przebieg sygnału w czasie
  • subplot(2,2,2): moduł jego transformaty Fouriera (narysować za pomocą funkcji py.stem, skalę na osi częstości można uzyskać wywołując funkcję F = FFT.fftfreq(len(s), 1/Fs) gdzie s to nasz sygnał),
  • subplot(2,2,3): Proszę wykreślić trzykrotnie periodycznie powielony oryginalny sygnał. Można go skonstruować wywołując funkcję: s_period = np.concatenate((s,s,s)).
  • subplot(2,2,4): moduł transformaty Fouriera s_period (narysować za pomocą funkcji py.stem, skalę na osi częstości można uzyskać wywołując funkcję F = FFT.fftfreq(3*len(s), 1/Fs)

• Powtórz te same kroki dla sinusa o częstości 10.3 Hz.

Pytania:

1. Czym różnią się przedłużenia sinusoidy 10 Hz od sinusoidy 10.3 Hz? Proszę zwrócić uwagę na miejsca sklejania sygnałów.
2. Skąd bierze się widoczna różnica w widmie sinusoidy 10 Hz i 10.3 Hz?

*

Długość sygnału a rozdzielczość widma FFT

Z dotychczasowych rozważań o transformacie Fouriera ograniczonych w czasie sygnałów dyskretnych wynika, że w widmie reprezentowane są częstości od [math]-F_N[/math] do [math]F_N[/math] gdzie [math]F_N[/math] to częstości Nyquista. Dostępnych binów częstości jest N - tyle samo ile obserwowanych punktów sygnału. Zatem zwiększenie długości sygnału w czasie poprawia "rozdzielczość" reprezentacji częstotliwościowej sygnału.

Załóżmy, że dysponujemy jedynie sekwencją N próbek pewnego sygnału. Rozważymy teraz jakie można przyjąć strategie przedłużania tego sygnału w celu zwiększenia gęstości binów częstotliwościowych i jakie te strategie mają konsekwencje.

Przedłużanie sygnału zerami

Inną popularną metodą na zwiększanie ilości binów w transformacie Fouriera jest przedłużanie sygnału zerami (zero-padding). Jest to szczególny przypadek następującego podejścia: Nasz "prawdziwy" sygnał jest długi. Oglądamy go przez prostokątne okno, które ma wartość 1 na odcinku czasu dla którego próbki mamy dostępne i 0 dla pozostałego czasu (więcej o różnych oknach będzie na kolejnych zajęciach). W efekcie możemy myśleć, że oglądany przez nas sygnał to efekt przemnożenia "prawdziwego" sygnału przez okno. Efekty takiego przedłużania proszę zbadać przy użyciu poniższego kodu.

• Jak można zinterpretować wyniki tego eksperymentu w świetle twierdzenia o splocie?

# -*- coding: utf-8 -*-

import pylab as py
import numpy as np
import numpy.fft as FFT

def sin(f = 1, T = 1, Fs = 128, phi =0 ):
	'''sin o zadanej częstości (w Hz), długości, fazie i częstości próbkowania
	Domyślnie wytwarzany jest sygnał reprezentujący 
	1 sekundę sinusa o częstości 1Hz i zerowej fazie próbkowanego 128 Hz
	'''
 
	dt = 1.0/Fs
	t = np.arange(0,T,dt)
	s = np.sin(2*np.pi*f*t + phi)
	return (s,t)	
	
(s1,t) = sin(f = 15.0, T =0.1, Fs = 100, phi = 0)
(s2,t)= sin(f = 20.0, T =0.1, Fs = 100, phi = 0)
s=s1+s2
py.clf()
py.subplot(2,2,1)
py.plot(t,s)
py.subplot(2,2,2)
S = FFT.fft(s)
F = FFT.fftfreq(len(s),0.01)
py.stem(F,np.abs(S)/len(S))
py.xlim((-50,50))
py.ylim((0,0.7))
z= np.zeros(len(s))
py.subplot(2,2,3)
s_period = np.concatenate((s,z,z,z,z,z,z,z,z,z))
t_period = np.arange(0,len(s_period)/100.0,0.01)
py.plot(t_period,s_period)

py.subplot(2,2,4)
S_period = FFT.fft(s_period)
F_period = FFT.fftfreq(len(s_period),0.01)
py.stem(F_period,np.abs(S_period)/len(S))
py.stem(F,np.abs(S)/len(S),linefmt='r-', markerfmt='ro')
py.xlim((-50,50))
py.ylim((0,0.7))
py.show()