图像和视频处理手册

图像和视频处理手册

英文原文

4.13 Gradient and Laplacian Edge Detection

1 Introduction

2 Gradient-based Methods

2.1 Continuous Gradient

2.2 Discrete Gradient Operators

3 Laplacian-based Methods

3.1 Continuous Laplacian

3.2 Discrete Laplacian Operators

3.3 The Laplacian of Guassian (Marr-Hildreth Operator)

3.4 Difference of Gaussian

4 Canny’s Method

5 Approaches for Color and Multispectral Images

6 Summary

References

1 Introduction

One of the most fundamental image analysis operations is edge detection. Edges are ofen vital clues toward the analysis and interpretation of image information, both in biological vision and computer image analysis. Some sort of edge detection capability is present in the visual systems of a wide variety of creatures, so it is obviously useful in their abilities to perceive their surroundings.

For this discussion, it is important to define what is and is not meant by the term “edge”. The everyday notion of an edge is usually a physical one, caused by either the shapes of physical objects in three dimensions or by their inherent material properties. Described in geometric terms, there are two types of physical edges:（1）the set of points along which there is an abrupt change in local orientation of a physical surface, and (2) the set of points describing the boundary between two or more materially distinct regions of a physical surface. Most of our perceptual senses, including vision, operate at a distance and gather information using receptors that work in, at most, two dimensions. Only the sense of touch, which requires direct contact to stimulate the skin’s pressure sensors, is capable of direct perception of objects in three-dimensional (3D) space. However, some physical edges of the second type may not be perceptible by touch because material differences-for instance different colors of paint-do not always produce distinct tactile sensations. Everyone first develops a working understanding of physical edges in early childhood by touching and handling every object within reach.

The imaging process inherently performs a projection from a 3D scene to a two-dimensional (2D) representation of that scene, according to the viewpoint of the imaging device. Because of this projection process, edges in images have a somewhat different meaning than physical edges. Although the precise definition depends on the

application context, an edge can generally be defined as a boundary or contour that separates adjacent image regions having relatively distinct characteristics according to some feature of interest. Most often this feature is gray level or luminance, but others, such as reflectance, color, or texture, are sometimes used. In the most common situation where luminance is of primary interest, edge pixels are those at the locations of abrupt gray-level change. To eliminate single-point impulses from consideration as edge pixels, one usually requires that edges be sustained along a contour; i.e., an edge point must be part of an edge structure having some minimum extent appropriate for the scale of interest. Edge detection is the process of determining which pixels are the edge pixels. The result of the edge detection process is typically an edge map, a new image that describes each original pixel’s edge classification and perhaps additional edge attributes, such as magnitude and orientation.

There is usually a strong correspondence between the physical edges of a set of objects and the edges in images containing views of those objects. Infants and young children learn this as they develop hand-eye coordination, gradually associating visual patterns with touch sensations as they feel and handle items in their vicinity. There are many situations, however, in which edges in an image do not correspond to physical edges. Illumination differences are usually responsible for this effect-for example, the boundary of a shadow cast across an otherwise uniform surface.

Conversely, physical edges do not always give rise to edges in images. This can also be caused by certain cases of lighting and surface properties. Consider what happens when one wishes to photograph a scene rich with physical edges-for example, a craggy mountain face consisting of a single type of rock. When this scene is imaged while the sun is directly behind the camera, no shadows are visible in the scene and hence shadow-dependent edges are nonexistent in the photo. The only edges in such a photo are produced by the differences in material reflectance, texture, or color. Since our rocky subject material has little variation of these types, the result is a rather dull photograph, because of the lack of apparent depth caused by the missing edges. Thus, images can exhibit edges having no physical counterpart, and they can also miss capturing edges that do. Although edge information can be very useful in the initial stages of such image processing and analysis tasks as segmentation, registration, and object recognition, edges are not completely reliable for these purposes.

If one defines an edge as an abrupt gray-level change, then the derivative, or gradient, is a natural basis for an edge detector. Figure 1 illustrates the idea with a continuous, one-dimensional (1D) example of a bright central region against a dark

background. The left-hand portion of the gray-level function f c(x)shows a smooth transition from dark to bright as x increases. There must be a point x0that marks the transition from the low-amplitude region on the left to the adjacent high-amplitude region in the center. The gradient approach to detecting this edge is to locate x0

where f c′(x)reaches a local maximum or, equivalently, f c′(x)reaches a local extremum, as shown in the second plot of Fig. 1. The second derivative, or Laplacianapproach, locates x0where a zero-crossing of f c"(x)occurs, as in the third plot of Fig. 1. The right-hand side of Fig. 1 illustrates the case for a falling edge located at x1.

To use the gradient or the Laplacian approaches as the basis for practical image edge detectors, one must extend the process to two dimensions, adapt to the discrete case, and somehow deal with the difficulties presented by real images. Relative to the 1D edges shown in Fig. 1, edges in 2D images have the additional quality of direction. One usually wishes to find edges regardless of direction, but a directionally sensitive edge detector can be useful at times. Also, the discrete nature of digital images requires the use of an approximation to the derivative. Finally, there are a number of problems that can confound the edge detection process in real images. These include noise, crosstalk or interference between nearby edges, and inaccuracies resulting from the use of a discrete grid. False edges, missing edges, and errors in edge location and orientation are often the result.

Because the derivative operator acts as a high-pass filter, edgedetectors based on it are sensitive to noise. It is easy for noiseinherent in an image to corrupt the real edges by shifting their apparent locations and by adding many false edge pixels. Unless care is taken, seemingly moderate amountsof noise are capable of overwhelming the edge detection process, rendering the results virtually useless. The wide variety of edge detection algorithmsdeveloped over the past three decades exists, in large part, because of the many ways proposed for dealing with noise andits effects. Most algorithms employ noise-suppression filteringof some kind before applying the edge detector itself. Some decomposethe image into a set of low-pass or bandpass versions,apply the edge detector to each, and merge the results. Still othersuse adaptive methods, modifying the edge detector’s parametersand behavior according to the noise characteristics of the imagedata.

An important tradeoff exists between correct detection of theactual edges and precise location of their positions. Edge detectionerrors can occur in two forms: false positives, in whichnonedge pixels are misclassified as edge pixels, and false

negatives,which are the reverse. Detection errors of both types tendto increase with noise, making good noise suppression very importantin achieving a high detection accuracy. In general, the potential for noise suppression improves with the spatial extent of the edge detection filter. Hence, the goal of maximumdetection accuracy calls for a large-sized filter. Errors in edge localizationalso increase with noise. To achieve good localization,however, the filter should generally be of small spatial extent.The goals of detection accuracy and location accuracy are thusput into direct conflict, creating a kind of uncertainty principlefor edge detection.

In this chapter, we cover the basics of gradient and Laplacianedge detection methods in some detail. Following each, we alsodescribe several of the more important and useful edge detectionalgorithms based on that approach. While the primary focus ison gray-level edge detectors, some discussion of edge detectionin color and multispectral images is included.

2 Gradient-Based Methods

2.1 Continuous Gradient

The core of gradient edge detection is, of course, the gradientoperator,?In continuous form, applied to a continuous-spaceimage, f c(x,y), the gradient is defined as

?f c x,y=ef c(x,y)

ex

i x+

ef c(x,y)

ey

i y（1）

Where i x, and i y are the unit vectors in the x and y directions.Notice that the gradient is a vector, having both magnitude anddirection. Its magnitude,?f c(x0,y0), measures the maximumrate of change in the intensity at the location(x0,y0). Its directionis that of the greatest increase in intensity; i.e., it points“uphill.”

To produce an edge detector, one may simply extend the 1-Dcase described earlier. Consider the effect of finding the localextrema of ?f c x,y or the local maxima of

?f c x,y=ef c(x,y)

ex

2

+

ef c(x,y)

ey

2

（2）

The precise meaning of “local” is very important here. If themaximaof Eq. (2) are foundover a2-D neighborhood, the resultis a set of isolated points rather than the desired edge contours.The problem stems from the fact that the gradient magnitude

isseldom constant along a given edge, so finding the 2-D local maximayields only the locally strongest of the edge contour points.To fully construct edge contours, it is better to apply Eq. (2) to a1-D local neighborhood, namely a line segment, whose directionis chosen to cross the edge. The situation is then similar to that ofFig. 1, where the point of locally maximum gradient magnitudeis the edge point. Now the issue becomes how to select the bestdirection for the line segment used for the search.

The most commonly used method of producing edge segmentsor contours from Eq. (2) consists of two stages: thresholdingand thinning. In the thresholding stage, the gradient magnitudeat every point is compared to a predefined threshold value,T. All points satisfying the following criterion are classified ascandidate edge points:

?f c(x,y)≥T (3)

The set of candidate edge points tends to form strips, whichhave positive width. Since the desire is usually for zero-widthboundary segments or contours to describe the edges, a subsequentprocessing stage is needed to thin the strips to the finaledge contours. Edge contours derived from continuous-spaceimages should have zero width because any local maxima of ?f c(x,y), along a line segment that crosses the edge, cannotbe adjacent points. For the case of discrete-space images, thenonzero pixel size imposes a minimum practical edge width.

Edge thinning canbe accomplished in a number of ways, depending on the application, but thinning by nonmaximum suppression is usually the best choice. Generally speaking, we wish to suppress any point that is not, in a 1-D sense, a local maximum in gradient magnitude. Since a 1-D local neighborhood search typically produces a single maximum, those points that are local maxima will form edge segments only one point wide. One approach classifies an edge-strip point asan edge point if its gradient magnitude is a local maximum in at least one direction. However, this thinning method sometimes has the side effect of creating false edges near strong edge lines [12]. It isalso somewhat inefficient because of the computation requiredto check along a number of different directions. Abetter, moreefficient thinning approach checks only a single direction, thegradient direction, to test whether a given point is a local maximumin gradient magnitude. The points that pass this scrutinyare classified as edge points. Looking in the gradient directionessentially searches perpendicular to the edge itself, producinga scenario similar to the 1-D case shown in Fig. 1. The methodis efficient because it is not necessary to search in multiple directions.It also tends to produces edge segments having goodlocalization accuracy. These characteristics make the gradientdirection, local extremum method quite

popular. The following

steps summarize its implementation.

1. Using one of the techniques described in the next section,compute ?f for all pixels.

2. Determine candidate edge pixels bythresholding all pixels’gradient magnitudes by T.

3. Thin by checking whether each candidate edge pixel’s gradientmagnitude is a local maximum along its gradientdirection. If so, classify it as an edge pixel.

Consider the effect of performing the thresholding and thinningoperations in isolation. If thresholding alone were done,the computational cost of thinning would be saved and the edgeswould show as strips or patches instead of thin segments. Ifthinning were done without thresholding, that is, if edge pointswere simply those having locally maximum gradient magnitude,many false edge points would likely result because of noise. Noisetends to create false edge points because some points in edge-freeareas happen to have locally maximum gradient magnitudes.The thresholding step of Eq. (3) is often useful to reducenoise prior to thinning. A variety of adaptive methods havebeen developed that adjust the threshold according to certainimage characteristics, such as an estimate of local signal-to-noiseratio. Adaptive thresholding can often do a better job ofnoise suppression while reducing the amount of edge fragmentation.

The edge maps in Fig. 3, computed from the original imagein Fig. 2, illustrate the effect of the thresholding and subsequentthinning steps.

The selection of the threshold value T is a tradeoff betweenthe wish to fully capture the actual edges in the image and thedesire to reject noise. Increasing T decreases sensitivity to noise atthe cost of rejecting the weakest edges, forcing the edge segmentsto become more broken and fragmented. By decreasing T, onecan obtain more connected and richer edge contours, but thegreater noise sensitivity is likely to produce more false edges. Ifonly thresholding is used, as in Eq. (3) and Fig. 3(a), the edgestrips tend to narrow as T increases and widen as it decreases.Figure 4 compares edge maps obtained from several differentthreshold values.

Sometimes a directional edge detector is useful. One can beobtained by decomposing the gradient into horizontal and verticalcomponents and applying them separately. Expressed in thecontinuous domain, the operators become:

ef c(x,y)

≥T for edges in the y direction

ex

ef c(x,y)

ey

≥T for edges in the x direction

An example of directional edge detection is illustrated in Fig. 5.

Adirectional edge detector can be constructed for any desireddirection by using the directional derivative along a unit vector，

?f c

=?f c x,y?n

=ef c(x,y)

cosθ+

ef c(x,y)

sinθ (4)

whereθis the angle of n relative to the positive x axis. The directionalderivative is most sensitive to edges perpendicular to n .

The continuous-space gradient magnitude produces an isotropicor rotationally symmetric edge detector, equally sensitiveto edges in any direction. It is easy to show why ?f isisotropic. In addition to the original X-Y coordinate system, letus introduce a new system, X'-Y', which is rotated by an angleof φ, relative to X-Y. Let n x′and n y′be the unit vectors in the x'and y' directions, respectively. For the gradient magnitude to beisotropic, the same result must be produced in both coordinatesystems, regardless ofφ. Using Eq. (4) along with abbreviatednotation, we find the partial derivatives with respect to the newcoordinate axes are

f x′=?f?n x′=f x cosφ+f y sinφ

f y′=?f?n y′=f x sinφ+f y cosφ

Now let us examine the gradient magnitude in the new coordinatesystem:

?f=f x′2+f y′2

=(f x cosφ+f y sinφ)2+(?f x sinφ+f y cosφ)2

=f x2cos2φ+2f x f y cosφsinφ+f y2sin2φ+f x2sin2φ?2f x f y sinφcosφ+f y2cos2φ=f x2cos2φ+sin2φ+f y2cos2φ+sin2φ

f x2+f y2

So the gradient magnitude in the new coordinate system matches

that in the original system, regardless of the rotation angle,φ.

Occasionally, one may wish to reduce the computation load ofEq. (2) by approximating the square root with a computationallysimpler function. Three possibilities are

?f c(x,y)≈max ?f c(x,y)

?x

,

?f c(x,y)

?y

(5)

≈?f c x,y

+

?f c x,y

(6)

≈max ?f c x,y

?x

,

?f c x,y

?y

+1

4

min

?f c x,y

?x

,

?f c x,y

?y

(7)

One should be aware that approximations of this type may alter the properties of the gradient somewhat. For instance, the approximated gradient magnitudes of Eqs.(5), (6),and (7) are not isotropic and produce their greatest errors for purely diagonally oriented edges. All three estimates are correct only for the pure horizontal and vertical cases. Otherwise, Eq. (5) consistently underestimates the true gradient magnitude while Eq. (6)overestimates it. This makes Eq. (5) biased against diagonal edges and Eq. (6)biased toward them. The estimate of Eq. (7) is by far the most accurate of the three.

2.2 Discrete Gradient Operators

In the continuous-space image, f c(x,y), let xand y representthe horizontal and vertical axes, respectively. Let the discretespacerepresentation of f c(x,y)be f(n1,n2), with n1describingthe horizontal position and n2describing the vertical. For useon discrete-space images, the continuous gradient's derivativeoperators must be approximated in discrete form. The approximationtakes the form of a pair of orthogonally oriented filters,h1(n1,n2)and h2(n1,n2), which must be separately convolvedwith the image. Based on Eq. (l), the gradient estimate is

?f n1,n2=f1n1,n2i n

1+f2n1,n2i n

2

where

f1n1,n2=f n1,n2?h1n1,n2

f2n1,n2=f n1,n2?h2n1,n2

Two filters are necessary because the gradient requires the computationof an orthogonal pair of directional derivatives. Thegradient magnitude and direction estimates can then be computedas follows:

?f n1,n2=f12n1,n2+f22n1,n2

∠f n1,n2=tan?1f2(n1,n2)

112

(8)

Each of the filters implements a derivative and should not respondto a constant, so the sum of its coefficients must alwaysbe zero. Amore general statement of this property is describedlater in this chapter by Eq. (10).

There are many possible derivative-approximation filters foruse in gradient estimation. Let us start with the simplest case.Two simple approximation schemes for the horizontal derivative are, for the first and central differences, respectively,

f1n1,n2=f n1,n2?f n1?1,n2

f1n1,n2=1

f n1+1,n2?f(n1?1,n2)

The scaling factor of 1/2 for the central difference is causedby the two-pixel distance between the nonzero samples. Theorigin positions for both filters are usually

set at （n1,n2）. Thegradient magnitude threshold value can be easily adjusted tocompensate for any scaling, so we omit the scale factor fromhere on. Both of these differences respond most strongly to verticallyoriented edges and do not respond to purely horizontaledges. The case for the vertical direction is similar, producinga derivative approximation that responds most strongly tohorizontal edges. These derivative approximations can be expressedas filter kernels, whose impulse responses, h1(n1,n2)and h2(n1,n2), are as follows for the first and central differences,respectively:

h1n1,n2=00

?11, h2n1,n2=10

?10

h1n1,n2=

000

?101

000

, h2n1,n2=

010

000

0?10

Boldface elements indicate the origin position.

If used to detect edges, the pair of first difference filters abovepresents the problem that the zero crossings of its two [-1,1] derivative kernels lie at different positions. This prevents the twofilters from measuring horizontal and vertical edge characteristicsat the same location, causing error in the estimated gradient.The central difference, caused by the common center of itshorizontal and vertical differencing kernels, avoids this positionmismatch problem. This benefit comes at the costs of larger filtersize and the fact that the measured gradient at a pixel (n1,n2)does not

actually consider the value of that pixel.

Rotating the first difference kernels by an angle of π4andstretching the grid a bit produces the h1(n1,n2)and h2(n1,n2)kernels for the Roberts operator:

01?1010 0?1

The Roberts operator's component filters are tuned for diagonaledges rather than vertical and horizontal ones. For use in an edgedetector based on the gradient magnitude, it is important onlythat the two filters be orthogonal. They need not be aligned withthe n1 and n2axes. The pair of Roberts filters have a commonzero-crossing point for their differencing kernels. This commoncenter eliminates the position mismatch error exhibited by thehorizontal-vertical first difference pair, as described earlier. Ifthe origins of the Roberts kernels are positioned on the +1 samples,as is sometimes found in the literature, then no commonenter point exists for their first differences.

The Roberts operator, like any simple first-difference gradientoperator, has two undesirable characteristics. First, the zerocrossing of its [-1,1] diagonal kernel lies off grid, but the edgelocation must be assigned to an actual pixel location, namelythe one at the filter's origin. This can create edge location biasthat may lead to location errors approaching the interpixel distance.If we could use the central difference instead of the firstdifference, this problem would be reduced because the centraldifference operator inherently constrains its zero crossing to anexact pixel location.

The other difficulty caused by the first difference is its noisesensitivity. In fact, both the first- and central-difference derivativeestimators are quite sensitive to noise. The noise problemcan be reduced somewhat by incorporating smoothing into eachfilter in the direction normal to that of the difference. Consideran example based on the central difference in one direction forwhich we wish to smooth along the orthogonal direction witha simple three-sample average. To that end, let us define theimpulse responses of two filters:

?a n1= 1 1 1,?b n2=?1 0 1

Since h a is a function only of n1and h b depends only on n2,one can simply multiply them as an outer product to form aseparable derivative filter that incorporates smoothing:

h a n1h b n2=h1(n1,n2)

1 1 1?101=

?101

?101

?101

Repeating this process for the orthogonal case produces thePrewitt operator:

?101?101?101

111 000?1?1?1

The Prewitt edge gradient operator simultaneously accomplishesdifferentiation in one coordinate direction, using the centraldifference, and noise reduction in the orthogonal direction, bymeans of local averaging. Because it uses the central differenceinstead of the first difference, there is less edge-location bias.

In general, the smoothing characteristics can be adjusted bychoosing an appropriate low-pass filter kernel in place of thePrewitt's three-sample average. One such variation is the Sobeloperator, one of the most widely used gradient edge detectors:

?101?202?101

121 000?1?2?1

Sobel's operator is often a better choice than Prewitt's becausethe low-pass filter produced by the [1 2 1] kernel results in asmoother frequency response compared to that of [1 1 1].

The Prewitt and Sobel operators respond differently to diagonaledges than to horizontal or vertical ones. This behavioris a consequence of the fact that their filter coefficients do notcompensate for the different grid spacings in the diagonal andthe horizontal directions. The Prewitt operator is less sensitiveto diagonal edges than to vertical or horizontal ones. The oppositeis true for the Sobeloperator . A variation designed forequal gradient magnitude response to diagonal, horizontal, andvertical edges is the Frei-Chen operator:

?101? 202?101

121 000?1? 2?1

However, even the Frei-Chen operator retains some directionalsensitivity in gradient magnitude, so it is not truly isotropic.

The residual anisotropy is caused by the fact that the differenceoperators used to approximate Eq. (1) are not rotationally symmetric.Merron and Brady [15] describe a simple method forgreatly reducing the residual directional bias by using a set offour difference operators instead of two. Their operators are orientedin increments of π

4radians, adding a pair of diagonalones to the original horizontal and vertical pair. Averaging thegradients produced by the diagonal operators with those of thenondiagonal ones allows their complementary directional biasesto reduce the

overall anisotropy. However, Ziou and Wang have described how an isotropic gradient applied to a discretegrid tends to introduce some anisotropy. They have also analyzedthe errors of gradient magnitude and direction as a function ofedge translation and orientation for several detectors. Figure 6shows the results of performing edge detection on an exampleimage by applying the discrete gradient operators discussedso far.

Haralick's facet model [8] provides another way of calculatingthe gradient in order to perform edge detection. In the slopedfacet model, a small neighborhood is parameterized byαn2+βn1+γdescribing the plane that best fits the gray levels in thatneighborhood. The plane parametersαand βcan be used tocompute the gradient magnitude:

?f c(n1,n2)=22

The facet model also provides means for computing directionalderivatives, zero crossings, and a variety of other useful operations.

Improved noise suppression is possible with increased kernelsize. The additional coefficients can be used to betterapproximate the desired continuous-space noise-suppression filter.Greater filter extent can also be used to reduce directionalsensitivity by more accurately modeling an ideal isotropic filter.However, increasing the kernel size will exacerbate edge localizationproblems and create interference between nearby edges.Noise suppression can be improved by other methods as well.Papers by Bovik[3] and Hardie and Boncelet [9] are just two thatdescribe the use of edge-enhancing prefilters, which simultaneouslysuppress noise and steepen edges prior to gradient edgedetection.

3 Laplacian-Based Methods

3.1 Continuous Laplacian

The Laplacian is defined as

?2f c x,y=???f c x,y=?2f c(x,y)

2

+

?2f c(x,y)

2

(9)

The zero crossings of ?2f c x,y occur at the edge points of f c x,y because of the second derivative action (see Fig. 1).Laplacian-based edge detection has the nice

property that it producesedges of zero thickness, making edge-thinning steps unnecessary.This is because the zero crossings themselves definethe edge locations.

The continuous Laplacian is isotropic, favoring no particularedge orientation. Consequently, its second partial terms inEq. (9) can be oriented in any direction as long as they remainperpendicular to each other. Consider an ideal, straight, andnoise-free edge oriented in an arbitrary direction. Let us realignthe first term of Eq.

(9) parallel to that edge and the second termperpendicular to it. The first term then generates no response atall because it acts only along the edge. The second term producesazero crossing at the edge position along its edge-crossing profile.

An edge detector based solely on the zero crossings of the continuouLaplacian produces closed edge contours if the image,f(x,y),meets certain smoothness constraints [20]. The contoursare closed because edge strength is not considered, so eventhe slightest, most gradual intensity transition produces a zerocrossing. In effect, the zero-crossing contours define the boundariesthat separate regions of nearly constant intensity in theoriginal image. The second derivative zero crossings occur at thelocal extrema of the first derivative (see Fig. l), but many zerocrossings are not local maxima of the gradient magnitude. Somelocal minima of the gradient magnitude give rise to phantomedges, which can be largely eliminated by appropriately thresholdingthe edge strength. Figure 7 illustrates a 1-D example of aphantom edge.

Noise presents a problem for the Laplacian edge detectorin several ways. First, the second-derivative action of Eq. (9)makes the Laplacian even more sensitive to noise than thefirst-derivative-based gradient. Second, noise produces manyfalse edge contours because it introduces variation to theconstant-intensity regions in the noise-free image. Third, noisealters the locations of the zero-crossing points, producing locationerrors along the edge contours. The problem of noise-inducedfalse edges can be addressed by applying an additionaltest to the zero-crossing points. Only the zero crossingsthat satisfy this new criterion are considered edge points. Onecommonly-used technique classifies a zero crossing as an edgepoint if the local gray-level variance exceeds a threshold amount.Another method is to select the strong edges by thresholding thegradient magnitude or the slope of the Laplacian output at thezero crossing. Both criteria serve to reject zero crossing pointswhich are more likely caused by noise than by a real edge in theoriginal scene. Of course, thresholding the zero crossings in thismanner tends to break up the closed contours.

Like any derivative filter, the continuous-space Laplacian filter,h c(x,y), has this important property:

h c x,y dxdy =0 10 +∞

?∞+∞?∞

In other words, h c (x,y)is a surface bounding equal volumesabove and below zero.Consequently,?2f c (x,y)will also haveequal volumes above and below zero. This property eliminatesany response that is due to the constant or DC bias contained in f c (x,y). Without DC bias rejection, the filter's edge detectionperformance would be compromised.

3.2 Discrete Laplacian Operators

It is useful to construct a filter to serve as the Laplacian operatorwhen applied to a discrete-space image. Recall that the gradient,which is a vector, required a pair of orthogonal filters. The Laplacianis a scalar. Therefore, a single filter,h(n 1,n 2), is sufficientfor realizing a Laplacian operator. The Laplacian estimate for animage,f(n 1,n 2),is then

?

2f n 1,n 2 =f(n 1,n 2)?h(n 1,n 2) One of the simplest Laplacian operators can be derived asfollows. First needed is an approximation to the derivative in x,so let us use a simple first difference.

?f c (x,y)?x

→f x n 1,n 2 =f n 1+1,n 2 ?f n 1,n 2 （11） The second derivative in x can be built by applying the firstdifference to Eq. (11). However, we discussed earlier how thefirst difference produces location errors because its zero crossinglies off grid. This second application of a first difference can beshifted to counteract the error introduced by the previous one:

?2f c x,y ?x 2

→f xx n 1,n 2 =f x n 1,n 2 ?f x n 1?1,n 2 12 Combining the two derivative-approximation stages fromEqs. (1 1) and (12) produces

?2f c x,y 2

→f xx n 1,n 2 =f n 1+1,n 2 ?2f n 1,n 2 +f n 1?1,n 2

= 1?2 1 (13)

Proceeding in an identical manner for y yields

?2f c x,y ?y 2

→f yy n 1,n 2 =f n 1,n 2+1 ?2f n 1,n 2 +f(n 1,n 2?1)

=

1

?2

1

(14)

Combining the x and y second partials of Eqs. (13) and (14)produces a filter,h(n1,n2), which estimates the Laplacian:

?2f c x,y→?2f n1,n2

=f xx n1,n2+f yy(n1,n2)

=f n1+1,n2+f n1?1,n2+f n1,n2+1

+f n1,n2?1?4f(n1,n2)

=1?2 1+

1

?2

1

=

010

1?41

010

Other Laplacian estimation filters can be constructed by usingthis method of designing a pair of appropriate 1-D secondderivative filters and combining them into a single 2-D filter.The results depend on the choice of derivative approximator,the size of the desired filter kernel, and the characteristics ofany noise-reduction filtering applied. Two other 3 x 3examplesare

111 1?81 111,

?12?1 2?42?12?1

In general, a discrete-space smoothed Laplacian filter can be easilyconstructed by sampling an appropriate continuous-spacefunction, such as the Laplacian of Gaussian. When constructinga Laplacian filter, make sure that the kernel's coefficients sum tozero in order to satisfy the discrete form of Eq. (10). Truncationeffects may upset this property and create bias. If so, the filter coefficientsshould be adjusted in a waythat restores proper balance.

Locating zero crossings in the discrete-space image,?2f n1,n2is fairly straightforward. Each pixel should be comparedto its eight immediate neighbors; a four-way neighborhoodcomparison, while faster, may yield broken contours. If apixel, p,differs in sign with its neighbor, q, an edge lies betweenthem. The pixel, p, is classified as a zero crossing if

?2f(p)≤?2f(q) (15)

3.3 The Laplacian of Gaussian(Marr-Hildreth Operator)

It is common for a single image to contain edges having widelydifferent sharpnesses and scales, from blurry and gradual to crispand abrupt. Edge scale

information is often useful as an aid towardimage understanding. For instance, edges at low resolutiontend to indicate gross shapes, whereas texture tends to becomeimportant at higher resolutions. An edge detected over a widerange of scale is more likely to be physically significant in thescene than an edge found only within a narrow range of scale.Furthermore, the effects of noise are usually most deleterious atthe finer scales.

Marr and Hildrethadvocated the need for an operatorthat can be tuned to detect edges at a particular scale. Theirmethod is based on filtering the image with a Gaussian kernelselected for a particular edge scale. The Gaussian smoothingoperation serves to band-limit the image to a small range of frequencies,reducing the noise sensitivity problem when detectingzero crossings. The image is filtered over a variety of scales andthe Laplacian zero crossings are computed at each. This producesa set of edge maps as a function of edge scale. Each edge point

can be considered to reside in a region of scale space, for which edge point location is a function of x, y, and σ. Scale space has been successfully used to refine and analyze edge maps [22].

The Gaussian has some very desirable properties that facilitatethis edge detection procedure. First, the Gaussian functionis smooth and localized in both the spatial and frequencydomains, providing a good compromise between the need foravoiding false edges and for minimizing errors in edge position.In fact, Torre and Poggio [20] describe the Gaussian asthe only real-valued function that minimizes the product ofspatial- and frequency-domain spreads. The Laplacian ofGaussian essentially acts as a bandpass filter because of its differentialand smoothing behavior. Second, the Gaussian is separable,which helps make computation very efficient.

Omitting the scaling factor, the Gaussian filter can be writtenas

g c x,y=exp ?x2+y2

2σ2

16

Its frequency response, GΩ

x ,Ω

y

,is also Gaussian:

GΩ

x ,Ω

y

=2πσ2exp ?

σ2

2

(Ω

x

2

+Ω

y

2

)

The σparameter is inversely related to the cutoff frequency.

Because the convolution and Laplacian operations are bothlinearand shift invariant, their computation order can be interchanged：

?2f c x,y?g c x,y=?2g c x,y?f c x,y (17)

Here we take advantage of the fact that the derivative is a linearoperator. Therefore, Gaussian filtering followed by differentiationis the same as filtering with the derivative of a Gaussian.The right-hand side of Eq. (17) usually provides for more efficientcomputation since ?2g c x,y can be prepared in advanceas a result of its image independence. The Laplacian of Gaussian(LOG) filter,h c x,y , therefore has the folIowing impulse

response:

h c x,y =?2g c x,y

=x 2+y 2?2σ2σ4exp x 2+y 2

2σ2

(18) To implement the LOG in discrete form, one may constructa filter,h(n 1,n 2)by sampling Eq. (18) after choosing a valueof σ, then convolving with the image. If the filter extent is notsmall, it is usuallymore efficient to workin the frequency domainby multiplying the discrete Fourier transforms of the filter andthe image, then inverse transforming the result. The fast Fouriertransform, or FFT, is the method of choice for computing thesetransforms.

Although the discrete form of Eq. (18) is a 2-D filter, Chenet al. [6] have shown that it is actually the sum of two separablefilters because the Gaussian itself is a separable function. By constructingand applying the appropriate 1-D filters successively tothe rows and columns of the image, the computational expenseof 2-D convolution becomes unnecessary. Separable convolutionto implement the LoG is roughly 1-2 orders of magnitudemore efficient than 2-D convolution. If an image is M x Min size, the number of operations at each pixel is M 2 for 2-Dconvolution and only 2 M if done in a separable, 1-D manner.

Figure 8 shows an example of applying the LOG using various σ values. Figure 8(d) includes a gradient magnitude threshold,which suppresses noise and breaks contours. Lim [12] describesan adaptive thresholding scheme that produces better results.

Equation (18) has the shape of a sombrero or “Me xi can hat.”Figure 9 shows a perspectiveplot of 2(,)c g x y ?and its frequencyresponse, {}2(,)c F g x y ?. This profile closely mimics the responseof the spatial receptive field found in biological vision.Biological receptive fields have been shown to have a circularlysymmetric impulse response, with a central excitory region surroundedby an inhibitory band.

When sampling the LOG to produce a discrete version, it isimportant to size the

filter large enough to avoid significant truncationeffects. A good rule of thumb is to make the filter at leastthree times the width of the LOG’s central excitory lobe

[16].Siohan [19] describes two approaches for the practical designof LOG filters. The errors in edge location produced by the LOGhave been analyzed in some detail by Berzins[2]

3.4 Difference of Gaussian

The Laplacian of Gaussian of Eq. (18) can be closely approximatedby the difference of two Gaussians having properly chosenscales. The difference of Gaussian (DOG) filter is

12(,)(,)(,),c c c h x y g x y g x y =- where

21

1.6σσ≈ and 12,c c g g are evaluated by using Eq. (16). However, the LOG isusually preferred because it is theoretically optimal and its separabilityallows for efficient computation [14]. For the same accuracyof results, the DOG requires a slightly larger filter size [10].

The technique of unsharp masking, used in photography, isbasically a difference of Gaussians operation done with lightand negatives. Unsharp masking involves making a somewhatblurry exposure of an original negative onto a new piece of film.When the film is developed, it contains a blurred and inverted-rightnessversion of the original negative. Finally, a print ismade from these two negatives sandwiched together, producinga sharpened image with the edges showing an increasedcontrast.

Nature uses the difference of Gaussians as a basis for the architectureof the retina’s visual receptive field. The spatial -domainimpulse response of a photoreceptor cell in the mammalianretina has a roughly Gaussian shape. The photoreceptor outputfeeds into horizontal cells in the adjacent layer of neurons. Eachhorizontal cell averages the responses of the receptors in its immediateneighborhood, producing a Gaussian-shaped impulseresponse with a higher σ than that of a single photoreceptor.Both layers send their outputs to the third layer, where bipolarneurons

subtract the high-σneighborhood averages from thecentral photoreceptors’ low-σresponses. This produces a biologicalrealization of the difference-of-Gaussian filter, approximatingthe behavior of the Laplacian of Gaussian. The retinaactually implements DOG bandpass filters at several spatial frequencies [13].

4 Canny’s Method

Canny’s method [4] uses the concepts of both the first andsecond derivatives in a very effective manner. His is a classicapplication of the gradient approach to edge detection in thepresence of additive white Gaussian noise, but it also incorporateselements of the Laplacian approach. The method has threesimultaneous goals: low rate of detection errors, good edge localization,and only a single detection response per edge. Cannyassumed that false-positive and false-negative detection errorsare equally undesirable and so gave them equal weight. He furtherassumed that each edge has nearly constant cross sectionand orientation, but his general method includes a way to effectivelydeal with the cases of curved edges and corners. Withthese constraints, Canny determined the optimal 1-D edge detectorfor the step edge and showed that its impulse response canbe approximated fairly well by the derivative of a Gaussian.

An important action of Canny’s edge detector is to preventmultiple responses per true edge. Without this criterion, theoptimal step-edge detector would have an impulse response inthe form of a truncated signum function. (The signum functionproduces + 1 for any positive argument and - 1 for any negativeargument.) But this type of filter has high bandwidth, allowingnoise or texture to produce several local maxima in the vicinityof the actual edge. The effect of the derivative of Gaussian is toprevent multiple responses by smoothing the truncated signumin order to permit only one response peak in the edge neighborhood.The choice of variance for the Gaussian kernel controlsthe filter width and the amount of smoothing. This defines thewidth of the neighborhood in which only a single peak is tobe allowed. The variance selected should be proportional to theamount of noise present. If the variance is chosen too low, the filtercan produce multiple detections for a single edge; if too high,edge localization suffers needlessly. Because the edges in a givenimage are likely to differ in signal-to-noise ratio, a single-filterimplementation is usually not best for detecting them. Hence, athorough edge detection procedure should operate at differentscales.