import java.awt.image.BufferedImage;
2import java.util.Arrays;
3import java.awt.*;
4import java.awt.image.*;
5import javax.swing.*;
6
7/**
8 * <p><em>This software has been released into the public domain.
9 * <strong>Please read the notes in this source file for additional information.
10 * </strong></em></p>
11 *
12 * <p>This class provides a configurable implementation of the Canny edge
13 * detection algorithm. This classic algorithm has a number of shortcomings,
14 * but remains an effective tool in many scenarios. <em>This class is designed
15 * for single threaded use only.</em></p>
16 *
17 * <p>Sample usage:</p>
18 *
19 * <pre><code>
20 * //create the detector
21 * CannyEdgeDetector detector = new CannyEdgeDetector();
22 * //adjust its parameters as desired
23 * detector.setLowThreshold(0.5f);
24 * detector.setHighThreshold(1f);
25 * //apply it to an image
26 * detector.setSourceImage(frame);
27 * detector.process();
28 * BufferedImage edges = detector.getEdgesImage();
29 * </code></pre>
30 *
31 * <p>For a more complete understanding of this edge detector's parameters
32 * consult an explanation of the algorithm.</p>
33 *
34 * @author Tom Gibara
35 *
36 */
37
38public class CannyEdgeDetector {
39
40 // statics
41
42 private final static float GAUSSIAN_CUT_OFF = 0.005f;
43 private final static float MAGNITUDE_SCALE = 100F;
44 private final static float MAGNITUDE_LIMIT = 1000F;
45 private final static int MAGNITUDE_MAX = (int) (MAGNITUDE_SCALE * MAGNITUDE_LIMIT);
46
47 // fields
48
49 private int height;
50 private int width;
51 private int picsize;
52 private int[] data;
53 private int[] magnitude;
54 private BufferedImage sourceImage;
55 private BufferedImage edgesImage;
56
57 private float gaussianKernelRadius;
58 private float lowThreshold;
59 private float highThreshold;
60 private int gaussianKernelWidth;
61 private boolean contrastNormalized;
62
63 private float[] xConv;
64 private float[] yConv;
65 private float[] xGradient;
66 private float[] yGradient;
67
68 // constructors
69
70 /**
71 * Constructs a new detector with default parameters.
72 */
73
74 public CannyEdgeDetector() {
75 lowThreshold = 2.5f;
76 highThreshold = 7.5f;
77 gaussianKernelRadius = 2f;
78 gaussianKernelWidth = 16;
79 contrastNormalized = false;
80 }
81
82 // accessors
83
84 /**
85 * The image that provides the luminance data used by this detector to
86 * generate edges.
87 *
88 * @return the source image, or null
89 */
90
91 public BufferedImage getSourceImage() {
92 return sourceImage;
93 }
94
95 /**
96 * Specifies the image that will provide the luminance data in which edges
97 * will be detected. A source image must be set before the process method
98 * is called.
99 *
100 * @param image a source of luminance data
101 */
102
103 public void setSourceImage(BufferedImage image) {
104 sourceImage = image;
105 }
106
107 /**
108 * Obtains an image containing the edges detected during the last call to
109 * the process method. The buffered image is an opaque image of type
110 * BufferedImage.TYPE_INT_ARGB in which edge pixels are white and all other
111 * pixels are black.
112 *
113 * @return an image containing the detected edges, or null if the process
114 * method has not yet been called.
115 */
116
117 public BufferedImage getEdgesImage() {
118 return edgesImage;
119 }
120
121 /**
122 * Sets the edges image. Calling this method will not change the operation
123 * of the edge detector in any way. It is intended to provide a means by
124 * which the memory referenced by the detector object may be reduced.
125 *
126 * @param edgesImage expected (though not required) to be null
127 */
128
129 public void setEdgesImage(BufferedImage edgesImage) {
130 this.edgesImage = edgesImage;
131 }
132
133 /**
134 * The low threshold for hysteresis. The default value is 2.5.
135 *
136 * @return the low hysteresis threshold
137 */
138
139 public float getLowThreshold() {
140 return lowThreshold;
141 }
142
143 /**
144 * Sets the low threshold for hysteresis. Suitable values for this parameter
145 * must be determined experimentally for each application. It is nonsensical
146 * (though not prohibited) for this value to exceed the high threshold value.
147 *
148 * @param threshold a low hysteresis threshold
149 */
150
151 public void setLowThreshold(float threshold) {
152 if (threshold < 0) throw new IllegalArgumentException();
153 lowThreshold = threshold;
154 }
155
156 /**
157 * The high threshold for hysteresis. The default value is 7.5.
158 *
159 * @return the high hysteresis threshold
160 */
161
162 public float getHighThreshold() {
163 return highThreshold;
164 }
165
166 /**
167 * Sets the high threshold for hysteresis. Suitable values for this
168 * parameter must be determined experimentally for each application. It is
169 * nonsensical (though not prohibited) for this value to be less than the
170 * low threshold value.
171 *
172 * @param threshold a high hysteresis threshold
173 */
174
175 public void setHighThreshold(float threshold) {
176 if (threshold < 0) throw new IllegalArgumentException();
177 highThreshold = threshold;
178 }
179
180 /**
181 * The number of pixels across which the Gaussian kernel is applied.
182 * The default value is 16.
183 *
184 * @return the radius of the convolution operation in pixels
185 */
186
187 public int getGaussianKernelWidth() {
188 return gaussianKernelWidth;
189 }
190
191 /**
192 * The number of pixels across which the Gaussian kernel is applied.
193 * This implementation will reduce the radius if the contribution of pixel
194 * values is deemed negligable, so this is actually a maximum radius.
195 *
196 * @param gaussianKernelWidth a radius for the convolution operation in
197 * pixels, at least 2.
198 */
199
200 public void setGaussianKernelWidth(int gaussianKernelWidth) {
201 if (gaussianKernelWidth < 2) throw new IllegalArgumentException();
202 this.gaussianKernelWidth = gaussianKernelWidth;
203 }
204
205 /**
206 * The radius of the Gaussian convolution kernel used to smooth the source
207 * image prior to gradient calculation. The default value is 16.
208 *
209 * @return the Gaussian kernel radius in pixels
210 */
211
212 public float getGaussianKernelRadius() {
213 return gaussianKernelRadius;
214 }
215
216 /**
217 * Sets the radius of the Gaussian convolution kernel used to smooth the
218 * source image prior to gradient calculation.
219 *
220 * @return a Gaussian kernel radius in pixels, must exceed 0.1f.
221 */
222
223 public void setGaussianKernelRadius(float gaussianKernelRadius) {
224 if (gaussianKernelRadius < 0.1f) throw new IllegalArgumentException();
225 this.gaussianKernelRadius = gaussianKernelRadius;
226 }
227
228 /**
229 * Whether the luminance data extracted from the source image is normalized
230 * by linearizing its histogram prior to edge extraction. The default value
231 * is false.
232 *
233 * @return whether the contrast is normalized
234 */
235
236 public boolean isContrastNormalized() {
237 return contrastNormalized;
238 }
239
240 /**
241 * Sets whether the contrast is normalized
242 * @param contrastNormalized true if the contrast should be normalized,
243 * false otherwise
244 */
245
246 public void setContrastNormalized(boolean contrastNormalized) {
247 this.contrastNormalized = contrastNormalized;
248 }
249
250 // methods
251
252 public void process() {
253 width = sourceImage.getWidth();
254 height = sourceImage.getHeight();
255 picsize = width * height;
256 initArrays();
257 readLuminance();
258 if (contrastNormalized) normalizeContrast();
259 computeGradients(gaussianKernelRadius, gaussianKernelWidth);
260 int low = Math.round(lowThreshold * MAGNITUDE_SCALE);
261 int high = Math.round( highThreshold * MAGNITUDE_SCALE);
262 performHysteresis(low, high);
263 thresholdEdges();
264 writeEdges(data);
265 }
266
267 // private utility methods
268
269 private void initArrays() {
270 if (data == null || picsize != data.length) {
271 data = new int[picsize];
272 magnitude = new int[picsize];
273
274 xConv = new float[picsize];
275 yConv = new float[picsize];
276 xGradient = new float[picsize];
277 yGradient = new float[picsize];
278 }
279 }
280
281 //NOTE: The elements of the method below (specifically the technique for
282 //non-maximal suppression and the technique for gradient computation)
283 //are derived from an implementation posted in the following forum (with the
284 //clear intent of others using the code):
285 // [url=http://forum.java.sun.com/thread.jspa?threadID=546211&start=45&tstart=0]Java 2D - canny algorithm for edge detection [Locked][/url]
286 //My code effectively mimics the algorithm exhibited above.
287 //Since I don't know the providence of the code that was posted it is a
288 //possibility (though I think a very remote one) that this code violates
289 //someone's intellectual property rights. If this concerns you feel free to
290 //contact me for an alternative, though less efficient, implementation.
291
292 private void computeGradients(float kernelRadius, int kernelWidth) {
293
294 //generate the gaussian convolution masks
295 float kernel[] = new float[kernelWidth];
296 float diffKernel[] = new float[kernelWidth];
297 int kwidth;
298 for (kwidth = 0; kwidth < kernelWidth; kwidth++) {
299 float g1 = gaussian(kwidth, kernelRadius);
300 if (g1 <= GAUSSIAN_CUT_OFF && kwidth >= 2) break;
301 float g2 = gaussian(kwidth - 0.5f, kernelRadius);
302 float g3 = gaussian(kwidth + 0.5f, kernelRadius);
303 kernel[kwidth] = (g1 + g2 + g3) / 3f / (2f * (float) Math.PI * kernelRadius * kernelRadius);
304 diffKernel[kwidth] = g3 - g2;
305 }
306
307 int initX = kwidth - 1;
308 int maxX = width - (kwidth - 1);
309 int initY = width * (kwidth - 1);
310 int maxY = width * (height - (kwidth - 1));
311
312 //perform convolution in x and y directions
313 for (int x = initX; x < maxX; x++) {
314 for (int y = initY; y < maxY; y += width) {
315 int index = x + y;
316 float sumX = data[index] * kernel[0];
317 float sumY = sumX;
318 int xOffset = 1;
319 int yOffset = width;
320 for(; xOffset < kwidth ;) {
321 sumY += kernel[xOffset] * (data[index - yOffset] + data[index + yOffset]);
322 sumX += kernel[xOffset] * (data[index - xOffset] + data[index + xOffset]);
323 yOffset += width;
324 xOffset++;
325 }
326
327 yConv[index] = sumY;
328 xConv[index] = sumX;
329 }
330
331 }
332
333 for (int x = initX; x < maxX; x++) {
334 for (int y = initY; y < maxY; y += width) {
335 float sum = 0f;
336 int index = x + y;
337 for (int i = 1; i < kwidth; i++)
338 sum += diffKernel[i] * (yConv[index - i] - yConv[index + i]);
339
340 xGradient[index] = sum;
341 }
342
343 }
344
345 for (int x = kwidth; x < width - kwidth; x++) {
346 for (int y = initY; y < maxY; y += width) {
347 float sum = 0.0f;
348 int index = x + y;
349 int yOffset = width;
350 for (int i = 1; i < kwidth; i++) {
351 sum += diffKernel[i] * (xConv[index - yOffset] - xConv[index + yOffset]);
352 yOffset += width;
353 }
354
355 yGradient[index] = sum;
356 }
357
358 }
359
360 initX = kwidth;
361 maxX = width - kwidth;
362 initY = width * kwidth;
363 maxY = width * (height - kwidth);
364 for (int x = initX; x < maxX; x++) {
365 for (int y = initY; y < maxY; y += width) {
366 int index = x + y;
367 int indexN = index - width;
368 int indexS = index + width;
369 int indexW = index - 1;
370 int indexE = index + 1;
371 int indexNW = indexN - 1;
372 int indexNE = indexN + 1;
373 int indexSW = indexS - 1;
374 int indexSE = indexS + 1;
375
376 float xGrad = xGradient[index];
377 float yGrad = yGradient[index];
378 float gradMag = hypot(xGrad, yGrad);
379
380 //perform non-maximal supression
381 float nMag = hypot(xGradient[indexN], yGradient[indexN]);
382 float sMag = hypot(xGradient[indexS], yGradient[indexS]);
383 float wMag = hypot(xGradient[indexW], yGradient[indexW]);
384 float eMag = hypot(xGradient[indexE], yGradient[indexE]);
385 float neMag = hypot(xGradient[indexNE], yGradient[indexNE]);
386 float seMag = hypot(xGradient[indexSE], yGradient[indexSE]);
387 float swMag = hypot(xGradient[indexSW], yGradient[indexSW]);
388 float nwMag = hypot(xGradient[indexNW], yGradient[indexNW]);
389 float tmp;
390 /*
391 * An explanation of what's happening here, for those who want
392 * to understand the source: This performs the "non-maximal
393 * supression" phase of the Canny edge detection in which we
394 * need to compare the gradient magnitude to that in the
395 * direction of the gradient; only if the value is a local
396 * maximum do we consider the point as an edge candidate.
397 *
398 * We need to break the comparison into a number of different
399 * cases depending on the gradient direction so that the
400 * appropriate values can be used. To avoid computing the
401 * gradient direction, we use two simple comparisons: first we
402 * check that the partial derivatives have the same sign (1)
403 * and then we check which is larger (2). As a consequence, we
404 * have reduced the problem to one of four identical cases that
405 * each test the central gradient magnitude against the values at
406 * two points with 'identical support'; what this means is that
407 * the geometry required to accurately interpolate the magnitude
408 * of gradient function at those points has an identical
409 * geometry (upto right-angled-rotation/reflection).
410 *
411 * When comparing the central gradient to the two interpolated
412 * values, we avoid performing any divisions by multiplying both
413 * sides of each inequality by the greater of the two partial
414 * derivatives. The common comparand is stored in a temporary
415 * variable (3) and reused in the mirror case (4).
416 *
417 */
418 if (xGrad * yGrad <= (float) 0 /*(1)*/
419 ? Math.abs(xGrad) >= Math.abs(yGrad) /*(2)*/
420 ? (tmp = Math.abs(xGrad * gradMag)) >= Math.abs(yGrad * neMag - (xGrad + yGrad) * eMag) /*(3)*/
421 && tmp > Math.abs(yGrad * swMag - (xGrad + yGrad) * wMag) /*(4)*/
422 : (tmp = Math.abs(yGrad * gradMag)) >= Math.abs(xGrad * neMag - (yGrad + xGrad) * nMag) /*(3)*/
423 && tmp > Math.abs(xGrad * swMag - (yGrad + xGrad) * sMag) /*(4)*/
424 : Math.abs(xGrad) >= Math.abs(yGrad) /*(2)*/
425 ? (tmp = Math.abs(xGrad * gradMag)) >= Math.abs(yGrad * seMag + (xGrad - yGrad) * eMag) /*(3)*/
426 && tmp > Math.abs(yGrad * nwMag + (xGrad - yGrad) * wMag) /*(4)*/
427 : (tmp = Math.abs(yGrad * gradMag)) >= Math.abs(xGrad * seMag + (yGrad - xGrad) * sMag) /*(3)*/
428 && tmp > Math.abs(xGrad * nwMag + (yGrad - xGrad) * nMag) /*(4)*/
429 ) {
430 magnitude[index] = gradMag >= MAGNITUDE_LIMIT ? MAGNITUDE_MAX : (int) (MAGNITUDE_SCALE * gradMag);
431 //NOTE: The orientation of the edge is not employed by this
432 //implementation. It is a simple matter to compute it at
433 //this point as: Math.atan2(yGrad, xGrad);
434 } else {
435 magnitude[index] = 0;
436 }
437 }
438 }
439 }
440
441 //NOTE: It is quite feasible to replace the implementation of this method
442 //with one which only loosely approximates the hypot function. I've tested
443 //simple approximations such as Math.abs(x) + Math.abs(y) and they work fine.
444 private float hypot(float x, float y) {
445 if (x == 0f) return y;
446 if (y == 0f) return x;
447 return (float) Math.sqrt(x * x + y * y);
448 }
449
450 private float gaussian(float x, float sigma) {
451 return (float) Math.exp(-(x * x) / (2f * sigma * sigma));
452 }
453
454 private void performHysteresis(int low, int high) {
455 //NOTE: this implementation reuses the data array to store both
456 //luminance data from the image, and edge intensity from the processing.
457 //This is done for memory efficiency, other implementations may wish
458 //to separate these functions.
459 Arrays.fill(data, 0);
460
461 int offset = 0;
462 for (int x = 0; x < width; x++) {
463 for (int y = 0; y < height; y++) {
464 if (data[offset] == 0 && magnitude[offset] >= high) {
465 follow(x, y, offset, low);
466 }
467 offset++;
468 }
469 }
470 }
471
472 private void follow(int x1, int y1, int i1, int threshold) {
473 int x0 = x1 == 0 ? x1 : x1 - 1;
474 int x2 = x1 == width - 1 ? x1 : x1 + 1;
475 int y0 = y1 == 0 ? y1 : y1 - 1;
476 int y2 = y1 == height -1 ? y1 : y1 + 1;
477
478 data[i1] = magnitude[i1];
479 for (int x = x0; x <= x2; x++) {
480 for (int y = y0; y <= y2; y++) {
481 int i2 = x + y * width;
482 if ((y != y1 || x != x1)
483 && data[i2] == 0
484 && magnitude[i2] >= threshold) {
485 follow(x, y, i2, threshold);
486 return;
487 }
488 }
489 }
490 }
491
492 private void thresholdEdges() {
493 for (int i = 0; i < picsize; i++) {
494 data[i] = data[i] > 0 ? -1 : 0xff000000;
495 }
496 }
497
498 private int luminance(float r, float g, float b) {
499 return Math.round(0.299f * r + 0.587f * g + 0.114f * b);
500 }
501
502 private void readLuminance() {
503 int type = sourceImage.getType();
504 if (type == BufferedImage.TYPE_INT_RGB || type == BufferedImage.TYPE_INT_ARGB) {
505 int[] pixels = (int[]) sourceImage.getData().getDataElements(0, 0, width, height, null);
506 for (int i = 0; i < picsize; i++) {
507 int p = pixels[i];
508 int r = (p & 0xff0000) >> 16;
509 int g = (p & 0xff00) >> 8;
510 int b = p & 0xff;
511 data[i] = luminance(r, g, b);
512 }
513 } else if (type == BufferedImage.TYPE_BYTE_GRAY) {
514 byte[] pixels = (byte[]) sourceImage.getData().getDataElements(0, 0, width, height, null);
515 for (int i = 0; i < picsize; i++) {
516 data[i] = (pixels[i] & 0xff);
517 }
518 } else if (type == BufferedImage.TYPE_USHORT_GRAY) {
519 short[] pixels = (short[]) sourceImage.getData().getDataElements(0, 0, width, height, null);
520 for (int i = 0; i < picsize; i++) {
521 data[i] = (pixels[i] & 0xffff) / 256;
522 }
523 } else if (type == BufferedImage.TYPE_3BYTE_BGR) {
524 byte[] pixels = (byte[]) sourceImage.getData().getDataElements(0, 0, width, height, null);
525 int offset = 0;
526 for (int i = 0; i < picsize; i++) {
527 int b = pixels[offset++] & 0xff;
528 int g = pixels[offset++] & 0xff;
529 int r = pixels[offset++] & 0xff;
530 data[i] = luminance(r, g, b);
531 }
532 } else {
533 throw new IllegalArgumentException("Unsupported image type: " + type);
534 }
535 }
536
537 private void normalizeContrast() {
538 int[] histogram = new int[256];
539 for (int i = 0; i < data.length; i++) {
540 histogram[data[i]]++;
541 }
542 int[] remap = new int[256];
543 int sum = 0;
544 int j = 0;
545 for (int i = 0; i < histogram.length; i++) {
546 sum += histogram[i];
547 int target = sum*255/picsize;
548 for (int k = j+1; k <=target; k++) {
549 remap[k] = i;
550 }
551 j = target;
552 }
553
554 for (int i = 0; i < data.length; i++) {
555 data[i] = remap[data[i]];
556 }
557 }
558
559 private void writeEdges(int pixels[]) {
560 //NOTE: There is currently no mechanism for obtaining the edge data
561 //in any other format other than an INT_ARGB type BufferedImage.
562 //This may be easily remedied by providing alternative accessors.
563 if (edgesImage == null) {
564 edgesImage = new BufferedImage(width, height, BufferedImage.TYPE_INT_ARGB);
565 }
566 edgesImage.getWritableTile(0, 0).setDataElements(0, 0, width, height, pixels);
567 }