From 5a988fdae79c66c0b8ffec6a1328be42bd8b0c32 Mon Sep 17 00:00:00 2001
From: Vighnesh Birodkar <vighneshbirodkar@gmail.com>
Date: Wed, 5 Mar 2014 00:46:58 +0530
Subject: [PATCH] Changed formula, this improves detected size accuracy

---
 skimage/feature/blob.py            | 26 +++++++++++++-------------
 skimage/feature/tests/test_blob.py | 28 +++++++++++++++++++---------
 2 files changed, 32 insertions(+), 22 deletions(-)

diff --git a/skimage/feature/blob.py b/skimage/feature/blob.py
index 11d4c717..4ba18b33 100644
--- a/skimage/feature/blob.py
+++ b/skimage/feature/blob.py
@@ -123,7 +123,7 @@ def _prune_blobs(blobs_array, overlap):
     return np.array([b for b in blobs_array if b[2] > 0])
 
 
-def blob_dog(image, min_sigma=1, max_sigma=25, sigma_ratio=1.6, threshold=2.0,
+def blob_dog(image, min_sigma=1, max_sigma=50, sigma_ratio=1.6, threshold=2.0,
              overlap=.5,):
     """Finds blobs in the given grayscale image.
 
@@ -165,30 +165,30 @@ def blob_dog(image, min_sigma=1, max_sigma=25, sigma_ratio=1.6, threshold=2.0,
     Examples
     --------
     >>> from skimage import data, feature
-    >>> feature.blob_dog(data.coins())
+    >>> feature.blob_dog(data.coins(),threshold=.5,max_sigma=40)
     array([[  45,  336, 1608],
-           [  51,  277, 1608],
            [  52,  155, 1608],
            [  52,  216, 1608],
            [  54,   42, 1608],
-           [  56,  101, 1608],
+           [  54,  276,  628],
+           [  58,  100,  628],
            [ 120,  272, 1608],
-           [ 124,  206, 1608],
-           [ 124,  339, 1608],
+           [ 124,  337,  628],
            [ 125,   45, 1608],
-           [ 125,  102, 1608],
-           [ 127,  154, 1608],
+           [ 125,  208,  628],
+           [ 127,  102,  628],
+           [ 128,  154,  628],
            [ 185,  347, 1608],
            [ 193,  213, 1608],
            [ 194,  277, 1608],
            [ 195,  102, 1608],
-           [ 196,   41, 1608],
-           [ 197,  154, 1608],
+           [ 196,   43,  628],
+           [ 198,  155,  628],
            [ 260,   46, 1608],
            [ 261,  173, 1608],
            [ 263,  245, 1608],
            [ 263,  302, 1608],
-           [ 267,  115, 1608],
+           [ 267,  115,  628],
            [ 267,  359, 1608]])
 
     """
@@ -208,9 +208,9 @@ def blob_dog(image, min_sigma=1, max_sigma=25, sigma_ratio=1.6, threshold=2.0,
     gaussian_images = [gaussian_filter(image, s) for s in sigma_list]
 
     # computing difference between two successive Gaussian blurred images
-    # multiplying with square of standard deviation provides scale invariance
+    # multiplying with standard deviation provides scale invariance
     dog_images = [(gaussian_images[i] - gaussian_images[i + 1])
-                  * sigma_list[i] ** 2 for i in range(k)]
+                  * sigma_list[i] for i in range(k)]
     image_cube = np.dstack(dog_images)
 
     local_maxima = _get_local_maxima_3d(image_cube, threshold)
diff --git a/skimage/feature/tests/test_blob.py b/skimage/feature/tests/test_blob.py
index e1d2baf9..54ed38ce 100644
--- a/skimage/feature/tests/test_blob.py
+++ b/skimage/feature/tests/test_blob.py
@@ -1,6 +1,7 @@
 import numpy as np
 from skimage.draw import circle
 from skimage.feature import blob_dog
+import math
 
 
 def test_blob_dog():
@@ -12,17 +13,26 @@ def test_blob_dog():
     xs, ys = circle(100, 300, 25)
     img[xs, ys] = 255
 
-    xs, ys = circle(200, 350, 30)
+    xs, ys = circle(200, 350, 45)
     img[xs, ys] = 255
 
-    blobs = blob_dog(img)
-    coords = blobs[:, 0:2]
+    blobs = blob_dog(img, min_sigma=5, max_sigma=50)
+    area = lambda x: x[2]
+    radius = lambda x: math.sqrt(x / math.pi)
+    s = sorted(blobs, key=area)
+    thresh = 5
 
-    if coords[coords == [400, 130]].shape != (2,):
-        assert False
+    b = s[0]
+    assert abs(b[0] - 400) <= thresh
+    assert abs(b[1] - 130) <= thresh
+    assert abs(radius(b[2]) - 5) <= thresh
 
-    if coords[coords == [100, 300]].shape != (2,):
-        assert False
+    b = s[1]
+    assert abs(b[0] - 100) <= thresh
+    assert abs(b[1] - 300) <= thresh
+    assert abs(radius(b[2]) - 25) <= thresh
 
-    if coords[coords == [200, 350]].shape != (2,):
-        assert False
+    b = s[2]
+    assert abs(b[0] - 200) <= thresh
+    assert abs(b[1] - 350) <= thresh
+    assert abs(radius(b[2]) - 45) <= thresh