# -*- 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]])])
# -----------------------------------------------------------------
# -----------------------------------------------------------------
# -----------------------------------------------------------------