Fun with sobel - Neon Effect - Thick Edge Detection

[fmw42]

see radon transform at viewtopic.php?f=2&t=14532&p=50617#p50617

[HugoRune]

I have been thinking about doing a hough transform with the tools already provided by imagemagick, and it seems relatively straightforward.

All that is missing is a way to simulate the accumulator array.

This is just a random idea, but is there some way to do a 2d histogram in imagemagick?
I have seen your method of generating a 1d histogram for a single channel image:
convert zelda3.jpg -colorspace gray histogram:- | convert - -scale 256x1! zelda3g_hist.png
I there a way to extend this to a 2d histogram for a 2 channels?

For example, given the red/green channel of an image, can I create a 2d histogram image with

• - 0-255 red values on the x axis
  - 0-255 green values on the y axis
  - the brightness of a pixel at x0,y0 determined by how often the color RGB(x0,y0,0) appears in the image

Such a 2d histogram could be used as an accumulator. With that and the fx operator, I think i could do a hough transform in one imagemagick command.
(of course, this migth fall under the category "just because you can does not mean you should")

EDIT: looking over this, I am actually not so sure anymore if it would work, I have to think about it again)

[fmw42]

If you decide to push forward, with regard to 2D plots, see my scripts, scatter and scatterchannels. They are scattergrams. Don't know if this helps or not.

[HugoRune]

Well, the idea was this:

- Generate an image with the theta values of all edge pixels (the orientation of the gradient)
- Generate an image with the rho values of all edge pixels (The shortest distance to the origin of a line throught that pixel with orientation theta)
- create a 2d histogram: count how often each possible (rho,theta) value pair occcurs, and plot this with rho on one axis, theta on the other axis, brightness=number of occurences.

This should be equivalent to the oriented hough transform described in the paper by Zhang and He. The classic hough transform would require a lot more effort to compute this way (repeat the whole process N times, for all possible theta values, then add the resulting 2d plots).

In theory it should work, but using the process shown below I did not get any good results, and the whole thing seems extremely impractical and inefficient to me anyway.

Code: Select all

convert rectrot2.png  -resize 50x50! -colorspace gray ^
    ( -clone 0 -bias 50% -convolve "-0.125,0,0.125,  -0.25,0,0.25,  -0.125,0,0.125" ) ^
    ( -clone 0 -bias 50% -convolve "-0.125,-0.25,-0.125,  0,0,0,  0.125,0.25,0.125" ) ^
    -delete 0 ^
    -fx "uu=(u-0.5)*2;vv=(v-0.5)*2; hypot(uu,vv)>0.05 ? atan2(vv,uu) /pi/2+0.5 : 0"  -write THETA.png ^
    -fx "uu=(u-0.5)*pi*2; u>0 ? ((i/w)*cos(uu) + (j/h)*sin(uu)) /sqrt(2)/2+0.5 : 0" RHO.png

scatter THETA.png RHO.png HOUGH.png

input ------------------------------ THETA --------------------------- RHO




There are 4 peaks which might respond to the right lines, but it is hard to tell with all that noise. The binary nature of the scatter script makes this difficult, since it is impossible to tell how many counts each pixel represents, but the high calulation time limits the useability of this approach anyway.

[fmw42]

This should be equivalent to the oriented hough transform described in the paper by Zhang and He. The classic hough transform would require a lot more effort to compute this way (repeat the whole process N times, for all possible theta values, then add the resulting 2d plots).

I agree and that was my first thought also

Does it help to quantize your theta values to the 8 primary compass directions? Also to threshold the magnitude to binary form partially to eliminate low priority edges (probably here it is not an issue as your image is not 'noisy'), but to make the edge strengths equal. This could help speed up processing the 2D histogram if you throw out zero magnitude edges (black pixels) or conversely just use edges with mag=1 (white pixels)

The binary nature of the scatter script makes this difficult

Yes, that is why I was not sure it would help other than give you an idea about how to produce a 2D histogram. But it needs extending so that it is grayscale and not binary.

So that was why I had hoped that the radon transform within -deskew might be immediately useable without too much effort. But it looks like it needs significant alteration and I have not yet replied to Magick with any pseudocode, but I have found several pieces of code and pseudocode to send him. However, he is always very busy and I don't expect it needs escalation to the top of his queue.

[HugoRune]

Some more experiments. I went back to python since I do not see a way to make the previous attempt work at acceptable speed and quality.

My conclusion is that the oriented hough transform is indeed much better suited than the classical for this particular problem.

Still not sure if this could or should be integrated with imagemagick. Seems like there are too many possible modifications for command line switches

original image:


oriented hough transform

This is simply using the orientation estimated from sobel. Very good results, much better peaks than the classical hough transform, less noise. Much faster to calculate too.
Only a few seconds for a medium-sized image, even in python.

classical hough transform

Very long calculation time with my unoptimized python implementation.
As expected, the borders are lost amid the many edges from the text.
That is why I needed the solbelseries script from the first post in this thread to filter all small structures.

A modification I tried:

A sort of combination of the two previous algorithms
Like the classical hough transfrom, I loop trough all possible theta values from -180° to +180°
But here I only add +1 to the accumulator bin when theta is within +-10° of the measured orientation. All other bins get -1

Seems very effective for certain noisy images. The peaks are also a bit sharper than with the oriented hough transform. As slow as the classical hough transform though.

Code: Select all

#!/usr/bin/python
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 3 of the License, or
# (at your option) any later version.

# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.

# You should have received a copy of the GNU General Public License
# along with this program.  If not, see <http://www.gnu.org/licenses/>.
#
################################################################
#
# This is a standalone program.  
# it requires 
# - python 2.6+ ( 2.5 or 2.4 works probably too)
# - openCV with PythonInterface installed
# 
################################################################
#
# USAGE: hough2.py SOURCEFILE [OUTFILE]
#

import sys
from math import sin,cos,sqrt,degrees,radians,atan2,hypot,pi
from opencv.cv import *
from opencv.highgui import *


# ORIENTED HOUGH TRANSFORM
def orientedHoughTransform(edges,sobelX,sobelY,accumulator):
	maxtheta = pi
	maxrho = hypot(edges.width, edges.height)
	maxacc = -1000000
	for x in range(0,edges.width):
		for y in range(0,edges.height):
			#print x,y, edges.width,edges.height
			if edges[y,x]>0:
				theta = atan2(sobelY[y,x],sobelX[y,x])
				rho = x*cos(theta) + y*sin(theta)
				t=int(round  ((theta/maxtheta/2  +0.5) * (accumulator.width-1) ))
				r=int(round  ((rho  /maxrho  /2  +0.5) * (accumulator.height-1)))
				acc=accumulator[r,t]
				acc+=1
				accumulator[r,t]=acc
				if acc>maxacc: maxacc=acc
				#print theta,sobelY[y,x],sobelX[y,x]
	cvConvertScale(accumulator,accumulator,100000/maxacc)

# ORIENTED HOUGH TRANSFORM WITH MODIFICATIONS
def orientedHoughTransform2(edges,sobelX,sobelY,accumulator):
	maxtheta = pi
	maxrho = hypot(edges.width, edges.height)
	maxacc = -1000000
	for x in range(0,edges.width):
		for y in range(0,edges.height):
			#print x,y, edges.width,edges.height
			if edges[y,x]>0:
				theta2 = atan2(sobelY[y,x],sobelX[y,x])
				for ii in range(-180,179):
					theta= radians(ii)
					rho = x*cos(theta) + y*sin(theta)
					t=int(  ((theta/maxtheta/2  +0.5) * (accumulator.width-1) ))
					r=int(  ((rho  /maxrho  /2  +0.5) * (accumulator.height-1)))

					acc=accumulator[r,t]
					if abs(theta-theta2)<radians(10): acc+=1
					else: acc-=1
					accumulator[r,t]=acc
					if acc>maxacc: maxacc=acc
				#print theta,sobelY[y,x],sobelX[y,x]
	print maxacc
	cvConvertScale(accumulator,accumulator,32767/maxacc)

				
if __name__ == "__main__":
	if len(sys.argv)>1:
		srcfile = sys.argv[1]
		
	srcImg=cvLoadImage(srcfile, 0);
	if not srcImg:
		print "Error opening image %s" % srcfile
		sys.exit(-1)
		
	edgeImg = cvCreateImage( cvGetSize(srcImg), 8, 1 );
	cvCanny( srcImg, edgeImg, 200, 100, 3 );

	sobelXImg = cvCreateImage( cvGetSize(srcImg), IPL_DEPTH_16S, 1 );
	cvSobel( srcImg, sobelXImg, 1, 0, 3 );

	sobelYImg = cvCreateImage( cvGetSize(srcImg), IPL_DEPTH_16S, 1 );
	cvSobel( srcImg, sobelYImg, 0, 1, 3 );

	accumulatorImg = cvCreateImage( cvSize(360,360), IPL_DEPTH_16S, 1 );
	cvZero( accumulatorImg );
	
	orientedHoughTransform(edgeImg, sobelXImg, sobelYImg, accumulatorImg)
	
	if len(sys.argv)>2:
		outfile = sys.argv[2]
		cvSaveImage(outfile,accumulatorImg)
	else:
		cvNamedWindow( "hough", 1 );
		cvShowImage( "hough", accumulatorImg );
		cvWaitKey(0);
		cvShowImage( "hough", edgeImg );
		cvWaitKey(0);
		# cvShowImage( "hough", sobelYImg );
		# cvWaitKey(0);

[fmw42]

very nice work.

why 5 (not 4) bright spots? Is that due to the open pages on the left side of the book?

I wonder if something like this could be done with AWK (Anthony?)


Just for curiosity and comparison, I tried the ImageJ radon transform on the original and an edge image (both cropped to square dimensions required by the radon plugin)

image:


radon:


simple edge:
convert pic4nb.png -monochrome -median 2 -edge 1 pic4nb_mono_med2_edge1.png


radon from edge image:

[HugoRune]

why 5 (not 4) bright spots? Is that due to the open pages on the left side of the book?

Not sure actually, it might be the middle of the book, or one of the other lines. I still have to finish the code that finds the local maxima in hough space and converts them back into lines.

I did not do any noise filtering on this image, just a canny filter with a low threshold to detect all edges. I wanted to see how it can deal with noisy images.
(the original image was resized before uploading)

I wonder if something like this could be done with AWK (Anthony?)

I never worked with AWK, but it seems like a bad match. Why AWK?

Just for curiosity and comparison, I tried the ImageJ radon transform on the original and an edge image (both cropped to square dimensions required by the radon plugin)

I never figured out what the advantages and disadvantages of the radon or the hough transform are for this sort of problem. Are there fundamental differences of the results (for the classical cases)?
Hough seems to be faster to calculate, if I only loop over the edges. and easier to thinker with. But that might be because the integrals confuse me.

simple edge:
convert pic4nb.png -monochrome -median 2 -edge 1 pic4nb_mono_med2_edge1.png

interestingly, the oriented hough transform does not work at all if it only has this sort of edge image as input.
To calculate precise gradient orientation, The sobel results of the original images are needed
If I run sobel over this edge image, the gradient varies to much to give any meaningful result.

[fmw42]

I never worked with AWK, but it seems like a bad match. Why AWK?

Just learning a bit of it and finding that for what it can do, it can do it so much faster than bash scripting with while or for loops.

I never figured out what the advantages and disadvantages of the radon or the hough transform are for this sort of problem. Are there fundamental differences of the results (for the classical cases)?
Hough seems to be faster to calculate, if I only loop over the edges. and easier to thinker with. But that might be because the integrals confuse me.

Not totally sure either. Seemed to be a bit less noise sensitive from what i have read. Also I knew that Magick had implemented some version of it and thought it might be something that could be extracted for direct use in IM easily (but not the case).

I made the examles, 1) knowing that you had not done much to reduce noise, 2) knowing that my edge result is different and not directly comparable to what you are doing, 3) it was available from ImageJ, and 4) just for curiosity to see the differences. Just learning more myself about both.

interestingly, the oriented hough transform does not work at all if it only has this sort of edge image as input.
To calculate precise gradient orientation, The sobel results of the original images are needed
If I run sobel over this edge image, the gradient varies to much to give any meaningful result.

I would not necessarily expect it to work as you are doing and edge extraction on an edge image. Perhaps better to try from the binary image I made before extracting the edge. Here it is:

[anthony]

Awk is suitable for handling TEXT data. Which is what we wer dealing with for histogram modification. Python and/or perl is still much better and faster than awk.

[fmw42]

What are the pros and cons of python vs perl? And where does Ruby fit in all this?

Just checked and my Mac seems to have come with at least one version of each. I read a book on Perl once, but never got around to trying it. Don't know much about python other than have heard about it.

Found these:

http://infohost.nmt.edu/tcc/help/lang/p ... sperl.html
http://xodian.net/serendipity/index.php ... -Ruby.html
http://mail.python.org/pipermail/python ... 95939.html
http://www.apple-group.com/forums/showt ... hp?p=37923
http://www.dmh2000.com/cjpr/

[HugoRune]

I would not necessarily expect it to work as you are doing and edge extraction on an edge image. Perhaps better to try from the binary image I made before extracting the edge. Here it is:

Not much better actually. As far as I can tell the problem is mainly the binary nature of the image. With a 3x3 sobel over binary data, there are only 8 possible values, so the hough space is strongly quantized in the theta direction.

What are the pros and cons of python vs perl? And where does Ruby fit in all this?

I am a big fan of python, but it is partialy a matter of taste. I think all of them are better suited than bash for complex algorithms.

I'd say the main advantage of perl are longer history so more available code, much of it rather old though. Also, several syntax shortcuts for handling things like regexp, text data and shell commands.

Main advantage of python would be much cleaner and more consistent syntax, which makes reading and writing code much easier. (My own python code above notwithstanding )

Ruby is probably similar to python, with some philosophical differences like how code blocks are marked. Did not really try it out thought. Less widley used than python, so lacking somewhat in available libraries and example code (unless it is for web applications with ruby on rails)

[fmw42]

Not much better actually. As far as I can tell the problem is mainly the binary nature of the image. With a 3x3 sobel over binary data, there are only 8 possible values, so the hough space is strongly quantized in the theta direction.

Had not considered that. I don't understand why it is quantized to 8 directions only if you are using arctan2(DY,DX) for getting theta. Should you not get the actual angle of the lines? I need to think on this more. Nevertheless, why are you getting doubling of the dots? Is that because each opposite side is at a different distance from the origin.

So at least this shows the difference between the radon and the hough transform. In the radon transform, one got 4 dots. You get 8 in your Hough.

Thanks for your experimentation and for your recommendation about Python. Guess I need to buy and read another book.

[HugoRune]

Had not considered that. I don't understand why it is quantized to 8 directions only if you are using arctan2(DY,DX) for getting theta. Should you not get the actual angle of the lines? I need to think on this more.

with binary data the 3x3 sobelX and the sobelY can have any integer value from -4 to +4, so atan2(sobelY,sobelX) has x<8^2 possible values? not sure about that.

Nevertheless, why are you getting doubling of the dots? Is that because each opposite side is at a different distance from the origin.

So at least this shows the difference between the radon and the hough transform. In the radon transform, one got 4 dots. You get 8 in your Hough.

Depends on what you mean.
The doubling of peaks in the last picture is only because the actual peak is inbetween these values, I think. However, since I am doing an oriented hough transform, it is not exactly comparable to the radon transform. In my case, a vertical border from white to black is distinct from a vertical border from black to white, while the classical hough and the radon transform of an edge image would treat both the same.
here is pic4nb_mono_med2.png with the combination of classical and oriented hough approach: all bins with theta whithin 20° of the estimated orientation value get +1, all others -1

In a classical hough or radon transform, the black peaks would be white too, so there is a doubling of peaks.
I think to get from my pictures to the output of a classical hough transform, you would need to cut them in half vertically, then flip the right half horizontally and overlay it over the left half.

[anthony]

Thanks for your experimentation and for your recommendation about Python. Guess I need to buy and read another book.

I know perl myself, I know it well enough that I have written CPAN modules for it. Specifically to do vector mathematics ("VectorReal" module) including the use of 'x' and '.' for cross and dot products, and tie it into another module "MatrixReal". As well as interupted multi-tasking IO handlers!

I tried learning python but the book was not very clear on variable indirection, and its complete dependance on external libraries through me. Also in my reading of python code I found it chaotic and hard to follow.

I have not tried Ruby yet.

In terms of ImageMagick. PerlMagick support is a standard part of the distribution, but always seemed to be a little forced in its image handling, especially with error checking.

I don't know about IM in Python

But their is a very large support group for RubyMagick with Ruby on Rails.

PHP with Imagemagick seems to have devolved to making Command Line or shell API calls for image processing. Basically as there has been at least three attempts to incorporate IM into a PHP module, but never with any real follow up support. As such they never seem to keep up the the new additions such as Image distortions, within ImageMagick.


My own preference is to use command line IM, but with appropriate scripting via whatever program or technique works best for whatever I am doing. I usally start with shell, then is complexity warrants shift to awk or calculations, perl if database store, matrix, file or IO handling, is needed need. Only once have I gone further to C when speed and direct image access became necessary. Usually though shell does the task needed.

Shells are a tool, and we are its carpenters.
Sometimes we need a hammer and sometimes a screwdriver. --- John Navarra