Analog & digital signals Analog Digital of

Analog & digital signals
Analog
Digital
Discrete function Vk of
discrete sampling variable tk,
with k = integer: Vk = V(tk).
0.3
0.3
0.2
0.2
Voltage [V]
Voltage [V]
Continuous function V of
continuous variable t (time,
space etc) : V(t).
0.1
0
-0.1
-0.2
0.1
0
ts ts
-0.1
-0.2
0
2
4
6
time [ms]
8
10
0
6
8
2
4
sampling time, tk [ms]
Uniform (periodic) sampling.
Sampling frequency fS = 1/ tS
Slides adapted from ME Angoletta, CERN
10
Digital vs analog proc’ing
Digital Signal Processing (DSPing)
Advantages
Limitations
• Often easier system upgrade.
• A/D & signal processors speed:
wide-band signals still difficult to
treat (real-time systems).
• Data easily stored.
• Finite word-length effect.
• Better control over accuracy
requirements.
• Obsolescence (analog
electronics has it, too!).
• More flexible.
• Reproducibility.
Slides adapted from ME Angoletta, CERN
Digital system example
General scheme
ms
Antialiasing
V
Sometimes steps missing
ms
A
(ex: economics);
k
- D/A + filter
(ex: digital output wanted).
A
k
V
Digital
Processing
D/A
Filter
Reconstruction
ms
Slides adapted from ME Angoletta, CERN
ANALOG
DOMAIN
V
ms
A/D
DIGITAL
DOMAIN
- Filter + A/D
Filter
ANALOG
DOMAIN
V
Digital system implementation
ANALOG INPUT
Antialiasing
Filter
KEY DECISION POINTS:
Analysis bandwidth, Dynamic range
• Pass / stop bands.
1
• Sampling rate.
A/D
• No. of bits. Parameters.
Digital
Processing
DIGITAL OUTPUT
• Digital format.
What to use for processing?
See slide “DSPing aim & tools”
Slides adapted from ME Angoletta, CERN
2
3
Sampling
1
How fast must we sample * a continuous
signal to preserve its info content?
Ex: train wheels in a movie.
25 frames (=samples) per second.
Train starts
wheels ‘go’ clockwise.
Train accelerates
wheels ‘go’ counter-clockwise.
Why?
Frequency misidentification due to low sampling frequency.
* Sampling: independent variable (ex: time) continuous →
discrete.
Quantisation: dependent variable (ex: voltage) continuous → discrete.
Here we’ll talk about uniform sampling.
Slides adapted from ME Angoletta, CERN
1
Sampling - 2
1.2
__ s(t) = sin(2πf t)
0
1
0.8
0.6
s(t) @ fS
0.4
0.2
f0 = 1 Hz, fS = 3 Hz
0
tt
-0.2
-0.4
-0.6
__ s (t) = sin(8πf t)
0
1
-0.8
-1
__ s (t) = sin(14πf t)
0
2
-1.2
s(t) @ fS represents exactly all sine-waves sk(t) defined by:
sk (t) = sin( 2π (f0 + k fS) t ) , ⏐k ⏐∈
Slides adapted from ME Angoletta, CERN
1
The sampling theorem
A signal s(t) with maximum frequency fMAX can be
Theo* recovered if sampled at frequency f > 2 f
S
MAX .
* Multiple proposers: Whittaker(s), Nyquist, Shannon, Kotel’nikov.
Naming gets
confusing !
Nyquist frequency (rate) fN = 2 fMAX or fMAX or fS,MIN or fS,MIN/2
Example
s(t) = 3 ⋅ cos(50 π t) + 10 ⋅ sin(300 π t) − cos(100π t)
F1
F2
F1=25 Hz, F2 = 150 Hz, F3 = 50 Hz
Condition on fS?
F3
fS > 300 Hz
fMAX
Slides adapted from ME Angoletta, CERN
1
Frequency domain (hints)
• Time & frequency:
frequency two complementary signal descriptions.
Signals seen as “projected’ onto time or frequency domains.
Example
Ear + brain act as frequency analyser: audio spectrum
split into many narrow bands
low-power sounds
detected out of loud background.
• Bandwidth:
Bandwidth indicates rate of change of a signal.
High bandwidth
signal changes fast.
Slides adapted from ME Angoletta, CERN
1
Sampling low-pass signals
Continuous spectrum
(a)
(a) Band-limited signal:
frequencies in [-B, B] (fMAX = B).
-B
0
(b)
B
f
Discrete spectrum
No aliasing
(b) Time sampling
frequency
repetition.
fS > 2 B
-B
0
B fS/2
f
Discrete spectrum
Aliasing & corruption
(c)
0
fS/2
no aliasing.
(c) fS
f
2B
aliasing !
Aliasing: signal ambiguity
in frequency domain
Slides adapted from ME Angoletta, CERN
Antialiasing filter
1
(a)
(a),(b) Out-of-band noise can alias
Signal of interest
Out of band
noise
Out of band
noise
-B
0
B
f
(b)
into band of interest. Filter it before!
(c) Antialiasing filter
Passband: depends on bandwidth of
interest.
Attenuation AMIN : depends on
(c)
-B
0
f
• ADC resolution ( number of bits N).
B fS/2
AMIN, dB ~ 6.02 N + 1.76
Antialiasing
filter
Passband
frequency
• Out-of-band noise magnitude.
Other parameters: ripple, stopband
frequency...
-B
0
B
f
Slides adapted from ME Angoletta, CERN
ADC - Number of bits N
2
Continuous input signal digitized into 2N levels.
Uniform, bipolar transfer function (N=3)
1113
2
Quantization step q =
1
0
-4
-3
-2
-1
0
1
2
3
4
Ex: VFSR = 1V , N = 12
V
-1
V FSR
2N
q = 244.1 µV
010
-2
001
-3
Voltage ( = q)
000
-4
LSB
VFSR
Scale factor (= 1 / 2N )
1
Percentage (= 100 / 2N )
q/2
0.5
0
-4
-3
-2
-1
0
1
-0.5
-q/2
2
3
4
Quantisation error
-1
Slides adapted from ME Angoletta, CERN
ADC - Quantisation error
2
0.3
• Quantisation Error eq in
[-0.5 q, +0.5 q].
Voltage [V]
0.2
0.1
0
0
2
4
6
-0.1
-0.2
time [ms]
8
10
• eq limits ability to resolve
small signal.
• Higher resolution means
lower eq.
Slides adapted from ME Angoletta, CERN
Frequency analysis: why?
• Fast & efficient insight on signal’s building blocks.
• Simplifies original problem - ex.: solving Part. Diff. Eqns. (PDE).
• Powerful & complementary to time domain analysis techniques.
• The brain does it?
time, t
analysis
General Transform as
problem-solving tool
frequency, f
F
S(f) = F[s(t)]
s(t)
s(t), S(f) :
Transform Pair
synthesis
Slides adapted from ME Angoletta, CERN
Fourier analysis - tools
Input Time Signal
Frequency spectrum
2.5
2
1.5
1
Periodic
0.5
0
0
1
2
3
4
time, t
5
6
7
8
Continuous
2.5
(period T)
Aperiodic
2
1.5
1
FS
Discrete
FT
Continuous
T
1
c k = ⋅ ∫ s(t) ⋅ e − j k ω t dt
T
0
− j2 π f t
+∞
S(f) = ∫ s(t) ⋅ e
dt
−∞
0.5
0
0
2
4
6
8
time, t
10
12
2.5
2
Periodic
1.5
1
0.5
(period T)
0
0
1
2
3
4
time, tk
5
6
7
8
Discrete
2.5
Aperiodic
2
1.5
1
0.5
0
0
2
4
time, tk
6
8
10
12
2πkn
−
N
1
j
−
1
~
N
ck = ∑ s[n] ⋅ e
N
n =0
DFS** Discrete
DTFT
Continuous
DFT** Discrete
Note: j =√-1, ω = 2π/T, s[n]=s(tn), N = No. of samples
S(f) =
+∞
∑ s[n] ⋅ e− j 2 π f n
n= −∞
−j
1 N−1
~
ck = ∑ s[n] ⋅ e
N
n =0
**
2πkn
N
Calculated via FFT
Slides adapted from ME Angoletta, CERN
A little history
¾ Astronomic predictions by Babylonians/Egyptians likely via trigonometric sums.
¾ 1669:
1669 Newton stumbles upon light spectra (specter = ghost) but fails to
recognise “frequency” concept (corpuscular theory of light, & no waves).
¾ 18th century:
century two outstanding problems
→ celestial bodies orbits: Lagrange, Euler & Clairaut approximate observation data
with linear combination of periodic functions; Clairaut,1754(!) first DFT formula.
→ vibrating strings: Euler describes vibrating string motion by sinusoids (wave
equation).
¾ 1807:
1807 Fourier presents his work on heat conduction ⇒ Fourier analysis born.
→ Diffusion equation ⇔ series (infinite) of sines & cosines. Strong criticism by peers
blocks publication. Work published, 1822 (“Theorie Analytique de la chaleur”).
Slides adapted from ME Angoletta, CERN
A little history -2
¾ 19th / 20th century:
century two paths for Fourier analysis - Continuous & Discrete.
CONTINUOUS
→
Fourier extends the analysis to arbitrary function (Fourier Transform).
→
Dirichlet, Poisson, Riemann, Lebesgue address FS convergence.
→
Other FT variants born from varied needs (ex.: Short Time FT - speech analysis).
DISCRETE: Fast calculation methods (FFT)
→
1805 - Gauss, first usage of FFT (manuscript in Latin went unnoticed!!!
Published 1866).
→
1965 - IBM’s Cooley & Tukey “rediscover” FFT algorithm (“An algorithm for
the machine calculation of complex Fourier series”).
→
Other DFT variants for different applications (ex.: Warped DFT - filter design &
signal compression).
→
FFT algorithm refined & modified for most computer platforms.
Slides adapted from ME Angoletta, CERN
Fourier Series (FS)
A periodic function s(t) satisfying Dirichlet’s conditions * can be expressed
as a Fourier series, with harmonically related sine/cosine terms.
is
s
+∞
he
t
n s(t) = a0 + ∑ [ak ⋅ cos (k ω t) − bk ⋅ sin (k ω t)]
sy
k =1 For all t but discontinuities
a0, ak, bk : Fourier coefficients.
k: harmonic number,
T: period, ω = 2π/T
s
si
y
T
al
n
1
(signal average over a period, i.e. DC term &
a a = ⋅ s(t)dt
0
zero-frequency component.)
T ∫
0
T
2
ak = ⋅ ∫ s(t) ⋅ cos(k ω t) dt
Note: {cos(kωt), sin(kωt) }k
T
0
form orthogonal base of
T
function space.
2
- bk = ⋅ ∫ s(t) ⋅ sin(k ω t) dt
T
0
* see next slide
Slides adapted from ME Angoletta, CERN
FS convergence
Dirichlet conditions
(a) s(t) piecewise-continuous;
(b) s(t) piecewise-monotonic;
In any period:
(c) s(t) absolutely integrable ,
T
∫
s(t) dt < ∞
0
Example:
square wave
Rate of convergence
T
if s(t) discontinuous then
|ak|<M/k for large k (M>0)
s(t)
T
(a)
(b)
Slides adapted from ME Angoletta, CERN
(c)
FS analysis - 1
T = 2π ⇒ ω = 1
2π
⎧π
⎫
1 ⎪
⎪
a0 =
⋅ ⎨ ∫ dt + ∫ ( −1)dt ⎬ = 0
2π ⎪
⎪⎭
π
⎩0
2π
⎧π
⎫
1 ⎪
⎪
ak = ⋅ ⎨ ∫ cos kt dt − ∫ cos kt dt ⎬ = 0
π ⎪
⎪⎭
π
⎩0
(zero average)
(odd function)
2π
⎫
⎧π
2
1 ⎪
⎪
⋅ { 1− cos kπ } =
- bk = ⋅ ⎨ ∫ sin kt dt − ∫ sin kt dt ⎬ = ... =
k
⋅
π
π ⎪
⎪⎭
π
⎩0
⎧ 4
⎪ k ⋅ π , k odd
⎪
= ⎨
⎪ 0 , k even
⎪
⎩
4
4
4
sw(t) = ⋅ sin t +
⋅ sin 3 ⋅ t +
⋅ sin 5 ⋅ t + ...
π
3⋅π
5⋅π
2π
1.5
square signal, sw(t)
FS of odd* function: square wave.
1
0.5
0
-0.5
0
2
4
6
8
10
t
-1
-1.5
* Even & Odd functions
s(x)
Even :
s(-x) = s(x)
x
s(x)
Odd :
s(-x) = -s(x)
Slides adapted from ME Angoletta, CERN
x
FS synthesis
Square wave reconstruction
from spectral terms
1.5
7
3
15
911
sw1
(t)
sin(kt)
(t)===∑
⋅sin(kt)
sin(kt)
]]]]
∑∑[[--[b-bkbkk⋅⋅sin(kt)
7
3
5
11
9(t)
kkk==1=11
square signal, sw(t)
1
0.5
0
-0.5
-1
-1.5
0
2
4
t
6
8
Convergence may be slow (~1/k) - ideally need infinite terms.
Practically, series truncated when remainder below computer tolerance
(⇒ error).
error BUT … Gibbs’ Phenomenon.
Slides adapted from ME Angoletta, CERN
10
Gibbs phenomenon
1.5
sw 79 (t) =
79
∑ [- bk ⋅ sin(kt)]
k =1
1
square signal, sw(t)
Overshoot exist @
each discontinuity
0.5
0
-0.5
-1
-1.5
0
2
4
t
6
• First observed by Michelson, 1898. Explained by Gibbs.
• Max overshoot pk-to-pk = 8.95% of discontinuity magnitude.
Just a minor annoyance.
• FS converges to (-1+1)/2 = 0 @ discontinuities, in this case.
case
Slides adapted from ME Angoletta, CERN
8
10
FS time shifting
⎧ 4
⎪ k ⋅ π , k odd, k = 1, 5, 9...
⎪
⎪
ak = ⎨ − 4
, k odd, k = 3, 7, 11...
⎪ k⋅π
⎪
⎪
0
, k even.
⎩
- bk = 0
(even function)
Note: amplitudes unchanged BUT
phases advance by k⋅π/2.
square signal, sw(t)
(zero average)
2π
1
0.5
0
-0.5
0
2
4
6
8
10
t
-1
-1.5
rk
4/π
4/3π
θk
f1
3f1
5f1
7f1
f
f1
3f1
5f1
7f1
f
π
ph
as
e
a 0= 0
1.5
am
pl
it
ud
e
FS of even function:
π/2-advanced square-wave
Slides adapted from ME Angoletta, CERN
Complex FS
Euler’s notation:
e-jt = (ejt)* = cos(t) - j·sin(t)
“phasor”
e jt + e − jt
cos(t) =
2
e jt − e − jt
sin(t) =
2⋅ j
s
T
si
1
y
l
a c k = ⋅ s(t) ⋅ e - j k ω t dt
n
a
T
∫
0
is
s
he
t
n
sy s(t) =
Complex form of FS (Laplace 1782). Harmonics
ck separated by ∆f = 1/T on frequency plot.
∞
jk ω t
c
⋅
e
∑ k
k = −∞
z=re
Note:
Note c-k = (ck)*
Link to FS real coeffs.
c 0 = a0
ck =
b
r
θ
a
1
1
⋅ (ak + j ⋅ bk ) = ⋅ (a −k − j ⋅ b −k )
2
2
Slides adapted from ME Angoletta, CERN
jθ
r = a2 + b2
θ = arctan(b/a)
FS properties
Time
Homogeneity
a·s(t)
Additivity
s(t) + u(t)
Linearity
a·s(t) + b·u(t)
Time reversal
Multiplication *
Convolution *
Time shifting
Frequency
a·S(k)
S(k)+U(k)
a·S(k)+b·U(k)
s(-t)
S(-k)
∞
s(t)·u(t)
T
1
⋅ ∫ s(t − t ) ⋅ u( t ) dt
T
0
s(t − t )
Frequency shifting e
+j
2π m t
T ⋅ s(t)
∑ S(k − m)U(m)
m = −∞
S(k)·U(k)
e
−j
2π k ⋅t
T
⋅ S(k)
S(k - m)
Slides adapted from ME Angoletta, CERN
*
FS - “oddities”
Orthonormal base
Fourier components {uk} form orthonormal base of signal space:
T
*
uk = (1/√T) exp(jkωt) (|k| = 0,1 2, …+∞) Def.: Internal product ⊗: uk ⊗ um = uk ⋅ um
dt
∫
uk ⊗ um = δk,m (1 if k = m, 0 otherwise).
(Remember (ejt)* = e-jt )
o
Then ck = (1/√T) s(t) ⊗ uk i.e. (1/√T) times projection of signal s(t) on component uk
Negative frequencies & time reversal
k = - ∞, … -2,-1,0,1,2, …+ ∞,
ωk = kω, φk = ωkt, phasor turns anti-clockwise.
Negative k ⇒ phasor turns clockwise (negative phase φk ), equivalent to negative time t,
⇒ time reversal.
Careful:
Careful phases important when combining several signals!
Slides adapted from ME Angoletta, CERN
FS - power
Average power W :
1
W =
T
Parseval’s Theorem
W=
∞
∑ ck
2
k = −∞
1
= a0 2 +
2
T
∫
s(t) 2 dt ≡ s(t) ⊗ s(t)
o
• FS convergence ~1/k
∞
Example
Pulse train, duty cycle δ = 2 τ / T
s(t)
2τ
T
t
bk = 0
a0 = δ sMAX
ak = 2δsMAX sync(k δ)
⇒ lower frequency terms
∑ ⎛⎜⎝ ak 2 + bk 2 ⎞⎟⎠
k =1
Wk = |ck|2 carry most power.
• Wk vs. ωk: Power density spectrum.
spectrum
2
1
10
-1
10
-2
10
-3
Wk/W0
Wk = 2 W0 sync2(k δ)
kf
0
50
W0 = (δ sMAX)2
sync(u) = sin(π u)/(π u)
100
150
200
∞ W ⎫
⎧⎪
⎪
W = W0 ⋅ ⎨1+ ∑ k ⎬
⎪⎩ k =1 W0 ⎪⎭
Slides adapted from ME Angoletta, CERN
FS of main waveforms
Slides adapted from ME Angoletta, CERN
Discrete Fourier Series (DFS)
Band-limited signal s[n], period = N.
DFS defined as:
is
2π k n
s
y
N
−
1
l
j
−
a
1
N
an ~
ck =
s[n] ⋅ e
N
∑
n =0
~
~
Note: ck+N = ck ⇔ same period N
i.e. time periodicity propagates to frequencies!
s
si
e
2π k n
th
N
−
1
j
n
~
sy s[n] =
ck ⋅ e N
∑
k =0
DFS generate periodic ck
with same signal period
Orthogonality in DFS:
2π n(k -m)
N
−
1
j
1
N
e
= δ k,m
∑
N
n =0
Kronecker’s delta
N consecutive samples of s[n]
completely describe s in time
or frequency domains.
Synthesis: finite sum ⇐ band-limited s[n]
Slides adapted from ME Angoletta, CERN
DFS analysis
DFS of periodic discrete
1-Volt square-wave
s[n]: period N, duty factor L/N
0 1 2 3 4 5 6 7 8 9 10
0
L
N
am
pl
it
ud
e
-5
1
ck
0.24
0.24
0.2
n
1
~
0.6
0.6
0.6
0.6
0.24
0.24
0 1 2 3 4 5 6 7 8 9 10
ph
as
e
L
⎧
,
k = 0, + N, ± 2N,...
⎪
N
⎪
⎪
⎪
~
ck = ⎨
π k (L −1)
⎛ π kL ⎞
⎪ −j
sin
⎜
⎟
N
N ⎠
⎪e
⎝
, otherwise
⋅
⎪
N
⎛π k⎞
sin ⎜
⎪
⎟
⎝ N ⎠
⎩
s[n]
1
θk
0.4π
0.2π
Discrete signals ⇒ periodic frequency spectra.
Compare to continuous rectangular function
(slide # 10, “FS analysis - 1”)
k
0.4π
0.2π
0
2
4 5 6 7 8 9 10
-0.2π
-0.4π
Slides adapted from ME Angoletta, CERN
n
-0.2π
-0.4π
DFS properties
Time
Homogeneity
a·s[n]
Additivity
s[n] + u[n]
Linearity
a·s[n] + b·u[n]
Multiplication *
s[n] ·u[n]
N−1
Frequency
a·S(k)
S(k)+U(k)
a·S(k)+b·U(k)
1 N−1
⋅ ∑ S(h)U(k - h)
N h=0
∑ s[m] ⋅ u[n − m]
Convolution *
S(k)·U(k)
m =0
Time shifting
s[n - m]
Frequency shifting
e
+j
2π h t
T ⋅ s[n]
e
−j
2π k ⋅m
T
⋅ S(k)
S(k - h)
Slides adapted from ME Angoletta, CERN
DFT – Window characteristics
•
Finite discrete sequence ⇒ spectrum convoluted with rectangular window spectrum.
•
Leakage amount depends on chosen window & on how signal fits into the window.
(1) Resolution: capability to distinguish different tones. Inversely proportional to mainlobe width. Wish: as high as possible.
(2) Peak-sidelobe level: maximum response outside the main lobe.
Determines if small signals are hidden by nearby stronger ones.
Wish: as low as possible.
(3) Sidelobe roll-off: sidelobe decay
(2)
(1)
(3)
per decade. Trade-off with (2).
Several windows used (applicationdependent): Hamming, Hanning,
Blackman, Kaiser ...
Rectangular window
Slides adapted from ME Angoletta, CERN
DFT of main windows
Windowing reduces leakage by
minimising sidelobes magnitude.
In time it reduces endpoints discontinuities.
Sampled sequence
Non
windowed
Windowed
Some window functions
Slides adapted from ME Angoletta, CERN
DFT - Window choice
Common windows characteristics
Window type
-3 dB Mainlobe width
-6 dB Mainlobe width
[bins]
1.21
Max sidelobe
level
Sidelobe roll-off
[dB/decade]
Rectangular
[bins]
0.89
Hamming
1.3
1.81
- 41.9
20
Hanning
1.44
2
- 31.6
60
Blackman
1.68
2.35
-58
60
Observed signal
[dB]
-13.2
20
Window wish list
Far & strong interfering components
⇒
Near & strong interfering components
⇒
Accuracy measure of single tone
⇒
high roll-off rate.
small max sidelobe level.
wide main-lobe
NB: Strong DC component can shadow nearby small signals. Remove it!
Slides adapted from ME Angoletta, CERN
DFT - Window loss remedial
Smooth data-tapering windows cause information loss near edges.
Solution:
sliding (overlapping) DFTs.
• Attenuated inputs get next
window’s full gain & leakage
reduced.
2 x N samples (input signal)
DFT #1
DFT #2
DFT #3
• Usually 50% or 75% overlap
(depends on main lobe width).
Drawback: increased
total processing time.
DFT AVERAGING
Slides adapted from ME Angoletta, CERN
DFT - parabolic interpolation
Rectangular window
Hanning window
1.968
0.977
1.967
1.966
0.976
1.965
1.964
0.975
1.963
1.962
198
199
200
201
202
203
0.974
199
200
201
202
¾ Parabolic interpolation often enough to find position of
peak (i.e. frequency).
¾ Other algorithms available depending on data.
Slides adapted from ME Angoletta, CERN
203
204
Systems spectral analysis (hints)
System analysis: measure input-output relationship.
Linear Time Invariant
x[n]
DIGITAL LTI
SYSTEM
δ[n]
y[n]
h[n]
x[n]
h[n]
1
0
y[n] = x[n] ∗ h[n] =
n
DIGITAL
LTI
SYSTEM
h[n]
0
h[t] = impulse response
∞
∑ x[n − m] ⋅ h[m]
y[n] predicted from { x[n], h[t] }
m =0
X(f)
H(f)
Y(f) = X(f) · H(f)
n
H(f) : LTI transfer function
Transfer function can be estimated by Y(f) / X(f)
Slides adapted from ME Angoletta, CERN
Estimating H(f)
G xx (f) = X(f) ⋅ X* (f)
G yx (f) = Y(f) ⋅ X* (f)
(hints)
Power Spectral Density of x[t]
(FT of autocorrelation).
Cross Power Spectrum of x[t] & y[t]
(FT of cross-correlation).
Y(f) Y(f) ⋅ X* (f) G yx
H(f) =
=
=
*
X(f) X(f) ⋅ X (f) G xx
Transfer Function
(ex: beam !)
It is a check on
H(f) validity!
Slides adapted from ME Angoletta, CERN