Opencv Threshold with Mask

OpenCV threshold with mask

In general, you can simply compute the threshold using cv::threshold, and then copy the src image on dst using the inverted mask.

// Apply cv::threshold on all image
thresh = cv::threshold(src, dst, thresh, maxval, type);

// Copy original image on inverted mask
src.copyTo(dst, ~mask);

With THRESH_OTSU, however, you also need to compute the threshold value only on the masked image. The following code is a modified version of static double getThreshVal_Otsu_8u(const Mat& _src) in thresh.cpp:

double otsu_8u_with_mask(const Mat1b src, const Mat1b& mask)
{
    const int N = 256;
    int M = 0;
    int i, j, h[N] = { 0 };
    for (i = 0; i < src.rows; i++)
    {
        const uchar* psrc = src.ptr(i);
        const uchar* pmask = mask.ptr(i);
        for (j = 0; j < src.cols; j++)
        {
            if (pmask[j])
            {
                h[psrc[j]]++;
                ++M;
            }
        }
    }

    double mu = 0, scale = 1. / (M);
    for (i = 0; i < N; i++)
        mu += i*(double)h[i];

    mu *= scale;
    double mu1 = 0, q1 = 0;
    double max_sigma = 0, max_val = 0;

    for (i = 0; i < N; i++)
    {
        double p_i, q2, mu2, sigma;

        p_i = h[i] * scale;
        mu1 *= q1;
        q1 += p_i;
        q2 = 1. - q1;

        if (std::min(q1, q2) < FLT_EPSILON || std::max(q1, q2) > 1. - FLT_EPSILON)
            continue;

        mu1 = (mu1 + i*p_i) / q1;
        mu2 = (mu - q1*mu1) / q2;
        sigma = q1*q2*(mu1 - mu2)*(mu1 - mu2);
        if (sigma > max_sigma)
        {
            max_sigma = sigma;
            max_val = i;
        }
    }
    return max_val;
}

You then can wrap all in a function, here called threshold_with_mask, that wraps all different cases for you. If there is no mask, or the mask is all-white, then use cv::threshold. Otherwise, use one of the above mentioned approaches. Note that this wrapper works only for CV_8UC1 images (for simplicity sake, you can easily expand it to work with other types, if needed), and accepts all THRESH_XXX combinations as original cv::threshold.

double threshold_with_mask(Mat1b& src, Mat1b& dst, double thresh, double maxval, int type, const Mat1b& mask = Mat1b())
{
    if (mask.empty() || (mask.rows == src.rows && mask.cols == src.cols && countNonZero(mask) == src.rows * src.cols))
    {
        // If empty mask, or all-white mask, use cv::threshold
        thresh = cv::threshold(src, dst, thresh, maxval, type);
    }
    else
    {
        // Use mask
        bool use_otsu = (type & THRESH_OTSU) != 0;
        if (use_otsu)
        {
            // If OTSU, get thresh value on mask only
            thresh = otsu_8u_with_mask(src, mask);
            // Remove THRESH_OTSU from type
            type &= THRESH_MASK;
        }

        // Apply cv::threshold on all image
        thresh = cv::threshold(src, dst, thresh, maxval, type);

        // Copy original image on inverted mask
        src.copyTo(dst, ~mask);
    }
    return thresh;
}

Here is the full code for reference:

#include <opencv2/opencv.hpp>
#include <iostream>
using namespace std;
using namespace cv;

// Modified from thresh.cpp
// static double getThreshVal_Otsu_8u(const Mat& _src)

double otsu_8u_with_mask(const Mat1b src, const Mat1b& mask)
{
    const int N = 256;
    int M = 0;
    int i, j, h[N] = { 0 };
    for (i = 0; i < src.rows; i++)
    {
        const uchar* psrc = src.ptr(i);
        const uchar* pmask = mask.ptr(i);
        for (j = 0; j < src.cols; j++)
        {
            if (pmask[j])
            {
                h[psrc[j]]++;
                ++M;
            }
        }
    }

    double mu = 0, scale = 1. / (M);
    for (i = 0; i < N; i++)
        mu += i*(double)h[i];

    mu *= scale;
    double mu1 = 0, q1 = 0;
    double max_sigma = 0, max_val = 0;

    for (i = 0; i < N; i++)
    {
        double p_i, q2, mu2, sigma;

        p_i = h[i] * scale;
        mu1 *= q1;
        q1 += p_i;
        q2 = 1. - q1;

        if (std::min(q1, q2) < FLT_EPSILON || std::max(q1, q2) > 1. - FLT_EPSILON)
            continue;

        mu1 = (mu1 + i*p_i) / q1;
        mu2 = (mu - q1*mu1) / q2;
        sigma = q1*q2*(mu1 - mu2)*(mu1 - mu2);
        if (sigma > max_sigma)
        {
            max_sigma = sigma;
            max_val = i;
        }
    }

    return max_val;
}

double threshold_with_mask(Mat1b& src, Mat1b& dst, double thresh, double maxval, int type, const Mat1b& mask = Mat1b())
{
    if (mask.empty() || (mask.rows == src.rows && mask.cols == src.cols && countNonZero(mask) == src.rows * src.cols))
    {
        // If empty mask, or all-white mask, use cv::threshold
        thresh = cv::threshold(src, dst, thresh, maxval, type);
    }
    else
    {
        // Use mask
        bool use_otsu = (type & THRESH_OTSU) != 0;
        if (use_otsu)
        {
            // If OTSU, get thresh value on mask only
            thresh = otsu_8u_with_mask(src, mask);
            // Remove THRESH_OTSU from type
            type &= THRESH_MASK;
        }

        // Apply cv::threshold on all image
        thresh = cv::threshold(src, dst, thresh, maxval, type);

        // Copy original image on inverted mask
        src.copyTo(dst, ~mask);
    }
    return thresh;
}

int main()
{
    // Load an image
    Mat1b img = imread("D:\\SO\\img\\nice.jpg", IMREAD_GRAYSCALE);

    // Apply OpenCV version
    Mat1b cvth;
    double cvth_value = threshold(img, cvth, 100, 255, THRESH_OTSU);

    // Create a binary mask
    Mat1b mask(img.rows, img.cols, uchar(0));
    rectangle(mask, Rect(100, 100, 200, 200), Scalar(255), CV_FILLED);

    // Apply threshold with a mask
    Mat1b th;
    double th_value = threshold_with_mask(img, th, 100, 255, THRESH_OTSU, mask);

    // Show results
    imshow("cv::threshod", cvth);
    imshow("threshold_with_balue", th);
    waitKey();

    return 0;
}

Mask Image region for Otsu Threshold with OpenCV

After following the answer Miki linked, I realized that one can index with conditions in Python. The explicit loop takes a second, the indexing is miliseconds.

     tempThresImg = img[img !=  255]

Using a mask with an adaptive threshold?

I've written some Python (sorry not c++) code that will allow for masked adaptive thresholding. Its not very fast, but it does what you want, and you may be able to use it as a basis for C++ code. It works as follows:

Sets masked pixels in the image to zero.
Determines the number of unmasked neighbours within the convolution block for each pixel.
Performs a convolution, and averages it by the number of unmasked neighbours within the block. This yields the average value within a pixels neighbourhood block.
Thresholds, by comparing the image to the mean neighbourhood values, mean_conv
Adds the masked off (non-thresholded) part of the image back on.

Sample Image

The images show, the initial image, the mask, the final processed image.

Here's the code:

import cv
import numpy
from scipy import signal

def thresh(a, b, max_value, C):
    return max_value if a > b - C else 0

def mask(a,b):
    return a if b > 100 else 0

def unmask(a,b,c):
    return b if c > 100 else a

v_unmask = numpy.vectorize(unmask)
v_mask = numpy.vectorize(mask)
v_thresh = numpy.vectorize(thresh)

def block_size(size):
    block = numpy.ones((size, size), dtype='d')
    block[(size - 1 ) / 2, (size - 1 ) / 2] = 0
    return block

def get_number_neighbours(mask,block):
    '''returns number of unmasked neighbours of every element within block'''
    mask = mask / 255.0
    return signal.convolve2d(mask, block, mode='same', boundary='symm')

def masked_adaptive_threshold(image,mask,max_value,size,C):
    '''thresholds only using the unmasked elements'''
    block = block_size(size)
    conv = signal.convolve2d(image, block, mode='same', boundary='symm')
    mean_conv = conv / get_number_neighbours(mask,block)
    return v_thresh(image, mean_conv, max_value,C)

image = cv.LoadImageM("image.png", cv.CV_LOAD_IMAGE_GRAYSCALE)
mask = cv.LoadImageM("mask.png", cv.CV_LOAD_IMAGE_GRAYSCALE)

#change the images to numpy arrays
original_image = numpy.asarray(image)
mask = numpy.asarray(mask)
# Masks the image, by removing all masked pixels.
# Elements for mask > 100, will be processed
image = v_mask(original_image, mask)
# convolution parameters, size and C are crucial. See discussion in link below.
image = masked_adaptive_threshold(image,mask,max_value=255,size=7,C=5)
# puts the original masked off region of the image back
image = v_unmask(original_image, image, mask)
#change to suitable type for opencv
image = image.astype(numpy.uint8)
#convert back to cvmat
image = cv.fromarray(image)

cv.ShowImage('image', image)
#cv.SaveImage('final.png',image)
cv.WaitKey(0)

After writing this I found this great link that has a good explanation with plenty of image examples, I used their text image for the above example.

Note. Numpy masks do not seem to be respected by scipy signal.convolve2d(), so the above workarounds were necessary.

How to apply threshold counters to only specified/masking region in the image using opencv

You can apply functions to ROIs using Numpy indexing as explained in this tutoral:

top,left = 338,282
h,w = 216,106

gray_img = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)
crop_img = cv2.rectangle(img, (left,top), (left+w,top+h), (0,255,0), 3)

retval, thresh = cv2.threshold(gray_img, 127, 255, 0)
img_contours, _ = cv2.findContours(thresh[top:top+h,left:left+w], cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)
cv2.drawContours(img[top:top+h,left:left+w], img_contours, -1, (255, 0, 0))

Please note that with Numpy indexing y comes first then x, as opposed to point coordinates with (x,y).

With this image we get: Sample Image

Otsu thresholding inside mask

I had a try at this using the threshold_otsu() method from skimage and a Numpy masked array. I don't know if there are faster ways - the skimage is normally pretty well optimised. If anyone else wants to take my sample data and try other ideas on it, please feel free - although there is a service charge of one upvote ;-)

#!/usr/bin/env python3

import cv2
import numpy as np
import numpy.ma as ma
from skimage.filters import threshold_otsu

# Set up some repeatable test data, 4 blocks 100x100 pixels each of random normal np.uint8s centred on 32, 64, 160,192
np.random.seed(42)
a=np.random.normal(size = (100,100), loc = 32,scale=10).astype(np.uint8)
b=np.random.normal(size = (100,100), loc = 64,scale=10).astype(np.uint8)
c=np.random.normal(size = (100,100), loc = 160,scale=10).astype(np.uint8)
d=np.random.normal(size = (100,100), loc = 192,scale=10).astype(np.uint8)
# Stack (concatenate) the 4 squares horizontally across the page
im = np.hstack((a,b,c,d))
# Next line is just for debug
cv2.imwrite('start.png',im)

That gives us this:

Sample Image

# Now make a mask revealing only left half of image, centred on 32 and 64
mask=np.zeros((100,400))
mask[:,200:]=1
masked = ma.masked_array(im,mask)
print(threshold_otsu(masked.compressed()))       # Prints 47

# Now do same revealing only right half of image, centred on 160 and 192
masked = ma.masked_array(im,1-mask)
print(threshold_otsu(masked.compressed()))       # Prints 175

The histogram of the test data looks like this, x-axis is 0..255

Sample Image

Adapting to your own sample data, I get this:

#!/usr/bin/env python3

import cv2
import numpy as np
import numpy.ma as ma
from skimage.filters import threshold_otsu

# Load images
im   = cv2.imread('eye.tif', cv2.IMREAD_UNCHANGED)
mask = cv2.imread('mask.tif', cv2.IMREAD_UNCHANGED)

# Calculate Otsu threshold on entire image
print(threshold_otsu(im))                       # prints 130

# Now do same for masked image
masked = ma.masked_array(im,mask>0)
print(threshold_otsu(masked.compressed())).     # prints 124

OPENCV mean value of image using only values under a threshold

You can create a mask using cv2.threshold. You need to set all values above 200 to 0 and the ones below it to 255, hence, the flag needs to be set to cv2.THRESH_BINARY_INV.

import cv2

# assuming row has your values
ret, thresh = cv2.threshold(row, 200, 255, cv2.THRESH_BINARY_INV)
result = cv2.mean(row, thresh)

This will only calculate the mean of pixels less equal to 200.

Opencv Threshold with Mask