This question was inspired by this question where it is necessary to numerically compute the Fourier transform of a Gaussian optical pulse with a Gaussian chirp function.
$$E(t)=e^{-t^2} \cos(50 t - e^{-2 t^2} 8 π)$$.
I was thinking that the Fourier transform, while a complete picture of the pulse, wasn't as useful as I would want it to be. I think a spectrogram, showing the time and frequency information of the pulse, would be helpful as well. The only spectrogram I am familiar with is the Wigner function (useful reading here as well)
$$ W_x(t,f) = \int_{-\infty}^\infty E\left(t+\frac{\tau}{2}\right) E^*\left(t-\frac{\tau}{2}\right) e^{-i 2\pi\tau\,f} \, d\tau$$
What is the best way to implement this in Mathematica?
My first attempt was to do it directly,
ω0 = 50;
pulse[t_] := Exp[-t^2] Exp[-I (ω0 t - 8 π Exp[-2 t^2] )];
pulsecc[t_] := Exp[-t^2] Exp[I (ω0 t - 8 π Exp[-2 t^2] )];
wignerfunc[t_, w_] :=
NIntegrate[(pulse[t + τ/2] pulsecc[
t - τ/
2]) Exp[-I w τ], {τ, -∞, ∞}];
Plot[wignerfunc[0, w], {w, 30, 80}]
But I get convergence errors,
NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small. >>
NIntegrate::ncvb: NIntegrate failed to converge to prescribed accuracy after 9 recursive bisections in τ near {τ} = {-0.925724}. NIntegrate obtained 2.42861*10^-17+8.97719*10^-17 I and 3.6945011190775015`*^-15 for the integral and error estimates. >>
Is there an option to give NIntegrate
that makes this converge?
This paper describes a discrete Wigner transform, that I imagine would take a 2D array (t
along one dimension, τ
along the other) and spit out the spectrogram, but signal processing always seems opaque to me.
ContinuousWaveletTransform
. $\endgroup$Spectrogram
would compare to the result below. $\endgroup$