|
|
|
|
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Michael Raiford <mra### [at] hotmailcom> wrote:
> Heh, as more come to me. I'm tossing together a little C# app to play
> around with the FFT and be able to examine it in both rectangular and
> polar forms.
Aha, I knew I had that lying around somewhere. If you just want to play with
the behavior and not implement it yourself, there's always Matlab, or its free
equivalent, Octave. Here's the simplest low-pass filter:
n=101; tmin=0; tmax=10;
%calculate time & frequency
tbase = linspace(tmin,tmax,n);
dt = (tmax-tmin)/(n-1);
ni = floor((n-1)/2);
fbase = [0:floor(n/2),ni:-1:1] /ni/2.0/dt;
%calculate y and fft(y)
y = rand(size(tbase));
Y = fft(y);
%filter
Y_filt = Y;
Y_filt(find(fbase>2.0)) = 0;
y_filt = ifft(Y_filt);
%plot
figure(1); plot(tbase,y,tbase,y_filt);
figure(2); plot(fbase,abs(Y_filt));
- Ricky
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
> That much I got, and the order of the data seemed to make sense to me,
> what's racking my brain is the butterfly... That seems to be the "magic"
> piece that makes the whole thing work.
Well, you *could* take the time-domain signal and cut it into two
halves. But that's not very easy to undo in the frequency domain. So
what they usually do is take all the odd samples and all the even samples.
Taking all the even samples is just halving the sampling rate. [Assuming
0-based indexing!] Halving the time domain means doubling the frequency
domain. (Think about it: You take a sound, cut out half the samples,
that's like doubling the playback speed. All the pitches are twice as high!)
Taking all the odd samples is like shifting one sample to the side and
then halving the sampling rate as per above. In the frequency domain,
such a time shift is a multiplication by a sinewave.
So, to undo all this, double both signals, multiply one by a sinewave,
then add them back together. Hence, the butterfly.
>> Combining in the frequency domain must undo what was done
>> by splitting in the frequency domain...
Grr... "Combining in the *frequency domain* must ondo what was done by
splitting in the *time domain*." Get it right, boy! >_<
>> Multiplication in the frequency domain is convulotion in the time
>> domain, and vice versa. FFT convolution is mildly more tricky than it
>> looks at first. (You gotta chop the input into "blocks" to run it
>> through the FFT, and you need to overlap them after you're done or the
>> "blur" will stop at the block edges - which is wrong.)
>
> I'm supposing this is where Hamming and Blackmann windows come into play?
No. Those come into play if you're trying to desing a digital filter.
(E.g., a "perfect" lowpass filter has an infinite frequency response,
but we can't do that. If you just chop the ends off it screws up the
frequency response. Using (say) a Hamming window reduces the effect -
but doesn't eliminate it.)
Recall that the convolution of a signal by some kernel is *longer* than
the original signal. This is the effect you have to repeat when doing
FFT convolution. (The details escape me at this exact moment... Suffice
it to say that your inverse-FFT results have to be overlapped to get a
smooth result.)
> I knew you'd have something to add to this, since you seem to be quite
> advanced in mathematics.
Meh. I wanted to build a digital filter, so I read a rather excellent
book on the subject. ;-)
Wanna see the Z-transform?
--
http://blog.orphi.me.uk/
http://www.zazzle.com/MathematicalOrchid*
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Orchid XP v8 wrote:
>
> Meh. I wanted to build a digital filter, so I read a rather excellent
> book on the subject. ;-)
But, you see ... you have an interest in such things.
> Wanna see the Z-transform?
Z transform?
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
>> Meh. I wanted to build a digital filter, so I read a rather excellent
>> book on the subject. ;-)
>
> But, you see ... you have an interest in such things.
I read Wikipedia - and discovered that it is *completely useless* for
learning this kind of thing! But then I found a book online, and read
that. ;-)
>> Wanna see the Z-transform?
>
> Z transform?
Weeell... The Fourier transform comes in several different "flavours",
but essentially you have
- The continuous Fourier transform takes a formula and turns it into a
different formula.
- The discrete Fourier transform takes some numbers and turns them into
some other numbers.
There is a generalisation of the [continuous] Fourier transform called
the Laplace transform. The correspondin discrete version is called the
Z-transform.
The Laplace transform is the trippy mumma you use for designing Infinite
Inpulse Response filters.
See, Finite Impulse Response (FIR) filters are very easy to design, and
very flexible, but they take quite a bit of compute power. If you want
precise filtering (e.g., for scientific purposes) then you're going to
use FIR.
On the other hand, Infinite Impulse Response filters use feedback to
generate an impulse response which is effectively "infinitely long". The
downside is that only certain impulse responses can be created by
feedback. Oh, and that controlling feedback is damned tricky.
Enter the Laplace transform. This makes it relatively easy to figure out
how much feedback to apply.
It just so happens that the impulse response of a perfect lowpass filter
"nearly" matches the kinds of impulse responses that an IIR can produce.
The net result of this is that with only a tiny amount of calculation,
you can get pretty good results.
In other words, IIR is massively, massively faster than FIR. (It's also
far less precise, nowhere near as flexible, and way harder to design.)
While we're on the subject, the Laplace transform turns differential
equations into algebraic equations, making them drastically easier to
solve. As I understand it, that's why mathematicians came up with it in
the first place.
Jesus Christ, SOMEBODY HIRE ME!! >_<
--
http://blog.orphi.me.uk/
http://www.zazzle.com/MathematicalOrchid*
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Orchid XP v8 wrote:
> So, to undo all this, double both signals, multiply one by a sinewave,
> then add them back together. Hence, the butterfly.
I want you to know that's the first explanation of it I've ever been able to
follow. I'm not a very mathematical person.
You definitely should write some Haskell documentation. :-)
--
Darren New, San Diego CA, USA (PST)
Why is there a chainsaw in DOOM?
There aren't any trees on Mars.
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
>> So, to undo all this, double both signals, multiply one by a sinewave,
>> then add them back together. Hence, the butterfly.
>
> I want you to know that's the first explanation of it I've ever been
> able to follow.
Shamelessly paraphrased from The DSP Guide. ;-) But thank you.
It seems to be that (in general) there are two ways to explain any
mathematical result:
- As a sequence of mechanical transformations of symbols.
- As a vague but intuitive examination of *why* this (or at least
something roughly this shape) should work.
Non-mathematicians tend to respond better to the latter, and it's
terribly hard to glean from the former.
I still vividly remember the day I figured out why that famous formula
actually solves quadratic equations. Maybe I'll share it with you?
> I'm not a very mathematical person.
...says the guy who actually knows WTF a nondeterministic Turing machine
*is*! :-P
> You definitely should write some Haskell documentation. :-)
Heh, yeah.
Sadly, before you can document something, you must comprehend it. Also,
the Haskell library documentation uses literate proramming - i.e., the
documentation is inside the library source code, so to update it you'd
have to submit a patch to the library maintainer, yadda yadda yackt...
--
http://blog.orphi.me.uk/
http://www.zazzle.com/MathematicalOrchid*
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Orchid XP v8 wrote:
> I still vividly remember the day I figured out why that famous formula
> actually solves quadratic equations. Maybe I'll share it with you?
OK, so "everybody knows" (?!) that ax^2 + bx + c = 0 can be transformed
into the magical expression
x = (-b +/- Sqrt(b^2 - 4ac)) / 2a
But why is that, exactly? Well, looking in my sister's A-level maths
book [unlike me, my sister was *taught* mathematics] I discover a long
sequence of algebraic manipulations that slowly transforms the one into
the other.
Well, that proves it *works*, but *why* does it work??
One day, I figured it out. Watch this:
The (linear) polynomial (x - A) has one solution: A. And the quadratic
polynomial (x - A)(x - B) has two solutions: A and B. So here we just
*constructed* a quadratic with *known* solutions!
(This may seem completely unremarkable, but remember I have no
mathematical training. It seemed pretty impressive to me!)
Now, if we expand that out, we get
x^2 - Ax - Bx + AB = 0
and with some rearranging,
x^2 + (-A-B)x + (AB) = 0
If we match that up against
x^2 + bx + c = 0
then we have
b = -A-B
c = AB
Or, put another way, -b = A+B. In other words, if we take
ax^2 + bx + c = 0
and divide through a, we now know the sum and product of the solutions!
But how do we figure out the solutions themselves given only their sum
and product?
Well, given that -b = A+B, there are an infinite number of possible
solutions. Similarly, c = AB on its own has infinitely many solutions.
Hopefully there is only two solutions to *both* equations. But what are
they? Hmm.
The next revelation: Since -b = A+B, then -b/2 = (A+B)/2. In other
words, the arithmetic mean of A and B. So -b/2 is a point exactly half
way between the two solutions. If I can just figure out how far apart
they are, I can compute them!
Clearly you can't determine the seperation from the sum. It's not
obvious that you can determine it from the product either though. But
that product must come into the equation somehow. Hmm.
Of course, I already know what the final solution is, just not why. That
gave me a clue.
So what we have is A+B, but what we *want* is A-B. Well now, what
happens if we square these?
(A+B)^2 = A^2 + 2AB + B^2
(A-B)^2 = A^2 - 2AB - B^2 (A mathematician would know this off by
heart, but it took me a while to work out.)
Ooo, that's *interesting*... The difference between the two is only in
the product in the middle! And we *have* that product!
From here, it's just a matter of tying up loose ends. b^2 = (A+B)^2.
The difference between sum^2 and difference^2 is exactly 4AB, so we have
(A-B)^2 = (A+B)^2 - 4AB = b^2 - 4c
Since we want to add/subtract half the difference to half the sum, we have
x = -b/2 +/- Sqrt(b^2 - 4c)/2
= (-b +/- Sqrt(b^2 - 4c))/2
Ah, but we already eliminated "a", didn't we? Can we put that right back
into the main formula? Well, comparing with the "standard" formula and
shifting a little algebra, it becomes clear that we can.
x = (-b +/- Sqrt(b^2 - 4ac))/2a
QED.
Bouyed up by this, I immediately set about performing the same trick for
the cubic. Obviously, this did not end well.
However figured out that you can compute the difference from the sum and
the product... GENIUS!! o_O I would never have known in a million years...
--
http://blog.orphi.me.uk/
http://www.zazzle.com/MathematicalOrchid*
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Orchid XP v8 wrote:
> - As a sequence of mechanical transformations of symbols.
That's the "formal" mechanism I was talking about in the programming
language thread.
> - As a vague but intuitive examination of *why* this (or at least
> something roughly this shape) should work.
Yeah. Math tends to go for the "do all the gruntwork first, building bottom
up to the final statement, and finish", usually without explaining why you
would care about the final statement. :-)
> I still vividly remember the day I figured out why that famous formula
> actually solves quadratic equations. Maybe I'll share it with you?
Sure.
>> I'm not a very mathematical person.
>
> ...says the guy who actually knows WTF a nondeterministic Turing machine
> *is*! :-P
That's pretty trivial mathematics compared to tensor calculus type stuff.
:-) I'm not too bad at discrete math.
> documentation is inside the library source code, so to update it you'd
> have to submit a patch to the library maintainer, yadda yadda yackt...
Write it, patch it, suck out the docs, say "isn't this great? Here's the patch."
--
Darren New, San Diego CA, USA (PST)
Why is there a chainsaw in DOOM?
There aren't any trees on Mars.
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
Darren New <dne### [at] sanrrcom> wrote:
> That's pretty trivial mathematics compared to tensor calculus type stuff.
> :-) I'm not too bad at discrete math.
Tensor calculus isn't so bad. Among other things, I gave my nephew my tensor
calculus book this past summer for his baby shower. Gotta get kids started
early these days, y'know.
- Ricky
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
| |
|
|
>> That's pretty trivial mathematics compared to tensor calculus type stuff.
>> :-) I'm not too bad at discrete math.
>
> Tensor calculus isn't so bad. Among other things, I gave my nephew my tensor
> calculus book this past summer for his baby shower. Gotta get kids started
> early these days, y'know.
I still can't figure out what a "tensor" is. :-P
--
http://blog.orphi.me.uk/
http://www.zazzle.com/MathematicalOrchid*
Post a reply to this message
|
|
| |
| |
|
|
|
|
| |
|
|