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
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
|
#!/usr/bin/env python3
"""
Fitting example: running same fit using various minimizer and their settings.
"""
import bornagain as ba
from bornagain import deg, angstrom, nm
def get_sample(P):
"""
A sample with uncorrelated cylinders and prisms on a substrate.
"""
cylinder_height = P["cylinder_height"]
cylinder_radius = P["cylinder_radius"]
prism_height = P["prism_height"]
prism_base_edge = P["prism_base_edge"]
# defining materials
vacuum = ba.RefractiveMaterial("Vacuum", 0, 0)
material_substrate = ba.RefractiveMaterial("Substrate", 6e-6, 2e-8)
material_particle = ba.RefractiveMaterial("Particle", 6e-4, 2e-8)
# collection of particles
cylinder_ff = ba.Cylinder(cylinder_radius, cylinder_height)
cylinder = ba.Particle(material_particle, cylinder_ff)
prism_ff = ba.Prism3(prism_base_edge, prism_height)
prism = ba.Particle(material_particle, prism_ff)
layout = ba.ParticleLayout()
layout.addParticle(cylinder, 0.5)
layout.addParticle(prism, 0.5)
interference = ba.InterferenceNone()
layout.setInterference(interference)
# vacuum layer with particles and substrate form multilayer
vacuum_layer = ba.Layer(vacuum)
vacuum_layer.addLayout(layout)
substrate_layer = ba.Layer(material_substrate, 0)
sample = ba.Sample()
sample.addLayer(vacuum_layer)
sample.addLayer(substrate_layer)
return sample
def get_simulation(P):
"""
A GISAXS simulation with beam and detector defined
"""
beam = ba.Beam(1e8, 1*angstrom, 0.2*deg)
n = 100
detector = ba.SphericalDetector(n, -1*deg, 1*deg, n, 0, 2*deg)
return ba.ScatteringSimulation(beam, get_sample(P), detector)
def fake_data():
"""
Generating "real" data from simulated image with default parameters.
"""
P = {
'cylinder_height': 5*nm,
'cylinder_radius': 5*nm,
'prism_height': 5*nm,
'prism_base_edge': 5*nm
}
simulation = get_simulation(P)
result = simulation.simulate()
return result
def run_fitting():
"""
main function to run fitting
"""
data = fake_data()
# prints info about available minimizers
print(ba.MinimizerFactory().catalogToString())
# prints detailed info about available minimizers and their options
print(ba.MinimizerFactory().catalogDetailsToString())
fit_objective = ba.FitObjective()
fit_objective.addFitPair(get_simulation, data, 1)
fit_objective.initPrint(10)
P = ba.Parameters()
P.add("cylinder_height", 4.*nm, min=0.01)
P.add("cylinder_radius", 6.*nm, min=0.01)
P.add("prism_height", 4.*nm, min=0.01)
P.add("prism_base_edge", 12.*nm, min=0.01)
minimizer = ba.Minimizer()
minimizer.setMinimizer("Minuit2", "Migrad", "MaxFunctionCalls=5") # fast test
# Uncomment one of the line below to adjust minimizer settings
"""
Setting Minuit2 minimizer with Migrad algorithm, limiting number of iterations.
Minimization will try to respect MaxFunctionCalls value.
"""
# minimizer.setMinimizer("Minuit2", "Migrad", "MaxFunctionCalls=50")
"""
Setting two options at once.
Strategy=2 promises more accurate fit.
"""
# minimizer.setMinimizer("Minuit2", "Migrad", "MaxFunctionCalls=500;Strategy=2")
"""
Setting Minuit2 minimizer with Fumili algorithm.
"""
# minimizer.setMinimizer("Minuit2", "Fumili")
"""
Setting Levenberg-Marquardt algorithm.
"""
# minimizer.setMinimizer("GSLLMA")
result = minimizer.minimize(fit_objective.evaluate_residuals, P)
fit_objective.finalize(result)
print("Fitting completed.")
print("chi2:", result.minValue())
for fitPar in result.parameters():
print(fitPar.name(), fitPar.value, fitPar.error)
if __name__ == '__main__':
run_fitting()
|
