Wednesday, December 16, 2009

Band Pass Filtering: SCORE SQUARED!

In the last post, I showed how a truncated cosine wave did a better job than an exponential function at smoothing data. The goal was to produce a band pass filter, which filtered fluctations occuring at both high and low frequencies, leaving only those in a desired band. Differentiation was used for the high-pass part. But the differentiation canceled out the low-pass effect of the exponential smoothing. The truncated cosine smoothing was strong enough to retain its low-pass character even after the differentiation step.

Here I consider another form of smoothing: a cosine squared. The exponential weighting function jumps discontinuously to zero. The truncated cosine has a continuous value to zero but a discontinuous slope when it is clamped at zero. The cosine squared, on the other hand, has a continuous slope as well when it is held at zero.

The three filter shapes are shown in the following plot. They are shown with the time shift I mentioned last time so to be causal. Only past-time data are used.

So the exponential and cosine filters have already been analyzed here. Time to look at the cosine squared filter. Similarly to what we did with the cosine filter:

f*(t) = (2 / πτ) ∫ d Δt · cos²(Δt / τ) f(t ‒ Δt),

where here the limits of integration are again from ‒π/2 τ to +π/2 τ. Note the coefficient is designed to cancel the integral of cos²: that has an average value of 1/2 over the range of πτ (the area of the convolution function should be one). And once again I will assert a time delay of π/2 τ is needed to make the filter causal.

The complex representation in this case becomes:

cos²(Δt / τ) = [2 + exp(2 i Δt / τ) + exp(‒2 i Δt / τ)] / 4

The model for the function remains:

f(t) = f0 cos ωt,

for which I will again use the complex representation:

f(t) = f0 Re{ exp iωt },

where Re{} denotes the real part of the complex number. .

Following the approach used last time:

f*(t) = (2 / πτ) ∫ d Δt · cos²(Δt / τ) f(t ‒ Δt),

= (f0 / πτ) Re { ∫ d Δt · [2 + exp(i 2Δt / τ) + exp(‒i 2Δt / τ)] exp iωt / 2 },

f*(t) = (f0 / πτ) Re { ∫ d Δt · [2 exp iω(t ‒ Δt) + exp(i 2Δt / τ + iω (t ‒ Δt) ) + exp(‒i Δt / 2τ + iω(t ‒ Δt))] / 2}.

Back to bookkeeping for the integration:

f*(t) = (f0 / π) Re { [‒2 exp iω(t ‒ Δt) / ωτ + exp(i 2Δt / τ + iω (t ‒ Δt)) / (2 ‒ ωτ) ‒ exp(‒i 2Δt / τ + iω (t ‒ Δt)) / (2 + ωτ) ] / 2i }.

Consolidate the denominator:

f*(t) = (f0 / π) Re { [ ([‒2 (4 ‒ ω²τ²) exp iω(t ‒ Δt) + ωτ (2 + ωτ) exp(i 2Δt / τ + iω (t ‒ Δt)) ‒ ωτ (2 ‒ ωτ) exp(‒i 2Δt / τ + iω (t ‒ Δt)) ] / [2i ωτ (4 - ω²τ²)] },

This is a bit more complex than the truncated cosine was! Still, it's not too hard to unwrap into its trigonometric components:

f*(t) = (f0 / π) Re { ((4 ‒ ω²τ²) i + 2ωτ sin(2Δt / τ) ‒ iω²τ² cos(2Δt / τ)) exp(iω (t ‒ Δt)) / ωτ(4 ‒ ω²τ²) }.

Now clean up this mess by converting it to a definite integral. Integration limits are Δt = ‒π/2 τ and +π/2 τ:

f*(t) = (f0 / π) Re {(4 ‒ ω²τ²) i + 2ω sin(π) ‒ iω²τ² cos(π)) exp(iω (t ‒ π/2 τ)) ‒ ((4 ‒ ω²τ²) i + 2ω sin(‒π) ‒ iω²τ² cos(‒π)) exp(iω (t + π/2 τ)) } / ωτ (4 ‒ ω²τ²),

The sine terms go to 0, the cosines to ‒1, and pull out the common exp iωt, then cancel iω²τ² terms:

f*(t) = (f0 / π) Re {[4i exp(‒π/2 iωτ) ‒ 4i exp(π/2 iωτ)] exp iωt} / ωτ (4 ‒ ω²τ²),

That's starting to look more civil! Here it is, the final result for the low-pass filter:

f*(t) = (f0 / π) 8 sin(π/2 ωτ) cos ωt / ωτ (4 ‒ ω²τ²).

Okay.... so the big check here is whether it does the right thing for small ωτ. At small values of ωτ the sin function is proportional to its argument, then I can clean out ωτ's and 4π's from the numerator and denominator:

f*(t) ≈ f0 cos ωt,

Nice! It passed that test. For |ωτ| « 1, it gives me back the original function.

Then there's the funny business near |ωτ| = 2. Here the argument of the sin is near π, so for ωτ = 2 + ε for small ε, we have:

sin(π/2 ωτ) = sin(π + (π/2) ε) ≈ ‒(π/2) ε.

The denominator ωτ (4 ‒ ω²τ²) then becomes (2 + ε) ( 4 - 4 + 4 ε + ε² ) which is to first order 8 ε. So for ωτ = 2 + ε for small ε,

f*(t) ≈ (f0 / π) 8 (‒ π/2) ε cos ωt / 8ε = f0 cos ωt / 2.

So no problem at all: the function has simply fallen to half its maximum amplitude at this frequency.

At large frequency (|ωτ| » 1) it's easy to see the function falls off proportional to ω³τ³. This is really exciting: with the exponential, it fell off at high frequency proportional to ωτ, with the truncated cosine it fell off proportional to ω²τ², while with this function it falls off proportional to ω³τ³. This means even after differentiation the function will still fall off proportional to ω²τ², which is very nice low-pass filtering. We just need to make sure we don't set the cut-off frequency too low such that we lose signals of interest.

Now if we want a band-pass filter, we differentiate and multiply by τ:

τ ∂/∂t f*(t) = ‒(f0 / π) 8 sin(π/2 ωτ) sin ωt / (4 ‒ ω²τ²).

It's easy to see this falls off at large ωτ due to the ω²τ² in the denominator. It also falls off at small ωτ due to the sin ωt in the numerator. So we're good.