from math import acos, sqrt, pi
import numpy as np
import matplotlib.pyplot as plt
import random


## Enveloppe convexe

def gauche(p0, p1, p2):
    """
    Entrées: trois points p0, p1, p2
    Sortie: True si (p0, p1, p2) tourne à gauche, False sinon
    """
    x0, y0 = p0
    x1, y1 = p1
    x2, y2 = p2
    # signe du produit vectoriel
    return ((y1 -y2)*(x0-x1) - (x1-x2)*(y0-y1)) > 0

def graham(pts):
    """
    Entrée: pts une liste de points
    Sortie: liste des points de l'enveloppe convexe de pts
    """
    if not pts:
        return []
    # trouver et enlever le point d'ordonnée minimale
    p0 = pts[0]
    for p in pts:
        if p[1] < p0[1]:
            p0 = p
        elif p[1] == p0[1] and p[0] > p0[0]:
            p0 = p
    pts.remove(p0)

    # trier les points par angle avec p0 et l'horizontale
    def angle(p):
        """
        renvoie l'angle de (p0, p) par rapport à la droite horizontale
        passant par p0.
        hypothèse: p0 est plus bas que p
        """
        x1, y1 = p0
        x2, y2 = p

        # cos = adjacent / hypothénuse
        cos_theta = (x2-x1) / sqrt((x2-x1)**2 + (y2-y1)**2)
        return acos(cos_theta)

    # on pouvait implémenter un des algos de tri vus en cours
    # comme le tri sélection ou le tri insertion, mais on peut
    # aussi utiliser la fonction sort pré-existante de python.
    # On peut donner en argument à cette fonction une fonction
    # utilisée pour la comparaison: en spécifiant key=angle, la
    # fonction va comparer deux points q et p selon les valeurs
    # de angle(p) et angle(q).
    pts.sort(key=angle)
    # une remarque: on aurait pu trier par cosinus décroissant plutôt que par
    # angle croissant, ce qui aurait donné le même résultat mais aurait permis
    # d'économiser le calcul de acos qui est coûteux.

    # itérer sur tous les points. On simule le passage pour p0
    # en le mettant directement sur la pile
    S = [p0] # pile des points
    for p in pts:
        while len(S)>1 and not gauche(S[-2], S[-1], p):
            S.pop()
        S.append(p)
    return S


def plot_pts(pts, link, color):
    """
    Entrées: pts une liste de points, link un booléen
    Sortie: Trace (sans afficher) les points de pts,
            en les reliant si link = True.
    """
    x = [p[0] for p in pts]
    y = [p[1] for p in pts]

    plt.plot(x, y, 'o', color=color)
    if link and x:
        # Rajouter le premier point à la fin pour faire une boucle
        plt.plot(x+[x[0]], y+[y[0]], color=color, linestyle='-')



def test_graham(A=50, N=100):
    pts = [(random.randint(-A, A), random.randint(-A, A)) for i in range(N)]
    # ne garder que les points dans le cercle
    pts = [(x, y) for (x, y) in pts if x*x+y*y < A*A]

    enveloppe = graham(pts)
    plot_pts(pts, False)
    plot_pts(enveloppe, True)
    plt.show()

## Théorème de Gauss-Lucas
"""
numpy possède une classe (c'est à dire un type) représentant les polynomes,
que l'on importe avec la commande suivante:
"""
from numpy.polynomial import Polynomial
"""
A partir d'un polynome p, on peut obtenir ses racines avec p.roots(),
obtenir le polynome dérivé avec p.deriv().
On peut créer un polynome avec
    p = Polynomial(l)
où l est la liste des coefficients. Par exemple,
p = Polynomial([1, 2, 3]) sera le polynôme X² + 2X + 3

De plus, on peut manipuler des complexes en python: l'unité imaginaire
s'écrit 1j. Par exemple:
    z = a + 1j * b
"""

def rand_poly(n):
    """
    Renvoie un polynome de degré n-1, dont chaque coefficient est choisi
    uniformément dans [-3, 3].
    """
    l = []
    for i in range(n):
        re = random.random()*3-6
        im = random.random()*3-6
        l.append(re + 1j * im)
    return Polynomial(l)

"""
Mise en place d'un système de couleurs:
Les fonctions de matplotlib.pyplot utilisent un format pour les couleurs
particulier: une couleur est représentée par une chaîne de caractères
de la forme "#rrggbb" où rr, gg, et bb sont les montants de rouge, vert et
bleu, écrits en hexadécimal.
Par exemple, la couleur #00ac1e a
- 0/255 de rouge
- 172/255 de vert car 0xac = 10*16 + 12 = 172
- 30/255 de bleu  car 0x1e = 30
"""

def pad(c):
    """
    Entrée: c un entier entre 0 et 255
    Renvoie une chaîne de caractère représentant c en hexadécimal,
    toujours sur 2 caractères.
    ex:
    - 2   -> 0x2  -> 02
    - 255 -> 0xff -> ff
    - 0   -> 0x0  -> 00
    """
    s = hex(c)
    return s[2:] + "0"*(4-len(s))

def to_hex(r, g, b):
    """
    Entrée: r, g, b des flottants entre 0 et 255
    Renvoie la couleur (r, g, b) sous format HTML
    #rrggbb
    """
    r, g, b = int(r), int(g), int (b)
    return f"#{pad(r)}{pad(g)}{pad(b)}"

def show_roots(p, color):
    """
    Entrée: p un polynome, color une couleur
    Place les racines de p dans le plan, relie
    les points de l'enveloppe convexe.
    """
    rts = p.roots()
    pts = [(a.real, a.imag) for a in rts]
    hull = graham(pts)
    plot_pts(pts, False, color)
    plot_pts(hull, True, color)


def gauss_lucas(n):
    """
    Trace les racines d'un polynome P aléatoire à n coefficients,
    ainsi que l'enveloppe convexe de ces racines.
    Répète l'opération pour P', P'', etc...
    """
    p = rand_poly(n)
    k = 0
    while(len(p) >= 3):
        k += 1
        # couleur de la k-ème dérivée
        color = to_hex(255-k*(255/n), 40 + k*(180/n), 50+k*(200/n))
        show_roots(p, color)
        p = p.deriv()
    plt.show()