# -*- coding: utf-8 -*-
#      diese Zeile legt die Codierung von Umlauten fest
##################################################################

# script PlotResonanz.py
''' 
   Funktionsdarstellung mit matplotlib
   Illustriert die Anwendung von Stil-Optionen im Modul `myPlotStyle`

.. author:: Günter Quast <g.quast@kit.edu> 
     für den Kurs Computergestützte Datenauswertung (CgDA)
'''
#--------------------------------------------------------------

# pyhton2 - python3 compatibility
from __future__ import division
from __future__ import absolute_import
from __future__ import print_function
from __future__ import unicode_literals
#

import numpy as np
import matplotlib.pyplot as plt
from myPltStyle import *           # collection of options


# Beispiele für Funktionsdefinitionen
#   Funktionen aus numpy verwenden, da sie auf array angewendet werden können

# 1. Resonanzkurve eines harmonisch getriebenen, gedaempften Oszillators
def resonanz(eta, D, A=1.):
  return A/np.sqrt( (1-eta*eta)**2 + (2.*eta*D)**2 )

# 2. Fermifunktion, e Energie, t Temperatur, mu Fermi-Energie
def fermi(e, t, EF=1.):
  return ( 1.0/(np.exp((e-EF)/t) + 1.0))

##### ---- main Program starts here -----

# create figure object
fig = plt.figure('MainFigure', figsize=(7.5, 10.))

## numpy-array(s) fuer Stuetzstellen und axis-Objekt(e) erzeugen
etamin,etamax=0.,3.
npoints=256
X1=np.linspace( etamin, etamax, npoints)
ax1 = fig.add_subplot(2,1,1)

emin,emax = 0., 5.
npoints=256
X2=np.linspace( emin, emax, npoints)
ax2 = fig.add_subplot(2,1,2)

# zunaechst eigene Stil-Optionen setzen
my_axStyle(fig, ax1)
my_axStyle(fig, ax2)

## Plot(s) der Funktion(en) erzeugen
for D in (0.05, 0.1, 0.2, 0.3, 0.5, 1., 2., 3.0):
  ax1.plot(X1, resonanz(X1, D), linewidth=2, label="D="+str(D))

for E in (0.01, 0.05, 0.1, 0.2, 0.3, 0.6, 1.0, 2.5, 5.0, 10.):
  ax2.plot(X2, fermi(X2, E), linewidth=2, label="$k_B T$="+str(E))
 

# und Grafiken noch beschriften
ax1.set_title('Resonanzkurve')
ax1.set_xlabel(r'$\omega/\omega_0$', size='x-large')
ax1.set_ylabel(r'$Resonanzamplidude$', size='x-large')
ax1.set_xlim(X1.min()*1.1, X1.max()*1.1)
ax1.set_ylim(0., 10.)
ax1.legend(loc='upper right')
ax1.grid(True)
# add grapical elements
ax1.axvline(x=1., ymin=0., ymax=10, color='r', lw=2, linestyle='--' )

ax2.set_title('Fermi-Dirac-Verteilung')
ax2.set_xlabel(r'$Energie$', size='x-large')
ax2.set_ylabel(r'$Wahrscheinlichkeit$', size='x-large')
ax2.set_xlim(X2.min()*1.1, X2.max()*1.1)
ax2.set_ylim(-0.1, 1.1)
ax2.legend(fontsize=12, loc='upper right')
ax2.grid(True)
# add grapical element(s)
ax2.axvline(x=1., ymin=-0.1, ymax=1.1, color='r', lw=2, linestyle='--' )
ax2.text(1.05 ,-0.09, r'$E_F=1.$', color='r', size='x-large')

plt.show() # show on screen
#fig.savefig('MyFig.pdf')