# DFT of Transient vs. Analytical Solution

I'm trying to compute the Fourier Transform of a lightning impulse waveform. The python code I'm using is below. The problem that I am running into is that the magnitude of the DFT and the magnitude of the analytically derived Fourier Transform are off by a factor equal to 1 / time period, where the time period is the length of time that I am using to model the signal for the DFT.

from math import *
import numpy as np
import matplotlib.pyplot as plt

def Impulse(t, V_0=1.04, alpha=0.0146e6, beta=2.467e6, t_chop=float("inf")):
if t<=t_chop:
return V_0 * (exp(-alpha*t)-exp(-beta*t))
else:
return 0

npImpulse = np.vectorize(Impulse, excluded=['V_0', 'alpha', 'beta', 't_chop'])

# Sampling rate and time vector
total_time = 500e-6 # seconds
n_samples = 2**24
sampling_rate = n_samples/total_time # Hz (samples/s)

t_values = np.linspace(0, total_time, n_samples)

# Time constants
tau1 = 68.2e-6
tau2 = 0.405e-6

y_values = npImpulse(t_values, V_0=1.0, alpha=1/tau1, beta=1/tau2)

# plot with various axes scales
plt.figure(1)
plt.plot(t_values, y_values)
plt.grid(True)
plt.show()

# FFT calculation
sp = np.fft.rfft(y_values) # real FFT, + freq only
freq = np.fft.rfftfreq(y_values.size, t_values[1]-t_values[0])

# Normalization
# This is where things are not right.  Should be dividing by n_samples, not samples/sec
sp = sp / n_samples

# Calculation of Fourier Transform analytically
sp_analytic_real = tau1/(1+(2*pi*freq*tau1)**2)-tau2/(1+(2*pi*freq*tau2)**2)
sp_analytic_imag = -(2*pi*freq*(tau1**2)/(1+(2*pi*freq*tau1)**2)-2*pi*freq*(tau2**2)/(1+(2*pi*freq*tau2)**2))
sp_analytic_mag = np.sqrt(sp_analytic_real**2 + sp_analytic_imag**2)

plt.figure(2)
plt.plot(freq, np.abs(sp), label='FFT')
plt.plot(freq, sp_analytic_mag, label='Analytic')
plt.xscale("log")
plt.yscale("log")
plt.legend()
plt.grid(True)
plt.show()


If I set my calculation time period to 1s, the two curves match. If I divide the magnitude of the DFT return by rfft by the sampling rate (or alternately multiply by the 1 / time period), they match. I can't for the life of me figure out why that should be. Can someone clue me in to the obvious that I am missing?

That's the case. Assuming you have a band-limited continuous-time signal $$x_c(t)$$, then the relation between its (analytically derived) CTFT, $$X_c(j\Omega)$$, and the (numerically computed) DTFT, $$X_d(e^{j\omega})$$, (or actually its DFT $$X_d[k]$$ ) of the samples, $$x_d[n] = x_c(n T_s)$$, of it taken at a period $$T_s$$ will be :
$$X_d(e^{j\omega}) = \frac{1}{T_s} X_c( \frac{\omega}{T_s}) ~~~, |\omega| < \pi$$
In practice you are computing the N-point DFT $$X[k]$$ of $$x_d[n]$$ which is the uniform samples of $$X_d(e^{j \omega})|_{\omega = \frac{2\pi}{N} k}$$, for $$k=0,1,...,N-1$$. And they have the same magnitude relation too.
$$X[k] =X_d(e^{j\frac{2\pi}{N} k}) = \frac{1}{T_s} X_c( \frac{\omega}{T_s}) ~~~, |\omega| < \pi$$
In other words, you will have a weight of $$K = \frac{1}{T_s}$$