Introduction to Fourier Transforms in Oceanography Data Analysis
E N D
Presentation Transcript
Laboratory in Oceanography: Data and Methods Intro to the Signal Processing Toolbox MAR550, Spring 2013 Miles A. Sundermeyer
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Fourier Transform • Suppose we have time or space series data ... • wish to quantify information content of signal • wish to separate periodic component from random component
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Fourier Transform (cont’d) • Basic assumptions • x(t) is one realization from an ensemble of realizations • x(t) has a mean and correlation function, • x(t) is stationary • mean and correlation function are independent of t • (i.e., “weakly” stationary) • make ergodic assumption – can replace an ensemble average with average over time of single realization (in general, don’t have multiple realizations)
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Fourier Transform (cont’d) • Define a finite Fourier transform as: • Define “Power Spectrum” as: • where * denotes the complex conjugate • The power spectrum quantifies the amount of energy contained in different frequencies in the time series. • The “theoretical” power spectrum has the property: where k denotes realizations within an ensemble
0 T ... ... -T 0 T 2T Intro to Signal Processing Toolbox Basics of Fourier Transforms • Fourier Transform (cont’d) • Problems with this: • have discrete data (digitized) • not infinite time series • only have one realization • In practice, we thus perform Fourier analysis on our single realization: • By doing this, implicitly assume our finite interval time series is periodic.
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Fourier Transform (cont’d) • Matlab uses Fourier transform equivalent to continuous integral transform on infinite domain: • Discrete transform on finite domain:
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Example: simple fft • >> x = 5 + 3*cos(2*pi*[0:7]/8) • >> X = fft(x); % forward fft • >> xnew = ifft(X); % inverse fft • >> [x' fft(x)' xnew'] • ans = • 8.0000 40.0000 8.0000 • 5.1213 12.0000 - 0.0000i 5.1213 • 5.0000 0.0000 + 8.0000i 5.0000 • 4.8787 0 + 0.0000i 4.8787 • 2.0000 0 2.0000 • 0.8787 0 - 0.0000i 0.8787 • 5.0000 0.0000 - 8.0000i 5.0000 • 9.1213 12.0000 + 0.0000i 9.1213 • Note: • Imaginary parts are all zero - no sine component • First fft value is freq (k-1) = 0, cos(0) = 1, => fft = (npts)*(mean(x)) = 8x5 = 40 • 2nd & 8th fft values are same & real, represent cosine variability with 8 points, • i.e., freq of 2p/8. Amp of cosine variability in orig signal = 2*X2/N • Other terms are zero since zero energy at other freqs.
Intro to Signal Processing Toolbox Basics of Fourier Transforms • Example: simple fft (cont’d) • Add a sine component and repeat • >> x = 5 + 3*cos(2*pi*[0:7]/8) -2*sin(4*pi*[0:7]/8) • >> X = fft(x); % forward fft • >> xnew = ifft(X); % inverse fft • >> [x' fft(x)' xnew'] • ans = • 8.0000 40.0000 8.0000 • 5.1213 12.0000 - 0.0000i 5.1213 • 5.0000 0.0000 + 8.0000i 5.0000 • 4.8787 0 + 0.0000i 4.8787 • 2.0000 0 2.0000 • 0.8787 0 - 0.0000i 0.8787 • 5.0000 0.0000 - 8.0000i 5.0000 • 9.1213 12.0000 + 0.0000i 9.1213 • Note: • X3 = 8i, X7 = -8i ... Xn and XN+2-n are complex conjugates • Imag parts of X2 and X7 => sine w/ freq 2*2p/N has amp 2*X3/8 = 2. • In General, frequencies represented by fft are: 2*pi(k-1)/N, k = 0:(N/2) • zero freq (mean), • 2*pi*(1/N) (lowest) ... 2*pi*((N/2 - 1)/N) (highest = Nyquist freq)
Intro to Signal Processing Toolbox Frequency Spectra • Example: Muddy Creek, Chatham, MA • stage data – fft/spectrum via 4 methods: • Harmonic analysis • 1/N X*X • Matlab’s ‘spectrum’ • Matlab’s ‘periodogram’
Intro to Signal Processing Toolbox Frequency Spectra • Variance Preserving Form • Variance preserving form: • f· Pxx plotted on a semilogx axis
Intro to Signal Processing Toolbox Cautions for Fourier Space – Gibbs Phenomenon
Intro to Signal Processing Toolbox Cautions for Fourier Space - Aliasing
Intro to Signal Processing Toolbox Cautions for Fourier Space - Aliasing signal freq Nyquist freq
Intro to Signal Processing Toolbox Signal Processing Toolbox • Convolution and filters • The convolution of two functions is defined as: • where ∗ denotes the convolution operation. • In Fourier space, the convolution is the product of the Fourier transforms of the functions:
Intro to Signal Processing Toolbox Signal Processing Toolbox • Convolution and filters (cont’d) • Matlab’s ‘fdesign’ function for filter building
Intro to Signal Processing Toolbox Signal Processing Toolbox Example: Low-Pass Filter
Intro to Signal Processing Toolbox Signal Processing Toolbox Example: Low-Pass Filter (cont’d)
Intro to Signal Processing Toolbox Signal Processing Toolbox Example: Windowing
Intro to Signal Processing Toolbox Signal Processing Toolbox Example: Windowing
Intro to Signal Processing Toolbox Signal Processing Toolbox • Spectral Estimators in Matlab • Spectral analysis includes three types of spectral estimators — power spectral density (PSD), mean-square spectrum (MSS) and pseudo spectrum. • Power spectral density (psd) measures power per unit of frequency and has power/frequency units. • Mean-square (power) spectrum (msspectrum) measures power at a specific frequency. • Pseudospectrum (pseudospectrum) returns a pseudo spectrum that does not have any units.
Intro to Signal Processing Toolbox Signal Processing Toolbox • Useful Tidbits: • fft, ifft - compute forward and inverse fft • spectrum - for computing various types of spectra • spectrum.welch - for computing windowed spectra • butter - for computing Butterworth filters • freqz - for computing Fourier representations of filters • filter, filtfilt - for time domain filtering • Some References: • Bendat, J. S., and A. G. Piersol: Random Data: Analysis and Measurement Procedures (1st Ed. 1971) • Priestly, M. B.: Spectral Analysis and Time Series. 1983.