Calculating the first derivative of an image using DFT

Question

I need to calculate the first derivative of a greyscale image (a 2D array) using a DFT function I built (which works). Unfortunately, the results don't seem to be correct.

In the fourier domain, the derivative d/dx is given as F(u,v)* 2*pi*i/N * u, where u is the x-axis transformed, N is the size of one of the matrix's dimensions, v being the other one.

Attached is the code. What bothers me is that I'm not getting the same results as I would differentiating by convolution with (1,-1) or (1,-1) as a column vector.

def derivative(fourier_signal):
"""
Derivative in fourier domain is multiplying by u or v, and 2pi*i/N
:param fourier_signal:
:return:
"""
N = np.shape(fourier_signal)[ZERO]
M = np.shape(fourier_signal)[ONE]
u = np.arange(N)
v = np.arange(M)
du = fourier_signal * (u*TWO_PI*1j)/N
dv = fourier_signal * (v*TWO_PI*1j)/M
return du, dv


def fourier_der(im):
    # Calculate DFT2
    dft_image = DFT2(im)
    # Function that Multiply by rows by u, columns by y
    du, dv = derivative(dft_image)
    shifted_du, shifted_dv = np.fft.fftshift(du), np.fft.fftshift(dv)
    dx, dy = IDFT2(shifted_du), IDFT2(shifted_dv)

I'm not looking for easy answers on how to do it, but rather a direction to as to why my output is incorrect.

Answer 1

They are not the same. Using 1D notation,

the discrete-time (backward) first difference is $x[n] - x[n-1]$ whose frequency domain DTFT equivalent is

$$ x[n]-x[n-1] \leftrightarrow X(e^{j\omega}) - e^{-j \omega} X(e^{j\omega}) =X(e^{j\omega})(1- e^{-j \omega}) $$

which becomes $$ x[n]-x[n-1] \longleftrightarrow X[k](1 - e^{-j \frac{2 \pi}{N} k})$$ using the DFT to implement it.

The FIR impulse response of the discrete-time system that implements the first difference is therefore,

$$ h[n] = \delta[n] - \delta[n-1]$$

The first derivative of a continuous-variable function $x(t)$ is $x'( t)$ and in CTFT domain it becomes : $$ x'(t) \longleftrightarrow j\Omega X(\Omega) $$

where the analog system frequency response is

$$H_d(\Omega) = j \Omega $$

which is not implementable in digital form, but a bandlimited approximation to it is attained under a sampling period of $T$ that yields an equivalent discrete-time frequency response of a (bandlimited) differentiator as

$$ H_d(e^{j \omega}) = j \frac{\omega}{T} ~~~, ~~~~ |\omega| < \pi$$

The associated discrete-time (IIR) impulse response is $$ h_d[n] = \text{IDTFT} \{ H_d(e^{j \omega}) \} = \text{IDTFT} \{ j \frac{\omega}{T} \} $$

Practically you would truncate and window $h_d[n]$ before using.

So they are not the same.