Commit a2d3f259 authored by Wuttke, Joachim's avatar Wuttke, Joachim
Browse files

gisas model with cylinders only

parent 8efbc0c2
......@@ -15,10 +15,10 @@ def fake_data():
"""
params = {
'cylinder_height': 5*nm,
'lg(intensity)': 6,
'lg(background)': -99, # background is not part of theory, but added below
'cylinder_height': 4*nm,
'cylinder_radius': 5*nm,
'prism_height': 5*nm,
'prism_base_edge': 5*nm
}
# Compute model distribution
......@@ -28,7 +28,7 @@ def fake_data():
# Draw noisy data
np.random.seed(0)
data = np.random.normal(theory, np.sqrt(theory)) + np.random.poisson(3.57, theory.shape)
data = np.random.normal(theory, np.sqrt(theory)) + np.random.poisson(10**0.5, theory.shape)
# Save to numpy
np.savetxt("gisas-model1.txt.gz", data)
......
......@@ -16,8 +16,6 @@ def get_sample(params):
"""
cylinder_height = params["cylinder_height"]
cylinder_radius = params["cylinder_radius"]
prism_height = params["prism_height"]
prism_base_edge = params["prism_base_edge"]
# defining materials
m_vacuum = ba.HomogeneousMaterial("Vacuum", 0, 0)
......@@ -25,13 +23,10 @@ def get_sample(params):
m_particle = ba.HomogeneousMaterial("Particle", 6e-4, 2e-8)
# collection of particles
cylinder_ff = ba.FormFactorCylinder(cylinder_radius, cylinder_height)
cylinder = ba.Particle(m_particle, cylinder_ff)
prism_ff = ba.FormFactorPrism3(prism_base_edge, prism_height)
prism = ba.Particle(m_particle, prism_ff)
ff = ba.FormFactorCylinder(cylinder_radius, cylinder_height)
cylinder = ba.Particle(m_particle, ff)
layout = ba.ParticleLayout()
layout.addParticle(cylinder, 0.5)
layout.addParticle(prism, 0.5)
layout.addParticle(cylinder)
# vacuum layer with particles and substrate form multi layer
vacuum_layer = ba.Layer(m_vacuum)
......@@ -47,19 +42,14 @@ def get_simulation(params):
"""
Returns a GISAXS simulation with beam and detector defined
"""
intensity = 1e6
beam = ba.Beam(intensity, 0.1*nm, ba.Direction(0.2*deg, 0))
det = ba.SphericalDetector(100, -1*deg, 1*deg, 100, 0, 2*deg)
beam = ba.Beam(10**params['lg(intensity)'], 0.1*nm, ba.Direction(0.2*deg, 0))
det = ba.SphericalDetector(100, -1.5*deg, 1.5*deg, 100, 0, 3*deg)
sample = get_sample(params)
return ba.GISASSimulation(beam, sample, det)
simulation = ba.GISASSimulation(beam, sample, det)
simulation.setBackground(ba.ConstantBackground(10**params['lg(background)']))
def load_data():
"""
Loads experimental data from file
"""
return np.loadtxt("gisas-model1.txt.gz", dtype=float)
return simulation
def run_fitting():
......@@ -67,7 +57,7 @@ def run_fitting():
Setup simulation and fit
"""
real_data = load_data()
real_data = np.loadtxt("gisas-model1.txt.gz", dtype=float)
fit_objective = ba.FitObjective()
fit_objective.addSimulationAndData(get_simulation, real_data)
......@@ -79,10 +69,10 @@ def run_fitting():
fit_objective.initPlot(10)
params = ba.Parameters()
params.add("cylinder_height", 4.*nm, min=0.01)
params.add("lg(intensity)", 5)
params.add("lg(background)", 1)
params.add("cylinder_height", 6.*nm, min=0.01)
params.add("cylinder_radius", 6.*nm, min=0.01)
params.add("prism_height", 4.*nm, min=0.01)
params.add("prism_base_edge", 6.*nm, min=0.01)
minimizer = ba.Minimizer()
result = minimizer.minimize(fit_objective.evaluate, params)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment