# Add chirp to time domain pulse

I have a Gaussian pulse $$E(\omega)$$ centered at 500 PHz with a bandwidth of 10 PHz. The pulse is transform limited initially.

I am trying to add chirp, through multiplying the FT of $$E(\omega)$$ by a quadratic phase factor given as

$$E(t)_{chirp} = E(t)\cdot\text{exp}(-ibt^2)$$

however, the exponential factor never seems to change the electric field?

def gaussian(x, mu, sig):
return np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.)))

# frequency in Hz
freq = np.linspace(0,2000*10**12,1000)
dt = np.diff(freq)[0]
wl = (3*10**8)/freq;
band = len(freq)

I_stokes = gaussian(freq,500*10**12,10*10**12)

phase = freq*0

E_stokes = np.sqrt(I_stokes)*np.exp(+1j*phase)
#FT:
n = len(E_stokes)
Fs = np.diff(freq)[0]
t = np.arange(0, n/Fs, 1/Fs)

Et = fftshift(ifft(E_stokes))

Phi = 1

E2t = Et**2

E2 = (fft(E2t))

I = E2t*np.exp((-1j*((t)**2)/2*Phi))

p = np.angle(I)

#plt.plot(p)

plt.plot(np.unwrap(np.angle(I))- np.unwrap(np.angle(E2t)))


Here is the difference between the angles of I and E2t:

Couple points that may be helpful:

The command plt.plot(I) is only plotting the real component.

The command phase = freq**0 results in a constant phase of 1 for all values of frequency, rather than $$bt^2$$ (was this intended)?

That said, E_stokes ends up being E_stokes = np.sqrt(I_stokes)

To see the difference in the waveforms, I recommend looking at the unwrapped phase as determined between I and E2t using

np.unwrap(np.angle(I))- np.unwrap(np.angle(E2t))

Below is the plot of the result in doing that.

Other differences may be in the specific fft ifft used (unknown origin in the OP's code as the imports are not listed). Below is the complete code I used to generate the plot, resulting in the expected quadrature phase at twice the frequency (due to the product done). As far as I am aware, I don't think I changed anything other than fft references and plotting the unwrapped phase:

import numpy as np
import scipy.fftpack as fft
import matplotlib.pyplot as plt
def gaussian(x, mu, sig):
return np.exp(-np.power(x - mu, 2.) / (2 * np.power(sig, 2.)))

# frequency in Hz
freq = np.linspace(0,2000*10**12,1000)
dt = np.diff(freq)[0]
wl = (3*10**8)/freq;
band = len(freq)
I_stokes = gaussian(freq,500*10**12,10*10**12)
phase = freq**0
E_stokes = np.sqrt(I_stokes)*np.exp(+1j*phase)
#FT:
n = len(E_stokes)
Fs = np.diff(freq)[0]
t = np.arange(0, n/Fs, 1/Fs)
Et = fft.fftshift(fft.ifft(E_stokes))
Phi = 1
E2t = Et**2
E2 = (fft.fft(E2t))
I = E2t*np.exp((-1j*((t*10**10)**2)/2*Phi))
plt.figure()
plt.plot(freq, np.unwrap(np.angle(I))- np.unwrap(np.angle(E2t)))

• Thanks for the feedback, so yes I should be plotting the magnitude squared. the constant phase in the pulse creation was intended, as i wanted to add the chirp in the time domain after the fact.
– Ryan
Jan 13, 2022 at 19:49
• I see. Well if you plot the absolute value then you will just get the envelope (amplitude) of your analytic signal. Compare the two (I and E2t) using np.angle() for each and you will see the impact of your added chirp. The instantaneous frequency is the derivative of your time domain phase. Jan 13, 2022 at 19:55
• I am doing that but the phase, or absolute intensity is not changing at all as per the original question.
– Ryan
Jan 13, 2022 at 20:26
• It changed for me: specifically comparing plt.plot(np.angle(I)) and plt.plot(np.angle(E2t)). Am I still missing what you are trying to do? Jan 13, 2022 at 20:31
• Also should E2t = Et**2 be E2t = np.abs(Et)**2 ? Or similarly Et*np.conj(Et)? If you do that (which is how I would typically compute a power qty if that is the purpose of that line) it will be a real result and the phase will only be the phase as given by the multiplied exponential. If not, what is the purpose of squaring Et (which would cause its phase to double)? Jan 13, 2022 at 20:33