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
|
#!/usr/bin/env python3
"""
Fitting example: fit with masks
"""
from matplotlib import pyplot as plt
import bornagain as ba
from bornagain import ba_fit, ba_fitmonitor, deg, nm, nm2
import lmfit
import numpy as np
def get_sample(P):
"""
Uncorrelated cylinders on a substrate, parameterized for fitting.
"""
substrate_mat = ba.RefractiveMaterial(
"Substrate", (0.28, 0.57, 0.82), 6e-6, 2e-8)
particle_mat = ba.RefractiveMaterial(
"Particle", (0.86, 0.24, 0.18), 6e-4, 2e-8)
particle = ba.Particle(
particle_mat, ba.Cylinder(P["radius"], P["height"]))
particle_layer = ba.Layer(ba.Vacuum())
particle_layer.deposit2D(ba.Dilute2D(0.01/nm2, particle))
sample = ba.Sample()
sample.addLayer(particle_layer)
sample.addLayer(ba.Layer(substrate_mat))
return sample
def get_simulation(P):
"""
GISAS simulation for the parameterized cylinder sample.
"""
n = 100
beam = ba.Beam(1e8, 0.1*nm, 0.2*deg)
detector = ba.SphericalDetector(n, -1*deg, 1*deg, n, 0, 2*deg)
simulation = ba.ScatteringSimulation(beam, get_sample(P), detector)
simulation.options().setUseAvgMaterials(False)
return simulation
def fake_data():
"""
Generate noisy synthetic data for a known cylinder sample.
"""
P = {"radius": 5*nm, "height": 10*nm}
return get_simulation(P).simulate().noisy(0.1, 0.1)
def get_masked_simulation(P):
"""
GISAS simulation with a Python-generated detector mask.
"""
n = 100
beam = ba.Beam(1e8, 0.1*nm, 0.2*deg)
sample = get_sample(P)
detector = ba.SphericalDetector(n, -1*deg, 1*deg, n, 0, 2*deg)
add_mask_to_detector(detector)
simulation = ba.ScatteringSimulation(beam, sample, detector)
simulation.options().setUseAvgMaterials(False)
return simulation
def add_mask_to_detector(detector):
"""
Adds a Python-generated detector bitmap mask to the detector.
"""
n_phi = detector.axis(0).size()
n_alpha = detector.axis(1).size()
y, x = np.ogrid[:n_alpha, :n_phi]
x0 = 0.5*n_phi
y0 = 0.5*n_alpha
r = 0.32*min(n_phi, n_alpha)
mask = np.ones((n_alpha, n_phi), dtype=bool)
head = (x - x0)**2 + (y - y0)**2 <= r**2
eye = (x - (x0 + 0.22*r))**2 + (y - (y0 + 0.35*r))**2 <= (0.11*r)**2
mouth = (x > x0) & (np.abs(y - y0) < 0.35*(x - x0))
mask[head] = False
mask[eye] = True
mask[mouth] = True
for center in (x0 + 1.15*r, x0 + 1.45*r, x0 + 1.75*r):
snack = (x - center)**2 + (y - y0)**2 <= (0.09*r)**2
mask[snack] = False
detector.setMask(mask)
if __name__ == '__main__':
data = fake_data()
exp_values = data.intensities()
sim_result = None # latest simulation, shared with the plot callback
def residuals(P):
"""
Runs a simulation for given parameters P, and returns residuals
vector; masked (NaN) pixels contribute zero.
"""
global sim_result
sim_result = get_masked_simulation(P.valuesdict()).simulate()
return ba_fit.valid_pixel_residual(exp_values,
sim_result.intensities())
plotter = ba_fitmonitor.PlotterGISAS()
def plot_iteration(P, iteration, resid):
if iteration % 10 == 0 and sim_result is not None:
plotter.plot(data, sim_result, P, float(np.sum(resid*resid)))
P = lmfit.Parameters()
P.add("radius", 6.*nm, min=4, max=8)
P.add("height", 9.*nm, min=8, max=12)
result = lmfit.minimize(residuals, P, iter_cb=plot_iteration)
print(lmfit.fit_report(result))
finalP = result.params.valuesdict()
ba.showSample3D(get_sample(finalP), sample_size=120*nm, seed=0)
plt.show()
|
