10. Exemples🔗

10.1. Mean power (fonction, argparse)🔗

 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
#!/usr/bin/env python3
# coding: utf-8

"""
Exemple de script (shebang, docstring, etc.) permettant une
utilisation en module (`import mean_power`) et en exécutable (`python
mean_power.py -h`);
"""


def mean_power(alist, power=1):
    r"""
    Retourne la racine `power` de la moyenne des éléments de `alist` à
    la puissance `power`:

    .. math:: \mu = (\frac{1}{N}\sum_{i=0}^{N-1} x_i^p)^{1/p}

    `power=1` correspond à la moyenne arithmétique, `power=2` au *Root
    Mean Squared*, etc.

    Exemples:
    >>> mean_power([1, 2, 3])
    2.0
    >>> mean_power([1, 2, 3], power=2)
    2.160246899469287
    """

    # *mean* = (somme valeurs**power / nb valeurs)**(1/power)
    mean = (sum( val ** power for val in alist ) / len(alist)) ** (1 / power)

    return mean


if __name__ == '__main__':

    # start-argparse
    import argparse

    parser = argparse.ArgumentParser()
    parser.add_argument('list', nargs='*', type=float, metavar='nombres',
                        help="Liste de nombres à moyenner")
    parser.add_argument('-i', '--input', nargs='?', type=argparse.FileType('r'),
                        help="Fichier contenant les nombres à moyenner")
    parser.add_argument('-p', '--power', type=float, default=1.,
                        help="'Puissance' de la moyenne (par défaut: %(default)s)")

    args = parser.parse_args()
    # end-argparse

    if args.input:       # Lecture des coordonnées du fichier d'entrée
        # Le fichier a déjà été ouvert en lecture par argparse (type=file)
        try:
            args.list = [float(x) for x in args.input
                         if not x.strip().startswith('#')]
        except ValueError:
            parser.error("Impossible de déchiffrer la ligne "
                         f"{x!r} du fichier {args.input!r}")

    # Vérifie qu'il y a au moins un nombre dans la liste
    if not args.list:
        parser.error("La liste doit contenir au moins un nombre")

    # Calcul
    moyenne = mean_power(args.list, args.power)

    # Affichage du résultat
    print("La moyenne puissance 1/{0} des {1} nombres à la puissance {0}"
          " est {2}.".format(args.power, len(args.list), moyenne))

Source: mean_power.py

10.2. Animal (POO)🔗

  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
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
#!/usr/bin/env python3
# coding: utf-8

"""
Exemple (tragique) de Programmation Orientée Objet.
"""


# Définition d'une classe ==============================

class Animal:
    """
    Un animal, défini par sa masse.
    """

    def __init__(self, masse):
        """
        Initialisation d'un Animal, a priori vivant.

        :param float masse: masse en kg (> 0)
        :raise ValueError: masse non réelle ou négative
        """

        self.estVivant = True

        self.masse = float(masse)
        if self.masse < 0:
            raise ValueError("La masse ne peut pas être négative.")


    def __str__(self):
        """
        Surcharge de la fonction `str()`.

        L'affichage *informel* de l'objet dans l'interpréteur, p.ex. `print(a)`
        sera résolu comme `a.__str__()`

        :return: une chaîne de caractères
        """

        return f"Animal {'vivant' if self.estVivant else 'mort'}, " \
            f"{self.masse:.0f} kg"


    def meurt(self):
        """
        L'animal meurt.
        """

        self.estVivant = False


    def grossit(self, masse):
        """
        L'animal grossit (ou maigrit) d'une certaine masse (valeur algébrique).

        :param float masse: prise (>0) ou perte (<0) de masse.
        :raise ValueError: masse non réelle.
        """

        self.masse += float(masse)


# Définition d'une classe héritée ==============================

class AnimalFeroce(Animal):
    """
    Un animal féroce est un animal qui peut dévorer d'autres animaux.

    La classe-fille hérite des attributs et méthodes de la
    classe-mère, mais peut les surcharger (i.e. en changer la
    définition), ou en ajouter de nouveaux:

    - la méthode `AnimalFeroce.__init__()` dérive directement de
      `Animal.__init__()` (même méthode d'initialisation);
    - `AnimalFeroce.__str__()` surcharge `Animal.__str__()`;
    - `AnimalFeroce.devorer()` est une nouvelle méthode propre à
      `AnimalFeroce`.
    """

    def __str__(self):
        """
        Surcharge de la fonction `str()`.
        """

        return "Animal féroce " \
            f"{'bien vivant' if self.estVivant else 'mais mort'}, " \
            f"{self.masse:.0f} kg"

    def devore(self, other):
        """
        L'animal (self) devore un autre animal (other).

        * Si other est également un animal féroce, il faut que self soit plus
          gros que other pour le dévorer. Sinon, other se défend et self meurt.
        * Si self dévore other, other meurt, self grossit de la masse de other
          (jusqu'à 10% de sa propre masse) et other maigrit d'autant.

        :param Animal other: animal à dévorer
        :return: prise de masse (0 si self meurt)
        """

        if isinstance(other, AnimalFeroce) and (other.masse > self.masse):
            # Pas de chance...
            self.meurt()
            prise = 0.
        else:
            other.meurt()             # Other meurt
            prise = min(other.masse, self.masse * 0.1)
            self.grossit(prise)       # Self grossit
            other.grossit(-prise)     # Other maigrit

        return prise


# Définition d'une autre classe héritée ==============================

class AnimalGentil(Animal):
    """
    Un animal gentil est un animal avec un petit nom.

    La classe-fille hérite des attributs et méthodes de la
    classe-mère, mais peut les surcharger (i.e. en changer la
    définition), ou en ajouter de nouveaux:

    - la méthode `AnimalGentil.__init__()` surcharge l'initialisation originale
      `Animal.__init__()`;
    - `AnimalGentil.__str__()` surcharge `Animal.__str__()`;
    """

    def __init__(self, masse, nom='Youki'):
        """
        Initialisation d'un animal gentil, avec son masse et son nom.
        """

        # Initialisation de la classe parente (nécessaire pour assurer
        # l'héritage)
        Animal.__init__(self, masse)

        # Attributs propres à la classe AnimalGentil
        self.nom = nom

    def __str__(self):
        """
        Surcharge de la fonction `str()`.
        """

        return f"{self.nom}, un animal gentil " \
            f"{'bien vivant' if self.estVivant else 'mais mort'}, " \
            f"{self.masse:.0f} kg"

    def meurt(self):
        """
        L'animal gentil meurt, avec un éloge funéraire.
        """

        Animal.meurt(self)
        print(f"Pauvre {self.nom} meurt, paix à son âme...")


if __name__ == '__main__':

    # Exemple d'utilisation des classes définies ci-dessus

    print("Une tragédie en trois actes".center(70, '='))

    print("Acte I: la vache prend 10 kg.".center(70, '-'))
    vache = Animal(500.)        # Instantiation d'un animal de 500 kg
    vache.grossit(10)           # La vache grossit de 10 kg
    print(vache)

    print("Acte II: Dumbo l'éléphant".center(70, '-'))
    elephant = AnimalGentil(1000., "Dumbo")  # Instantiation d'un animal gentil
    print(elephant)

    print("Acte III: le féroce lion".center(70, '-'))
    lion = AnimalFeroce(200)    # Instantiation d'un animal féroce
    print(lion)

    print("Scène tragique: le lion dévore l'éléphant...".center(70, '-'))
    lion.devore(elephant)       # Le lion dévore l'éléphant

    print(elephant)
    print(lion)

Exécution du code:

$ ./animal.py
=====================Une tragédie en trois actes======================
--------------------Acte I: la vache prend 10 kg.---------------------
Animal vivant, 510 kg
----------------------Acte II: Dumbo l'éléphant-----------------------
Dumbo, un animal gentil bien vivant, 1000 kg
-----------------------Acte III: le féroce lion-----------------------
Animal féroce bien vivant, 200 kg
-------------Scène tragique: le lion dévore l'éléphant...-------------
Pauvre Dumbo, paix à son âme...
Dumbo, un animal gentil mais mort, 980 kg
Animal féroce bien vivant, 220 kg

Source: animal.py

10.3. Cercle circonscrit (POO, argparse)🔗

  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
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
#!/usr/bin/env python3
# coding: utf-8

"""
Calcule le cercle circonscrit à 3 points du plan.

Ce script sert d'illustration à plusieurs concepts indépendants:

- un exemple de script (shebang, docstring, etc.) permettant une
  utilisation en module (`import circonscrit`) et en exécutable
  (`python circonscrit.py -h`);
- des exemples de Programmation Orientée Objet: classe `Point` et la
  classe héritière `Vector`;
- un exemple d'utilisation du module `argparse` de la bibliothèque
  standard, permettant la gestion des arguments de la ligne de
  commande;
- l'utilisation de tests unitaires sous la forme de `doctest` (tests
  inclus dans les *docstrings* des éléments à tester).

  Pour exécuter les tests unitaires du module:

  - avec doctest: `python -m doctest -v circonscrit.py`;
  - avec pytests: `py.test --doctest-modules -v circonscrit.py`;
  - avec nose:    `nosetests --with-doctest -v circonscrit.py`.
"""

__author__ = "Yannick Copin <y.copin@ipnl.in2p3.fr>"
__version__ = "Time-stamp: <2014-01-12 22:19 ycopin@lyonovae03.in2p3.fr>"

# Définition d'une classe ==============================


class Point:

    """
    Classe définissant un `Point` du plan, caractérisé par ses
    coordonnées `x`,`y`.
    """

    def __init__(self, x, y):
        """
        Méthode d'instanciation à partir de deux coordonnées réelles.

        >>> Point(0, 1)          # doctest: +ELLIPSIS
        <....Point object at 0x...>
        >>> Point(1 + 3j)
        Traceback (most recent call last):
        ...
        TypeError: __init__() missing 1 required positional argument: 'y'
        """

        try:  # Convertit les coords en `float`
            self.x = float(x)
            self.y = float(y)
        except (ValueError, TypeError):
            raise TypeError(f"Invalid input coordinates ({x}, {y})")

    def __str__(self):
        """
        Surcharge de la fonction `str()`: l'affichage *informel* de l'objet
        dans l'interpréteur, p.ex. `str(self)` sera résolu comme
        `self.__str__()`

        Retourne une chaîne de caractères.

        >>> str(Point(1,2))
        'Point (x=1.0, y=2.0)'
        """

        return f"Point (x={self.x}, y={self.y})"

    def isOrigin(self):
        """
        Teste si le point est à l'origine en testant la nullité des deux
        coordonnées.

        Attention aux éventuelles erreurs d'arrondis: il faut tester
        la nullité à la précision numérique près.

        >>> Point(1,2).isOrigin()
        False
        >>> Point(0,0).isOrigin()
        True
        """

        import sys

        eps = sys.float_info.epsilon  # Le plus petit float non nul

        return ((abs(self.x) <= eps) and (abs(self.y) <= eps))

    def distance(self, other):
        """
        Méthode de calcul de la distance du point (`self`) à un autre point
        (`other`).

        >>> A = Point(1,0); B = Point(1,1); A.distance(B)
        1.0
        """

        # hypot(dx, dy) = sqrt(dx**2 + dy**2)
        return ((self.x - other.x)**2 + (self.y - other.y)**2)**0.5


# Définition du point origine O
O = Point(0, 0)


# Héritage de classe ==============================


class Vector(Point):

    """
    Un `Vector` hérite de `Point` avec des méthodes additionnelles
    (p.ex. la négation d'un vecteur, l'addition de deux vecteurs, ou
    la rotation d'un vecteur).
    """

    def __init__(self, A, B):
        """
        Définit le vecteur `AB` à partir des deux points `A` et `B`.

        >>> Vector(Point(1,0), Point(1,1)) # doctest: +ELLIPSIS
        <....Vector object at 0x...>
        >>> Vector(0, 1)
        Traceback (most recent call last):
        ...
        AttributeError: 'int' object has no attribute 'x'
        """

        # Initialisation de la classe parente
        Point.__init__(self, B.x - A.x, B.y - A.y)

        # Attribut propre à la classe dérivée
        self.sqnorm = self.x ** 2 + self.y ** 2  # Norme du vecteur au carré

    def __str__(self):
        """
        Surcharge de la fonction `str()`: `str(self)` sera résolu comme
        `Vector.__str__(self)` (et non pas comme `Point.__str__(self)`)

        >>> A = Point(1, 0); B = Point(1, 1); str(Vector(A, B))
        'Vector (x=0.0, y=1.0)'
        """

        return f"Vector (x={self.x}, y={self.y})"

    def __add__(self, other):
        """
        Surcharge de l'opérateur binaire `{self} + {other}`: l'instruction
        sera résolue comme `self.__add__(other)`.

        On construit une nouvelle instance de `Vector` à partir des
        coordonnées propres à l'objet `self` et à l'autre opérande
        `other`.

        >>> A = Point(1, 0); B = Point(1, 1)
        >>> str(Vector(A, B) + Vector(B, O))  # = Vector(A, O)
        'Vector (x=-1.0, y=0.0)'
        """

        return Vector(O, Point(self.x + other.x, self.y + other.y))

    def __sub__(self, other):
        """
        Surcharge de l'opérateur binaire `{self} - {other}`: l'instruction
        sera résolue comme `self.__sub__(other)`.

        Attention: ne surcharge pas l'opérateur unaire `-{self}`, géré
        par `__neg__`.

        >>> A = Point(1, 0); B = Point(1, 1)
        >>> str(Vector(A, B) - Vector(A, B))  # Différence
        'Vector (x=0.0, y=0.0)'
        >>> -Vector(A, B)                     # Négation
        Traceback (most recent call last):
        ...
        TypeError: bad operand type for unary -: 'Vector'
        """

        return Vector(O, Point(self.x - other.x, self.y - other.y))

    def __eq__(self, other):
        """
        Surcharge du test d'égalité `{self}=={other}`: l'instruction sera
        résolue comme `self.__eq__(other)`.

        >>> Vector(O, Point(0, 1)) == Vector(Point(1, 0), Point(1, 1))
        True
        """

        # On teste ici la nullité de la différence des 2
        # vecteurs. D'autres tests auraient été possibles -- égalité
        # des coordonnées, nullité de la norme de la différence,
        # etc. -- mais on tire profit de la méthode héritée
        # `Point.isOrigin()` testant la nullité des coordonnées (à la
        # précision numérique près).
        return (self - other).isOrigin()

    def __abs__(self):
        """
        Surcharge la fonction `abs()` pour retourner la norme du vecteur.

        >>> abs(Vector(Point(1, 0), Point(1, 1)))
        1.0
        """

        # On pourrait utiliser sqrt(self.sqnorm), mais c'est pour
        # illustrer l'utilisation de la méthode héritée
        # `Point.distance`...
        return Point.distance(self, O)

    def rotate(self, angle, deg=False):
        """
        Rotation (dans le sens trigonométrique) du vecteur par un `angle`,
        exprimé en radians ou en degrés.

        >>> Vector(Point(1, 0), Point(1, 1)).rotate(90, deg=True) == Vector(O, Point(-1, 0))
        True
        """

        from cmath import rect  # Bibliothèque de fonctions complexes

        # On calcule la rotation en passant dans le plan complexe
        z = complex(self.x, self.y)
        phase = angle if not deg else angle / 57.29577951308232  # [rad]
        u = rect(1., phase)  # exp(i*phase)
        zu = z * u           # Rotation complexe

        return Vector(O, Point(zu.real, zu.imag))


def circumscribedCircle(M, N, P):
    """
    Calcule le centre et le rayon du cercle circonscrit aux points
    M, N, P.

    Retourne: (centre [Point], rayon [float])

    Lève une exception `ValueError` si le rayon ou le centre du cercle
    circonscrit n'est pas défini.

    >>> M = Point(-1, 0); N = Point(1, 0); P = Point(0, 1)
    >>> C, r = circumscribedCircle(M, N, P)  # Centre O, rayon 1
    >>> C.distance(O), round(r, 6)
    (0.0, 1.0)
    >>> circumscribedCircle(M, O, N)         # Indéfini
    Traceback (most recent call last):
    ...
    ValueError: Undefined circumscribed circle radius.
    """

    MN = Vector(M, N)
    NP = Vector(N, P)
    PM = Vector(P, M)

    # Rayon du cercle circonscrit
    m = abs(NP)  # |NP|
    n = abs(PM)  # |PM|
    p = abs(MN)  # |MN|

    d = (m + n + p) * (-m + n + p) * (m - n + p) * (m + n - p)
    if d > 0:
        rad = m * n * p / d**0.5
    else:
        raise ValueError("Undefined circumscribed circle radius.")

    # Centre du cercle circonscrit
    d = -2 * (M.x * NP.y + N.x * PM.y + P.x * MN.y)
    if d == 0:
        raise ValueError("Undefined circumscribed circle center.")

    om2 = Vector(O, M).sqnorm  # |OM|**2
    on2 = Vector(O, N).sqnorm  # |ON|**2
    op2 = Vector(O, P).sqnorm  # |OP|**2

    x0 = -(om2 * NP.y + on2 * PM.y + op2 * MN.y) / d
    y0 = (om2 * NP.x + on2 * PM.x + op2 * MN.x) / d

    return (Point(x0, y0), rad)  # (centre [Point], R [float])


if __name__ == '__main__':

    # start-argparse
    import argparse

    parser = argparse.ArgumentParser(
        usage="%(prog)s [-p/--plot] [-i/--input coordfile | x1,y1 x2,y2 x3,y3]",
        description="Compute the circumscribed circle to 3 points in the plan.")
    parser.add_argument('coords', nargs='*', type=str, metavar='x,y',
                        help="Coordinates of point")
    parser.add_argument('-i', '--input', nargs='?', type=argparse.FileType('r'),
                        help="Coordinate file (one 'x,y' per line)")
    parser.add_argument('-P', '--plot', action="store_true", default=False,
                        help="Draw the circumscribed circle")
    parser.add_argument('-T', '--tests', action="store_true", default=False,
                        help="Run doc tests")
    parser.add_argument('--version', action='version', version=__version__)

    args = parser.parse_args()
    # end-argparse

    if args.tests:                                    # Auto-test mode
        import sys, doctest

        fails, tests = doctest.testmod(verbose=True)  # Run doc tests
        sys.exit(fails > 0)

    if args.input:  # Lecture des coordonnées du fichier d'entrée
        # Le fichier a déjà été ouvert en lecture par argparse (type=file)
        args.coords = [coords for coords in args.input
                       if not coords.strip().startswith('#')]

    if len(args.coords) != 3:  # Vérifie le nb de points
        parser.error("Specify 3 points by their coordinates 'x,y' "
                     f"(got {len(args.coords)})")

    points = []  # Liste des points
    for i, arg in enumerate(args.coords, start=1):
        try:  # Déchiffrage de l'argument 'x,y'
            x, y = (float(t) for t in arg.split(','))
        except ValueError:
            parser.error(f"Cannot decipher coordinates #{i}: {arg!r}")

        points.append(Point(x, y))      # Création du point et ajout à la liste
        print(f"#{i:d}: {points[-1]}")  # Affichage du dernier point

    # Calcul du cercle cisconscrit (lève une ValueError en cas de problème)
    center, radius = circumscribedCircle(*points)  # Délistage
    print(f"Circumscribed circle: {center}, radius: {radius}")

    if args.plot:                   # Figure
        import matplotlib.pyplot as P

        fig = P.figure()
        ax = fig.add_subplot(1, 1, 1, aspect='equal')
        # Points
        ax.plot([p.x for p in points], [p.y for p in points], 'ko')
        for i, p in enumerate(points, start=1):
            ax.annotate(f"#{i}", (p.x, p.y),
                        xytext=(5, 5), textcoords='offset points')
        # Cercle circonscrit
        c = P.matplotlib.patches.Circle((center.x, center.y), radius=radius,
                                        fc='none', ec='k')
        ax.add_patch(c)                    # Cercle
        ax.plot(center.x, center.y, 'r+')  # Centre

        P.show()
$ ./circonscrit.py -h
usage: circonscrit.py [-p/--plot] [-i/--input coordfile | x1,y1 x2,y2 x3,y3]

positional arguments:
  x,y                   Coordinates of point

optional arguments:
  -h, --help            show this help message and exit
  -i [INPUT], --input [INPUT]
                        Coordinate file (one 'x,y' per line)
  -p, --plot            Draw the circumscribed circle
  --version             show program's version number and exit

$ ./circonscrit.py -p 0,0 1,0 1,1
Cercle circonscrit.

Source: circonscrit.py

10.4. Matplotlib🔗

10.4.1. Figure (relativement) simple🔗

Exemple de figure simple.

Fig. 10.1 Figure: exemple de figure (charte graphique: seaborn)🔗

 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
#!/usr/bin/env python3
# Time-stamp: <2018-01-09 15:03:52 ycopin>

"""
Exemple un peu plus complexe de figure, incluant 2 axes, légendes, axes, etc.
"""

import numpy as N
import matplotlib.pyplot as P

x = N.linspace(-N.pi, 3*N.pi, 2*360)

# Signal carré
y = N.sign(N.sin(x))            # = ± 1

# 3 premiers termes de la décomposition en série de Fourier
y1 = 4/N.pi * N.sin(x)          # Fondamentale
y2 = 4/N.pi * N.sin(3*x) / 3    # 1re harmonique
y3 = 4/N.pi * N.sin(5*x) / 5    # 2de harmonique
# Somme des 3 premières composantes
ytot = y1 + y2 + y3

# Figure
fig = P.figure()                # Création de la Figure

# 1er axe: composantes
ax1 = fig.add_subplot(2, 1, 1,  # 1er axe d'une série de 2 × 1
                      ylabel="Composantes",
                      title="Décomposition en série de Fourier")
ax1.plot(x, y1, label="Fondamental")
ax1.plot(x, y2, label=u"1re harmonique")
ax1.plot(x, y3, label=u"2de harmonique")
ax1.legend(loc="upper left", fontsize="x-small")

# 2nd axe: décomposition
ax2 = fig.add_subplot(2, 1, 2,  # 2d axe d'une série de 2 × 1
                      ylabel=u"Décomposition",
                      xlabel="x [rad]")
ax2.plot(x, y, lw=2, color='k', label="Signal carré")
ax2.plot(x, ytot, lw=2, ls=':', color='k', label="Somme des composantes")
ax2.legend(loc="upper left", fontsize="x-small")

# Sauvegarde de la figure (pas d'affichage intéractif)
fig.savefig("figure.png")

Source: figure.py

10.4.2. Filtres du 2nd ordre🔗

Filtre passe-haut du 2nd ordre

Fig. 10.2 Figure: Filtre passe-haut du 2nd ordre.🔗

  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
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
#!/usr/bin/env python3

import numpy as N
import matplotlib.pyplot as P


def passeBas(x, Q=1):
    """
    Filtre passe-bas en pulsation réduite *x* = omega/omega0, facteur de
    qualité *Q*.
    """

    return 1 / (1 - x ** 2 + x / Q * 1j)


def passeHaut(x, Q=1):

    return -x ** 2 / (1 - x ** 2 + x / Q * 1j)


def passeBande(x, Q=1):

    return 1 / (1 + Q * (x - 1 / x) * 1j)


def coupeBande(x, Q=1):

    return (1 - x ** 2) / (1 - x ** 2 + x / Q * 1j)


def gainNphase(f, dB=True):
    """
    Retourne le gain (éventuellement en dB) et la phase [rad] d'un
    filtre de fonction de transfert complexe *f*.
    """

    g = N.abs(f)                        # Gain
    if dB:                              # [dB]
        g = 20 * N.log10(g)
    p = N.angle(f)                      # [rad]

    return g, p


def asympGain(x, pentes=(0, -40)):

    lx = N.log10(x)
    return N.where(lx < 0, pentes[0] * lx, pentes[1] * lx)


def asympPhase(x, phases=(0, -N.pi)):

    return N.where(x < 1, phases[0], phases[1])


def diagBode(x, filtres, labels,
             title='', plim=None, gAsymp=None, pAsymp=None):
    """
    Trace le diagramme de Bode -- gain [dB] et phase [rad] -- des filtres
    de fonction de transfert complexe *filtres* en fonction de la pulsation
    réduite *x*.
    """

    fig = P.figure()
    axg = fig.add_subplot(2, 1, 1,        # Axe des gains
                          xscale='log',
                          ylabel='Gain [dB]')
    axp = fig.add_subplot(2, 1, 2,        # Axe des phases
                          sharex=axg,
                          xlabel=r'x = $\omega$/$\omega_0$', xscale='log',
                          ylabel='Phase [rad]')

    lstyles = ['--', '-', '-.', ':']
    for f, label, ls in zip(filtres, labels, lstyles):  # Tracé des courbes
        g, p = gainNphase(f, dB=True)       # Calcul du gain et de la phase
        axg.plot(x, g, lw=2, ls=ls, label="Q=" + str(label))  # Gain
        axp.plot(x, p, lw=2, ls=ls)                           # Phase

    # Asymptotes
    if gAsymp is not None:              # Gain
        axg.plot(x, asympGain(x, gAsymp), 'k:', lw=2, label='_')
    if pAsymp is not None:              # Phase
        # axp.plot(x, asympPhase(x,pAsymp), 'k:')
        pass

    axg.legend(loc='best', prop=dict(size='small'))

    # Labels des phases
    axp.set_yticks(N.arange(-2, 2.1) * N.pi / 2)
    axp.set_yticks(N.arange(-4, 4.1) * N.pi / 4, minor=True)
    axp.set_yticklabels([r'$-\pi$', r'$-\pi/2$', r'$0$', r'$\pi/2$', r'$\pi$'])
    # Domaine des phases
    if plim is not None:
        axp.set_ylim(plim)

    # Ajouter les grilles
    for ax in (axg, axp):
        ax.grid()                           # x et y, majors
        ax.grid(which='minor')              # x et y, minors

    # Ajustements fins
    gmin, gmax = axg.get_ylim()
    axg.set_ylim(gmin, max(gmax, 3))

    fig.subplots_adjust(hspace=0.1)
    axg.xaxis.set_major_formatter(P.matplotlib.ticker.ScalarFormatter())
    P.setp(axg.get_xticklabels(), visible=False)

    if title:
        axg.set_title(title)

    return fig


if __name__ == '__main__':

    x = N.logspace(-1, 1, 1000)              # de 0.1 à 10 en 1000 pas

    # Facteurs de qualité
    qs = [0.25, 1 / N.sqrt(2), 5]            # Valeurs numériques
    labels = [0.25, r'$1/\sqrt{2}$', 5]      # Labels

    # Calcul des fonctions de transfert complexes
    pbs = [ passeBas(x, Q=q) for q in qs ]
    phs = [ passeHaut(x, Q=q) for q in qs ]
    pcs = [ passeBande(x, Q=q) for q in qs ]
    cbs = [ coupeBande(x, Q=q) for q in qs ]

    # Création des 4 diagrammes de Bode
    figPB = diagBode(x, pbs, labels, title='Filtre passe-bas',
                     plim=(-N.pi, 0),
                     gAsymp=(0, -40), pAsymp=(0, -N.pi))
    figPH = diagBode(x, phs, labels, title='Filtre passe-haut',
                     plim=(0, N.pi),
                     gAsymp=(40, 0), pAsymp=(N.pi, 0))
    figPC = diagBode(x, pcs, labels, title='Filtre passe-bande',
                     plim=(-N.pi / 2, N.pi / 2),
                     gAsymp=(20, -20), pAsymp=(N.pi / 2, -N.pi / 2))
    figCB = diagBode(x, cbs, labels, title='Filtre coupe-bande',
                     plim=(-N.pi / 2, N.pi / 2),
                     gAsymp=(0, 0), pAsymp=(0, 0))

    P.show()

Source: filtres2ndOrdre.py