Add Triangle thresholding algorithm
authorMiroslav Benes <mbenes@suse.cz>
Wed, 22 Oct 2014 14:54:39 +0000 (16:54 +0200)
committerMiroslav Benes <mbenes@suse.cz>
Wed, 22 Oct 2014 14:54:39 +0000 (16:54 +0200)
Add Triangle method for automatic threshold computation next to the existing
Otsu's method. Triangle deals better with images whose histogram does not
contain dominant peak.

See paper Zack GW, Rogers WE, Latt SA.: Automatic measurement of sister
chromatid exchange frequency. J Histochem Cytochem. 1977 Jul;25(7):741-53.

modules/cudaarithm/doc/element_operations.rst
modules/imgproc/doc/miscellaneous_transformations.rst
modules/imgproc/include/opencv2/imgproc.hpp
modules/imgproc/include/opencv2/imgproc/types_c.h
modules/imgproc/src/thresh.cpp

index d3dba67..5ccca78 100644 (file)
@@ -430,7 +430,7 @@ Applies a fixed-level threshold to each array element.
 
     :param maxval: Maximum value to use with  ``THRESH_BINARY`` and  ``THRESH_BINARY_INV`` threshold types.
 
-    :param type: Threshold type. For details, see  :ocv:func:`threshold` . The ``THRESH_OTSU`` threshold type is not supported.
+    :param type: Threshold type. For details, see  :ocv:func:`threshold` . The ``THRESH_OTSU`` and ``THRESH_TRIANGLE`` threshold types are not supported.
 
     :param stream: Stream for the asynchronous version.
 
index bbd692e..b0dedfa 100644 (file)
@@ -712,11 +712,11 @@ types of thresholding supported by the function. They are determined by ``type``
 
               \texttt{dst} (x,y) =  \fork{0}{if $\texttt{src}(x,y) > \texttt{thresh}$}{\texttt{src}(x,y)}{otherwise}
 
-Also, the special value ``THRESH_OTSU`` may be combined with
-one of the above values. In this case, the function determines the optimal threshold
-value using the Otsu's algorithm and uses it instead of the specified ``thresh`` .
+Also, the special values ``THRESH_OTSU`` or ``THRESH_TRIANGLE`` may be combined with
+one of the above values. In these cases, the function determines the optimal threshold
+value using the Otsu's or Triangle algorithm and uses it instead of the specified ``thresh`` .
 The function returns the computed threshold value.
-Currently, the Otsu's method is implemented only for 8-bit images.
+Currently, the Otsu's and Triangle methods are implemented only for 8-bit images.
 
 
 .. image:: pics/threshold.png
index f23080a..470e28f 100644 (file)
@@ -109,7 +109,8 @@ enum { THRESH_BINARY     = 0, // value = value > threshold ? max_value : 0
        THRESH_TOZERO     = 3, // value = value > threshold ? value : 0
        THRESH_TOZERO_INV = 4, // value = value > threshold ? 0 : value
        THRESH_MASK       = 7,
-       THRESH_OTSU       = 8  // use Otsu algorithm to choose the optimal threshold value
+       THRESH_OTSU       = 8, // use Otsu algorithm to choose the optimal threshold value
+       THRESH_TRIANGLE   = 16  // use Triangle algorithm to choose the optimal threshold value
      };
 
 //! adaptive threshold algorithm
index de8fb62..544dee0 100644 (file)
@@ -551,8 +551,11 @@ enum
     CV_THRESH_TOZERO      =3,  /* value = value > threshold ? value : 0           */
     CV_THRESH_TOZERO_INV  =4,  /* value = value > threshold ? 0 : value           */
     CV_THRESH_MASK        =7,
-    CV_THRESH_OTSU        =8  /* use Otsu algorithm to choose the optimal threshold value;
+    CV_THRESH_OTSU        =8, /* use Otsu algorithm to choose the optimal threshold value;
                                  combine the flag with one of the above CV_THRESH_* values */
+    CV_THRESH_TRIANGLE    =16  /* use Triangle algorithm to choose the optimal threshold value;
+                                 combine the flag with one of the above CV_THRESH_* values, but not
+                                 with CV_THRESH_OTSU */
 };
 
 /* Adaptive threshold methods */
index 6e0639c..08feaba 100644 (file)
@@ -986,6 +986,110 @@ getThreshVal_Otsu_8u( const Mat& _src )
     return max_val;
 }
 
+static double
+getThreshVal_Triangle_8u( const Mat& _src )
+{
+    Size size = _src.size();
+    int step = (int) _src.step;
+    if( _src.isContinuous() )
+    {
+        size.width *= size.height;
+        size.height = 1;
+        step = size.width;
+    }
+
+    const int N = 256;
+    int i, j, h[N] = {0};
+    for( i = 0; i < size.height; i++ )
+    {
+        const uchar* src = _src.ptr() + step*i;
+        j = 0;
+        #if CV_ENABLE_UNROLLED
+        for( ; j <= size.width - 4; j += 4 )
+        {
+            int v0 = src[j], v1 = src[j+1];
+            h[v0]++; h[v1]++;
+            v0 = src[j+2]; v1 = src[j+3];
+            h[v0]++; h[v1]++;
+        }
+        #endif
+        for( ; j < size.width; j++ )
+            h[src[j]]++;
+    }
+
+    int left_bound = 0, right_bound = 0, max_ind = 0, max = 0;
+    int temp;
+    bool isflipped = false;
+
+    for( i = 0; i < N; i++ )
+    {
+        if( h[i] > 0 )
+        {
+            left_bound = i;
+            break;
+        }
+    }
+    if( left_bound > 0 )
+        left_bound--;
+
+    for( i = N-1; i > 0; i-- )
+    {
+        if( h[i] > 0 )
+        {
+            right_bound = i;
+            break;
+        }
+    }
+    if( right_bound < N-1 )
+        right_bound++;
+
+    for( i = 0; i < N; i++ )
+    {
+        if( h[i] > max)
+        {
+            max = h[i];
+            max_ind = i;
+        }
+    }
+
+    if( max_ind-left_bound < right_bound-max_ind)
+    {
+        isflipped = true;
+        i = 0, j = N-1;
+        while( i < j )
+        {
+            temp = h[i]; h[i] = h[j]; h[j] = temp;
+            i++; j--;
+        }
+        left_bound = N-1-right_bound;
+        max_ind = N-1-max_ind;
+    }
+
+    double thresh = left_bound;
+    double a, b, dist = 0, tempdist;
+
+    /*
+     * We do not need to compute precise distance here. Distance is maximized, so some constants can
+     * be omitted. This speeds up a computation a bit.
+     */
+    a = max; b = left_bound-max_ind;
+    for( i = left_bound+1; i <= max_ind; i++ )
+    {
+        tempdist = a*i + b*h[i];
+        if( tempdist > dist)
+        {
+            dist = tempdist;
+            thresh = i;
+        }
+    }
+    thresh--;
+
+    if( isflipped )
+        thresh = N-1-thresh;
+
+    return thresh;
+}
+
 class ThresholdRunner : public ParallelLoopBody
 {
 public:
@@ -1086,6 +1190,7 @@ double cv::threshold( InputArray _src, OutputArray _dst, double thresh, double m
 
     Mat src = _src.getMat();
     bool use_otsu = (type & THRESH_OTSU) != 0;
+    bool use_triangle = (type & THRESH_TRIANGLE) != 0;
     type &= THRESH_MASK;
 
     if( use_otsu )
@@ -1094,6 +1199,12 @@ double cv::threshold( InputArray _src, OutputArray _dst, double thresh, double m
         thresh = getThreshVal_Otsu_8u(src);
     }
 
+    if( use_triangle )
+    {
+        CV_Assert( src.type() == CV_8UC1 );
+        thresh = getThreshVal_Triangle_8u(src);
+    }
+
     _dst.create( src.size(), src.type() );
     Mat dst = _dst.getMat();