Fourier Series and Fourier Transform
------------------------------------

Any stationary signal (mean-variance do not change over time, e.g., periodic signal) 
may be represented as integral of sinusoidal functions of different frequencies (e^(i*2Pi/n)), i is squroot of -1, n is continuous frequency variable):  Refrence: "Fourier Series Wiki"

Discrete signal may be represented by summation over discrete frequencies n

Co-efficients inside the integrtal or summation are the components for each frequency

The coefficients can be extracted by Fourier transform, integrate with the respective 
sinsoidal function for the corresponding frequency: this is Fourier Transform

Discrete Fourier Transform: Consider the time series as a polynomial and evaluate it at sinosoidal functions for all frequencies n = 2*pi / S, for S samples in the input time-series

"DFT" is a divide and conquer algorithm, 
which is essentially a fast polynomial multiplication algorithm



Polynomials
-----------

2 Representations:
Co-efficients (A0, A1, A2, . . ., An-1)  
Evaluated at n points y(x0), y(x1), ..., y(xn-1)

Exchangable: Given A's evaluate at n points,   
or  Given n points solve for n co-efficients
Traditionally: O(n^2) algorithms for each.
D&Q  DFT:  O(n log n)

Given A[0..n-1] DFT evaluates==> on wn^0, wn^1, wn^2, . . ., wn^n-1,
where wn is n-th root of the equation  x^n -1 = 0
Takes O(n) for each eval, O(n^2) for all

Recursive divide conquer DFT takes O(n log n) time.

Why on those complex points?
Any set of n points would, but, those points let us 
develop O(nlogn) DFT, 
and convert it to fast iterative and parallelizable FFT

Another advantage:
wn = e^(i*2pi(theta)/n) = cos(theta) + i*sin(theta)
This allows us to do Fourier analysis of a signal

A[] => y() by DFT (or FFT) transforms frequency spectrum to time series, and
y() => A[] by inverseDFT (or IFFT) transforms from time series to its frequency spectrum



Divide and Conquer-based Efficient Polynomial Multiplication
------------------------------------------------------------

A(x) = Sum [ai *x^i] over i=0->n-1, a's are n coefficients
B(x) = Sum [bi *x^i] over i=0->n-1, b's are n coefficients

Standard way of multiplying A and B:
Two loops over i=0->n-1, j=0->n-1
C(x) = Sum over i=0->2n-2, total 2n-1 coefficients, each a combination of ai's times bj's
With carry-forwards added appropriately
Complexity O(n^2)

Efficient way, p828
(0)Pad each polynomial to 2n-2 -th order by adding an=0, an+1=0, ..., a2n-1=0
(1)Evaluate each polynomial at the 2n-th complex roots of 1
w2n^0=1, w2n^1=wn, w2n^2, w2n^3, ..., w2n^(2n-1)
(2)Then, element-to-element multiply values at each evaluated pts above: A(wn^k)*b(wn^k) =: C(wn^k), for k=0->2n-1
(3)Find coefficients of polynomial C(x) by using w0 through w2n

Steps 1 & 3 each takes O(2nlog2n), & Step-2 takes O(2n)
Total O(nlogn)

------------------------------------------------------------
Three advantages why we evaluate at complex roots of 1:
1) Fourier transform happens as a side effect (main usage of the algorithm)
2) Half of the loop need be run, but more importantly, algorithm becomes highly parallelizable (FFT)
3) Inverse can use the same algorithm with O(nlogn) complexity
Note: Forward DFT is O(nlogn) for evaluating polynomial at any n points, not just on w_n's
------------------------------------------------------------

FFT: Further improvement by iterative butterfly technique, instead of recursive algorithm.
Each butterfly may be a thread, &
butterflies at each level may run parallely on different processors,
thus, reducing time.

===============================================

Understand,
Note: 2n --> n after padding

Complex roots of unity, Fig 30.2, p831

Corollary 30.4 of cancellation lemma, p832

Multiply and example, two polynomial each of order 2, to verify Step 2 above

Matrix mult, p836 O(n^2) traditional way

Recurrence eqation for D&Q, Eq 30.9 p834

Recursive-FFT, p835

Why does it work? yk bottom-of p835, yk+(n/2) top-of p836 computations, & y0 below-alg p835

Complexity? Set up recurrence and solve.

Inverse FFT, eqn 30.11 p837, needed for Step 3

By-The-Way: Poly-mult = Convolution, as in signal processing;
& inverse-FFT is co-eff extraction, as in Fourier analysis for finding frequency components of a time-series
===============================================


Evaluating a Frequency Component of a Time Series
-------------------------------------------------
Given time series y0, y1, y2, ..., yn-1,
at equal intervals (sampling), 
as if at time = 0,1,2,3,...,n-1
or, at   wn^0, wn^2, ..., wn^(n-1)

We want to find amplitude of the frequency component at fk.

Eavaluate (inverse DFT or inverseFFT),
Sum over i=0->n, yi* e^i*(j2Pi/n)*e^(-j2Pifk), 
where j is the sqr-root of -1
Real part is amplitude, imaginery part gives phase at fk,
power is squareroots(Re^2 +Im^2)

Same as doing inverse-FFT above.

Complexity, same as Step 3 above, O(nlogn)