diff --git a/Examples/fit55_SpecularIntro/FitWithUncertainties.py b/Examples/fit55_SpecularIntro/FitWithUncertainties.py
index f0f20cd8ed48e52cf534cd92e7f3fa0fadca5cc0..f554228c0c52f29ac575bcbfbcbee90296c8b93f 100755
--- a/Examples/fit55_SpecularIntro/FitWithUncertainties.py
+++ b/Examples/fit55_SpecularIntro/FitWithUncertainties.py
@@ -18,105 +18,20 @@ from matplotlib import pyplot as plt
 import bornagain as ba
 from bornagain import ba_fitmonitor
 
-
-def get_sample(params):
-    """
-    Creates a sample and returns it
-    :param params: a dictionary of optimization parameters
-    :return: the sample defined
-    """
-
-    # substrate (Si)
-    si_sld_real = 2.0704e-06  # \AA^{-2}
-
-    # layers' parameters
-    n_repetitions = 10
-    # Ni
-    ni_sld_real = 9.4245e-06  # \AA^{-2}
-    ni_thickness = 70*ba.angstrom
-    # Ti
-    ti_sld_real = -1.9493e-06  # \AA^{-2}
-    ti_thickness = params["ti_thickness"]
-
-    # defining materials
-    m_vacuum = ba.MaterialBySLD()
-    m_ni = ba.MaterialBySLD("Ni", ni_sld_real, 0)
-    m_ti = ba.MaterialBySLD("Ti", ti_sld_real, 0)
-    m_substrate = ba.MaterialBySLD("SiSubstrate", si_sld_real, 0)
-
-    # vacuum layer and substrate form multi layer
-    vacuum_layer = ba.Layer(m_vacuum)
-    ni_layer = ba.Layer(m_ni, ni_thickness)
-    ti_layer = ba.Layer(m_ti, ti_thickness)
-    substrate_layer = ba.Layer(m_substrate)
-    multi_layer = ba.MultiLayer()
-    multi_layer.addLayer(vacuum_layer)
-    for _ in range(n_repetitions):
-        multi_layer.addLayer(ti_layer)
-        multi_layer.addLayer(ni_layer)
-    multi_layer.addLayer(substrate_layer)
-    return multi_layer
-
-
-def get_real_data(filename):
-    """
-    Loading data from genx_interchanging_layers.dat
-    Returns a Nx2 array (N - the number of experimental data entries)
-    with first column being coordinates,
-    second one being values.
-    """
-
-    real_data = np.loadtxt(filename, usecols=(0, 1), skiprows=3)
-
-    # translating axis values from double incident angle (degs)
-    # to incident angle (radians)
-    real_data[:, 0] *= np.pi/360
-
-    global expdata
-    expdata = real_data.copy()
-
-
-def get_real_data_axis():
-    """
-    Get axis coordinates of the experimental data
-    :return: 1D array with axis coordinates
-    """
-    return expdata[:, 0]
-
-
-def get_real_data_values():
-    """
-    Get experimental data values as a 1D array
-    :return: 1D array with experimental data values
-    """
-    return expdata[:, 1]
-
-
-def get_simulation(params):
-    """
-    Create and return specular simulation with its instrument defined
-    """
-    wavelength = 1.54*ba.angstrom  # beam wavelength
-
-    simulation = ba.SpecularSimulation()
-    scan = ba.AlphaScan(wavelength, get_real_data_axis())
-    simulation.setScan(scan)
-    simulation.setSample(get_sample(params))
-    return simulation
-
+import FitSpecularBasics as fsb
 
 def run_fitting():
     """
     Setup simulation and fit
     """
 
-    real_data = get_real_data_values()
+    real_data = fsb.get_real_data_values()
     # setting artificial uncertainties (uncertainty sigma equals a half
     # of experimental data value)
     uncertainties = real_data*0.5
 
     fit_objective = ba.FitObjective()
-    fit_objective.addSimulationAndData(get_simulation, real_data,
+    fit_objective.addSimulationAndData(fsb.get_simulation, real_data,
                                        uncertainties)
 
     plot_observer = ba_fitmonitor.PlotterSpecular(units=ba.Axes.RQ4)
@@ -140,6 +55,6 @@ def run_fitting():
 if __name__ == '__main__':
     datadir = os.getenv('BORNAGAIN_EXAMPLE_DATA_DIR', '')
     data_fname = os.path.join(datadir, "genx_interchanging_layers.dat.gz")
-    get_real_data(data_fname)
+    fsb.get_real_data(data_fname)
     run_fitting()
     plt.show()