Experimental weighting

examples/joint-fit_single-dataset.py

@@ -0,0 +1,77 @@
+"""
+Joint Refinement Example (Advanced API)
+
+This example demonstrates a more flexible and advanced usage of the EasyDiffraction
+library by explicitly creating and configuring some objects. It is more suitable for
+users comfortable with Python programming and those interested in custom workflows.
+"""
+from numpy.testing import assert_almost_equal
+
+from easydiffraction import (
+    Project,
+    SampleModel,
+    Experiment
+)
+
+# Create and configure sample model
+
+model = SampleModel("pbso4")
+model.space_group.name.value = "P n m a"
+model.cell.length_a.value = 8.4693
+model.cell.length_b.value = 5.3910
+model.cell.length_c.value = 6.9506
+model.atom_sites.add("Pb", "Pb", 0.1876, 0.25, 0.167, b_iso=1.37)
+model.atom_sites.add("S", "S", 0.0654, 0.25, 0.684, b_iso=0.3777)
+model.atom_sites.add("O1", "O", 0.9082, 0.25, 0.5954, b_iso=1.9764)
+model.atom_sites.add("O2", "O", 0.1935, 0.25, 0.5432, b_iso=1.4456)
+model.atom_sites.add("O3", "O", 0.0811, 0.0272, 0.8086, b_iso=1.2822)
+
+# Create and configure experiments
+
+# Experiment: Neutron powder diffraction (full dataset)
+expt = Experiment(id="npd", radiation_probe="neutron", data_path="examples/data/pbso4_powder_neutron_cw_full.dat")
+expt.instrument.setup_wavelength = 1.91
+expt.instrument.calib_twotheta_offset = -0.1406
+expt.peak.broad_gauss_u = 0.139
+expt.peak.broad_gauss_v = -0.4124
+expt.peak.broad_gauss_w = 0.386
+expt.peak.broad_lorentz_x = 0
+expt.peak.broad_lorentz_y = 0.0878
+expt.linked_phases.add("pbso4", scale=1.46)
+expt.background_type = "line-segment"
+for x, y in [
+    (11.0, 206.1624),
+    (15.0, 194.75),
+    (20.0, 194.505),
+    (30.0, 188.4375),
+    (50.0, 207.7633),
+    (70.0, 201.7002),
+    (120.0, 244.4525),
+    (153.0, 226.0595),
+]:
+    expt.background.add(x, y)
+
+# Create project and add sample model and experiments
+project = Project()
+project.sample_models.add(model)
+project.experiments.add(expt)
+
+# Set calculator, minimizer and refinement strategy
+project.analysis.current_calculator = "cryspy"
+project.analysis.current_minimizer = "lmfit (leastsq)"
+project.analysis.fit_mode = 'joint'
+# project.analysis.joint_fit.add("expt1", weight=0.4)  # Default weight could be 0.5
+# project.analysis.joint_fit.add("expt2", weight=0.6)  # Default weight could be 0.5
+project.analysis.joint_fit.setdefault("npd1", 0.5)  # Default weight could be 0.5
+project.analysis.joint_fit.setdefault("npd2", 0.5)  # Default weight could be 0.5
+
+# Define free parameters
+model.cell.length_a.free = True
+model.cell.length_b.free = True
+model.cell.length_c.free = True
+
+# Run refinement
+project.analysis.fit()
+
+# Assert results
+assert_almost_equal(project.analysis.fit_results.reduced_chi_square, 4.66, decimal=1)

examples/joint-fit_split-single-dataset.py

@@ -0,0 +1,105 @@
+"""
+Joint Refinement Example (Advanced API)
+
+This example demonstrates a more flexible and advanced usage of the EasyDiffraction
+library by explicitly creating and configuring some objects. It is more suitable for
+users comfortable with Python programming and those interested in custom workflows.
+"""
+from numpy.testing import assert_almost_equal
+
+from easydiffraction import (
+    Project,
+    SampleModel,
+    Experiment
+)
+
+# Create and configure sample model
+
+model = SampleModel("pbso4")
+model.space_group.name.value = "P n m a"
+model.cell.length_a.value = 8.4693
+model.cell.length_b.value = 5.3910
+model.cell.length_c.value = 6.9506
+model.atom_sites.add("Pb", "Pb", 0.1876, 0.25, 0.167, b_iso=1.37)
+model.atom_sites.add("S", "S", 0.0654, 0.25, 0.684, b_iso=0.3777)
+model.atom_sites.add("O1", "O", 0.9082, 0.25, 0.5954, b_iso=1.9764)
+model.atom_sites.add("O2", "O", 0.1935, 0.25, 0.5432, b_iso=1.4456)
+model.atom_sites.add("O3", "O", 0.0811, 0.0272, 0.8086, b_iso=1.2822)
+
+# Create and configure experiments
+
+# Experiment 1: Neutron powder diffraction (first half of the dataset)
+expt1 = Experiment(id="npd1", radiation_probe="neutron", data_path="examples/data/pbso4_powder_neutron_cw_first-half.dat")
+expt1.instrument.setup_wavelength = 1.91
+expt1.instrument.calib_twotheta_offset = -0.1406
+expt1.peak.broad_gauss_u = 0.139
+expt1.peak.broad_gauss_v = -0.4124
+expt1.peak.broad_gauss_w = 0.386
+expt1.peak.broad_lorentz_x = 0
+expt1.peak.broad_lorentz_y = 0.0878
+expt1.linked_phases.add("pbso4", scale=1.46)
+expt1.background_type = "line-segment"
+for x, y in [
+    (11.0, 206.1624),
+    (15.0, 194.75),
+    (20.0, 194.505),
+    (30.0, 188.4375),
+    (50.0, 207.7633),
+    (70.0, 201.7002),
+    (120.0, 244.4525),
+    (153.0, 226.0595),
+]:
+    expt1.background.add(x, y)
+
+# Experiment 2: Neutron powder diffraction (second half of the dataset)
+expt2 = Experiment(id="npd2", radiation_probe="neutron", data_path="examples/data/pbso4_powder_neutron_cw_second-half.dat")
+expt2.instrument.setup_wavelength = 1.91
+expt2.instrument.calib_twotheta_offset = -0.1406
+expt2.peak.broad_gauss_u = 0.139
+expt2.peak.broad_gauss_v = -0.4124
+expt2.peak.broad_gauss_w = 0.386
+expt2.peak.broad_lorentz_x = 0
+expt2.peak.broad_lorentz_y = 0.0878
+expt2.linked_phases.add("pbso4", scale=1.46)
+expt2.background_type = "line-segment"
+for x, y in [
+    (11.0, 206.1624),
+    (15.0, 194.75),
+    (20.0, 194.505),
+    (30.0, 188.4375),
+    (50.0, 207.7633),
+    (70.0, 201.7002),
+    (120.0, 244.4525),
+    (153.0, 226.0595),
+]:
+    expt2.background.add(x, y)
+
+# Create project and add sample model and experiments
+project = Project()
+project.sample_models.add(model)
+project.experiments.add(expt1)
+project.experiments.add(expt2)
+
+# Set calculator, minimizer and refinement strategy
+project.analysis.current_calculator = "cryspy"
+project.analysis.current_minimizer = "lmfit (leastsq)"
+project.analysis.fit_mode = 'joint'
+# project.analysis.joint_fit.add("expt1", weight=0.4)  # Default weight could be 0.5
+# project.analysis.joint_fit.add("expt2", weight=0.6)  # Default weight could be 0.5
+print(project.analysis.joint_fit_experiments["npd1"])
+project.analysis.joint_fit_experiments["npd1"] = 0.5  # Default weight could be 0.5
+project.analysis.joint_fit_experiments["npd2"] = 0.5  # Default weight could be 0.5
+
+# Define free parameters
+model.cell.length_a.free = True
+model.cell.length_b.free = True
+model.cell.length_c.free = True
+#expt1.linked_phases["pbso4"].scale.free = True
+#expt2.linked_phases["pbso4"].scale.free = True
+
+# Run refinement
+project.analysis.fit()
+
+# Assert results
+assert_almost_equal(project.analysis.fit_results.reduced_chi_square, 4.66, decimal=1)
+

examples/joint-refinement_weighting.py

@@ -0,0 +1,98 @@
+"""
+Joint Refinement Example (Advanced API)
+
+This example demonstrates a more flexible and advanced usage of the EasyDiffraction
+library by explicitly creating and configuring some objects. It is more suitable for
+users comfortable with Python programming and those interested in custom workflows.
+"""
+
+from easydiffraction import (
+    Project,
+    SampleModel,
+    Experiment
+)
+
+# Create and configure sample model
+
+model = SampleModel("pbso4")
+model.space_group.name.value = "P n m a"
+model.cell.length_a.value = 8.5
+model.cell.length_b.value = 5.35
+model.cell.length_c.value = 6.9
+model.atom_sites.add("Pb", "Pb", 0.1876, 0.25, 0.167, b_iso=1.37)
+model.atom_sites.add("S", "S", 0.0654, 0.25, 0.684, b_iso=0.3777)
+model.atom_sites.add("O1", "O", 0.9082, 0.25, 0.5954, b_iso=1.9764)
+model.atom_sites.add("O2", "O", 0.1935, 0.25, 0.5432, b_iso=1.4456)
+model.atom_sites.add("O3", "O", 0.0811, 0.0272, 0.8086, b_iso=1.2822)
+
+# Create and configure experiments
+
+# Experiment 1: Neutron powder diffraction
+expt1 = Experiment(id="npd", radiation_probe="neutron", data_path="examples/data/d1a_pbso4.dat")
+expt1.instrument.setup_wavelength = 1.91
+expt1.instrument.calib_twotheta_offset = -0.1406
+expt1.peak.broad_gauss_u = 0.139
+expt1.peak.broad_gauss_v = -0.412
+expt1.peak.broad_gauss_w = 0.386
+expt1.peak.broad_lorentz_x = 0
+expt1.peak.broad_lorentz_y = 0.088
+expt1.linked_phases.add("pbso4", scale=1.0)
+expt1.background_type = "line-segment"
+for x, y in [
+    (11.0, 206.1624),
+    (15.0, 194.75),
+    (20.0, 194.505),
+    (30.0, 188.4375),
+    (50.0, 207.7633),
+    (70.0, 201.7002),
+    (120.0, 244.4525),
+    (153.0, 226.0595),
+]:
+    expt1.background.add(x, y)
+
+# Experiment 2: X-ray powder diffraction
+expt2 = Experiment(id="xrd", radiation_probe="xray", data_path="examples/data/lab_pbso4.dat")
+expt2.instrument.setup_wavelength = 1.540567
+expt2.instrument.calib_twotheta_offset = -0.05181
+expt2.peak.broad_gauss_u = 0.304138
+expt2.peak.broad_gauss_v = -0.112622
+expt2.peak.broad_gauss_w = 0.021272
+expt2.peak.broad_lorentz_x = 0
+expt2.peak.broad_lorentz_y = 0.057691
+expt2.linked_phases.add("pbso4", scale=0.005)
+expt2.background_type = "chebyshev polynomial"
+for x, y in [
+    (0, 119.195),
+    (1, 6.221),
+    (2, -45.725),
+    (3, 8.119),
+    (4, 54.552),
+    (5, -20.661),
+]:
+    expt2.background.add(x, y)
+
+# Create project and add sample model and experiments
+project = Project()
+project.sample_models.add(model)
+project.experiments.add(expt1)
+project.experiments.add(expt2)
+
+# Set calculator, minimizer and refinement strategy
+project.analysis.current_calculator = "cryspy"
+project.analysis.current_minimizer = "lmfit (leastsq)"
+project.analysis.fit_mode = 'joint'
+# project.analysis.joint_fit.add("expt1", weight=0.4)  # Default weight could be 0.5
+# project.analysis.joint_fit.add("expt2", weight=0.6)  # Default weight could be 0.5
+project.analysis.joint_fit.setdefault("xrd",0.4)  # Default weight could be 0.5
+project.analysis.joint_fit.setdefault("npd",0.6)  # Default weight could be 0.5
+
+# Define free parameters
+model.cell.length_a.free = True
+model.cell.length_b.free = True
+model.cell.length_c.free = True
+expt1.linked_phases["pbso4"].scale.free = True
+expt2.linked_phases["pbso4"].scale.free = True
+
+# Run refinement
+project.analysis.fit()
+

src/easydiffraction/analysis/analysis.py

@@ -8,6 +8,7 @@
 from .calculators.calculator_factory import CalculatorFactory
 from .minimization import DiffractionMinimizer
 from .minimizers.minimizer_factory import MinimizerFactory
+from .joint_fit_experiments import JointFitExperiments
 
 
 class Analysis:
@@ -118,7 +119,13 @@ def fit_mode(self, strategy):
         if strategy not in ['single', 'joint']:
             raise ValueError("Fit mode must be either 'single' or 'joint'")
         self._fit_mode = strategy
-        print(paragraph("Current ffit mode changed to"))
+        if strategy == 'joint':
+            if not hasattr(self, 'joint_fit_experiments'):
+                # Pre-populate all experiments with weight 0.5
+                self.joint_fit_experiments = JointFitExperiments()
+                for id in self.project.experiments.ids:
+                    self.joint_fit_experiments.add(id, weight=0.5)
+        print(paragraph("Current fit mode changed to"))
         print(self._fit_mode)
 
     def show_available_fit_modes(self):
@@ -212,11 +219,11 @@ def fit(self):
             return
 
         # Run the fitting process
-        experiment_ids = list(experiments._items.keys())
+        experiment_ids = experiments.ids
 
         if self.fit_mode == 'joint':
             print(paragraph(f"Using all experiments 🔬 {experiment_ids} for '{self.fit_mode}' fitting"))
-            self.fitter.fit(sample_models, experiments, calculator)
+            self.fitter.fit(sample_models, experiments, calculator, weights=self.joint_fit_experiments)
         elif self.fit_mode == 'single':
             for expt_id in list(experiments._items.keys()):
                 print(paragraph(f"Using experiment 🔬 '{expt_id}' for '{self.fit_mode}' fitting"))
@@ -249,4 +256,4 @@ def show_as_cif(self):
         print(paragraph(f"Analysis 🧮 info as cif"))
         print(top)
         print("\n".join(padded_lines))
-        print(bottom)
+        print(bottom)

src/easydiffraction/analysis/joint_fit_experiments.py

@@ -0,0 +1,29 @@
+from easydiffraction.core.collection import Collection
+from easydiffraction.core.parameter import Descriptor
+
+class JointFitExperiments(Collection):
+    """
+    Collection manager for experiments that are fitted together
+    in a `joint` fit.
+    """
+    def __init__(self):
+        super().__init__()
+
+    def add(self, id: str, weight: float):
+        """Add an experiment with it's associated weight"""
+        # Save both id and weight as immutable Descriptors
+        self.id = Descriptor(
+            value=id,
+            cif_name="id"
+        )
+        self.weight = Descriptor(
+            value=weight,
+            cif_name="weight"
+        )
+        self._items[id] = weight
+
+    def __getitem__(self, key):
+        return self._items[key]
+
+    def __setitem__(self, key, value):
+        self._items[key] = value