obspy icon indicating copy to clipboard operation
obspy copied to clipboard

Zerophase filtering not zerophase for numpy matrices

Open ftilmann opened this issue 3 years ago • 2 comments

With zerophase=True the low pass filter is applied twice in opposite directions to give a zerophase output, as in this example with two delta spikes showing the expected behaviour:

import numpy as np
import matplotlib.pyplot as plt
from obspy.signal.filter import lowpass

ts=np.zeros([100])
ts[50]=1
ts[25]=0.5

ts_1=lowpass(ts, 0.05,1)      # single phase filter for comparison
ts_2=lowpass(ts, 0.05,1,zerophase=True)

plt.figure()
plt.plot(ts)
plt.plot(ts_1)
plt.plot(ts_2)

However, when instead a matrix presented, the single phase filter works as expected on the last dimension, but with zerophase=True instead of working in two directions, the zerophase filter is applied twice in the same direction, see this example

md=np.zeros([3,100])
md[:,50]=1
md[:,25]=0.5

md_1=lowpass(md, 0.05,1)
md_2=lowpass(md, 0.05,1,zerophase=True)

plt.figure()
plt.plot(md[0,:])
plt.plot(md_1[0,:])
plt.plot(md_2[0,:])

The following fix seems to work for this little test (although I have not subjected it to formal testing) also for the one-dimensional array. In the definition of lowpass, replace the line

        return sosfilt(sos, firstpass[::-1])[::-1]

with

        return np.flip(sosfilt(sos, np.flip(firstpass)))

Of course equivalent fixes should be implemented for all filter functions. I don't have a good overview and am not familiar with the testing framework, and feel therefore not qualified to make a pull request. I hope someone can implement and test this adaption.

ObsPy version: 1.2.2 (installation via pip)

ftilmann avatar Mar 08 '21 13:03 ftilmann

Hi, thanks for the short report. Probably nobody tried to use the filter functions with matrices up til now. Seems like a nice addition. We could also add the axis kwarg and pass it to the relevant functions.

trichter avatar Mar 10 '21 11:03 trichter

Yes, adding an axis kwarg seems to make perfect sense.

ftilmann avatar Mar 11 '21 14:03 ftilmann

This issue has been mentioned on ObsPy Forum. There might be relevant details there:

https://discourse.obspy.org/t/unexpected-behaviour-of-bandpass-with-zero-phase-option/1692/3

obspy-bot avatar Apr 13 '23 14:04 obspy-bot

Hello, I've recently come across this issue as well using the filter functions on 2-D data (ambient noise cross-correlation time-distance gathers). Just want to confirm @ftilmann 's solution worked for me! :-) See before/after below:

image

savardge avatar Apr 13 '23 14:04 savardge