Reply by Robert Orban January 28, 20092009-01-28
In article <gl6mfd$8uq$1@ork-un.noris.net>, thomas.weber@logomotive.eu says...
> > >I've expected here some more useful information, like read this book, look >at this page or some example code. > >Things like > > "Working in frequency domain is no big deal", > > "Learning something new will take time" or > > "there are a couple of key details you need to know to get things right" > >are not really helpful to start with the topic. > > > >I operate offline but if I apply one of the filter in our signal processing >software (NI Diadem) than it takes some seconds and in the existing Fortran >code some minutes (I have to apply this to many files and channels, so it is >hard to life with a factor in computation time of about 60), therefore I >look for a better filter code. > >BTW It don't have to be the best.
If the filter is FIR, then you can exploit the DOT_PRODUCT intrinsic in modern Fortran (i.e., Fortran 90 or higher). In the time domain, an FIR filter is a convolution of the filter coefficients and the input data, so in Fortran, computing the filter's output sample is as simple as writing output_sample=DOT_PRODUCT(coefficients(1:n),data(j:(j+n-1)) where coefficients is a real array containing n filter coefficients and data is a real array containing the last n input samples. In Fortran, the cleanest way of updating the data array with each new sample is by using the EOSHIFT intrinsic, although it may be more computationally efficient to write a circular buffer explicitly. The advantage of using the DOT_PRODUCT and EOSHIFT intrinsics is that they are probably highly optimized by the compiler.
Reply by Robert Orban January 28, 20092009-01-28
In article <gl1k28$ivp$1@ork-un.noris.net>, thomas.weber@logomotive.eu 
says...
> > >Hi, > > > >I've to apply some special shaped filter to long signals, the filter >coefficients I know already but the existing filter routine is probably to >slow. I want to update the filter routine without becoming a signal >processing expert ;-) > > > >What is the right way to apply a filter > >a) in time domain or > >b) in frequency domain > > > >It is possible to use some existing library code in FORTRAN ?
S/FILSYN will not only design just about any filter you can imagine but will automatically create a Fortran subroutine that implements the filter. (You call the subroutine once per sample.) http://www.alkeng.com/sfilsyn.html (Warning: It ain't cheap.)
Reply by Thomas January 23, 20092009-01-23
Fred Marshall schrieb:
> I once was able to speed up a "filter" by better than a factor of 10 by > going to the frequency domain. It's what was predicted and it's what we > got. > > The key is that the filter and the data are relatively long sequences. The > longer the sequences, the more computational saving. > > Here's how it's done: > > 1) Determine a FIR filter. This is a step unto itself. Call the length of > the filter "N". > > 2) Determine the length of data sequence that you will process at any one > time. Call the length of the data sequence "M". > > 3) Calculate M + N -1. This, or a greater value, will be the length of the > FFT / IFFT sequences. Decide on this value "R" - maybe a power of 2 if your > FFT is finicky about that. > > 3) Pad the filter with zeros to length "R". > > 4) Pad the data with zeros to length "R". > > 5) FFT the padded filter sequence of length "R" and save this complex > sequence. > > 6) FFT each padded data sequence of length "R". > > 7) Multiply the filter sequence by the data sequence (length "R"). > > 8) IFFT the resulting product (length "R"). > > The result will be a real time sequence which is nonzero for M + N -1 > samples. > > Note that N-1 nonzero samples at each end are filter transients. This is > the same as if you applied the filter in time to a data sequence of length M > plus ending with enough zeros to see the tail end. > > Fred > >
thanks a lot for your introduction, depending on length of signal I reached factors between 40 and 110, as I guessed the old routine was not optimal for this problem. Thomas
Reply by Fred Marshall January 21, 20092009-01-21
I once was able to speed up a "filter" by better than a factor of 10 by 
going to the frequency domain.  It's what was predicted and it's what we 
got.

The key is that the filter and the data are relatively long sequences.  The 
longer the sequences, the more computational saving.

Here's how it's done:

1) Determine a FIR filter.  This is a step unto itself.  Call the length of 
the filter "N".

2) Determine the length of data sequence that you will process at any one 
time.  Call the length of the data sequence "M".

3) Calculate M + N -1.  This, or a greater value, will be the length of the 
FFT / IFFT sequences.  Decide on this value "R" - maybe a power of 2 if your 
FFT is finicky about that.

3) Pad the filter with zeros to length "R".

4) Pad the data with zeros to length "R".

5) FFT the padded filter sequence of length "R" and save this complex 
sequence.

6) FFT each padded data sequence of length "R".

7) Multiply the filter sequence by the data sequence (length "R").

8) IFFT the resulting product (length "R").

The result will be a real time sequence which is nonzero for M + N -1 
samples.

Note that N-1 nonzero samples at each end are filter transients.  This is 
the same as if you applied the filter in time to a data sequence of length M 
plus ending with enough zeros to see the tail end.

Fred 


Reply by Thomas January 21, 20092009-01-21
Randy Yates schrieb:
> "Thomas" <thomas.weber@logomotive.eu> writes: > >> Hi, >> >> >> >> I've to apply some special shaped filter to long signals, the filter >> coefficients I know already but the existing filter routine is probably to >> slow. I want to update the filter routine without becoming a signal >> processing expert ;-) >> >> >> >> What is the right way to apply a filter >> >> a) in time domain or >> >> b) in frequency domain >> >> >> >> It is possible to use some existing library code in FORTRAN ? >> >> >> >> Thank you very much. > > Thomas, > > Why do you have to use fortran? Have you considered using GNU Octave > instead? > > h = coef; % your own coefficients here > x = wavread('mywave.wav'); > y = filter(h,1,x); > > Done. > > I presume the filter routine in octave is relatively efficient.
This sounds good, I'll try this soon. Till now I didn't know this software. As I pointed out that I do not need the very best routine only a better one, octave could be a solution. Thanks Thomas
Reply by Randy Yates January 21, 20092009-01-21
"Thomas" <thomas.weber@logomotive.eu> writes:

> Hi, > > > > I've to apply some special shaped filter to long signals, the filter > coefficients I know already but the existing filter routine is probably to > slow. I want to update the filter routine without becoming a signal > processing expert ;-) > > > > What is the right way to apply a filter > > a) in time domain or > > b) in frequency domain > > > > It is possible to use some existing library code in FORTRAN ? > > > > Thank you very much.
Thomas, Why do you have to use fortran? Have you considered using GNU Octave instead? h = coef; % your own coefficients here x = wavread('mywave.wav'); y = filter(h,1,x); Done. I presume the filter routine in octave is relatively efficient. -- % Randy Yates % "She's sweet on Wagner-I think she'd die for Beethoven. %% Fuquay-Varina, NC % She love the way Puccini lays down a tune, and %%% 919-577-9882 % Verdi's always creepin' from her room." %%%% <yates@ieee.org> % "Rockaria", *A New World Record*, ELO http://www.digitalsignallabs.com
Reply by Thomas January 21, 20092009-01-21
"Rune Allnor" <allnor@tele.ntnu.no> schrieb im Newsbeitrag 
news:19032fc6-7c27-40e7-95ea-e277df6dc4d9@z27g2000prd.googlegroups.com...
On 21 Jan, 09:32, "Thomas" <thomas.we...@logomotive.eu> wrote:
> I've expected here some more useful information, like read this book, look > at this page or some example code. > > Things like > > "Working in frequency domain is no big deal", > > "Learning something new will take time" or > > "there are a couple of key details you need to know to get things right" > > are not really helpful to start with the topic.
You didn't ask for help to learn a new topic. You asked which of the time or frequency domains would be 'better' in your case.
> I operate offline but if I apply one of the filter in our signal > processing > software (NI Diadem)
'Your' DSP software? Some of your colleagues ought to be able to answer your questions, then. Either as users or designers of that software. Why not ask them?
> than it takes some seconds and in the existing Fortran > code some minutes (I have to apply this to many files and channels, so it > is > hard to life with a factor in computation time of about 60), therefore I > look for a better filter code.
Then learn the FFT. The Numerical Recipies books are as good as any to start learning from. Or pick up a text on DSP, like anything by Oppenheim and Schafer as (co)authors or the book by Proakis and Manolakis. If you want an introductory book, try Lyons' "Understanding DSP", although this is more an introduction explanation that paints the big picture, rather than a detailed ecplanation of technical concepts. I am not sure how helpful it will be if you want hands-on results quickly. Again, talk with your colleagues. There will be lots of tiny details you will need to wrap your mind around to get things right. Rune I'm working in a smal company and of course I spoke with my colleagues, but they work with the software to measure the data and apply the installed filters only. (There is no DSP expert required for this work.) We need a differnt filter characteristic and I should implement it into a existing Fortran code (written by anybody else years ago). This is already done, but it is to slow. Hence I tryed to find out something to speed it up. Thomas
Reply by Rune Allnor January 21, 20092009-01-21
On 21 Jan, 09:32, "Thomas" <thomas.we...@logomotive.eu> wrote:
> I've expected here some more useful information, like read this book, look > at this page or some example code. > > Things like > > &#4294967295; &#4294967295; "Working in frequency domain is no big deal", > > &#4294967295; &#4294967295; "Learning something new will take time" or > > &#4294967295; &#4294967295; "there are a couple of key details you need to know to get things right" > > are not really helpful to start with the topic.
You didn't ask for help to learn a new topic. You asked which of the time or frequency domains would be 'better' in your case.
> I operate offline but if I apply one of the filter in our signal processing > software (NI Diadem)
'Your' DSP software? Some of your colleagues ought to be able to answer your questions, then. Either as users or designers of that software. Why not ask them?
> than it takes some seconds and in the existing Fortran > code some minutes (I have to apply this to many files and channels, so it is > hard to life with a factor in computation time of about 60), therefore I > look for a better filter code.
Then learn the FFT. The Numerical Recipies books are as good as any to start learning from. Or pick up a text on DSP, like anything by Oppenheim and Schafer as (co)authors or the book by Proakis and Manolakis. If you want an introductory book, try Lyons' "Understanding DSP", although this is more an introduction explanation that paints the big picture, rather than a detailed ecplanation of technical concepts. I am not sure how helpful it will be if you want hands-on results quickly. Again, talk with your colleagues. There will be lots of tiny details you will need to wrap your mind around to get things right. Rune
Reply by Thomas January 21, 20092009-01-21
I've expected here some more useful information, like read this book, look 
at this page or some example code.

Things like

    "Working in frequency domain is no big deal",

    "Learning something new will take time" or

    "there are a couple of key details you need to know to get things right"

are not really helpful to start with the topic.



I operate offline but if I apply one of the filter in our signal processing 
software (NI Diadem) than it takes some seconds and in the existing Fortran 
code some minutes (I have to apply this to many files and channels, so it is 
hard to life with a factor in computation time of about 60), therefore I 
look for a better filter code.

BTW It don't have to be the best.



Thomas



"Rune Allnor" <allnor@tele.ntnu.no> schrieb im Newsbeitrag 
news:665ea1b7-48c9-4108-b9f0-181a0eb689b9@w1g2000prk.googlegroups.com...
> On 19 Jan, 12:31, "Thomas" <thomas.we...@logomotive.eu> wrote: >> "Rune Allnor" <all...@tele.ntnu.no> schrieb im >> Newsbeitragnews:65688a96-a7f1-4bb8-b7c7-e9f80783fb89@v5g2000prm.googlegroups.com... >> On 19 Jan, 11:20, "Thomas" <thomas.we...@logomotive.eu> wrote: >> >> > Hi, >> >> > I've to apply some special shaped filter to long signals, the filter >> > coefficients I know already but the existing filter routine is probably >> > to >> > slow. I want to update the filter routine without becoming a signal >> > processing expert ;-) > ... >> The filters are band limited comfort filter till order 8 and the signals >> are >> typical measured accelerations from riding tests. >> >> I work with Intel Visual Fortran 10. and was able to survive till now >> with >> the "Numerical Recipes"-Library only but I hope I may download something. > > These data will be processed off-line? If so, processing time will > almost certainly not be an issue. In that case use what method you > are comfortable with. Learning something new will take time and > effort (and pain) anyway, so just use what you already know - it will > work. > >> The time domain operation would be my favourite version. > > Then stick with it. > > Working in frequency domain is no big deal once you know your way > around the concepts. However, there are a couple of key details you > need to know to get things right, and also a few numerical quirks > that can become awkward if you are not aware what you are doing. > Based on several posts here in the past, it seems that people who > work on their own have a hard time getting such details right, and > tend to become confused. No need to go there if you don't have to. > > Rune
Reply by Rune Allnor January 19, 20092009-01-19
On 19 Jan, 12:31, "Thomas" <thomas.we...@logomotive.eu> wrote:
> "Rune Allnor" <all...@tele.ntnu.no> schrieb im Newsbeitragnews:65688a96-a7f1-4bb8-b7c7-e9f80783fb89@v5g2000prm.googlegroups.com... > On 19 Jan, 11:20, "Thomas" <thomas.we...@logomotive.eu> wrote: > > > Hi, > > > I've to apply some special shaped filter to long signals, the filter > > coefficients I know already but the existing filter routine is probably to > > slow. I want to update the filter routine without becoming a signal > > processing expert ;-)
...
> The filters are band limited comfort filter till order 8 and the signals are > typical measured accelerations from riding tests. > > I work with Intel Visual Fortran 10. and was able to survive till now with > the "Numerical Recipes"-Library only but I hope I may download something.
These data will be processed off-line? If so, processing time will almost certainly not be an issue. In that case use what method you are comfortable with. Learning something new will take time and effort (and pain) anyway, so just use what you already know - it will work.
> The time domain operation would be my favourite version.
Then stick with it. Working in frequency domain is no big deal once you know your way around the concepts. However, there are a couple of key details you need to know to get things right, and also a few numerical quirks that can become awkward if you are not aware what you are doing. Based on several posts here in the past, it seems that people who work on their own have a hard time getting such details right, and tend to become confused. No need to go there if you don't have to. Rune