bugGNU Octave - Bugs: bug #47038, [octave forge] (signal) Ringing,...

 
 

bug #47038: [octave forge] (signal) Ringing, possibly unwanted, with filtfilt

Submitted by:  Alexander Klein <matalex>
Submitted on:  Mon 01 Feb 2016 12:23:16 PM UTC  
 
Category: Octave Forge PackageSeverity: 3 - Normal
Priority: 5 - NormalItem Group: Inaccurate Result
Status: Need InfoAssigned to: Mike Miller <mtmiller>
Originator Name: Open/Closed: Open
Release: otherOperating System: Any

Add a New Comment(Rich Markup)
   

You are not logged in

Please log in, so followups can be emailed to you.

 

Wed 10 Feb 2016 09:49:30 AM UTC, comment #2:

I'd say this definitely needs to be looked at. On my system I had a filtfilt.m version in my path from probably 5+ years ago and it is extremely simple, i.e., essentially calling filter in both directions after tagging some zeros at the end of the input signal to capture the filter transient and use that as the start of the reverse filter. Here is the comment from that old version:

## TODO: In Matlab filtfilt `reduces filter startup transients by carefully
## TODO: choosing initial conditions, and by prepending onto the input
## TODO: sequence a short, reflected piece of the input sequence'.
## TODO: Once filtic is written, use that here.

There is a filtic() function, but I don't see it called anywhere in the latest filtfilt.m. Not sure it applies.

I experimented a bit, and I found that the main contributing factor is this line:

If I change that line to the following:

the transient at the front end goes away. I guess that is no surprise, but why the 2 of 2*x(1,c)? Consider if filter length is 2

x(1,1) = 3
x(2,1) = 1.5
x(3,1) = -0.5
...

Then

v(1) = 3 + 3 - (-0.5) = 6.5
v(2) = 3 + 3 - (1.5) = 4.5
v(3) = 3
v(4) = 1.5
v(5) = -0.5

The point is that there is always some additional bias in this padding meant to prime the filter.

I'm not following the logic of using si*v(1) for the initial state either. I took out the padding at the front and still see a small transient that isn't at the level of the random signal. Guess I'd have to read the paper referenced within, but I don't see how this initialization scheme makes sense.

Further note, your example runs very slow because filtfilt.m is using a for-loop whereas filter might not be. It is clear that the F = flipud ( filter ( b, a, flipud ( filter ( b, a, M ) ) ) ); line is faster than filtfilt.

Dan Sebald <sebald>
Mon 01 Feb 2016 07:24:57 PM UTC, comment #1:

I observe what you describe with your given example. Can someone test this example in Matlab to see what the comparison looks like there?

Mike Miller <mtmiller>
Project AdministratorIn charge of this item.
Mon 01 Feb 2016 12:23:16 PM UTC, original submission:

Filtfilt may introduce strange artefacts into the spectra and time series of filtered signals.

This seems to result from the code computing the initial filter state, and manifests itself as a ringing artefact at both ends of the filtered signal, which does not appear when simply applying "filter" in both directions.

See the attached file for sample code plotting a short demo.

While being far from catastrophic for most applications, the effect is certainly strange enough to cause a certain amount of headache to unprepeared users, as it does not seem to be documented anywhere, and it might break assumptions regarding expectation values towards both ends of time series.

Is this really how we want in this function to work?

Alexander Klein <matalex>

 

(Note: upload size limit is set to 16384 kB, after insertion of the required escape characters.)

Attach File(s):
   
   
Comment:
   

Attached Files
file #36229:  source_me added by matalex (777B - application/octet-stream)

 

Depends on the following items: None found

Items that depend on this one: None found

 

Carbon-Copy List
  • -unavailable- added by jwe (Updated the item)
  • -unavailable- added by sebald (Posted a comment)
  • -unavailable- added by mtmiller (Posted a comment)
  • -unavailable- added by matalex (Submitted the item)
  •  

    Do you think this task is very important?
    If so, you can click here to add your encouragement to it.
    This task has 0 encouragements so far.

    Only project members can vote.

     

    Please enter the title of George Orwell's famous dystopian book (it's a date):

     

     

    Follow 7 latest changes.

    Date Changed By Updated Field Previous Value => Replaced By
    Sun 13 Aug 2017 12:58:27 PM UTCjweSummarysignal package: Ringing, possibly unwanted, with filtfilt=>[octave forge] (signal) Ringing, possibly unwanted, with filtfilt
    Tue 29 Mar 2016 05:59:35 PM UTCmtmillerAssigned toNone=>mtmiller
    Mon 01 Feb 2016 07:24:57 PM UTCmtmillerStatusNone=>Need Info
      Release3.8.2=>other
      Operating SystemBSD=>Any
      SummaryRinging, possibly unwanted, with filtfilt=>signal package: Ringing, possibly unwanted, with filtfilt
    Mon 01 Feb 2016 12:23:16 PM UTCmatalexAttached File-=>Added source_me, #36229

    Back to the top


    Powered by Savane 3.1-cleanup1