0
$\begingroup$

I'm trying to illustrate the principle of constructing a periodic signal using an elementary pattern. This is the code I use (python 2.7 + numpy + scipy.signal) :

import numpy as np 
import matplotlib.pyplot as plt
import scipy.fftpack as sfft
import scipy.signal as spsig

plt.close('all')

npts = 2**14 #number of points for both time-domain signals and FFT's
tsim = 100e-9 #ending time (in s)
fs = npts/tsim #sampling frequency

t = np.linspace(0,tsim,npts) #time bins vector
t1 = 30e-12 #pulse rise time
t2 = 80e-12 #pulse fall time
delta_t = 200e-12 #pulse width (including t1 and t2)
fsig = 1e9 #pulse repetition frequency
A = 1. #pulse amplitude

# here I construct the elementary pulse signal, it works fine
sig_elem = np.zeros(t.size)
sig_elem[(t>=0)&(t<t1)] = (A/t1)*t[(t>=0)&(t<t1)]
sig_elem[(t>=t1)&(t<delta_t-t2)] = A 
sig_elem[(t>=delta_t-t2)&(t<delta_t)] = (-A/t2)*t[(t>=delta_t-t2)&(t<delta_t)] + (A/t2)*delta_t

periode = 1/fsig #this is the pulse repetition period
comb = np.reshape(np.array([np.round(t*1e9, 2)%1==0]), (npts,)) #this is where I construct the comb, I use numerical tricks to make it work 
sig = spsig.fftconvolve(sig_elem, comb, 'full') #and the convolution itself, which should give a periodization of my elementary signal

##### FIGURE ####
plt.figure(1)
plt.plot(t*1e9, sig_elem, 'r',lw = 0.8)
plt.xlabel('t (ns)')
plt.xlim(-10,100)
plt.ylim(-1,5)
plt.grid(True)

# single sided FFT of the elementary signal. I restrict the plot frequency at 50 GHz since the sampling frequency is insanely high
sig_elem_fft = np.abs(sfft.fftshift(sfft.fft(sig_elem, npts)))
plt.figure(2)
plt.plot(sfft.fftshift(sfft.fftfreq(npts, 1/fs))*1e-9, sig_elem_fft, 'b',lw = 0.8)
plt.xlabel('f (GHz)')
plt.xlim(0,50)
plt.grid(True)

#plot of the dirac "comb". It seems good
plt.figure(3)
plt.scatter(t*1e9, comb)
plt.grid(True)
plt.xlim(-10,110)
plt.ylim(-0.5,1.5)

# plot of the convolve result, this is where I don't understand what's happening
plt.figure(4)
plt.plot(sig, 'r',lw = 0.8)
plt.xlabel('t (ns)')
#plt.xlim(-10,100)
#plt.ylim(-1,5)
plt.grid(True)

What I'm doing is defining a trapezoidal pulse :

$$ \begin{cases} s(t) = \frac{At}{t_1} & \text{ if }0 \leq t < t_1\\ A & \text{ if } t_1\leq t < \Delta t-t_2\\ \frac{A}{t_2}\left(\Delta t - t\right) & \text{ if } \Delta t-t_2 \leq t < \Delta t \end{cases} $$

$t_1$, $t_2$ and $\Delta t$ being respectively rise time, fall time and total pulse width. A is an arbitrary amplitude. I implement it on $2^{14}$ (16384) points, making the total duration equal to 100 ns. Which makes the sampling frequency around 164 GHz.

This is what I get, the rest of the signal being equal to zero. enter image description here

After creating this elementary pulse, I create a "dirac comb", using the same number of points, and using a 3 GHz repetition frequency. Then, both signals are convolved using an FFT method (meaning both signals FFT's are calculated, multiplied element-wise and the result IFFT is computed afterwards).

What I should get is a periodic signal consisting of $s(t)$ defined above being repeated each 1 ns. However, I get a strange signal, where the amplitude is sometimes A, sometimes 2*A. The period is however consistent.

![enter image description here

$\endgroup$

1 Answer 1

1
$\begingroup$

Your "numeric trick" is not working. After running above code, have a look at :

plt.figure(); plt.plot(comb); plt.xlim((324,336)); plt.ylim((-0.1, 1.1))

enter image description here

Your "Dirac" train has peaks with width 2, hence you have the double-height amplitude pulses.

You can fix this by e.g. doing

comb2 = np.zeros(npts)
step = 1e-9*fs  # distance between peaks in terms of samples (floating point)
peakPos = step*np.arange(int(npts / step))  # integer position of the peaks
comb2[np.round(peakPos).astype(int)] = 1  # set the peaks
plt.figure()
plt.plot(comb)
plt.plot(comb2); plt.ylim((-0.1, 1.1)); plt.xlim((324, 336))
$\endgroup$

Your Answer

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

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