1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
| #!/usr/env python
# -*- coding: utf-8 -*-
##############################################################################
import imgen
imgen.setup()
##############################################################################
import numpy as np
######
# Paramètres du problème
######
t_min = 0.
t_max = 20.
c = 0.3 # Gm/s, vitesse de la lumière
wc1 = 20. # GHz, pulsation de coupure 1
WC = [50, 35, 20]
# Spectre du signal
w0 = 80. # GHz, pulsation centrale
dw = 10. # largeur de la gaussienne en pulsation
Nw = 10000
w_min = max(0, w0 - 10*dw)
w_max = w0 + 10*dw
W = np.linspace(w_min, w_max, Nw) # domaine de pulsations
def A(w):
return np.exp(- ((w - w0)**2) / (2 * dw**2)) # spectre gaussien
# FenĂŞtre spatiale
x_min = -0. # m
x_max = c * t_max # m
Nx = 1000
X = np.linspace(x_min, x_max, Nx)
######
# Calcul de la relaiton de dispersion
######
def K(w, w_p):
k1 = np.sqrt((0.j + w**2 - w_p**2)) / c
return np.conjugate(k1) # On doit prendre les parties imaginaires négatives
######
# Calcul du signal S(X, t)
######
def S(t, wc1):
p_max = int(w_max/wc1 + 1)
S = np.zeros_like(X, dtype=np.complex)
for i in range(len(W)):
for p in range(1, p_max+1):
w = W[i]
k = K(w, p*wc1)
a = A(w)
S += a * np.exp(1.j * (w*t - k*X))
return S
######
# Calcul de c*t
######
def C(t):
return c*t
##############################################################################
import matplotlib.pyplot as plt
# Affichage
T = np.linspace(t_min, t_max, 5)
fig, axes = plt.subplots(len(T) + 1, len(WC))
for i, ax in enumerate(axes[0]):
wc = WC[i]
for n in range(len(T)):
axes[0,i].set_title(u"$\omega_c = %i$ GHz" % (wc))
t = T[n]
# Onde
axes[n, i].plot(X, S(t, wc).real)
# Vitesse c
axes[n, i].plot(C(t), 0, marker='o', color='r')
# Labels
axes[n, i].set_yticks([])
axes[n, 0].set_ylabel("t = %3.f" % t)
# Limites
axes[n, i].set_xlim(x_min, x_max)
# Labels
axes[-2, i].set_xlabel(u"x / Gm")
axes[-1, i].set_xlabel(u"$\omega$ / GHz")
# Relation de dispersion et TF du signal
axes[-1, i].plot(W, K(W, wc).real/(np.max(K(W, 1*wc).real)), color='orange')
axes[-1, i].plot(W, A(W).real, color='olive')
for p in range(1, int(w_max/wc + 1)+1):
axes[-1, i].plot(W, K(W, p*wc).real/(np.max(K(W, p*wc).real)), color='orange')
axes[-1, i].legend([u"Relation de dispersion", u"TF du signal"])
axes[-1, i].autoscale(enable=True)
##############################################################################
# Save
imgen.done(__file__)
|