This repository has been archived by the owner on Apr 13, 2024. It is now read-only.
-
Notifications
You must be signed in to change notification settings - Fork 1
/
simulated_patterns.py
executable file
·98 lines (71 loc) · 3.57 KB
/
simulated_patterns.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
import glob
import numpy as np
from pymatgen.core.lattice import Lattice
from pymatgen.analysis.diffraction.xrd import XRDCalculator
from pymatgen.io.cif import CifParser
from model_functions import GaussianSum
from random import uniform, choice
class simulated_pattern(object):
def __init__(self, cif, max_2theta=40.0):
cif = CifParser(cif)
self.struct = cif.get_structures(primitive=False)[0]
self.max_2theta = max_2theta
def change_lattice_constants(self, uc_params):
a,b,c,alpha,beta,gamma = uc_params
pi = np.pi
ax = a
ay = 0.0
az = 0.0
bx = b * np.cos(gamma * pi / 180.0)
by = b * np.sin(gamma * pi / 180.0)
bz = 0.0
cx = c * np.cos(beta * pi / 180.0)
cy = (c * b * np.cos(alpha * pi / 180.0) - bx * cx) / by
cz = (c ** 2.0 - cx ** 2.0 - cy ** 2.0) ** 0.5
unit_cell = np.array([[ax,ay,az],[bx,by,bz],[cx,cy,cz]])
new_lattice = Lattice(unit_cell)
self.struct.lattice = new_lattice
def randomly_perturb_lattice_constants(self, pmin, pmax):
uc_params0 = self.struct.lattice.abc + self.struct.lattice.angles
perturbed_params_len = [p + choice((-1,1))*uniform(pmin,pmax)*p for p in uc_params0[0:3]]
#perturbed_params_ang = [p + choice((-1,1))*uniform(pmin/2.0,pmax/2.0)*p for p in uc_params0[3:]]
perturbed_params_ang = list(uc_params0[3:])
perturbed_params = perturbed_params_len + perturbed_params_ang
self.change_lattice_constants(perturbed_params)
def find_peaks(self, wavelength='CuKa'):
XRD = XRDCalculator(symprec=0, wavelength=wavelength).get_pattern(self.struct, scaled=True, two_theta_range=(0,self.max_2theta))
return np.array([XRD.x, XRD.y]).T
def fit_GaussianModel(self, uc_params, mode='aniso', sd=0.01, wavelength='CuKa'):
if mode == 'iso':
uc_params = np.array([uc_params[0],uc_params[0],uc_params[0],uc_params[1],uc_params[2],uc_params[3]])
elif mode == 'ab':
uc_params = np.array([uc_params[0],uc_params[0],uc_params[1],uc_params[2],uc_params[3],uc_params[4]])
elif mode == 'ac':
uc_params = np.array([uc_params[0],uc_params[1],uc_params[0],uc_params[2],uc_params[3],uc_params[4]])
elif mode == 'bc':
uc_params = np.array([uc_params[0],uc_params[1],uc_params[1],uc_params[2],uc_params[3],uc_params[4]])
self.change_lattice_constants(uc_params)
peaks = self.find_peaks(wavelength=wavelength)
data = np.c_[peaks, np.full((len(peaks),), fill_value=sd)]
model = lambda x: GaussianSum(x, data)
return model
def unchanged_GaussianModel(self, xvals, sd=0.01, wavelength='CuKa'):
peaks = self.find_peaks(wavelength=wavelength)
data = np.c_[peaks, np.full((len(peaks),), fill_value=sd)]
model = lambda x: GaussianSum(x, data)
yvals = model(xvals)
return np.array([xvals,yvals]).T
def model_pattern_predictions(self, xvals, uc_params, mode='aniso', wavelength='CuKa'):
model = self.fit_GaussianModel(uc_params, mode=mode, wavelength=wavelength)
yvals = model(xvals)
yvals /= max(yvals)
return np.array([xvals,yvals]).T
if __name__ =='__main__':
cifs = glob.glob('*.cif')
for cif in cifs:
name = cif.split('.')[0]
topo = name.split('_')[1]
SIMP = simulated_pattern(cif, max_2theta=10.0)
xvals = np.linspace(0.5,10.0,2500)
pattern = SIMP.unchanged_GaussianModel(xvals, wavelength=0.45256)
np.savetxt('sim_' + topo + '.txt', pattern, delimiter=' ')