11
$\begingroup$

Code: https://github.com/echometerain/sound-derivative/blob/master/sd.py

I wrote a program which would subtract the current audio frame by the former one. Why does this act like a high-pass filter?

Spectrogram view

Waveform view (slightly offset)

$\endgroup$
0

5 Answers 5

13
$\begingroup$

Yes both the time derivative $d/dt$ and the discrete time difference $x[n]-x[n-1]$ are indeed "high pass filters": consider the extreme case of the lowest frequency, which is DC or a constant that doesn’t change with time. The derivative of a constant is zero (it cannot pass through) and we see simply that a moving difference of two values that are the same will be zero. As the frequency approaches zero, the maginitude will reduce by a factor of 20 dB/decade consistent with a first order filter response. We also see from the Laplace Transform of a time domain derivative of a function results in a multiplication by $s$ in Laplace. The frequency response for a transfer function $H(s)$ is found by setting $s=j\omega$, and here we again see how the result will be a high pass.

A difference equation $y[n] = x[n] - x[n-1]$ is one of many discrete time approximations or mappings of the operation $d/dt$ (continuous time-domain derivative), in this case using the backward Euler method. For those familiar with the s and z planes, we see how a zero at the origin in s maps to $z=1$ in the z plane, and how both of those locations are consistent with "DC" on the frequency axis! Over the frequencies of interest within the first Nyquist zone (out to Nyquist or $f_s/2$ where $f_s$ is the sampling rate, the single sample difference is a high pass filter. This is the complement of a single sample moving summation of two samples which is a low pass filter. We may also note that the combination of a single sample summation and single sample difference is a 2 point DFT (each resulting sample in the two point DFT is the result of a complementary low pass and high pass function). We also note that beyond Nyquist the response will return to zero in magnitude; which is exactly what would occur for any discrete-time "high-pass filter" with real coefficients (from periodicity in the frequency domain for discrete-time signals).

Similarly, a time domain integration which is $1/s$ in Laplace is a low pass filter. An accumulator $y[n] = x[n]+ y[n-1]$ is a discrete time approximation or mapping of a time domain integration, and we see here how a pole at the origin in s maps to $z=1$ in the z plane.

Below is a log log plot of the frequency magnitude response comparing the continuous time derivative in orange to the discrete time difference filter in blue out to the Nyquist frequency of $f_s/2$. (The frequency axis was normalized to the sampling rate in radians per second as $\omega = 2\pi f/f_s$ with $f$ in Hz, and the time domain derivative here is $ \omega$). From this we can clearly observe the high pass function as a first order high-pass with a slope of 20 dB/decade, and we see for this simple mapping how the more we oversample a waveform, the closer the difference between adjacent samples is to a true derivative in both magnitude and phase (which is consistent with the derivative as that difference in the limit as the difference approaches zero!). The closer we get to Nyquist the larger the deviation, with effects that can be consistently explained with frequency aliasing.

frequency response

As an interesting comparison, below is the result of mapping the differentiator using the Bilinear Transform, otherwise called Tustin's method. This results in an exact match in the phase response and distortion in the magnitude response as this particular mapping compresses (warps) the frequency axis spanning from DC to infinity in continuous-time to the range from DC to Nyquist (fs/2) in discrete time.

Tustin's method results in the difference equation given by the following (the factor of 2 is to match the scaling by two for the result $x[n]-x[n-1]$):

$$y[n] = 2(x[n] - x[n-1] - y[n-1])$$

So in addition to what the OP has done with a simple difference of the input, we also subtract the previous result.

The magnitude and phase response for this result is shown below:

Tustin's Method

Note that the continuous time derivative $d/dt$ can NOT be implemented as an exact match in both amplitude and phase in discrete form but only approximated such as described in this post. For other example approximations of the derivative for discrete time implementation see these articles on dsprelated.com by Rick Lyons and Neil Robertson:

https://www.dsprelated.com/showarticle/1447.php#tabs1-newcomment

https://www.dsprelated.com/showarticle/814.php#tabs1-newcomment

$\endgroup$
12
$\begingroup$

$\frac{d}{dt}\sin(\omega t) = \omega\cos(\omega t)$ as a consequence of the chain rule. $\sin$ and $\cos$ are the same apart from the phase shift, so the derivative is a filter with a response that's directly proportional to frequency.

What you're doing is a finite difference, which is only a discrete approximation to the derivative, but what holds true for the derivative is approximately true for the finite difference. The response starts at zero and increases up to the Nyquist frequency, although it falls below linear with increasing frequency as shown in this answer.

$\endgroup$
10
  • 3
    $\begingroup$ @OverLordGoldDragon derivative is linear, meaning deriv(a+b) = deriv(a)+deriv(b), meaning that deriv(sum of sine) = deriv(sine 1)+deriv(sine 2)+... and since we can decompose any signal into a sum of sines, yes, the derivative of a general signal IS a sum of sine derivatives. $\endgroup$
    – user253751
    Jun 13, 2022 at 10:12
  • $\begingroup$ @OverLordGoldDragon check your linked graphic for x(n)-x(n+1). It is a high-pass filter for frequencies below half the sampling rate! ... which we know are the only ones that exist in the file. $\endgroup$
    – user253751
    Jun 13, 2022 at 10:12
  • $\begingroup$ @user253751 It's not that simple. My answer shows the exact frequency response, and the derivative is not directly proportional to frequency. $\endgroup$ Jun 13, 2022 at 10:15
  • 1
    $\begingroup$ @OverLordGoldDragon the derivative is directly proportional to frequency, but this isn't a derivative. $\endgroup$
    – user253751
    Jun 13, 2022 at 10:24
  • 3
    $\begingroup$ It seemed to me that the OP knows that the finite difference is only an approximation to the derivative, and that linearity applies, and they asked about the derivative, so I wrote an answer that gets to the heart of that. But I added a little bit more detail for clarity. $\endgroup$
    – hobbs
    Jun 13, 2022 at 13:58
6
$\begingroup$

@OverLordGoldDragon Your x(n)−x(n−1) expression is called a "first-difference discrete differentiator". Its frequency magnitude response is the blue dashed curve below:

enter image description here

$\endgroup$
2
  • $\begingroup$ Yes, it is what OP is doing. I did forget that my other answer normalized the plot for the animation; the peak is indeed at 2. $\endgroup$ Jun 13, 2022 at 13:27
  • 1
    $\begingroup$ @OverLOrdGoldDragon Yes, the peak magnitude response is two, which is also the area underneath a half cycle of a sine wave. I find that fact to be interesting. $\endgroup$ Jun 14, 2022 at 9:30
4
$\begingroup$

You're doing $x(n) - x(n - 1)$, the finite difference / "discrete derivative", not to be confused with $\frac{d}{dt}x(t)$. The frequency response is

where $H(\omega) = 1 - e^{-j 1\omega}$ following the time-shift property $x(t - t_0) \Leftrightarrow e^{-j t_0 \omega} X(\omega)$.

Note this is a very poor high-pass filter. It peaks at Nyquist and vanishes at DC, but still passes a significant amount of low frequencies. In case of $\frac{d}{dt}$ we get, within discretization error, a straight line from DC to Nyquist, which is also a poor high-pass (but a better one for some purposes, and can be implemented discrete).

$d/dt$ vs $\text{diff}$

Conflating continuous and discrete derivatives can be a grave error, especially for OP's context (STFT). It's not only about magnitude response, and the derivative is a terrible approximation of finite difference's phase response. To call either an approximation of the other in the first place involves important caveats that none of the other answers discussed, that can create unnecessary confusion at especially implementation stage.

Explanation

$x(n) - x(n - 1)$ is simply a transformation of $x$. Whether it's "lowpass" or "highpass" is a frequency-domain description, so we describe the transform directly in frequency domain.

Here we're simply adding time-shifted and sign-flipped version of $x$ to itself. In frequency domain, this is described by the time-shift property (shown earlier). So we have:

$$ \begin{align} x(n) & \Leftrightarrow X(\omega) \\ x(n) - x(n - 1) & \Leftrightarrow X(\omega) - e^{-j 1 \omega}X(\omega) \\ \end{align} $$

To see the effect this has on any $x$, independent of $x$, we simply divide it out - which gives the frequency response:

$$ \begin{align} H(\omega) &= Y(\omega) / X(\omega) \\ &= \left( X(\omega) \cdot (1 - e^{-j1\omega})\right) / X(\omega) \\ &= \boxed{1 - e^{-j1\omega}} \end{align} $$

and thus, for any $x$ with frequencies $X(\omega)$, the effect of the transform is $X(\omega) \rightarrow H(\omega)X(\omega)$.

For the more general case of $x(n) - x(n - 1) + x(n - 2) - ...$, and code, see this answer.

$\endgroup$
3
  • 1
    $\begingroup$ So this could be thought of a comb filter (like in flanging), with only first half of the first tooth effectively used? $\endgroup$ Jun 14, 2022 at 3:47
  • 2
    $\begingroup$ @user 1079505 Yes yes. It's a comb filter whose delay line contains only one delay element. $\endgroup$ Jun 14, 2022 at 9:32
  • $\begingroup$ @user27621 The question only describes finite difference. That $d/dt$ approximates it in some sense is only a convenience: a continuous-time result does not warrant a similar discrete-time result. If OP processes STFT with diff where $d/dt$ is needed, or vice versa, the result can be completely wrong. -- @ implemented discrete, yes depending on application; that's a separate topic. $\endgroup$ Jun 15, 2022 at 19:27
2
$\begingroup$

In addition to the excellent existing answers we can also consider the continuous domain: If we have a signal as a fourier series (here in sine-cosine form)

$$f(t) = \frac{a_0}{2} + \sum_{n=1}^\infty \left( a_n \cos(nt) + b_n \sin(n t)\right)$$

And compare it with the derivative

$$f'(t) = \sum_{n=1}^\infty n\left( b_n \cos(nt) -a_n \sin(n t)\right)$$

we can immediately observe a few things: First of all the constant offset $a_0$ disappeared which is already consistent with your observation of having a high-pass filter. Then we see that for each frequency $n$ we get the expected 90° phase shift but in addition we get a factor $n$.

This means that higher frequencies get a greater weight than lower frequencies, which (after normalization) means we do indeed get a high-pass filter.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.