IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8, NO. 6,NOVEMBER 1986
A Computational Approach to Edge Detection JOHN CANNY, MEMBER, IEEE
679
Abstract-Thispaperdescribesacomputationalapproachtoedge
detection. The success of the approach depends on the definition of a
comprehensivesetofgoalsforthecomputationofedgepoints.These
goalsmustbepreciseenoughtodelimitthedesiredbehaviorofthe
detector while making minimal assumptions about the form of the so-
lution.Wedefinedetectionandlocalizationcriteriaforaclassofedges,
andpresentmathematicalformsforthesecriteriaasfunctionalsonthe
operatorimpulseresponse.A thirdcriterionisthenaddedtoensure
thatthedetectorhasonlyoneresponseto-asingleedge.Weusethe
criteriainnumericaloptimizationtoderivedetectorsforseveralcom- cessedbyanedgedetectorbeforematchingisdone[19], mon imagefeatures,includingstepedges.On specializingtheanalysis
tostepedges,wefindthatthereisanaturaluncertaintyprinciplebe-
tweendetectionandlocalizationperformance,whicharethetwomain
goals.Withthisprinciplewederiveasingleoperatorshapewhichis
optimalatanyscale.Theoptimaldetectorhasasimpleapproximate
implementationinwhichedgesaremarkedatmaximaingradientmag-
nitudeofaGaussian-smoothedimage.We extendthissimpledetector
usingoperatorsofseveralwidthstocopewithdifferentsignal-to-noise
ratiosintheimage.Wepresentageneralmethod,calledfeaturesyn-
thesis,forthefine-to-coarseintegrationofinformationfromoperators
atdifferentscales.Finallyweshowthatstepedgedetectorperfor- occurintheimageshouldnotbemissedandthattherebe manceimprovesconsiderablyastheoperatorpointspreadfunctionis nospuriousresponses.Inaltheabovecases,systemper- extended along the edge. This detection scheme uses several elongated
operatorsateachpoint,andthedirectionaloperatoroutputsarein- tegratedwiththegradientmaximum detector.
Index Terms-Edge detection, feature extraction, image processing, machine vision, multiscale image analysis.
I. INTRODUCTION
EDGE detectorsofsomekind,particularlystepedge detectors, have been an essential part of many com- putervisionsystems.Theedgedetectionprocessserves tosimplifytheanalysisofimagesbydrasticallyreducing
formancewillbehamperedbyedgedetectorerrors.The second criterion is that the edge points be well localized. That is, the distance between the points marked by the detector and the “center” of the true edge should be min- imized. This is particularly true of stereo and shape from motion,wheresmalldisparitiesaremeasuredbetweenleft and right images or between images produced at slightly differenttimes.
In this paper we will develop a mathematical form for thesetwocriteriawhichcanbeusedtodesigndetectors forarbitraryedges.Wewillalsodiscoverthatthefirsttwo criteriaarenot”tight”enough,andthatitisnecessary
theamountofdatatobeprocessed,whileatthesametime
preservingusefulstructuralinformationaboutobject toaddathirdcriteriontocircumventthepossibilityof boundaries.Thereiscertainlyagreatdealofdiversityin multipleresponsestoasingleedge.Usingnumericalop- theapplicationsofedgedetection,butitisfeltthatmany timization,wederiveoptimaloperatorsforridgeandroof
applicationsshareacommonsetofrequirements.These requirementsyieldanabstractedgedetectionproblem,the solutionofwhichcanbeappliedinanyoftheoriginal problemdomains.
Weshouldmentionsomespecificapplicationshere.The Binford-Hornlinefinder[14]usedtheoutputofanedge
Manuscript received December 10, 1984; revised November 27, 1985. Recommended for acceptance by S. L. Tanimoto. This work was supported in part by the System Development Foundation, in part by the Office of NavalResearchunderContractN00014-81-K-0494,andinpartbytheAd- vanced Research Projects Agency under Office of Naval Research Con- tracts N00014-80-C-0505 and N00014-82-K-0334.
The author is with the Artificial Intelligence Laboratory, Massachusetts InstituteofTechnology, Cambridge, MA 02139.
edges.Wewillthenspecializethecriteriaforstepedges andgiveaparametricclosedformforthesolution.Inthe processwewilldiscoverthatthereisanuncertaintyprin- ciple relating detection and localization of noisy step edges,andthatthereisadirecttradeoffbetweenthetwo. Oneconsequenceofthisrelationshipisthatthereisasin- gle unique “shape” of impulse response for an optimal stepedgedetector,andthatthetradeoffbetweendetection and localization can be varied by changing the spatial width of the detector. Several examples of the detector performanceonrealimageswillbegiven.
I. ONE-DIMENSIONAL FORMULATION
To facilitate the analysis we first consider one-dimen-
IEEE Log Number 8610412.
sional edge profiles. That is, we will assume that 0162-8828/86/1100-0679$01.00ý 1986IEEE
two-
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
detectorasinputtoaprogramwhichcouldisolatesimple geometric solids. More recently the model-based vision systemACRONYM[3]usedanedgedetectorasthefront endtoasophisticatedrecognitionprogram.Shapefrom motion[29],[13]canbeusedtoinferthestructureof three-dimensionalobjectsfromthemotionofedgecon- tours or edge points in the image plane. Several modem theoriesofstereopsisassumethatimagesareprepro-
[20].Beattie[1]describesanedge-basedlabelingscheme forlow-levelimageunderstanding.Finally,somenovel methodshavebeensuggestedfortheextractionofthree- dimensionalinformationfromimagecontours,namely shapefromcontour[27]andshapefromtexture[31].
Inaloftheseexamplestherearecommon criteriarel- evanttoedgedetectorperformance.Thefirstandmost obviousislowerrorrate.Itisimportantthatedgesthat
680
IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8,NO. 6.NOVEMBER 1986
(a)
(b)
(c)
(d)
(e)
be
~.0 – 0
problemthenbecomesoneoffindingthefilterwhichgives thebestperformancewithrespecttothecriteriagivenbe- low. For example, the filter in Fig. l(d) performs much better than Fig. l(b) on this example, because the re- sponse ofthelatterexhibitsseverallocalmaxima inthe regionoftheedge.
Insummary, thethreeperformancecriteriaare as fol- lows:
1)Gooddetection.Thereshouldbealowprobability
criteriondidnotcapturethemultipleresponserequire- mentandithadtobemadeexplicit.
12 Az
Fig.1.(a)A noisystepedge.(b)Differenceofboxesoperator. (c)Dif- ference of boxes operator applied to the edge. (d) First derivative of Gaussian operator. (e) First derivative of Gaussian applied to the edge.
dimensionaledgeslocallyhaveaconstantcross-section insome direction.Thiswouldbetrueforexample,of smoothedgecontoursorofridges,butnottrueofcorners. Wewillassumethattheimageconsistsoftheedgeand additive white Gaussian noise.
Thedetectionproblemisformulatedas follows:We be-
ginwithan edgeofknowncross-sectionbathedinwhite
GaussiannoiseasinFig.l(a),whichshowsastepedge.
We convolve this with a filter whose impulse response
couldbeillustratedbyeitherFig.1(b)or(d).Theoutputs itlycapturedinthefirstcriterionsincewhentherearetwo oftheconvolutionsare shown,respectively,inFig.l(c) responses tothesame edge,one ofthemmustbeconsid- and(e).Wewillmarkthecenterofanedgeatalocal eredfalse.However,themathematicalformofthefirst maximumintheoutputoftheconvolution.Thedesign
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
offailingtomarkrealedgepoints,andlowprobabilityof falselymarkingnonedgepoints.Sinceboththeseproba- bilitiesaremonotonicallydecreasingfunctionsoftheout- putsignal-to-noiseratio,thiscriterioncorrespondsto maximizing signal-to-noise ratio.
2) Good localization. The points marked as edge points bytheoperatorshouldbeas closeas possibletothecenter ofthetrueedge.
3)Onlyone responsetoa singleedge.Thisisimplic-
A.DetectionandLocalizationCriteria
A crucialstepinour methodistocapturetheintuitive criteriagivenaboveinamathematicalformwhichisread- ilysolvable.We dealfirstwithsignal-to-noiseratioand localization.Lettheimpulseresponseofthefilterbef(x), anddenotetheedgeitselfbyG(x).We willassume that theedgeiscenteredatx = 0.Thentheresponseofthe
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION 681
filtertothisedgeatitscenterHGisgivenbyaconvolution ric,andthatitsderivativesofoddorders[whichappear integral: in the coefficients of even order in (5)] are zero at the +w origin.Equations(4)and(5)give
HG=J G(-x)f(x)dx (1) -w
HG(O)x0 = -H(XO) (6)
assumingthefilterhasafiniteimpulseresponsebounded Now Hx(xo)isaGaussianrandomquantitywhosevari- by[-W,W].Theroot-mean-squaredresponsetothe anceisthemean-squaredvalueofHn(xo),andisgiven
noisen(x)only,willbe by
H,=noLw f2(x)dx] (2)
E[H,(XO)2]= no2
+ Ww -w
f ‘2(X)dx
(7)
wheren2isthemean-squarednoiseamplitudeperunit whereE[y]istheexpectationvalueofy.Combiningthis length.Wedefineourfirstcriterion,theoutputsignal-to- resultwith(6)andsubstitutingforHG(0)gives
noiseratio,asthequotientofthesetwo responses. | G(-x)f(x)dx
+w
n2 f,2(X)dX
of xo. The localization is defined as the reciprocal of 6xo. r+W
SNR =
nO, f2(x)dZx
For the localization criterion, we want some measure which increases as localization improves, and we will use the reciprocal of the root-mean-squared distance of the marked edge from the center of the true edge. Since we have decided to mark edges at local maxima in the re- sponse of the operatorf(x), the first derivative of the re- sponse will be zero at these points. Note also that since edgesarecenteredatx= 0,intheabsenceofnoisethere shouldbealocalmaximumintheresponseatx = O.
LetHn(x)betheresponseofthefiltertonoiseonly,and HG(x)beitsresponsetotheedge,andsupposethereisa localmaximuminthetotalresponseatthepointx=xO. Then we have
Hn(XO)+ HG(x0)= 0. (4) The Taylor expansion of H&(xo) about the origin gives
H&(xo)= HG(O)+ HG(0)x0+ O(x0). (5)
Localization
3 G'(-x)f'(x)dx
nO f”‘2(x)dx
(9)
I+W
| W Thedisplacementxooftheactualmaximumisassumed no f2(X)dx no f 2(X)dx
ByassumptionHG(0)= 0,i.e.,theresponseofthefil- ter in the absence of noise has a local maximum at the origin, so the first term in the expansion can be ignored.
G(-x)f(x)dx|
G'(-x)f'(x)dx
to be small so we will ignore quadratic and higher terms. In fact by a simple argument we can show that ifthe edge G(x)iseithersymmetricorantisymmetric,aleventerms in xo vanish. Suppose G(x) is antisymmetric, and express
f(x) as a sum of a symmetric component and an antisym- metriccomponent.Theconvolutionofthesymmetric componentwithG(x)contributesnothingtothenumerator oftheSNR,butitdoescontributetothenoisecom- ponentinthedenominator.Therefore,iff(x)hasany symmetriccomponent,itsSNRwillbeworsethana purelyantisymmetricfilter.A dualargumentholdsfor symmetricedges,sothatiftheedgeG(x)issymmetricor antisymmetric,thefilterf(x)willfollowsuit.Thenetre- sultofthisisthattheresponseHG(x)isalwayssymmet-
There may be some additional constraints on the solution, such as the multiple response constraint (12) described next.
B.EliminatingMultipleResponses
In our specification of the edge detection problem, we decidedthatedgeswouldbemarkedatlocalmaximain theresponseofalinearfilterappliedtotheimage.The detectioncriteriongiveninthelastsectionmeasuresthe effectivenessofthefilterindiscriminatingbetweensignal andnoiseatthecenterofanedge.Itdoesnottakeinto accountthebehaviorofthefilternearbytheedgecenter. Thefirsttwocriteriacanbetriviallymaximizedasfol-
E[X2]
(3)
where 6xo is an approximation to the standard deviation
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
L G'(-x)f'(x)dx
2 X2
Equations (3) and (9) are mathematical forms for the firsttwocriteria,andthedesignproblemreducestothe maximization of both of these simultaneously. In order to dothis,wemaximizetheproductof(3)and(9).Wecould conceivablyhavecombined(3)and(9)usinganyfunction thatismonotonicintwoarguments,buttheuseofthe product simplifies the analysis for step edges, as should becomeclearinSectionI.Forthepresentwewillmake use of the product of the criteria for arbitrary edges, i.e., we seektomaximize
+w +W –
(10)
-W2 (8)
682 IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8, NO. 6,NOVEMBER 1986
lows. From the Schwarz inequality for integrals we can showthatSNR (3)isboundedaboveby
n-I |t; w G2(x)dx andlocalization(9)by
I w
no1 81|G’2(x)dCX.
Both bounds are attained, and the product of SNR and localizationismaximizedwhenf(x)= G(-x)in[-
W].
Thus, according to the first two criteria, the optimal
detectorforstepedgesisatruncatedstep,ordifference ofboxesoperator.Thedifferenceofboxeswasusedby RosenfeldandThurston[25],andinconjunctionwithlat- eralinhibitionbyHerskovitsandBinford[11].However ithas a very high bandwidth and tends to exhibit many maxima in its response to noisy step edges, which is a seriousproblemwhentheimagingsystemaddsnoiseor whentheimageitselfcontainstexturedregions.Theseex- traedgesshouldbeconsiderederroneousaccordingtothe firstofourcriteria.However,theanalyticformofthis criterion was derived from the response at a single point (thecenteroftheedge)anddidnotconsidertheinterac- tionoftheresponsesatseveralnearbypoints.Ifweex- aminetheoutputofadifferenceofboxesedgedetector wefindthattheresponsetoanoisystepisaroughlytri- angularpeakwithnumeroussharpmaximainthevicinity of the edge (see Fig. 1).
Thesemaximaaresoclosetogetherthatitisnotpos- sibletoselectoneastheresponsetothestepwhileiden- tifyingtheothersasnoise.Weneedtoaddtoourcriteria the requirement that the function f will not have “too many” responses to a single step edge in the vicinity of the step. We need to limit the number of peaks in the response so that there will be a low probability of declar- ingmorethanoneedge.Ideally,wewouldliketomake thedistancebetweenpeaksinthenoiseresponseapprox- imatethewidthoftheresponseoftheoperatortoasingle step.Thiswidthwillbesomefractionoftheoperator width W.
R(O)= g2(x)dx and R”(0)= – g’2(x)dx -00 b zr 00
themean distancebetweenzero-crossingsoff’willbe +OD~ )1/2
Inordertoexpressthisasafunctionalconstraintonf, straightforward.Inparticular,ifthefunctionfisrepre- weneedtoobtainanexpressionforthedistancebetween sentedbyadiscretetimesequence,evaluationof(10) adjacentnoisepeaks.Wefirstnotethatthemeandistance requiresonlythecomputationoffourinnerproducts betweenadjacentmaximaintheoutputistwicethedis- betweensequences.Thissuggeststhatnumericaloptimi- tancebetweenadjacentzero-crossingsinthederivativeof zationcanbedonedirectlyonthesampledoperatorim- theoperatoroutput.Thenwemakeuseofaresultdueto pulseresponse.
Rice [24] that the average distance between zero-cross- The output will not be an analytic form for the operator, ingsoftheresponseofafunctiongtoGaussiannoiseis butanimplementationofadetectorfortheedgeofinter-
( R(O)12
est will require discrete point-spread functions anyway. It isalsopossibletoincludeadditionalconstraintsbyusing a penalty method [15]. In this scheme, the constrained
whereR(i)istheautocorrelationfunctionofg.Inour optimizationisreducedtoone,orpossiblyseveral,un- casewearelookingforthemeanzero-crossingspacing constrainedoptimizations.Foreachconstraintwedefine forthefunctionf’.Now since apenaltyfunctionwhichhasanonzerovaluewhenone
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
r x,,(f)= (+
\ oo
‘2(x)dx
tO ft2(x)dx
(12)
The distance between adjacent maxima in the noise re- sponseoff,denotedXmax,willbetwicexzc.We setthis distance to be some fraction k of the operator width.
Xmax(f)= 2x,,(f)= kW. (13)
Thisisanaturalformfortheconstraintbecausethere- sponseofthefilterwillbeconcentratedinaregionof width2W,andtheexpectednumberofnoisemaximain thisregionisNnwhere
2W 2
N=-m=k (14)
Fixingkfixesthenumberofnoisemaximathatcouldlead toafalseresponse.
We remark here that the intermaximum spacing (12) scaleswiththeoperatorwidth.Thatis,wefirstdefinean operatorf,whichistheresultofstretchingfbyafactor ofw,fw(x)=f(xlw).Thenaftersubstitutinginto(12)we findthattheintermaximumspacingforf,isx,(fj)= wxzc(f).Therefore,ifafunctionfsatisfiesthemultiple response constraint (13) for fixed k, then the function f, willalsosatisfyit,assumingW scaleswithw.Forany fixedk,themultipleresponsecriterionisinvariantwith respecttospatialscalingoff.
I. FINDING OPTIMAL DETECTORS BY NUMERICAL OPTIMIZATION
In general itwill be dificult (or impossible) to find a closedformforthefunctionfwhichmaximizes(10)sub- jecttothemultipleresponseconstraint.EvenwhenGhas aparticularlysimpleform(e.g.,itisastepedge),the formoffmaybecomplicated.However,ifwearegiven a candidate function f, evaluation of (10) and (12) is
Xmax k
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION 1.8
(a)
(b)
8.0U 4 8 1.
683
8’0 te 8 1. lb 28 Zq Z8 3Z 38 4 + 4 5Z 56 60 64 1. 8
-.e 4 77
(a)
03.0…ncn
(b)
18
-m-
oftheconstraintsisviolated.Wethenfindthefwhich maximizes
SNR(f)*Localization(f)-L piPi(f) (15)
plementationalreadyincludesdetectorsforthese.The only reason for using a ridge detector is that there are ridgesinimagesthataretoosmalltobedealtwitheffec- tively by the narrowest edge operator. These occur fre- quentlybecausetherearemany edges(e.g.,scratchesand cracksorprintedmatter)whichlieatorbeyondtheres- olutionofthecameraandresultincontoursonlyoneor twopixelswide.
duced, until the constraints are “almost” satisfied.
An example of the method applied to the problem of detecting”ridge”profilesisshowninFig.2.Foraridge,
thefunctionGisdefinedtobeaflatplateauofwidthw, withsteptransitionstozeroattheends.Theauxiliary constraintsare
* Themultipleresponseconstraint.Thisconstraintis takendirectlyfrom(12),anddoesnotdependontheform oftheedge.
* The operator should have zero dc component. That isitshouldhavezerooutputtoconstantinput.
Sincethewidthoftheoperatorisdependentonthe widthoftheridge,thereisasuggestionthatseveralwidths ofoperatorsshouldbeused.Thishasnotbeendonein thepresentimplementationhowever.Awideridgecanbe consideredtobetwocloselyspacededges,andtheim-
one for zero response to constant input.
A roofedgedetectorhasnotbeenincorporatedintothe
implementationoftheedgedetectorbecauseitwasfound thatidealroofedgeswererelativelyrare.Inanycasethe ridgedetectorisanapproximationtotheidealroofdetec- tor, and isadequate to cope with roofs. The situation may bedifferentinthecaseofanedgedetectordesignedex- plicitlytodealwithimagesofpolyhedra,liketheBin- ford-Horn line-finder [14].
Themethodjustdescribedhasbeenusedtofindoptimal operators for both ridge and roof profiles and in addition itsuccessfullyfindstheoptimalstepedgeoperatorde- rivedinSectionIV.Itshouldbepossibletouseittofind operatorsforarbitraryone-dimensionaledges,andit shouldbepossibletoapplythemethodintwodimensions tofindoptimaldetectorsforvarioustypesofcorner.
16
28/ 24 28 32 3b .¡ì
Fig.2. A ridgeprofileandtheoptimaloperatorforit.
85
Fig. 3. A roof profile and an optimal operator for roofs.
wherePiisafunctionwhichhasapositivevalueonly
whenaconstraintisviolated.Thelargerthevalueof,t
themorenearlytheconstraintswillbesatisfied,butatthe
same time the greater the likelihood that the problem will
beill-conditioned.Asequenceofvaluesof,uimayneed
tobeused,withthefinalformofffromeachoptimization concavejunctionsoftwoplanarfacesinpolyhedralob- usedasthestartingformforthenext.The1uiareincreased jects.TheresultsareshowninFig.3.Againthereare ateachiterationsothatthevalueofPi(f)willbere- twosubsidiaryconstraints,oneformultipleresponsesand
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
A similar procedure was used to find an optimal oper- atorforroofedges.Theseedgestypicallyoccuratthe
684 68 iFlEITRANSACTIONSON PAT1TRNANNAIYSISANDM1ACHINEINTELIGENCE,VOt.PAVMI-8N().6.NOVEMBER1986
IV.A D[ETECIORFORSTEPEDGES stepedgedetector.Throughspatialscalingoffwecan Wenowspecializetheresultsofthelastsectiontothe tradeoffdetectionperformanceagainstlocalization,but
casewheretheinputG(x)isstepedge.Specificallyweset G(x) Au (x)whereit,(Y)isthenthderivativeofadelta function, and A is the amplitude of the sterp That is,
It(X)-(0 fox< 0; (16)
(A, fro.-x-: 0 i
and substituting for G(x) in (3) and (9) gives A f(x)dx
wecannotimprovebothsimultaneously.Thissuggests thatanaturalchoiceforthecompositecriterionwouldbe the product of (19) and (20). since this product would be invariant under changes in scale.
0
1X f(x)
2+(w Lx +fw (22)
SNR
Localization
r- -W? F(+WK
no,01 2(X)dX
(17) The solutions to the maximization of this expression willbeaclassoffunctionsalrelatedby spatialscaling. In fact this result is independent of the method of com- bination of the criteria. To see this we assume that there
18 isafunctionfwhichgivesthebestlocalizationAfora
(18) narticularE.Thatiswefindfsuchthat F"ILI%IL41"I Ad. -tilLtL lo, vv%., illl%-Lj OL4%.,Il LII"L
2(f) clandA(f')ismaximized. (23)
and define two performance measures pendon thefilteronly:
A SNR - -2(f)
A A(fJ)
+W. ff(x)dx
Localization -
(19)
~f,2(x)
no
Aif'(O)i
no
Bothofthesecriteriaimprovedirectlywiththeratio theimage.We now remove thisdependenceon theimage
isfixedtoadifferentvalue,i.e.,
E(fv)= C2 while A(f,)ismaximized. (24)
Ifwenowdefinef1(x)intermsoffi(x)asf1(x)= fJ,(xw) where
}S-c2 lc
thentheconstraintonftin(24)translatestoaconstraint on f, which is identical to (23), and (24) can be rewritten as
E(f1)= c1 and w A(f1)ismaximized (25)
whichhasthesolutionf -f Soifwefindasinglesuch functionf,we canobtainmaximallocalizationforany fixed signal-to-noise ratio by scaling f. The design prob-
-4 ilz l W.
1I f'(x)di
Nowsupposeweseekasecondfunctionf,whichgives A/nowhichmightbetermedthesignal-to-noiseratioof thebestpossiblelocalizationwhileitssignal-to-noiseratio
and A which de- \10
I2(f)-
dX
(20) lemforstepedgedetectionhasasingleuniquie(uptospa- Supposenow thatwe forma spatiallyscaledfilterf, tialscaling)solutionregardlessoftheabsolutevaluesof
fromf,wherefj(x)-f(/w).RecallfromtheendofSec- signaltonoiseratioorlocalization. butwe
tion11 thatthemultipleresponsecriterionisunaffectedby Theoptimalfilterisimplicitlydefinedby(22),
spatialscaling.Whenwe substituteft into(19)and(20) musttransformtheproblemslightlybeforewecanapply we obtairnfortheperformanceofthescaledfilter: thecalculusofvariations.Specifically,wetransformthe
maximizationof(22)intoaconstrainedminimizationthat I involves only integral functionals. All but one of the in- 2(ff) wE2(f) and A(f't) A(f'). (21) tegralsin(22)aresettoundeterminedconstantvalues.
We thenfindtheextremevalueoftheremainingintegral Thefirstoftheseequationsisquiteintuitive,andim- (sinceitwillcorrespondtoanextremeinthetotalexpres- pliesthatafilterwithabroadimpulseresponsewillhave sion)asafunctionoftheundeterminedconstants.The bettersignal-to-noiseratiothananarrowfilterwhenap- valuesoftheconstantsarethenchosensoastomaximize
w
pliedtoastepedge.Thesecondislessobvious,andit theoriginalexpression,whichisnowafunctiononlyof implies that a narrow filter will give better localization these constants. Given the constants, we can uniquely thanabroadone. Whatissurprisingisthatthechanges specifythefunctionf(x)whichgivesamaximumofthe areinverselyrelated,thatis,bothcriteriaeitherincrease compositecriterion.
ordecreaseby UThereisanuncertaintyprinciplere- Asecondmodificationinvolvesthelimitsoftheinte- latingthedetectionandlocalizationperformanceofthe grals.Thetwointegralsinthedenominatorof(22)have
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
2(f) A(f') -
dx
f ()d
-W -
Wf2wd
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION
685
limits at + W and - W, while the integral in the numer- ator has one limit at 0 and the other at - W. Since the functionfshouldbeantisymmetric,we canusethelatter limits for al integrals. The denominator integrals will havehalfthevalueoverthissubrangethattheywould haveoverthefulrange.Also,thisenablesthevalueof f'(0)tobesetasaboundarycondition,ratherthanex- pressedasanintegraloff".Iftheintegraltobemini- mizedsharesthesamelimitsastheconstraintintegrals, itispossibletoexploittheisoperimetricconstraintcon- dition(see[6,p.216]).Whenthisconditionisfulfiled, the constrained optimization can be reduced to an uncon- strained optimization using Lagrange multipliers for the constraint functionals. The problem of finding the maxi- mumof(22)reducestotheminimizationoftheintegral
Substituting,
T(x,f, f") =f2 + XIf'2 + X2f + X3f (28)
It may be seen from the form of this equation that the choiceofwhichintegralisextremizedandwhicharecon- straintsisarbitrary,thesolutionwillbethesame.Thisis anexampleofthereciprocitythatwasmentionedearlier. Thechoiceofanintegralfromthedenominatorissimply convenientsincethestandardformofvariationalproblem isaminimizationproblem.TheEulerequationthatcor- respondstothefunctional4'is
whereIfdenotesthepartialderivativeof4withrespect inthedenominatoroftheSNRterm,subjecttothecon- tof,etc.Wesubstitutefor4from(28)intheEulerequa-
straintthattheotherintegralsremainconstant.Bythe principleofreciprocity,we couldhavechosentoextrem- izeanyoftheintegralswhilekeepingtheothersconstant, andthesolutionshouldbethesame.
Weseeksomefunctionfchosenfromaspaceofad- missible functions that minimizes the integral
tiongiving:
2f(x) 2Xlf"(x)+ 2X2f.""(X)+ X3= 0- (30)
The solution of this differential equation is the sum of aconstantandasetoffourexponentialsoftheformex where 7yderives from the solution of the corresponding homogeneousdifferentialequation.Nowaymustsatisfy
subjectto
f(x)dx=cl
f"2(x)dX= C3
so
-U1
Thisequationmay haverootsthatarepurelyimaginary, purelyreal,orcomplexdependingonthevaluesofX1and X2.Fromthecompositefunctional4'wecaninferthat 2 ispositive(sincetheintegraloff"2istobeminimized) butitisnotclearwhatthe or
o0
0 f2(x)dx -w
(26)
tf2(X)dX= C2
f'(0)=C4. (27)
2- 2XI'y2+2X2y4 0
2 x X1-4X2
e=2+ __ (31)
Thespace ofadmissiblefunctionsinthiscase willbe thespace ofalcontinuousfunctionsthatsatisfycertain
boundaryconditions,namelythatf(0)= 0 andf(
sign magnitudeofX1should be.TheEulerequationsuppliesanecessaryconditionfor
- =0.Theseboundaryconditionsarenecessarytoensure thattheintegralsevaluatedover finitelimitsaccurately
theexistenceofaminimum,butitisnotasufficientcon- dition.Byformulatingsuchaconditionwecanresolve theambiguityinthevalueofXl.Todothiswemustcon- sider
representtheinfiniteconvolutionintegrals.Thatis,ifthe nthderivativeoffappears insome integral,thefunction
mustbecontinuousinits(n
-
range (- 0o, + oo). This implies that the values off and
itsfirst(n -
integration,sincetheyare zero outsidethisrange.
l)st derivative
over the
the second variation of the functional. Let xo
1)derivativesmustbezero atthelimitsof
The functional to be minimized isofthe form ibF(x,f, f',f")andwe havea seriesofconstraintsthatcan be
Then by Taylor's theorem (see also [6, p. 214]),
J[f + eg] = J[f] + EJ1[f,g] + 26'2[f+ Pg,g]
writtenintheformXbGi(X,f,f',f") ci. =
Sincethecon- straintsare isoperimetric,i.e.,theysharethesame limits
W)
ofintegrationas theintegralbeingminimized,we can wherepissomenumberbetween0andc,andgischosen formacompositefunctionalusingLagrangemultipliers fromthespaceofadmissiblefunctions,andwhere [6].Thefunctionalisa linearcombinationofthefunc-
tionalsthatappear intheexpressiontobeminimizedand
intheconstraints.Findinga solutiontotheunconstrained maximizationof*(x,f,f',f")isequivalentto findingthe solution to the constrained problem. The composite func- tionalis
*(x,f,f',f")= F(x,f,f',f")+XIG1(x,f,f',f") + X2G2(x,f,f',f")+ ...
J1[f,g]= , 'fg+4'1g'+*f"gdx xi
J21f,g] = X 'fg2 + Tffg,2+ 'ff"g,2 xo
+ 24Nfg'+ 2*f,fg'g"+ 2'f"gfdx. (32)
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
-
*f- d
* (29)
d2
dx 4's + dX2 f=10
2X2- 2X2
686 I6E TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INT7ELLIGENCE, VOL. PAMI-8. NO. 6. NOVEMBER 1986
NotethatJ,isnothingmorethantheintegralofgtimes the Euler equation forf (transformed using integration by parts)andwillbezeroiffsatisfiestheEulerequation.We cannowdefinethesecondvariation62Jas
62J 6JIfg]-
Thenecessaryconditionforaminimumis62J> 0.We computethesecondpartialderivativesofT from(28)and we get
a. + a4 + c= 0 a,easinw + a,eacosw + a3e tsinX
+ a4e a cosw -+c 0 a,w+a2U+a3W-a4U-S
aIea(axsinw +wcosw)+a2eo(acosw
– w sin w) + a3e- (-a sin w + w cos w)
+ a4e-Uc-(ucosw – w sinc)- 0. (36)
Jil[f g-
2g2+ 2X1g2+ 2X2g’2dx. 0.
rlO
These equations are linear in the four unknowns a1, a2, (3) a3,a4andwhensolvedtheyyield
Usingthefactthatgisanadmissiblefunctionandthere- a,- c(oe(3-a)sin2w-acocos2w+(-2w2sinha
fore vanishes at the integration limits, we transform the above using integration by parts to
g2-X1gg,+X.g,2adx 0
+2a2e-)sinw +2awsinhacosw
+we2o(u+l3)-3)/4(w2sinh2o- a sin2w) a2=c(a(3-a)cos2w+awsin2w-2awcosha
*sinw-2w2sinhacosco+2w2e sinha
+ a(a- 3)/4(w2sinh2a _ ae2sin2w) a3=c(-au(3+ae)sin2w+aw cos2w+(2w2sinha
+ 2a2e’)sinw + 2aw)sinha cosw
+we2a(j- a)-f3w)14(wS2sinh2a-a2sin2w) a4=c(-a($+a)cos2w-ao sin2w+2awcosha
2
which can be written as
2(g2X if +(X2- )dx0.
The integralisguaranteedtobepositiveiftheexpression beingintegratedispositiveforalx,so if
4X,>XI
X) gi2
thentheintegralwillbepositiveforalx and forarbitrary
g,andtheextremum willcertainlybeaminimum. Ifwe
referbackto(31)we findthatthisconditionisprecisely thatwhichgivescomplexrootsfor-y,sowehaveboth guaranteedtheexistenceofaminimumandresolveda possibleambiguityintheformofthesolution.Wecan nowproceedwiththederivationandassumefourcomplex ilarlyfora4froma2.
rootsoftheform+=y+a¡Àiwwithoa,wreal.Nowy2 Thefunctionfisnowparametrizedintermsofthecon- a2 w2+2iawandequatingrealandimaginaryparts stantsa,w,(,andc.Wehavestiltofindthevaluesof
with(31)we obtain
theseparameterswhichmaximizethequotientofintegrals that forms our composite criterion. To do this we first express each of the integrals in terms of the constants. Sincetheseintegralsareverylonganduninteresting,they arenotgivenherebutmaybefoundin[4].We havere- duced the problem of optimizing over an infinite-dimen- sional space of functions to a nonlinear optimization in threevariablesa,w,and3(notsurprisingly,thecom- bined criterion does not depend on c). Unfortunately the resulting criterion, which must stil satisfy the multiple responseconstraint,isprobablytoocomplextobesolved analytically, and numerical methods must be used to pro- videthefinalsolution.
The shape of f will depend on the multiple response constraint,i.e.,itwilldependonhowfarapartweforce f(-x)= -f(x).Thefourboundaryconditionsenableus theadjacentresponses.Fig.5showstheoperatorsthat resultfromparticularchoicesofthisdistance.Recallthat therewasnosinglebestfunctionforarbitraryw,buta class of functions which were obtained by scaling a pro-
f(0)= 0
f(-W)=0 f'(0)
=
s
f'(-W)
=
(35)
0
a2
2>and4aeo= A2
X_X 2 ?2 2 4X2
2I
Thegeneralsolutionintherange [- W,0]may now be written
Ie’xsinwx +a2e”xcoscox+a3e-x sinwx +a4e-“coswx +c.
ft(x)=
This function is subject to the boundary conditions
a
wheresisanunknownconstantequaltotheslopeofthe functionfattheorigin.Sincef(x)isasymmetric,we can extendtheabovedefinitiontotherange[- W,W]using
tosolveforthequantitiesa1 througha4intermsofthe unknownconstantsa, co,c,ands.Theboundarycondi- tionsmay berewritten
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
sinw+2w2sinhacosw- 2w2easinha
+ a(a-_3))/4(w2sinh2a – c2sinwC) (37)
wheref3istheslopesattheorigindividedbytheconstant c.Oninspectionoftheseexpressionswecanseethata3 canbeobtainedfroma1byreplacingaby-a,andsim-
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION ,687
totypefunctionbyw.We willwanttoforcetheresponses
further apart as the signal-to-noise ratio in the image is
lowered, but it is not clear what the value of signal-to-
noise ratio will be for a single operator. In the context in
whichthisoperatorisused,severaloperatorwidthsare
available,andadecisionprocedureisappliedtoselectthe
smallestoperatorthathasanoutputsignal-to-noiseratio
aboveafixedthreshold.Withthisarrangementtheoper-
atorswillspendmuchofthetimeoperatingclosetotheir Fig.4.Filterparametersandperformancemeasuresforthefiltersilus-
outputEthresholds.We trytochooseaspacingforwhich tratedinFig.5. theprobabilityofamultipleresponseerroriscomparable approximatedthisfilterusingthefirstderivativeofa to the probability of an error due to thresholding. Gaussian as described in the next section.
A rough estimate for the probability of a spurious max- The first derivative of Gaussian operator, or even filter
imumintheneighborhoodofthetruemaximumcanbe 6itself,shouldnotbetakenasthefinalwordinedge
formedasfollows.Ifwelookattheresponseofftoan detectionfilters,evenwithrespecttothecriteriawehave idealstepwefindthatitssecondderivativehasmagnitude used.Ifwe arewillingtotolerateaslightreductionin
Af'(0)atx=0.Therewillbeonlyonemaximumnear multipleresponseperformancer,wecanobtainsignifi-
thecenteroftheedgeif Af'(0) isgreaterthanthesec- cantimprovementsintheothertwocriteria.Forexample, ondderivativeoftheresponsetonoiseonly.Thislater filters4and5bothhavesignificantlybetterEAproduct quantity,denoteds,nisaGaussianrandomvariablewith thanfilter6,andonlyslightlylowerr.FromFig.5we
standard deviation can see that these filters have steeper slope at the origin,
n(j ift2(x)dx)
zation,althoughthishasnotbeenverifiedexperimentally.
-w A thoroughempiricalcomparisonoftheseotheroperators TheprobabilityPmthatthenoiseslopeSnexceedsAf'(0) remainstobedone,andthetheoryinthiscaseisunclear
pm to the probability of falsely marking an edge pfwhich is
–
pf=I -)E(39)
G(x)= exp (2 (-N2)
The reason for doing this is that there are very efficient waystocomputethetwo-dimensionalextensionofthefil- ter ifitcan be represented as some derivative of a Gauss-
bysettingpm=pf.Thisisanaturalchoicesinceitmakes ian.Thisisdescribedindetailelsewhere[4],butforthe adetectionerrororamultipleresponseerrorequally presentwewillcomparethetheoreticalperformanceofa likely. Then from (38) and (39) we have firstderivative of a Gaussian filterto the optimal operator.
If'(¡ã) = E. ( 40)
In practice it was impossible to find filters which sat-
isfiedthisconstraint,soinsteadwesearchforafiltersat- andthetermsintheperformancecriteriahavethevalues isfying
If'(0)I= rE (41) If'(O)l= OS
OS 0 +0
f(x)dx= 1
whererisascloseaspossibleto1.Theperformancein- dexesandparametervaluesforseveralfiltersaregivenin Fig.4.Theaicoefficientsforalthesefilterscanbefound
from(37),byfixingcto,say,c= 1.Unfortunately,the
largest value of r that could be obtained using the con- strainednumericaloptimizationwasabout0.576forfilter number6inthetable.Inourimplementation,wehave Theoverallperformanceindexforthisoperatoris
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
FiltcrParameters nx,zE1Ar a w=_
1 0.15 4.21 0.215 21.59550 0.12250 63.97566 2 0.3 2.87 0.313 12.47120 0.38284 31.26860
3
4 0.8 1.57 0.515 5.06500 2.56770 11.06100 5 1.0 1.33 0.561 3.455800.07161 4.80684
,+W \1/2 suggestingthattheperformancegainismostlyinlocali-
noas =
isgiveninterms ofthenormaldistributionfunction4) on how besttomake thetradeoff.
PM 1- (A1f'(O)I) (38) V.ANEFFICIENTAPPROXIMATION noas
We canchooseavalueforthisprobabilityasanac- Theoperatorderivedinthelastsectionasfilternumber
ceptableerrorrateandthiswilldeterminetheratiooff'(0) 6,andillustratedinFig.6,canbeapproximatedbythe toas.Wecanrelatetheprobabilityofamultipleresponse firstderivativeofaGaussianG'(x),where
0.5 2.13 0.417 7.85869 2.62856 18.28800
1.2 1.12 0.576 2.05220 1.569:39 2.91540 7 141 0.75 0.484 0.00297 3.50350 7.47700
The impulse response of the first derivative filter is
Os x-x\-
-0
_00
f(x)= – ~exp(-2 X~~~2r
(42)
f’2(x)dx –3 4a
ft2(x)dx =
-0 8a5
f2(x)dx=VI -00 2a
(43)
688688 ~~~~IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI 8, NO. 6, NOVEMBER 1986 I
(a)
(b)
8 ze 40 60 ae
.6 z Z.a z. 2.8 Z.o
320 380
-1.3141194 1,28S2 13
a la le
1.1515[57
355 0.6200538
99 ze -0.62e 37
TI, Z2. 3″
zle Z’!e Z.q ZW 3″ 3ze
2
3
4 qoo
5
6
7
as
ze qo so 80
149
Z26 Z49 Zee
3e 3qo Z” 3zv 3qe
380
350 3ae
while the r value is, from (41),
r=
–zz0.5 1 15I
60 as IN
220
Z45 Z”
EA –0.92 3w-
The performnanceof the first derivative of Gaussian op- (44) eratoraboveisworse thantheoptimaloperatorbyabout 20 percentand itsmultipleresponse measure r, isworse by about 10 percent. It would probably be difficult to de- tect a difference of this magnitude by looking at the per- formanceofthetwooperatorsonrealimages,andbe- cause thefirstderivativeofGaussianoperatorcan be computedwithmuchlesseffortintwodimensions,ithas
Fig.5.Optimalstepedgeoperatorsforvariousvaluesofx19..Fromtop tobottom,theyarex,,,,=0.15,0.3,0.5.0.8.1.0.1.2.1.4.
Fig. 6. (a) The optimal step edge operator. (b) The first derivative of a Gaussian.
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION 689
beenusedexclusivelyinexperiments.Theimpulsere- thendenotetheautocorrelationfunctionofg,asRI(r) sponsesofthetwooperatorscanbecomparedinFig.6. andthatofg2asR22(T),andtheircross-correlationas
AcloseapproximationofthefirstderivativeofGauss- ianoperatorwassuggestedbyMacleod[16]forstepedge detection. Macleod’s operator is a difference of two dis- placed two-dimensional Gaussians. It was evaluated in Fram and Deutsch [7] and compared very favorably with several other schemes considered in that paper. There are alsostronglinkswiththeLaplacianofGaussianoperator suggestedbyMarrandHildreth[18].Infact,aone-di- mensional Marr-Hildreth edge detector is almost identi- cal with the operator we have derived because maxima in theoutputofafirstderivativeoperatorwillcorrespondto zero-crossingsintheLaplacianoperatorasusedbyMarr and Hildreth. In two dimensions however, the directional propertiesofourdetectorenhanceitsdetectionandlocal- izationperformancecomparedtotheLaplacian.Another importantdifferenceisthattheamplitudeoftheresponse at a maximum provides a good estimate of edge strength, becausetheSNR criterionistheratioofthisresponseto the noise response. The Marr-Hildreth operator does not useanyformofthresholding,butanadaptivethreshold- ing scheme can be used to advantage with our firstderiv- ative operator. In the next section we describe such a scheme, which includes noise estimation and a novel method for thresholding edge points along contours.
R12(T),wherethecorrelationoftworealfunctionsisde- finedasfollows:
r+
Rij(T)= gi(x)g1(x+r)dx.
We have derived our optimal operator to deal with
knownimagefeaturesinGaussiannoise.Edgedetection
betweentexturedregionsisanotherimportantproblem. GaussianandRI,isthesecondderivativeofaGaussian Thisisstraightforwardifthetexturecanbemodelledas
theresponseofsomefiltert(x)toGaussiannoise.Wecan
then treat the texture as a noise signal, and the response
ofthefilterf(x)tothetextureisthesameastheresponse
ofthefilter(f*t)(x)toGaussiannoise.Makingthis
replacementineachintegralintheperformancecriteria
thatcomputesanoiseresponsegivesusthetextureedge
designproblem.Thegeneralizationtoothertypesoftex-
tureisnotaseasy,andforgooddiscriminationbetween
knowntexturetypes,abetterapproachwouldinvolvea
Markov image model as in [5].
VI. NOISE ESTIMATION AND THRESHOLDING
Toestimatenoisefromanoperatoroutput,weneedto
The filter K above is convolved with the output of the edgedetectionoperatorandtheresultissquared.Thenext stepistheestimationofthemean-squarednoisefromthe localvalues.Herethereareseveralpossibilities.Thesim- plestistoaveragethesquaredvaluesoversomeneigh- borhood,eitherusingamovingaveragefilterorbytaking anaverageovertheentireimage.Unfortunately,experi- encehasshownthatthefilterKisverysensitivetostep edges, and that as a consequence the noise estimate from any form of averaging is heavily colored by the density andstrengthofedges.
beabletoseparateitsresponsetonoisefrom
duetostepedges.Sincetheperformanceofthesystem
willbecriticallydependentontheaccuracyofthisesti-
mate,itshouldalsobeformulatedasanoptimization.
Wienerfilteringisamethodforoptimallyestimatingone
componentofatwo-componentsignal,andcanbeused
toadvantageinthisapplication.Itrequiresknowledgeof
theautocorrelationfunctionsofthetwocomponentsand lessthan80percent)willbedeterminedmainlythenoise ofthecombinedsignal.Oncethenoisecomponenthas energy,andthattheyarethereforeusefulestimatorsfor
theresponse
beenoptimallyseparated,weformaglobalhistogramof noise.Aglobalhistogramestimateisactuallyusedinthe noiseamplitude,andestimatethenoisestrengthfrom currentimplementationofthealgorithm.
some fixed percentile of the noise signal. Even with noise estimation, the edge detector will be Letgl(x)bethesignalwearetryingtodetect(inthis susceptibletostreakingifitusesonlyasinglethreshold. casethenoiseoutput),andg2(x)besomedisturbance Streakingisthebreakingupofanedgecontourcausedby (paradoxicallythiswillbetheedgeresponseofourfilter), the operator output fluctuating above and below the
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
We assumeinthiscasethatthesignalanddisturbance areuncorrelated,soR12(T)= 0.TheoptimalfilterisK(x), whichisimplicitlydefinedasfollows[30]:
r+ R1(T)=J(R1I(T-x)+R2(T- x)K(x)dx.
Since the autocorrelation of the output of a filter in re- sponsetowhitenoiseisequaltotheautocorrelationofits impulseresponse,wehave
RI1(x)= k__- 1)exp(-4$2)
Ifg2istheresponseoftheoperatorderivedin(42)toa stepedgethenwewillhaveg2(x)=kexp(-x12o2)and
R22(x)= k2exp 2
In the case where the amplitude of the edge is large comparedtothenoise,R22+ RI,isapproximatelya
ofthesamea.ThentheoptimalformofKisthesecond derivativeofanimpulsefunction.
In order to gain better separation between signal and noisewecanmakeuseofthefactthattheamplitudedis- tributionofthefilterresponsetendstobedifferentfor edgesandnoise.Byourmodel,thenoiseresponseshould haveaGaussiandistribution,whilethestepedgeresponse willbecomposedoflargevaluesoccurringveryinfre- quently.Ifwetakeahistogramofthefiltervalues,we shouldfindthatthepositionsofthelowpercentiles(say
690 IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8, NO. 6, NOVEMBER 1986
(a) (c)
(b) (d)
Fig. 7. (a) Parts image, 576 by 454 pixels. (b) Image thesholded at T,. (c) Image thresholded at 2 T,. (d) Image thresholded with hysteresis using both the thresholds in (a) and (b).
thresholdalongthelengthofthecontour. Supposewe Inthecurrentalgorithm,no attemptismadetopreseg- haveasinglethresholdsetatT1,andthatthereisanedge mentcontours.Insteadthethresholdingisdonewithhys- intheimagesuchthattheresponse oftheoperatorhas teresis.Ifanypartofacontourisaboveahighthreshold, meanvalueT1.Therewillbesomefluctuationoftheout- thosepointsareimmediatelyoutput,asistheentirecon- putamplitudeduetonoise,evenifthenoiseisveryslight. nectedsegmentofcontourwhichcontainsthepointsand Weexpectthecontourtobeabovethresholdonlyabout whichliesabovea lowthreshold.Theprobabilityof halfthetime.Thisleadstoabrokenedgecontour.While streakingisgreatlyreducedbecauseforacontourtobe thisisa pathologicalcase, streakingisa very common brokenitmustnow fluctuateabovethehighthresholdand problemwithedgedetectorsthatemploythresholding.It belowthelowthreshold.Alsotheprobabilityofisolated isvery difficulttoseta thresholdso thatthereissmall falseedgepointsisreducedbecausethestrengthofsuch
pointsmustbeaboveahigherthreshold.Theratioofthe hightolowthresholdintheimplementationisintherange two or three to one.
bybreakingitatmaximaincurvature.Thissegmentation isnecessary inthecase ofzero-crossingssincethezero- crossingsalwaysformclosedcontours,whichobviously do not always correspond to contours in the image.
probabilityofmarkingnoiseedgeswhileretaininghigh sensitivity.Anexampleoftheeffectofstreakingisgiven in Fig. 7.
Onepossiblesolutiontothisproblem,usedbyPentland
[22]withMarr-Hildrethzero-crossings,istoaverage the
edgestrengthofa contourover partofitslength.Ifthe
averageisabovethethreshold,theentiresegmentis stepedgeinspacewithonepositioncoordinate.Intwo marked.Iftheaverage isbelowthreshold,no partofthe dimensionsanedgealsohasanorientation.Inthissection contourappearsintheoutput.Thecontourissegmented
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
VII. Two OR MORE DIMENSIONS
Inone dimensionwe can characterizethepositionofa
wewillusetheterm”edgedirection”tomeanthedirec- tionofthetangenttothecontourthattheedgedefinesin twodimensions.Supposewe wishtodetectedgesofa particularorientation.Wecreateatwo-dimensionalmask forthisorientationbyconvolvinga linearedgedetection
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION 691
functionalignednormaltotheedgedirectionwithapro- thedirectionaloperatorG,butwe neednotknowthe jectionfunctionparalleltotheedgedirection.Asubstan- directionnbeforeconvolution.
tialsavingsincomputational efortispossibleifthepro- Theformofnonlinearsecondderivativeoperatorin(47) jectionfunctionisaGaussianwiththesameaasthe(first hasalsobeenusedbyTorreandPoggio[28]andbyHar- derivativeofthe)Gaussianusedasthedetectionfunction. alick[10].ItalsoappearsinPrewitt[23]inthecontext Itispossibletocreatesuchmasksbyconvolvingtheim- ofedgeenhancement.Aratherdifferenttwo-dimensional agewithasymmetrictwo-dimensionalGaussianandthen extensionisproposedbySpacek[26]whousesone-di- differentiatingnormaltotheedgedirection.Infactwedo mensionalfiltersalignednormaltotheedgedirectionbut
nothavetodothisineverydirectionbecausetheslopeof withoutextendingthemalongtheedgedirection.Spacek asmoothsurfaceinanydirectioncanbedeterminedex- startswithaone-dimensionalformulationwhichmaxi- actlyfromitsslopeintwodirections.Thisformofdirec- mizestheproductofthethreeperformancecriteriade- tionaloperator,whilesimpleandinexpensivetocompute, finedinSectionI,andleadstoastepedgeoperatorwhich formstheheartofthemoreelaboratedetectorwhichwill differsslightlyfromtheone we derivedinSectionIV. be described in the next few sections. Gennert [8] addresses the two-dimensional edge detector
Supposewewishtoconvolvetheimagewithanoper- problemdirectly,andappliesasetofdirectionalfirstde- atorGnwhichisthefirstderivativeofatwo-dimensional rivativeoperatorsateachpointintheimage.Theopera-
Gaussian G in some direction n, i.e., G=exp X2+Y2)
and
tors have limited extent along the edge direction and pro- ducegoodresultsatsharpchangesinedgeorientationand corners.
The operator (47) actually locates either maxima or minima by locating the zero-crossings in the second de- rivativeintheedgedirection.Inprincipleitcouldbeused
aG
G =-= n*VG.
(45) toimplementanedgedetectorinanarbitrarynumberof dimensions, by first convolving the image with a sym- Ideally,nshouldbeorientednormaltothedirectionof metricn-dimensionalGaussian.Theconvolutionwithan
anedgetobedetected,andalthoughthisdirectionisnot n-dimensionalGaussianishighlyefficientbecausethe knownapriori,wecanformagoodestimateofitfrom Gaussianisseparableintonone-dimensionalfilters.
the smoothed gradient direction But there are other more pressing reasons for using a smooth projection function such as a Gaussian. When we apply a linear operator to a two-dimensional image, we n IV(G*)I (46) formateverypointintheoutputaweightedsumofsome
where*denotesconvolution.Thisturnsouttobeavery oftheinputvalues.Fortheedgedetectordescribedhere,
goodestimatorforedgenormaldirectionforsteps,since thissum willbeadifferencebetweenlocalaverageson asmoothedstephasstronggradientnormaltotheedge. differentsidesoftheedge.Thisoutput,beforenonmaxi-
an
Itisexactforstraightlineedgesintheabsenceofnoise, mumsuppression,representsakindofmovingaverageof
andtheGaussiansmoothingkeepsitrelativelyinsensitive theimage.Ideallywe wouldliketouseaninfinitepro- to noise. jection function, but real edges are of limited extent. It is
Anedgepointisdefinedtobealocalmaximum(inthe thereforenecessarytowindowtheprojectionfunction[9]. directionn)oftheoperatorGnappliedtotheimageI.At Ifthewindowfunctionisabruptlytruncated,e.g.,ifitis
alocalmaximum, we have
a-G,*i= O an
and substituting for Gn from (45) and associating Gauss- ian convolution, the above becomes
rectangular,thefilteredimagewillnotbesmoothbecause of the very high bandwidth of this window. This effect is relatedtotheGibbsphenomenoninFouriertheorywhich occurs when a signal is transformed over a finite window. When nonmaximum suppressionisappliedtothisrough signal we find that edge contours tend to “wander” or that in severe cases they are not even continuous.
a2 an2G*I= 0.
The solution is to use a smooth window function. In (47) statistics,theHammingandHanningwindowsaretypi-
At such an edge point, the edge strength will be the mag- nitudeof
cally used for moving averages. The Gaussian is a rea- sonable approximation to both of these, and it certainly has very low bandwidth for a given spatial width. (The
IGn*I= IV(G*I)I. (48) Gaussianistheuniquefunctionwithminimalproductof
Becauseoftheassociativityofconvolution,wecanfirst bandwidthandfrequency.)Theeffectofthewindow convolvewithasymmetricGaussianGandthencompute functionbecomesverymarkedforlargeoperatorsizesand
directionalsecondderivativezerostolocateedges(47), itisprobablythebiggestsinglereasonwhyoperatorswith andusethemagnitudeof(48)toestimateedgestrength. largesupportwerenotpracticaluntiltheworkofMarr Thisisequivalenttodetectingandlocatingtheedgeusing andHildrethontheLaplacianofGaussian.
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
It is worthwhile here to compare the performance of
692 IEEE I’RANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8, NO. 6.NOVEMBER 1986
thiskindofdirectionalsecondderivativeoperatorwith theLaplacian.Firstwenotethatthetwo-dimensionalLa- placiancanbedecomposedintocomponentsofsecond derivative in two arbitrary orthogonal directions. If we choosetotakeoneofthederivativesinthedirectionof principalgradient.wefindthattheoperatoroutputwill containonecontributionthatisessentiallythesameasthe operatordescribedabove,andalsoacontributionthatis alignedalongtheedgedirection.Thissecondcomponent contributesnothingtolocalizationordetection(thesur- faceisroughlyconstantinthisdirection),butincreases theoutputnoise.
signal-to-noiseratiothresholdleadstoalowererrorrate, butwilltendtogivepoorerlocalizationsincefeweredges willberecordedfromthesmalleroperators.
In summary then, the first heuristic for choosing be- tweenoperatoroutputsisthatsmalloperatorwidths shouldbeusedwhenevertheyhavesuficientE.Thisis similartotheselectioncriterionproposedbyMarrand Hildreth[18]forchoosingbetweendifferentLaplacianof Gaussianchannels.Intheircasetheargumentwasbased ontheobservationthatthesmallerchannelshavehigher resolution,i.e.,thereislesspossibilityofinterference from neighboring edges. That argument is also very rel- evantinthepresentcontext,astodatetherehasbeenno considerationofthepossibilityofmorethanoneedgein agivenoperatorsupport.Interestingly,Rosenfeldand Thurston[25]proposedexactlytheoppositecriterionin thechoiceofoperatorforedgedetectionintexture.The
Inlatersectionswewilldescribeanedgedetectorwhich
incorporatesoperatorsofvaryingorientationandaspect
ratio,buttheseareasupersetoftheoperatorsusedinthe
simpledetectordescribedabove.Intypicalimages,most
oftheedgesaremarkedbytheoperatorsofthesmallest
width,andmostofthesebynonelongatedoperators.The argumentgivenwasthatthelargeroperatorsgivebetter simpledetectorperformswellenoughinthesecases,and averagingandtherefore(presumably)bettersignal-to- asdetectorcomplexityincreases,performancegainstend noiseratios.
to diminish. However, as we shall see in the following Taking the fine-to-coarse heuristic as a starting point, sections,therearecaseswhenlargerormoredirectional weneedtoformalocaldecisionprocedurethatwillen- operatorsshouldbeused,andthattheydoimproveper- ableustodecidewhethertomarkoneormoreedgeswhen formancewhentheyareapplicable.Thekeytomaking severaloperatorsinaneighborhoodareresponding.Ifthe suchacomplicateddetectorproduceacoherentoutputis operatorwiththesmallestwidthrespondstoanedgeand todesigneffectivedecisionproceduresforchoosingbe- ifithasasignal-to-noiseratioabovethethreshold,we tween operator outputs at each point in the image. shouldimmediatelymarkanedgeatthatpoint.We now
facetheproblemthattherewillalmostcertainlybeedges
marked by the larger operators, but that these edges will Havingdeterminedtheoptimalshapefortheoperator, probablynotbeexactlycoincidentwiththefirstedge.A
we now facetheproblemofchoosingthewidthofthe possibleanswertothiswouldbetosuppresstheoutputs operatorsoastogivethebestdetection/localization ofalnearbyoperators.Thishastheundesirableeffectof tradeoffinaparticularapplication.Ingeneralthesignal- preventingthelargechannelsforrespondingto”fuzzy” to-noiseratiowillbedifferentforeachedgewithinan edgesthataresuperimposedonthesharpedge. image,andsoitwillbenecessarytoincorporateseveral Insteadweusea”featuresynthesis”approach.Webe- widthsofoperatorinthescheme.Thedecisionastowhich ginbymarkingaltheedgesfromthesmallestoperators. operatortousemustbemadedynamicallybythealgo- Fromtheseedges,wesynthesizethelargeoperatorout- rithmandthisrequiresalocalestimateofthenoiseenergy putsthanwouldhavebeenproducediftheseweretheonly
edgesintheimage.Wethencomparetheactualoperator outputstothesyntheticoutputs.Wemarkadditionaledges onlyifthelargeoperatorhassignificantlygreaterre- sponsethanwhatwewouldpredictfromthesyntheticout- put.Thesimplestwaytoproducethesyntheticoutputsis
VIII. THE NEED FOR MULTIPLE WIDTHS
intheregionsurroundingthecandidateedge.Oncethe
noiseenergyisknown,thesignal-to-noiseratiosofeach oftheoperatorswillbeknown.Ifwethenuseamodelof theprobabilitydistributionofthenoise,wecaneffec-
tivelycalculatetheprobabilityofacandidateedgebeing
afalseedge(foragivenedge,thisprobabilitywillbe totaketheedgesmarkedbyasmalloperatorinapartic-
differentfordifferentoperatorwidths). Ifweassumethattheaprioripenaltyassociatedwith
ular direction, and convolve with a Gaussian normal to theedgedirectionforthisoperator.TheaofthisGaussian shouldbethesameastheaofthelargechanneldetection filter.
afalselydetectededgeisindependentoftheedgestrength, itisappropriatetothresholdthedetectoroutputsonprob-
ability of error rather than on magnitude of response. Once
theprobabilitythresholdisset,theminimumacceptable
signal-to-noiseratioisdetermined.However,theremay markedbyatthefirst,andthentofindtheedgesfromthe
beseveraloperatorswithsignal-to-noiseratiosabovethe threshold,andinthiscasethesmallestoperatorshouldbe chosen,sinceitgivesthebestlocalization.Wecanafford tobeconservativeinthesettingofthethresholdsince edgesmissedbythesmallestoperatorsmaybepickedup bythelargerones.Effectivelytheglobaltradeoffbetween errorrateandlocalizationremains,sincechoosingahigh
thirdscalethatwerenotmarkedbyeitherofthefirsttwo, etc.Thuswebuildupacumulativeedgemapbyadding thoseedgesateachscalethatwerenotmarkedbysmaller scales.Itturnsoutthatinmanycasesthemajorityofedges arepickedupbythesmallestchannel,andthelaterchan- nelsmarkmostlyshadowandshadingedges,oredgesbe- tweentexturedregions.
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
This procedure can be applied repeatedly to first mark theedgesfromthesecondsmallestscalethatwere not
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION
693
(a)
(c)
(b) (d)
Fig.8.(a)Edgesfrompartsimageata = 1.0. (b)Edgesata = 2.0.(c) Superposition of the edges. (d) Edges combined using feature synthesis.
Someexamplesoffeaturesynthesisappliedtosome detectionfunction.Infactboththedetectionandlocali- sampleimagesare showninFigs.8and9.Noticethat zationoftheoperatorimproveasthelengthoftheprojec- mostoftheedgesinFig.8aremarkedbythesmaller tionfunctionincreases.Wenowprovethisfortheoper- scaleoperator,andonlyafewadditionaledges,mostly atorsignal-to-noiseratio.Theproofforlocalizationis shadows,arepickedupbythecoarserscale.However similar.Wewillconsiderastepedgeinthexdirection whenthetwosetsofedgesaresuperimposed,we notice whichpassesthroughtheorigin.Thisedgecanberepre- thatinmany casestheresponsesofthetwooperatorsto sentedbytheequation
thesame edgearenotspatiallycoincident.Whenfeature
synthesisisappliedwe findthatredundantresponses of
thelargeroperatorareeliminatedleadingtoasharpedge whereu_1istheunitstepfunction,andAistheamplitude map. of the edge as before. Suppose that there is additive
Bycontrast,inFig.9theedgesmarkedbythetwoop- Gaussiannoiseofmean squaredvaluen2 perunitarea. eratorsare essentiallyindependent,anddirectsuperposi- Ifwe convolvethissignalwitha filterwhoseimplusere- tionoftheedgesgivesausefuledgemap. Whenwe apply sponse isf(x,y),thentheresponse totheedge(atthe featuresynthesistothesesetsofedgeswe findthatmost origin)is
oftheedgesatthecoarser scaleremain.BothFigs.8and 9were producedbytheedgedetectorwithexactlythe same setofparameters(otherthanoperatorsize),andthey were chosen to represent opposing extremes of image content across scale.
J c
f(x,y)dxdy.
f2hqhY)sdXdyw aGaussianwiththesame a astheGaussianusedforthe Thesignal-to-noiseratioisthequotientofthesetwo
IX. THE NEED FOR DIRECTIONAL OPERATORS
noo00
Sofarwe haveassumedthattheprojectionfunctionis
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
I(x,y)= Au i(y)
Therootmean squaredresponse tothenoiseonlyis +oo fr+ 1/2
694
IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, VOL. PAMI-8, NO. 6,NOVEMBER 1986
(a)
(c)
(b)
(d)
(e)
Fig. 9. (a) Handywipe image 576 by 454 pixels. (b) Edges from handy- wipeimageata = 1.0.(c)a = 5.0.(d)Superpositionoftheedges.(e) Edgescombinedusingfeaturesynthesis.
Authorized licensed use limited to: Griffith University. Downloaded on August 04,2020 at 03:22:24 UTC from IEEE Xplore. Restrictions apply.
CANNY: COMPUTATIONAL APPROACH TO EDGE DETECTION 695
(a)
al,iothionJSAR-lLr.TISi0-