2828import itertools
2929
3030import numpy as np
31+ import scipy
3132import matplotlib .pylab as plt # basic plotting
3233
3334import ase
3435import ase .build
3536
3637import pyTEMlib
3738from pyTEMlib .crystal_tools import electronFF
39+ from ..utilities import get_wavelength , get_rotation_matrix
40+
3841
3942inputKeys = ['acceleration_voltage_V' , 'zone_hkl' , 'Sg_max' , 'hkl_max' ]
4043optional_input_keys = ['crystal' , 'lattice_parameter_nm' , 'convergence_angle_mrad' ,
41- 'mistilt' , 'thickness' , 'dynamic correction' , 'dynamic correction K0' ]
44+ 'mistilt' , 'thickness' , 'dynamic correction' , 'dynamic correction K0' ]
4245
4346
4447def read_poscar (filename ):
@@ -141,7 +144,6 @@ def zuo_fig_3_18(verbose=True):
141144
142145 return atoms
143146
144- get_rotation_matrix = pyTEMlib .utilities .get_rotation_matrix
145147
146148def zone_mistilt (zone , angles ):
147149 """ Rotation of zone axis by mistilt
@@ -173,28 +175,11 @@ def get_metric_tensor(matrix):
173175
174176def vector_norm (g ):
175177 """ Length of vector
176-
177178 depreciated - use np.linalg.norm
178179 """
179180 return np .linalg .norm (g )
180181
181182
182- def get_wavelength (acceleration_voltage ) -> float :
183- """
184- Calculates the relativistic corrected de Broglie wavelength of an electron in Angstrom
185-
186- Parameter:
187- ---------
188- acceleration_voltage: float or int
189- acceleration voltage in volt
190- Returns:
191- -------
192- wavelength: float
193- wave length in Angstrom (= meter *10**10)
194- """
195- return pyTEMlib .utilities .get_wavelength (acceleration_voltage ) * 1e10
196-
197-
198183def get_all_miller_indices (hkl_max ):
199184 """ Get all Miller indices up to hkl_max"""
200185 h = np .linspace (- hkl_max , hkl_max , 2 * hkl_max + 1 ) # all evaluated single Miller Indices
@@ -220,6 +205,27 @@ def get_structure_factors(atoms, g_hkl):
220205 return np .array (structure_factors )
221206
222207
208+ def get_inner_potential (atoms ):
209+ """ inner potential in Volts """
210+ u_0 = 0 # in (Ang)
211+ # atom form factor of zero reflection angle is the inner potential in 1/A
212+ for atom in atoms :
213+ u_0 += form_factor (atom .symbol , 0 )
214+ scattering_factor_to_volts = ((scipy .constants .h * 1e10 )** 2
215+ / (2 * np .pi * scipy .constants .m_e * scipy .constants .e )
216+ / atoms .cell .volume )
217+ return u_0 * scattering_factor_to_volts
218+
219+ def ewald_sphere_center (acceleration_voltage , atoms , zone_hkl ):
220+ """ Ewald sphere center in 1/Angstrom """
221+ wavelength = get_wavelength (acceleration_voltage , unit = 'Å' ) # in Angstrom
222+ u_0 = get_inner_potential (atoms )
223+ incident_wave_vector = np .sqrt (1 / wavelength ** 2 + u_0 )#1/Ang
224+
225+ center = np .dot (zone_hkl , atoms .cell .reciprocal ())
226+ center = center / np .linalg .norm (center ) * incident_wave_vector
227+ return center
228+
223229def find_nearest_zone_axis (tags ):
224230 """Test all zone axis up to a maximum of hkl_max"""
225231
@@ -436,6 +442,7 @@ def gaussian(xx, pp):
436442 s1 = pp [2 ] / 2.3548
437443 prefactor = 1.0 / np .sqrt (2 * np .pi * s1 ** 2 )
438444 y = (pp [1 ] * prefactor ) * np .exp (- (xx - pp [0 ]) ** 2 / (2 * s1 ** 2 ))
445+
439446 return y
440447
441448
@@ -453,28 +460,28 @@ def get_unit_cell(atoms, tags):
453460def output_verbose (atoms , tags ):
454461 """ Verbose output of experimental parameters"""
455462 print ('Experimental Parameters:' )
456- print (f"Acceleration Voltage: { tags [ 'acceleration_voltage_V' ] * 1000 :.1f} kV" )
457- print (f"Convergence Angle: { tags [ 'convergence_angle_mrad' ] :.2f} mrad" ,
458- f" = { tags [ 'convergence_angle_A-1' ] :.2f} 1/Ang" )
459- print (f"Wavelength: { tags [ 'wave_length' ] * 1000 :.3f} pm" )
460- print (f"Incident Wave Vector: { tags [ 'incident_wave_vector' ] * 10 } 1/nm in material " ,
461- f"; in vacumm: { 1 / tags [ 'wave_length' ] :.4f} 1/nm" )
463+ print (f"Acceleration Voltage: { tags . get ( 'acceleration_voltage_V' , 0 ) * 1000 :.1f} kV" )
464+ print (f"Convergence Angle: { tags . get ( 'convergence_angle_mrad' , None ) :.2f} mrad" ,
465+ f" = { tags . get ( 'convergence_angle_A-1' , None ) :.2f} 1/Ang" )
466+ print (f"Wavelength: { tags . get ( 'wave_length' , 0 ) * 1000 :.3f} pm" )
467+ print (f"Incident Wave Vector: { tags . get ( 'incident_wave_vector' , 0 ) * 10 } 1/nm in material " ,
468+ f"; in vacumm: { 1 / tags . get ( 'wave_length' , 0 ) :.4f} 1/nm" )
462469 print ("\n Rotation to tilt zone axis onto z-axis:" )
463- print (f"Rotation alpha { np .rad2deg (tags [ 'y-axis rotation alpha' ] ):.1f} degree, "
464- f" beta { np .rad2deg (tags [ 'x-axis rotation beta' ] ):.1f} degree" )
465- print (f"from zone axis { tags [ 'zone_hkl' ] } " )
466- print (f"Tilting { 1 } by { np .rad2deg (tags [ 'mistilt_alpha' ] ):.2f} "
467- f" in alpha and { np .rad2deg (tags [ 'mistilt_beta' ] ):.2f} " +
470+ print (f"Rotation alpha { np .rad2deg (tags . get ( 'y-axis rotation alpha' , None ) ):.1f} degree, "
471+ f" beta { np .rad2deg (tags . get ( 'x-axis rotation beta' , None ) ):.1f} degree" )
472+ print (f"from zone axis { tags . get ( 'zone_hkl' , None ) } " )
473+ print (f"Tilting { 1 } by { np .rad2deg (tags . get ( 'mistilt_alpha' , None ) ):.2f} "
474+ f" in alpha and { np .rad2deg (tags . get ( 'mistilt_beta' , None ) ):.2f} " +
468475 " in beta direction results in :" )
469- nearest = tags [ 'nearest_zone_axes' ]
476+ nearest = tags . get ( 'nearest_zone_axes' , {})
470477 print ('Next nearest zone axes are:' )
471- for i in range (1 , nearest [ 'amount' ] ):
478+ for i in range (1 , nearest . get ( 'amount' , 0 ) ):
472479 print (f"{ nearest [str (i )]['hkl' ]} : mistilt:" ,
473480 f" { np .rad2deg (nearest [str (i )]['mistilt_alpha' ]):6.2f} , "
474481 f"{ np .rad2deg (nearest [str (i )]['mistilt_beta' ]):6.2f} " )
475- print ('Center of Ewald sphere ' , tags [ 'k0_vector' ] )
476- dif = atoms .info [ 'diffraction' ]
477- print ('Center or Laue circle' , dif [ 'Laue_circle' ] )
482+ print ('Center of Ewald sphere ' , tags . get ( 'k0_vector' , None ) )
483+ dif = atoms .info . get ( 'diffraction' , {})
484+ print ('Center or Laue circle' , dif . get ( 'Laue_circle' , None ) )
478485
479486
480487def make_pretty_labels (hkls , hex_label = False ):
0 commit comments