#!/usr/bin env python3
from formulary import *

# recreating a plot from adv. sc physics ex7/2
# showing different scattering mechanisms

def fionized(T):
    return 10**4 * T**(3/2)

@np.vectorize
def fpolar_optic(T, TDebye):
    return np.exp(TDebye/T) if T < TDebye else np.nan

def fpiezoelectric(T):
    return 0.8 * 10**7 * T**(-1/2)

def facoustic(T):
    return 10**9 * T**(-3/2)


def scattering():
    fig, ax = plt.subplots(1, 1, figsize=size_formula_normal_default)
    Ts = np.power(10, np.linspace(0, 3, 100))

    TDebye = 10**3
    mu_ionized = fionized(Ts)
    mu_polar_optic = fpolar_optic(Ts, TDebye)
    mu_piezoelectric = fpiezoelectric(Ts)
    mu_acoustic = facoustic(Ts)

    mu_sum = 1/(1/mu_ionized + 1/mu_polar_optic + 1/mu_piezoelectric + 1/mu_acoustic)
    ax.plot(Ts, mu_ionized, label="Ionized")
    ax.plot(Ts, mu_polar_optic, label="Polar-optic")
    ax.plot(Ts, mu_piezoelectric, label="Piezoelectric")
    ax.plot(Ts, mu_acoustic, label="Acoustic")
    ax.plot(Ts, mu_sum, label="Total", color="black", linestyle="dashed")
    ax.set_xscale("log")
    ax.set_yscale("log")
    ax.set_xlabel("$T$")
    ax.set_ylabel("$\\mu$")
    ax.legend()
    ax.set_ylim(1e4, 1e7)
    return fig


if __name__ == '__main__':
    export(scattering(), "cm_scattering")