Example of extracting grid nodes in MATLAB

For last few posts, I am writing about morphological image analysis with python with open source package scikit-image. In this post I want to show benefit of morphological image analysis with some real application. I been working with road extraction from satellite images for past few weeks and  one of the common task is to extract the grid corner nodes of the road network. So let’s assume the road grid has been already extracted and we want to extract intersection of grid network. Here I am showing three approaches for that based on synthetic data. Obviously with real data, the task is more complex and might require many slight modifications.

Here is the grid. You can generate such grid with checkerboard function in MATLAB with minor processing (edge detection and hole filling). Intersections of of grid lines have a property that it has 4 white pixels in 4-N neighborhood.  We can exploit this property to extract it. I am going to show to ways by which it can be done and in addition, another method based on extraction of corners will also be shown.

%make a checkerboard
I=checkerboard(20, 10,10)>0.5;
imshow(I)
% detect edge with canny
BW = edge(I,'canny');
figure, imshow(BW)
%make square structural element to fill holes with closing
SE = strel('square', 3);
BW1 = imclose(BW, SE);
figure, imshow(BW1)
% make a skeleton of the to lines 1 pixel thick
BW2 = bwmorph(BW1,'skel',Inf);
figure, imshow(BW2)

 


In mathematics and, in particular, functional analysis, convolution is a mathematical operation on two functions f and g, producing a third function that is typically viewed as a modified version of one of the original function .In simple word, convolution of image is a process where you scan your image from left to right and top to bottom within local neighborhood defined by user, and do mathematical calculation within the local neighborhood. Here we are going to use it for simply sum the number of pixels which are which are white within 3x3 neighborhood. As shown by above figure, if the point is an intersection grid point then at that point sum should be 5. Any points with less than 5 white pixels are not intersection grid points.

% with block processing and inline function with convolution
tic;
kernel = [0 1 0; ... %# Convolution kernel
1 1 1; ...
0 1 0];
sumX = conv2(double(BW2),kernel,'same');
result=sumX;
% only consider pixels which are in grid image
result (BW2==0)= 0
result(sumX<5)=0;
[r,c] = find(result>0);
toc;
figure,imshow((BW2)), hold on, plot(r,c,'r*');title ('with Convolution')
hold off

Way Two: With Mathematical Morphology (MM)

If you are working in the field of Geo-information and you still saying what the heck this ‘MM’ then seriously you should be acquainted with this field of image processing. You can use MM for classification, edge detection, filtering, segmentation, building detection and many other things. Grab a cup of coffee and Google MM in remote sensing; there is tons of stuff to read.  There are courses in MM taught in different universities which are one semester long, so you got an idea how vast the subject is. If you want to seriously know MM, then there is no better book than ‘ Morphological Image Analysis’ by P Soillie. Here we use ‘Erosion’ method of MM.  Time: 0.011184 seconds.

% method two with morphological analysis
tic;
SE = strel ('disk', 1);
result2 = imerode(BW, SE);
[r,c] = find(result>0);
toc;
figure,imshow((BW2)), hold on, plot(r,c,'b*'), title ('With MM Erosion')
hold off

Way Two:  Mathematical Morphology (MM)

There are many corner detection techniques such as Moravec, Harris, SIFT etc. Here, I am going to use Harris corner detector for intersection grid points. Time: 0.094848 seconds.

% detect corner by harris corner detector
tic;
C = corner(BW2,500 );
timet2 = toc;
toc;
figure,imshow(BW2);
hold on
plot(C(:,1), C(:,2), 'g*'); title('with harris corner detector')
hold off

Way Two:  With Harris Corner Detector

So you have seen there are number of ways to solve a particular problem. All above methods detected intersection points successfully. Among three methods, convolution required longest time for processing and MM required shortest time. The gain in computation cost between MM and Corner detection is 9 folds.  I did all this in MATLAB but it can easily coded in python with scikit-image.