Fix bug where direct beam is added twice

diffsims/generators/simulation_generator.py

@@ -211,7 +211,7 @@ def calculate_diffraction2d(
                     wavelength,
                     max_excitation_error,
                     shape_factor_width=shape_factor_width,
-                    with_direct_beam=with_direct_beam,
+                    with_direct_beam=False,  # Already handled in DiffractingVector.from_min_dspacing
                 )
 
                 # Calculate diffracted intensities based on a kinematic model.
@@ -223,7 +223,6 @@ def calculate_diffraction2d(
                     scattering_params=self.scattering_params,
                     debye_waller_factors=debye_waller_factors,
                 )
-
                 # Threshold peaks included in simulation as factor of zero beam intensity.
                 peak_mask = intensities > np.max(intensities) * self.minimum_intensity
                 intensities = intensities[peak_mask]
@@ -375,9 +374,9 @@ def get_intersecting_reflections(
         o = recip.data
         if with_direct_beam:
             o = np.vstack([o, [0, 0, 0]])
-            recip_hkl = np.vstack([recip_hkl, [0, 0, 0]])
             # Modify the data directly, as `data` is a property
             recip._data = np.vstack([recip._data, [0, 0, 0]])
+            recip_hkl = np.vstack([recip_hkl, [0, 0, 0]])
 
         diff = o - c
         dot = np.dot(u, diff.T)
@@ -437,7 +436,6 @@ def get_d(_c):
             intersection = (d < (upper + max_excitation_error)) & (
                 d > (lower - max_excitation_error)
             )
-
         # select these reflections
         intersected_vectors = ~rot * (recip[intersection].to_miller())
         intersected_vectors = DiffractingVector(

diffsims/tests/generators/test_simulation_generator.py

@@ -133,13 +133,30 @@ def test_repr(self, diffraction_calculator):
             "approximate_precession=True)"
         )
 
+    @pytest.mark.parametrize("with_direct_beam", [False, True])
+    def test_single_direct_beam(
+        self,
+        diffraction_calculator: SimulationGenerator,
+        local_structure,
+        with_direct_beam,
+    ):
+        diffraction = diffraction_calculator.calculate_diffraction2d(
+            local_structure, reciprocal_radius=5.0, with_direct_beam=with_direct_beam
+        )
+        # Check that there is only one direct beam
+        num_direct_beam = 0
+        for hkl in diffraction.coordinates.hkl:
+            hkl = tuple(hkl.flatten().tolist())
+            num_direct_beam += hkl == (0, 0, 0)
+        assert num_direct_beam == with_direct_beam
+
     def test_matching_results(
         self, diffraction_calculator: SimulationGenerator, local_structure
     ):
         diffraction = diffraction_calculator.calculate_diffraction2d(
             local_structure, reciprocal_radius=5.0
         )
-        assert diffraction.coordinates.size == 70
+        assert diffraction.coordinates.size == 69
 
     def test_precession_simple(
         self, diffraction_calculator_precession_simple, local_structure
@@ -148,7 +165,7 @@ def test_precession_simple(
             local_structure,
             reciprocal_radius=5.0,
         )
-        assert diffraction.coordinates.size == 250
+        assert diffraction.coordinates.size == 249
 
     def test_precession_full(
         self, diffraction_calculator_precession_full, local_structure
@@ -157,13 +174,15 @@ def test_precession_full(
             local_structure,
             reciprocal_radius=5.0,
         )
-        assert diffraction.coordinates.size == 250
+        assert diffraction.coordinates.size == 249
 
     def test_custom_shape_func(self, diffraction_calculator_custom, local_structure):
         diffraction = diffraction_calculator_custom.calculate_diffraction2d(
             local_structure,
             reciprocal_radius=5.0,
         )
+        # This shape factor model yields 0 intensity for the direct beam,
+        # which is less than the intensity threshold. The direct beam is therefore removed
         assert diffraction.coordinates.size == 52
 
     def test_appropriate_scaling(self, diffraction_calculator: SimulationGenerator):