Source code for haslib.crystal_class

# -*- coding: utf-8 -*-
# crystal_class.py
#
# Class definition for an object which stores all the information
# already obtained on a certain material
#
# Patrick Kraus
# 24.02.2012 (10 Months to Cristmas!)
# 27.02.2012 Added crystal turning and Leedpicture

# import the used libraries
from .science_basics import (kb, mel, hb)
from numpy import (zeros, array, cos, sin, dot, cross, pi, arctan,
                   hstack, arange, sqrt, vstack)
from matplotlib.pyplot import (axis, figure, show)


[docs]class Crystal(object): # ----------------------------------------------------------------- # ---------------------- INITIALROUTINE --------------------------- # ----------------------------------------------------------------- def __init__(self, a1, a2, ang, name="set sample name in crystal class"): """ crystal object Start Parameters: a1 ... Lattice Vector 1 in meters a2 ... Lattice Vector 2 in meters ang ... Angle between the lattice vectors (degrees) name ... Sample name(e.g. Bi_2Te_3) The parameters can be altered later by using ch_a1, ch_a2 and ch_ang or calculated from the stored elastic positions using calc_params """ self.name = name # Setting the initial parameters self.ang = pi*ang/180.0 self.a1 = array([a1, 0, 0]) self.a2 = array([a2*cos(self.ang), a2*sin(self.ang), 0]) self.a3 = array([0, 0, 1]) # Calculating the reciprocal lattice vectors self.calc_reciprocals() self.b1 = self.b1_orig self.b2 = self.b2_orig # Rotating the reciprocals to <10> and saving as default self.b1 = self.rotate_to(self.b1, 1, 0) self.b2 = self.rotate_to(self.b2, 1, 0) self.b1_orig = self.b1 self.b2_orig = self.b2 # Calculating the default Grid in <10> orientation self.ran = 3 self.calc_grid(self.ran) # Defining known bound state energies # At first only the Threshold is for sure self.interactE = array([0.0]) # Defining important reciprocal Lattice vectors for resonances self.Gvecs = array([[1, 0], [1, 1]]) # Define phonon container self.phoncont = array([0.0, 0.0, 0.0, 0.0, 0.0]) # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ---------------- ACTUALIZATION - ROUTINES ----------------------- # -----------------------------------------------------------------
[docs] def calc_reciprocals(self): # Calculating the reciprocal lattice Vectors self.b1_orig = (2*pi*cross(self.a2, self.a3)) / \ (dot(self.a1, cross(self.a2, self.a3))) self.b2_orig = (2*pi*cross(self.a3, self.a1)) / \ (dot(self.a1, cross(self.a2, self.a3))) self.b3_orig = (2*pi*cross(self.a1, self.a2)) / \ (dot(self.a1, cross(self.a2, self.a3)))
# ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ---------------- MANIPULATION - ROUTINES ------------------------ # ----------------------------------------------------------------- # First the simple functions for manipulating the crystal parameters
[docs] def ch_a1(self, new_a1): "Function to set a1 to a new value" self.a1 = array([new_a1, 0, 0]) # Renewing self.calc_reciprocals()
[docs] def ch_a2(self, new_a2): "Function to set a2 to a new value" self.a2 = array([new_a2*cos(self.ang), new_a2*sin(self.ang), 0]) # Renewing self.calc_reciprocals()
[docs] def ch_ang(self, new_ang): "Function to set ang to a new value - in degrees" self.ang = pi*new_ang/180.0 self.a2 = array([self.a2*cos(self.ang), self.a2*sin(self.ang), 0]) # Renewing self.calc_reciprocals()
[docs] def calc_angle(self, G1, G2, Flag=False): G = self.b1_orig*G1+self.b2_orig*G2 if abs(G[0]) < 1E-6: angle = 0 else: angle = -arctan(G[1]/G[0]) if Flag: print("Rotate the Sample ", 180.0*angle/pi, " degrees from <10>") return angle
[docs] def rotate(self, ang, vec): R = array([[cos(ang), -sin(ang), 0], [sin(ang), cos(ang), 0], [0, 0, 1]]) erg = dot(R, vec) return erg
[docs] def rotate_to(self, points, G1, G2): """ Function to rotate the position of a single point to a certain lattice orientation """ angle = self.calc_angle(G1, G2) points = self.rotate(angle, points) return points
[docs] def rotate_crystal(self, G1, G2): """ This function rotates the crystal to the desired direction """ self.b1 = self.rotate_to(self.b1_orig, G1, G2) self.b2 = self.rotate_to(self.b2_orig, G1, G2) self.calc_grid(self.ran)
[docs] def add_BSEnergy(self, Enew): "Add a single measured Bound state energy" self.interactE = hstack([self.interactE, Enew])
# ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------- PREDICTION - ROUTINES ------------------------- # -----------------------------------------------------------------
[docs] def calc_grid(self, ran=3): if ran != 3: self.ran = ran siz = (2*ran+1)**2 self.grid = zeros([siz, 4]) run = 0 for G1 in arange(2*ran+1)-ran: for G2 in arange(2*ran+1)-ran: self.grid[run, 0] = self.b1[0]*G1+self.b2[0]*G2 self.grid[run, 1] = self.b1[1]*G1+self.b2[1]*G2 self.grid[run, 2] = G1 self.grid[run, 3] = G2 run += 1
[docs] def leedpicture(self): "Function which plots the expected LEED picture" leed_fig = figure() leed_ax = leed_fig.add_subplot(111, axis_bgcolor='# 122416') leed_ax.plot(self.grid[:, 0], self.grid[:, 1], marker='o', linestyle='None', color='# 218128') axis('equal') leed_ax.set_xlim([-1.1*self.b1_orig[0], 1.1*self.b1_orig[0]]) leed_ax.set_ylim([-1.1*self.b1_orig[0], 1.1*self.b1_orig[0]]) show()
[docs] def givevec(self, G1, G2): """ Function to return the x and y coordinates of a Gvec in the current orientation """ # indx = find((self.grid[:, 2]==G1)&(self.grid[:, 3]==G2))[0] # return array([self.grid[indx, 0], self.grid[indx, 1]]) return G1*self.b1+G2*self.b2
# ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------------------------------------------------------- # ----------------- MEASUREMENT - ROUTINES ------------------------ # -----------------------------------------------------------------
[docs] def log_phonon(self, Appobj, Tofobj, QB=9230476431.8783379): # 7993827079.0401821): """ Function to log the measured phonos from a TOF measurement in the crystal object """ # Collecting the neccesary data m = Appobj.m # ki = sqrt(5*m*kb*Tofobj.TN/hb**2) Ei = 5*kb*Tofobj.TN/2.0 TSD = Appobj.TSD Ti = -(Tofobj.theta_app-Appobj.specpos)+TSD/2.0 Tf = TSD-Ti cont = Tofobj.phoncont # Now define the specific scancurve (inverted) def findpoint(Etest): return (sqrt(2*m/hb**2 * (Ei+Etest*mel))*sin(Tf*pi/180.0) - sqrt(2*m/hb**2 * Ei)*sin(Ti*pi/180.0)) # Now, for every entry in cont, log a point for index in range(cont.shape[0]): if cont.ndim > 1: Q = findpoint(cont[index, 0]) else: Q = findpoint(cont[index]) # Now refold for the brillouin-zone Q = Q % (QB*2.0) if Q > QB: Q = abs(Q-2*QB) self.phoncont = vstack([self.phoncont, array([abs(cont[index, 0]), cont[index, 1], Q, Tofobj.G[0], Tofobj.G[1]])])
# ----------------------------------------------------------------- # ----------------------------------------------------------------- # -----------------------------------------------------------------