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

#!/usr/bin/env python3
"""
3 layers system (substrate, teflon, air).
Vacuum layer is populated with spheres with size distribution.
"""
import bornagain as ba


class SampleBuilder:
    """
    Builds complex sample from set of parameters.
    """
    def __init__(self):
        """
        Init SampleBuilder with default sample parameters.
        """
        self.radius = 5.75*ba.nm
        self.sigma = 0.4
        self.distance = 53.6*ba.nm
        self.disorder = 10.5*ba.nm
        self.kappa = 17.5
        self.ptfe_thickness = 22.1*ba.nm
        self.hmdso_thickness = 18.5*ba.nm

    def create_sample(self, params):
        """
        Create sample from given set of parameters.

        Args:
            params: A dictionary containing parameter map.

        Returns:
            A multi layer.
        """
        self.radius = params["radius"]
        self.sigma = params["sigma"]
        self.distance = params["distance"]
        return self.multilayer()

    def multilayer(self):
        """
        Constructs the sample from current parameter values.
        """

        # defining materials
        m_vacuum = ba.HomogeneousMaterial("Vacuum", 0, 0)
        m_Si = ba.HomogeneousMaterial("Si", 5.7816e-6, 1.0229e-7)
        m_Ag = ba.HomogeneousMaterial("Ag", 2.2475e-5, 1.6152e-6)
        m_PTFE = ba.HomogeneousMaterial("PTFE", 5.20509e-6, 1.9694e-8)
        m_HMDSO = ba.HomogeneousMaterial("HMDSO", 2.0888e-6, 1.3261e-8)

        # collection of particles with size distribution
        nparticles = 20
        nfwhm = 2.0
        sphere_ff = ba.FormFactorFullSphere(self.radius)

        sphere = ba.Particle(m_Ag, sphere_ff)
        position = ba.kvector_t(0, 0, -1*self.hmdso_thickness)
        sphere.setPosition(position)
        ln_distr = ba.DistributionLogNormal(self.radius, self.sigma)
        par_distr = ba.ParameterDistribution(
            "/Particle/FullSphere/Radius", ln_distr, nparticles, nfwhm,
            ba.RealLimits.limited(0, self.hmdso_thickness/2))
        part_coll = ba.ParticleDistribution(sphere, par_distr)

        # interference function
        interference = ba.InterferenceFunctionRadialParaCrystal(
            self.distance, 1e6*ba.nm)
        interference.setKappa(self.kappa)
        interference.setDomainSize(2e4*nm)
        pdf = ba.FTDistribution1DGauss(self.disorder)
        interference.setProbabilityDistribution(pdf)

        # assembling particle layout
        layout = ba.ParticleLayout()
        layout.addParticle(part_coll, 1)
        layout.setInterferenceFunction(interference)
        layout.setTotalParticleSurfaceDensity(1)

        # roughness
        r_ptfe = ba.LayerRoughness(2.3*ba.nm, 0.3, 5*ba.nm)
        r_hmdso = ba.LayerRoughness(1.1*ba.nm, 0.3, 5*ba.nm)

        # layers
        vacuum_layer = ba.Layer(m_vacuum)
        hmdso_layer = ba.Layer(m_HMDSO, self.hmdso_thickness)
        hmdso_layer.addLayout(layout)
        ptfe_layer = ba.Layer(m_PTFE, self.ptfe_thickness)
        substrate_layer = ba.Layer(m_Si)

        # assembling multilayer
        multi_layer = ba.MultiLayer()
        multi_layer.addLayer(vacuum_layer)
        multi_layer.addLayerWithTopRoughness(hmdso_layer, r_hmdso)
        multi_layer.addLayerWithTopRoughness(ptfe_layer, r_ptfe)
        multi_layer.addLayer(substrate_layer)

        return multi_layer

 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

#!/usr/bin/env python3
"""
Fitting experimental data: spherical nanoparticles with size distribution
in 3 layers system (experiment at GALAXI).
"""
import bornagain as ba
from bornagain import nm as nm
from sample_builder import SampleBuilder

wavelength = 1.34*ba.angstrom
alpha_i = 0.463*ba.deg

# detector setup as given from instrument responsible
pilatus_npx, pilatus_npy = 981, 1043
pilatus_pixel_size = 0.172  # in mm
detector_distance = 1730.0  # in mm
beam_xpos, beam_ypos = 597.1, 323.4  # in pixels


def create_detector():
    """
    Returns a model of the GALAXY detector
    """
    u0 = beam_xpos*pilatus_pixel_size  # in mm
    v0 = beam_ypos*pilatus_pixel_size  # in mm
    detector = ba.RectangularDetector(pilatus_npx,
                                      pilatus_npx*pilatus_pixel_size,
                                      pilatus_npy,
                                      pilatus_npy*pilatus_pixel_size)
    detector.setPerpendicularToDirectBeam(detector_distance, u0, v0)
    return detector


def create_simulation(params):
    """
    Creates and returns GISAS simulation with beam and detector defined
    """
    simulation = ba.GISASSimulation()
    simulation.setDetector(create_detector())
    simulation.setBeamParameters(wavelength, alpha_i, 0)
    simulation.beam().setIntensity(1.2e7)
    simulation.setRegionOfInterest(85, 70, 120, 92.)
    simulation.addMask(ba.Rectangle(101.9, 82.1, 103.7, 85.2),
                       True)  # beamstop

    sample_builder = SampleBuilder()
    sample = sample_builder.create_sample(params)
    simulation.setSample(sample)

    return simulation


def load_real_data(filename="galaxi_data.tif.gz"):
    """
    Loads experimental data and returns numpy array.
    """
    return ba.IntensityDataIOFactory.readIntensityData(filename).array()


def run_fitting():
    real_data = load_real_data()

    fit_objective = ba.FitObjective()
    fit_objective.addSimulationAndData(create_simulation, real_data, 1)
    fit_objective.initPrint(10)
    fit_objective.initPlot(10)

    params = ba.Parameters()
    params.add("radius", 5.*nm, min=4, max=6, step=0.1*nm)
    params.add("sigma", 0.55, min=0.2, max=0.8, step=0.01)
    params.add("distance", 27.*nm, min=20, max=70)

    minimizer = ba.Minimizer()
    result = minimizer.minimize(fit_objective.evaluate, params)
    fit_objective.finalize(result)


if __name__ == '__main__':
    run_fitting()