diff --git a/CHANGELOG.rst b/CHANGELOG.rst
index cbc49cfc..f6b763ac 100644
--- a/CHANGELOG.rst
+++ b/CHANGELOG.rst
@@ -1,6 +1,11 @@
 Changelog
 =========
 
+Release 2.1.17
+--------------
+
+* Added a new function to automatically determine a diffraction mask, :func:`auto_masking` (#46).
+
 Release 2.1.16
 --------------
 
diff --git a/skued/__init__.py b/skued/__init__.py
index fdcc97b5..0f4845a9 100644
--- a/skued/__init__.py
+++ b/skued/__init__.py
@@ -2,7 +2,7 @@
 __author__ = "Laurent P. René de Cotret"
 __email__ = "laurent.renedecotret@mail.mcgill.ca"
 __license__ = "GPLv3"
-__version__ = "2.1.16"
+__version__ = "2.1.17"
 
 from .affine import (
     affine_map,
diff --git a/skued/image/center.py b/skued/image/center.py
index 0db83606..38f115f2 100644
--- a/skued/image/center.py
+++ b/skued/image/center.py
@@ -117,36 +117,25 @@ def autocenter(im, mask=None, normalize_bg=True):
     return np.array([r_, c_]) + correction / 2
 
 
-def _center_of_intensity(im, mask=None):
-    weights = im * mask.astype(im.dtype)
-
-    rr, cc = np.indices(im.shape)
-    r_ = np.average(rr, weights=weights)
-    c_ = np.average(cc, weights=weights)
-    return int(r_), int(c_)
-
-
 def auto_masking(im, threshold=0.1):
     """
-    Generate a mask based on the darkest fraction of an image
+    Generate a mask based on the darkest fraction of an image.
+
+    .. versionadded:: 2.1.17
 
     Parameters
     ----------
-    im : floats, ndarrays of shape (N,M)
+    im : ndarray of shape (N,M)
         image used to generate a mask
     threshold: float, optional
-        fraction of the lowest values to be masked, default = 15%
+        fraction of the lowest values to be masked, default = 10%
 
-    Yields
-    ------
+    Returns
+    -------
     mask : boolean, ndarrays of shape (N,M)
         Mask that evaluates to True on valid pixels.
-
     """
     # Find the median of the highest intensity value of the image to avoid hot spots
-    max_median = np.median([max(x) for x in np.real(im)])
-    # Set the threshold value
-    lower_limit = threshold * max_median
+    lower_limit = threshold * np.median(np.maximum(im, 0))
     # generate a mask
-    mask = im >= lower_limit
-    return mask
+    return im > lower_limit
diff --git a/skued/image/tests/test_center.py b/skued/image/tests/test_center.py
index 716d2943..0925ed66 100644
--- a/skued/image/tests/test_center.py
+++ b/skued/image/tests/test_center.py
@@ -1,10 +1,9 @@
 # -*- coding: utf-8 -*-
 
-from skued import autocenter, gaussian, kinematicsim
+from skued import autocenter, gaussian, kinematicsim, auto_masking
 from crystals import Crystal
 from scipy.ndimage import gaussian_filter
 import numpy as np
-import itertools as it
 import pytest
 
 np.random.seed(23)
@@ -124,3 +123,14 @@ def test_autocenter_single_crystal_ewald_walkoff(rc, cc):
     I += 0.01 * I.max() * np.random.random(size=I.shape)
 
     assert np.allclose(autocenter(I, mask=mask), (rc, cc), atol=1)
+
+
+@pytest.mark.parametrize("seed", range(10))
+def test_auto_masking(seed: int):
+    im = 10 * np.ones(shape=(32, 32))
+    np.random.seed(seed)
+    reference_mask = np.random.randint(low=0, high=1 + 1, size=im.shape)
+
+    mask = auto_masking(im=im * reference_mask).astype(int)
+
+    assert np.allclose(reference_mask, mask)