Smoothing Edges of a Binary Image

Question

How to smooth the edges of this binary image of blood vessels obtained after thresholding.

I tried a method somewhat similar to this method but did not quite get the result I expected.

Here's the code:

import cv2
import numpy as np

INPUT = cv2.imread('so-br-in.png',0)
MASK = np.array(INPUT/255.0, dtype='float32')

MASK = cv2.GaussianBlur(MASK, (5,5), 11)
BG = np.ones([INPUT.shape[0], INPUT.shape[1], 1], dtype='uint8')*255

OUT_F = np.ones([INPUT.shape[0], INPUT.shape[1], 1],dtype='uint8')

for r in range(INPUT.shape[0]):
    for c in range(INPUT.shape[1]):
        OUT_F[r][c]  = int(BG[r][c]*(MASK[r][c]) + INPUT[r][c]*(1-MASK[r][c]))

cv2.imwrite('brain-out.png', OUT_F)  

What can be done to improve the smoothing of these harsh edges?

EDIT

I'd like to smoothen the edges something like http://pscs5.tumblr.com/post/60284570543. How to do this in OpenCV?

Outdoors

See if this is what you want. – Shirleenshirlene

try to use numpy operation on mat it is mutch faster than a pixel by pixel operation – Catechism

Maybe you could describe in more detail what you are expecting ... – Photoemission

@Shirleenshirlene That didn't work on my image. If you've got any other ideas, please share. :) – Outdoors

@Photoemission Updated the question. Please check if you've got any suggestions. :) – Outdoors

# EDIT with a cv2.adaptiveThreshold() it can be a lot better – Catechism

what do you want to use the result for? If you want to use the binary image as a mask, "pixel smoothing" doesn't work because you'll have to interpret the result binary again instead of grayscale values. Did you think about smoothing the contours? e.g. for each contour pixel: replace its position by the mean of the surrounding positions – Taps

@Taps I want to use the result as ground truth for nerve regions. I need it to be a binary image in the end. – Outdoors

It would give more chance to obtain a better result if you also provided the original image, before thresholding. – Litterbug

@Litterbug i.sstatic.net/HWOA9.png – Outdoors

@Litterbug Did you check? – Outdoors

@AbdulFatir I would like to know how you extracted the blood vessels? You can elaborate if you have time – Eu

Answer 1

Here is the result I obtained with your image:

My method is mostly based on several cv::medianBlurapplied on a scaled-up image.

Here is the code:

cv::Mat vesselImage = cv::imread(filename); //the original image
cv::threshold(vesselImage, vesselImage, 125, 255, THRESH_BINARY);
cv::Mat blurredImage; //output of the algorithm
cv::pyrUp(vesselImage, blurredImage);

for (int i = 0; i < 15; i++)
    cv::medianBlur(blurredImage, blurredImage, 7);

cv::pyrDown(blurredImage, blurredImage);
cv::threshold(blurredImage, blurredImage, 200, 255, THRESH_BINARY);

The jagged edges are due to the thresholding. If you are comfortable with an output image that is non-binary (i.e. with 256 shades of grAy), you can just remove it and you get this image:

Bushey

That's a better result. I've upvoted your answer, if I don't get better answers in the next few days, I'll accept this answer. Any ideas about pscs5.tumblr.com/post/60284570543 – Outdoors

In this post, the only real operation applied is a GaussianBlur. Then he adjusts the levels to get sharper edges. In your case, the defects in the image are too pronounced to be corrected with a simple Gaussian blur. You could try to apply a Gaussian blur to my result image but I suppose it would look like the second image I posted. – Bushey

Answer 2

You can dilate then erode the areas http://docs.opencv.org/2.4/doc/tutorials/imgproc/erosion_dilatation/erosion_dilatation.html.

import cv2
import numpy as np
blur=((3,3),1)
erode_=(5,5)
dilate_=(3, 3)
cv2.imwrite('imgBool_erode_dilated_blured.png',cv2.dilate(cv2.erode(cv2.GaussianBlur(cv2.imread('so-br-in.png',0)/255, blur[0], blur[1]), np.ones(erode_)), np.ones(dilate_))*255)  

EDIT whith a scale facor off 4 before the stuff

Catechism

Hi. Can you post some example code and results if possible? – Outdoors

there is what i mean ! – Catechism

due to your different kernel sizes for erosion and dilation you bloated his contours pretty much. :) Not sure if this is a desirable result. – Reservoir

The source image is pretty bad so ... that what I felt the best but I kept parameters clearly out to let him tweak them – Catechism

This is certainly better, but even in this image the edges are not very smooth. Could you tell how to smoothen those edges? – Outdoors

You'll have to increase image resolution using resizeprior to those operations – Photoemission

i am curently working on the new source image, we can't do mutch better with the prior one. – Catechism

Thanks for the code, but boy, that's some ugly code in the last line! :( – Murrhine

Answer 3

i made some modifications on @dhanushka 's answer for another question and get these images.

Sorry it is C++ code but maybe you will convert it to Python.

You can change the parameters below to get different results.

// contour smoothing parameters for gaussian filter
int filterRadius = 10; // you can try to change this value
int filterSize = 2 * filterRadius + 1;
double sigma = 20; // you can try to change this value

#include <opencv2/highgui.hpp>
#include <opencv2/imgproc.hpp>
#include <iostream>

using namespace cv;
using namespace std;

int main( int argc, const char** argv )
{
    Mat im = imread(argv[1], 0);

    Mat cont = ~im;
    Mat original = Mat::zeros(im.rows, im.cols, CV_8UC3);
    Mat smoothed = Mat(im.rows, im.cols, CV_8UC3, Scalar(255,255,255));

    // contour smoothing parameters for gaussian filter
    int filterRadius = 5;
    int filterSize = 2 * filterRadius + 1;
    double sigma = 10;

    vector<vector<Point> > contours;
    vector<Vec4i> hierarchy;
    // find contours and store all contour points
    findContours(cont, contours, hierarchy, CV_RETR_CCOMP, CV_CHAIN_APPROX_NONE, Point(0, 0));
    for(size_t j = 0; j < contours.size(); j++)
    {
        // extract x and y coordinates of points. we'll consider these as 1-D signals
        // add circular padding to 1-D signals
        size_t len = contours[j].size() + 2 * filterRadius;
        size_t idx = (contours[j].size() - filterRadius);
        vector<float> x, y;
        for (size_t i = 0; i < len; i++)
        {
            x.push_back(contours[j][(idx + i) % contours[j].size()].x);
            y.push_back(contours[j][(idx + i) % contours[j].size()].y);
        }
        // filter 1-D signals
        vector<float> xFilt, yFilt;
        GaussianBlur(x, xFilt, Size(filterSize, filterSize), sigma, sigma);
        GaussianBlur(y, yFilt, Size(filterSize, filterSize), sigma, sigma);
        // build smoothed contour
        vector<vector<Point> > smoothContours;
        vector<Point> smooth;
        for (size_t i = filterRadius; i < contours[j].size() + filterRadius; i++)
        {
            smooth.push_back(Point(xFilt[i], yFilt[i]));
        }
        smoothContours.push_back(smooth);

        Scalar color;

        if(hierarchy[j][3] < 0 )
        {
            color = Scalar(0,0,0);
        }
        else
        {
            color = Scalar(255,255,255);
        }
        drawContours(smoothed, smoothContours, 0, color, -1);
    }
    imshow( "result", smoothed );
    waitKey(0);
}

Cemetery

This is correct answer but needs to be converted in Python. – Salaried

Is this program written in python? – Bristow

Be careful - this degree of smoothing has hidden radiologically important vessel irregularity. It is 'pretty' but is no longer a clinically interpretable image. (I am a neuroradiologist) – Flotow

Answer 4

What you can do is increase the resolution of your image (e.g. double or triple it using resize). After that, erosion and dilation as described in the other answer above will lead to finer results.

Photoemission

are we brain connected ? i just finished that whith a scale factor of 4 and it's better, but whit the new source image in our hand we can do better :) – Catechism

Answer 5

You most probably got a gray scale image of the blood vessels first and then thresholded. It still looks non-smooth because the original gray scale image had noise inside. Asking now for a smoothing of the edges will result in lower resolution. For example the dilution and erosion proposed in another answer might fuse neighboring vessels in the dilution step which then cannot be separated again in the erosion step.

It might be preferable to remove the noise in the gray scale image first (aka do a smoothing there) and do the thresholding as the last step.

Because you did not deliver the gray scale image I performed a mild smoothing (about one pixel width) here on the binary image and performed a thresholding again.

I did a smoothing (with a Gaussian kernel of fixed size) and a thresholding (with a thresholding parameter). I suggest you do that on the grayscale image data and adjust the two parameters until you like the result.

Matlab code in case it is of interest:

% read
img = imread('YyNQV.png');
img = double(img(:, :, 1) ~= 255); % png is RGB -> binary

% smooth
kernel = fspecial('gaussian', 10, 1.5);
kernel = kernel / sum(kernel(:)); % normalize to 1
img_smooth = conv2(img, kernel, 'same');

% binarize again
threshold = 0.4; % experiment with values between 0 and 1
img_smooth_threshold = img_smooth > threshold;

% save (exchange black and white)
imwrite(~img_smooth_threshold, 'YyNQV_smooth.png');

Longshore

Answer 6

This is algorithm from sturkmen's post above converted to Python

import numpy as np
import cv2 as cv

def smooth_raster_lines(im, filterRadius, filterSize, sigma):
    smoothed = np.zeros_like(im)
    contours, hierarchy = cv.findContours(im, cv.RETR_CCOMP, cv.CHAIN_APPROX_NONE)
    hierarchy = hierarchy[0]
    for countur_idx, contour in enumerate(contours):
        len_ = len(contour) + 2 * filterRadius
        idx = len(contour) - filterRadius

        x = []
        y = []    
        for i in range(len_):
            x.append(contour[(idx + i) % len(contour)][0][0])
            y.append(contour[(idx + i) % len(contour)][0][1])

        x = np.asarray(x, dtype=np.float32)
        y = np.asarray(y, dtype=np.float32)

        xFilt = cv.GaussianBlur(x, (filterSize, filterSize), sigma, sigma)
        xFilt = [q[0] for q in xFilt]
        yFilt = cv.GaussianBlur(y, (filterSize, filterSize), sigma, sigma)
        yFilt = [q[0] for q in yFilt]


        smoothContours = []
        smooth = []
        for i in range(filterRadius, len(contour) + filterRadius):
            smooth.append([xFilt[i], yFilt[i]])

        smoothContours = np.asarray([smooth], dtype=np.int32)


        color = (0,0,0) if hierarchy[countur_idx][3] > 0 else (255,255,255)
        cv.drawContours(smoothed, smoothContours, 0, color, -1)
    
    return(smoothed)

Bertberta

Answer 7

This post is a little old now, however I'd like to raise a word of warning for others reading this question.

I am a neuroradiologist and interpret these angiograms in my day job. I am not sure what the use case is for smoothing these images as presented by the OP, but you would need to be careful when smoothing these images- some of the irregularity is actually due to pathology and you may (depending on use case) not want to remove this.

Depending on use case, it would be best to work with your clinical radiology side to agree on what is an appropriate degree of smoothing.

Flotow