C / C++ / MFC

19-Jun-08 20:36

Hello,

I have curve data into an array and I need to apply a low-pass filter (Bessel and Butterworth) on it.
The sampling rate of my data is not constant (for example, I get a value at 1s, then 1.2s, then 1.8s...).
As I said, the curve data are in an array, and the corresponding timestamps are in another :
double time[] = { 1.0, 1.2, 1.8, ... }
double value[] = { 3.1, 4.4, 5.0, ... }

By reading wikis and web sites, it seems that applying a low-pass filter will attenuate the high frequencies (i.e. higher that the cutoff frequency). The attenuation will depend of the type of the filter (Bessel, Butterworth, ...) and its order.
Please tell me if I am wrong!

In other words and in my case, if I have two timestamps which are "to close" (for example : 1/(time[1]-time[0]) > cutoff), I need to modify the corresponding value (value [1]). The new value will be given by the transfer function (which depends, again, of the type of the filter and its order).
Is that right or I am totaly lost ?

If all of this is right, then, how to get the new value ? Taking a look to wikipedia tells me a lot about bessel/butterworth transfer functions, but I do not understand how to apply them.
Could someone provide me an algorithm ?

If all of this is wrong, then, could you please correct me ?

Thanks in advance.

CPallini19-Jun-08 21:09

CPallini

19-Jun-08 21:09

Possibly this page [^] may help (see reference to Numerical Recipes in C).
Smile | :)

jeremie.delaitre19-Jun-08 21:55

19-Jun-08 21:55

Well, thanks for the answer, but Numerical Recipes does not seems to respond to my problem.

Nelek19-Jun-08 22:55

19-Jun-08 22:55

[ Quote ]
it seems that applying a low-pass filter will attenuate the high frequencies (i.e. higher that the cutoff frequency). The attenuation will depend of the type of the filter (Bessel, Butterworth, ...) and its order. Please tell me if I am wrong!
[ /Quote ]

No, you are right. As I told you in the other post (sorry for the delayment, I was away for some days)

[ Quote ]
if I have two timestamps which are "to close" (for example : 1/(time[1]-time[0]) > cutoff), I need to modify the corresponding value (value [1]).
[ /Quote ]

I think you are not having the idea. The filters are to be used with the frecuency of the signals you need to process, but have nothing to do with the timestamps where you are getting the values. Let's try with a example:

......./\............./\............../
......|..|............|..|............|
......|..|............|..|............|
...../....\........./....\........../
..../......\......./......\......../
.../........\...../........\....../
../..........\.../..........\..../
.|............|..|............|..|
.|............|..|............|..|
/..............\/.............\/

You have a signal like that for example. With the Fourier theoreme you can separate the signal into different ones (sinus-signal) which addition will be the original.
The result can be one sinus very thin but taller and other one with more width but shorter.

That means, one signal has more frecuency (the thin one) and more amplitude (vertical size) than the other. If you apply a low-pass filter, depending on which one, its cut-off frecuency and its order you will have a different result. Let's say that the "taller" signal is 150Hz and the shorter one is 50Hz. If your cuttoff frecuency is 60Hz then your result will be a signal quite similar to the second signal that compose the example drawn. It won't be exactly the same because the "taller" signal won't be totally filtered, but quite attenuated (so its vertical size will be much lower).

How much does a filter attenuate? If I remember ok it was a "fixed" order-related rate. For filters of first order the attenuation was 20 db, for second order was 40db and so on. This is easier to see in BODE-Diagramms (remember my last answer?)

first order

-----------
............\
..............\
................\
..................\
....................\

second order

-----------
............\
.............\
..............\
...............\
................\

Both have the same cuttof frecuency (the point where the ramp begins), but the order is the ramp inclination, for 1st order is 20db/dec and for 2nd order is 40db/dec (dec = decade???? "decada" on spanish)

So... turning back to your question...
"if I have two timestamps which are "to close" (for example : 1/(time[1]-time[0]) > cutoff)" this has nothing to do with the frecuency of the signal, this are the time between two measures. Frecuency is the time that the signal needs to be in the same value in the same conditions two times (red marked in my example function). That means, same value (for example 1.2) both going upwards. Then is valid what you say "frec = (1 / (time [1] - time [0]))" and that should be compared with the cutoff frecuency.

But another time a problem. This is only valid for pure signals "sinus-cosinus", "triangle", "square"... if you have signals with "strange" forms... you will need to apply fourier to separate it into different pure signals of different frecuencies. Then apply the filters and see how each component is going to be affected, and remount the end result with the addition of their curves.

About the algorithm... sorry but it is not so easy, it depends on many factors so there is no possible easy "general" application to be given. You will need to investigate and to make it fit your needs.

A tip... Matlab has already implemented some things that can be usefull for you. For example, giving timestamps arrays and values it draws the graphic. With the graphic (if the signal is not very complicated) you will easily get the transfer function. You can "easily" (when you know how) transform to bode-diagramms or apply fourier.

P.S. I am seeing it from the electrical perspective. I mean, I can be wrong in focusing the solution from the "informatic side" but concepts are concepts, and you seem to be a bit lost with them.

Hope it helps

Regards.
--------
M.D.V. Wink | ;)

If something has a solution... Why do we have to worry about?. If it has no solution... For what reason do we have to worry about?
Help me to understand what I'm saying, and I'll explain it better to you
“The First Rule of Program Optimization: Don't do it. The Second Rule of Program Optimization (for experts only!): Don't do it yet.” - Michael A. Jackson
Rating helpfull answers is nice, but saying thanks can be even nicer.

jeremie.delaitre19-Jun-08 23:34

19-Jun-08 23:34

Ok.. I think my problem is a little more complex than expected !

In fact, what you say remember me a lesson I had few years ago... I will dig into it and see what I can do.

Thanks again for your help.

Nelek19-Jun-08 23:43

19-Jun-08 23:43

You are wellcome. Good luck. Smile | :)

Regards.
--------
M.D.V. Wink | ;)

jeremie.delaitre19-Jun-08 23:22

19-Jun-08 23:22

I finally found a website which speaks the same langage than mine :

http://www.dsptutor.freeuk.com/dfilt1.htm[^]

Now that I am smarter, I will try to recapitulate Big Grin | :-D

The transfer function gives me the "filtered values". I must apply it on all my values to get the filtered ones.
To apply the transfer function, I need to compute the value of each coefficient which depend of the type of the filter (Bessel, Butterworth...) and the order.
Next I apply the function on my values using the coefficients in a simple algorithm.

The problem is that what I said here seems to work only for constant sampling. In my case, the elapsed time between two values is not always the same. How can I handle that ?

And another big question is : where the "low-pass thing" enter the game ? Does it mean that for low-pass Bessel filter I need to apply the Bessel transfer function on value which are to close ? or does it mean that the coefficient are different from a Bessel low-pass filter and a Bessel high-pass filter ?

Nelek20-Jun-08 0:29

20-Jun-08 0:29

"The transfer function gives me the "filtered values". I must apply it on all my values to get the filtered ones."

Yes.

"To apply the transfer function, I need to compute the value of each coefficient which depend of the type of the filter (Bessel, Butterworth...) and the order."

correct

"Next I apply the function on my values using the coefficients in a simple algorithm."

You can do it at the same time you are getting the value as well. Under certain conditions (see below)

"The problem is that what I said here seems to work only for constant sampling."

This will most probably be because the constant sampling is given by the frecuency of the signal. but...

"In my case, the elapsed time between two values is not always the same. How can I handle that ?"

The moment you are taking the values (sampling) is not so important. The important thing is that you have to get the transfer function of the signal you are meassuring. Filters can only be used in periodical signals, that means time-repetitive working patrons. Once you have the function that define your signal... A simple transfer function if the signal is pure or a fourier-series function if your signal is a mix of different signals (this is not a must, but having it separated makes it easier). Then you just apply the filter for the value at that moment. Frecuency of a signal should be constant, so you can use "t" (time) or "s" (time in laplaced functions) parameters insteads of using "f".

"Does it mean that for low-pass Bessel filter I need to apply the Bessel transfer function on value which are to close?"

If your signal is not periodical... you must not apply filters. If your signal is periodical... you should not use it that way. Example: You have a pure sinus signal centered in the zero axis with amplitude top at 1 and bottom is -1 (don't remember if it is amplitude 1 or 2). With frecuency 1 Hz, ok? Then ONCE a second you will be in the same point of your graphic, (i.e. every full second at the top of the function value 1).

then your sampling will have always the same pattern, time 1s -> value 1; time 1,25s -> value 0; time 1,5s -> value -1; time 1,75s -> value 0; time 2s -> value 1

and so on with the intermediate values...

Then if you apply a filter which cutoff frecuency is near 1Hz you will have an attenuation of the amplitude but the frecuency will be the same, let's say that instead of 1 and -1 now you have 0,7 and -0,7.

Now if we make the signal a bit more difficult, let's say... a combination of that simple one, and another sinus with frecuency 10Hz and amplitude 0.25, -0.25 then you will have a little and fast sinus wave centered in the first one going up and down.

If we complicate it a bit more adding a slower signal with a big amplitude you will have something more strange as result. And this is the principle of Fouriers theoreme. As long as a signal periodical is, there is always the possibility to get that signal as a combination of different independant signals with different amplitudes and frecuencies, that are called (armónicos on spanish) ???armonics???. So when you apply a filter the result of the original signal will be a bit weird to be seen or imagined, but if you apply the filter to the fourier serie it will be easily recognised how every one of them change and the result will be the same.

The biggest problem you are going to have is to delimitate your scenario, find out the function (or addition of functions) of your signal and choose the right filter and its best order to fit your problem. Once made that, to apply the algorithm of the filter into your programm should be so easy as apply a mathematical formula to a parameter.

Regards.
--------
M.D.V. Wink | ;)

jeremie.delaitre20-Jun-08 2:49

20-Jun-08 2:49

Ok, I think I am close to get it.

If the user gives me :
- the data with their corresponding time stamps (corresponding to values picked from a periodical signal)
- the frequency of the input signal
- the wanted type of low pass filter (Bessel or Butterworth)
- the wanted order of the filter he wants
- the wanted cutoff frequency

I do :
- compute the coefficient of the lowpass filter according to the type of the filter, the order and the cutoff frequency. These coefficients represent the transfer function of the filter.
- decompose my signal into pure ones (Fourrier series) (I do not search informations on this for now)
- apply the transfert function on each pure signal
- build the filtered global signal from each filtered pure signal
- get the filtered value for each of my time stamp from the filtered global signal

Is that right ? (please tell me it is !!)

Nelek20-Jun-08 3:09

20-Jun-08 3:09

More or less is right. The only problem is that if the user just gives you discrete meassures from the signal with the timestamps... it is going to be pretty difficult to imagine how the signal is. You will need A LOT of points to be able to imagine/extrapolate how is the function that fits that values in that time. Very hard, but needed (in my opinion).

If it is possible that the user gives the function that determine the curve of the original signal will be much much much easier.

Steps about the fourier serie can be omited, but in some scenarios will may make things easier.

More or less the main points are correct (or at least that would be my approach as well).

I am sorry that I don't have my notes from the university when I learned it, I had very good points and examples in them. 8 years and 2000km lapse are guilty. Sigh | :sigh:

Regards.
--------
M.D.V. Wink | ;)

jeremie.delaitre20-Jun-08 3:23

Re: http://www.musicdsp.org

20-Jun-08 3:23

Well, I can say now : it is pretty clear !

I understand the problem of the input signal but it is not possible to have the function. The data provides from experimental measurements.

As you say in your other post, I will ask a clear specification of what my client wants.

In any cases, I have learned a lot of things so it was not time loss.

I really want to thank you again...so..thank you !

http://www.musicdsp.org

div!sion20-Jun-08 0:24

div!sion

20-Jun-08 0:24

Have a look at http://www.musicdsp.org/archive.php?classid=3#237

jeremie.delaitre20-Jun-08 3:10

20-Jun-08 3:10

I will have a look, thanks.

cp987620-Jun-08 0:55

cp9876

20-Jun-08 0:55

jeremie.delaitre wrote:
The sampling rate of my data is not constant (for example, I get a value at 1s, then 1.2s, then 1.8s...).

In that case you have to think very carefully what your data actually represents, and it is not clear from your post. Applying any of the standard linear filtering techniques will cause the amplitude of the signal to increase in areas where there are more samples, this is a simple consequence of linearity. The standard filters, like the ones you mention, are defined in the frequency domain, so there is an inherent assumption of linearity.

It may help your thinking to consider what should happen if you add a second sample very close to one of your existing ones. If you would expect the effect of this to be almost identical to doubling the amplitude of the original sample, then proceeed with the linear filtering option.

If you expect the effect of a second sample to simply increase the certainty with which you know the sample value, then what you really want is some form of curve fitting, and for this numerical recipes would be a good starting point.

Peter
"Until the invention of the computer, the machine gun was the device that enabled humans to make the most mistakes in the smallest amount of time."

jeremie.delaitre20-Jun-08 3:01

20-Jun-08 3:01

Yes you are right, I must think about it.
I understand the problem now as I understand a little bit more linear filtering.
I do not know if this if the case of my client, so I will ask him Smile | :)

In fact, I think he wants really linear filtering, but I did not explain me the problem very precisely...

Thanks for pointing me that problem.

Nelek20-Jun-08 3:12