Linear Feature Detection Method and Apparatus

FIELD OF THE INVENTION

p-0002The present invention relates to methods of automated image analysis and in particular to computer methods for feature detection in an image.

p-0003The invention has been developed primarily for use as a method of linear feature detection for image analysis in computer vision and pattern recognition applications and will be described hereinafter with reference to this application. However it will be appreciated that the invention may not limited to this particular field of use. Neurite outgrowth detection is one of the application examples.

BACKGROUND OF THE INVENTION

p-0004Any discussion of the prior art throughout the specification should in no way be considered as an admission that such prior art is widely known or forms part of the common general knowledge in the field. In particular the references cited throughout the in specification should in no way be considered as an admission that such art is prior art, widely known or forms part of the common general knowledge in the field.

p-0005Linear or elongated feature detection is a very important task in the areas of image analysis, computer vision, and pattern recognition. It has a very wide range of applications ranging from retinal vessel extraction, skin hair removal for melanoma detection, and fingerprint analysis in the medical and biometrics area, neurite outgrowth detection and compartment assay analysis in the biotech area as described in Ramm et al., Van de Wouwer et al. and Meijering et al. It is also used for characterizing, tree branches, tree bark, plant roots, and leaf vein/skeleton detection; and in the infrastructure areas for road crack detection, roads and valleys detection in satellite images as described in Fischler et al.

p-0006There are a number of techniques in the literature for linear feature detection. Quite a few are aimed at detecting retinal vessels. A recent review of some of the available vessel extraction techniques and algorithms can be found in Kirbas and Quek

p-0007One type of extraction technique, such as that disclosed in Bamberger and Smith, requires a series of directional filters corresponding to the direction of the structures present in the image. Some of the methods employed there include steerable filters as described in Gauch and Pier, 2D matched filters as described in Chaudhuri et al., maximum gradient profiles as described in Colchester et al., and directional morphological filtering such as that used by Soille and Talbot. This type of techniques can be termed as either template or model-based and tend to be slow.

p-0008Another approach to linear feature detection uses the classical gradient/curvature or Hessian-based detectors. These techniques include the use of thin nets or crest lines as described in Monga et al., and ridges as described in Lang et al., Eberly, and Gauch and Pizer, and are also generally computationally expensive.

p-0009A further method of linear feature detection involves the use of tracking techniques such as stick growing as described in Nelson, and tracking as described in Can et al. and Tolias and Panas. The tracking based approach requires initial locations of linear features, which typically need user intervention. For example, methods using edge operators to detect pairs of edges and graph searching techniques to find the centrelines of vessel segments are presented in Fleagle et al. and Sonka et al. These methods require the user to identify the specific areas of interest. There are also other related techniques which are termed edge or “roof” based, such as those described in Zhou et al., Nevatia and Babu, and Canny. Other known feature detection algorithms include pixel classification using neural network scheme through supervised training as described in Staal, S-Gabor filter and deformable splines as described in Klein, and mathematical morphology as described in Zana and Klein.

REFERENCES

p-0010R. H. Bamberger and M. J. T. Smith. A filter bank for the directional decomposition of images: Theory and design. <i>IEEE Transactions on Signal Processing, </i>40(4):882-893, April 1992.

p-0011M. V. Boland and R. F. Murphy. A neural network classifier capable of recognizing the patterns of all major subcellular structures in fluorescence microscope images of HeLa cells. <i>Bioinformatics, </i>17(12):1213-1223, December 2001.

p-0012A. Can, H. Shen, J. N. Turner, H. L. Tanenbaum, and B. Roysam. Rapid automated tracing and feature extraction from retinal fundus images using direct exploratory algorithms. <i>IEEE Transactions on Information Technology in Biomedicine, </i>3(2):125-138, June 1999.

p-0013J. Canny. Finding edges and lines in images. Technical Report AITR-720, MIT, Artificial Intelligence Laboratory, Cambridge, USA, June 1983.

p-0014S. Chaudhuri, S. Chatterjee, N. Katz; M. Nelson, and M. Goldbaum. Detection of blood vessels in retinal images using two-dimensional matched filters. <i>IEEE Transactions on Medical Imaging, </i>8(3):263-269, September 1989.

p-0015A. C. F. Colchester, R. Ritchings, and N. D. Kodikara. Image segmentation using maximum gradient profiles orthogonal to edges. <i>Image and Vision Computing, </i>8(3):211-217, August 1990.

p-0016D. Eberly. <i>Ridges in Inage and Data Analysis. </i>Kluwer Academic Publishers, 1996.

p-0017M. A. Fischler, J. M. Tenenbaum, and H. C. Wolf. Detection of roads and linear structures in low-resolution aerial imagery using a multisource knowledge integration technique. <i>Computer Graphics and Image Processing, </i>15(3):201-223, March 1981.

p-0018S. R. Fleagle, M. R. Johnson, C. J. Wilbricht, D. J. Skorton, R. F. Wilson C. W. White, M. L. Marcus, and S. M. Collins. Automated analysis of coronary arterial morphology in cineangiograms: Geometric and physiologic validation in humans. <i>IEEE Transactions on Medical Imaging, </i>8(4):387-400, December 1989.

p-0019J. M. Gauch and S. M. Pizer. Multiresolution analysis of ridges and valleys in grey-scale images. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>15(6):635-646, June 1993.

p-0020M. Jacob and M. Unser. Design of steerable filters for feature detection using Canny-like criteria. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>26(8): 1007-1019, August 2004.

p-0021C. Kirbas and F. Quek. A review of vessel extraction techniques and algorithms. <i>ACM Computing Surveys, </i>36(2):81-121, June 2004.

p-0022A. K. Klein, F. Lee, and A. A. Amini. Quantitative coronary angiography with deformable spline models. <i>IEEE Transactions on Medical Imaging, </i>16(5):468-482, October 1997.

p-0023V. Lang, A. G. Belyaev, I. A. Bogaevsici and T. L. Kunii. Fast algorithms for ridge detection In <i>Proceedings of the International Conference on Shape Modeling and Applications, </i>pages 189-197, Aizu-Wakamatsu, Japan, 3-6 Mar. 1997.

p-0024E. Meijering, M. Jacob, J.-C. F. Sarria, and M. Unser. A novel approach to neurite tracing in fluorescence microscopy images. In M. H. Hamza, editor, <i>Proceedings of the Fifth LASTED International Conference on Signal and Image Processing, </i>pages 491-495, Honolulu, USA, 13-15 Aug. 2003.

p-0025O. Monga, N. Armande, and P. Montesinos. Thin nets and crest lines: Application to satellite data and medical images. <i>Computer Vision and Image Understanding, </i>67(3):285-295, September 1997.

p-0026R. C. Nelson. Finding line segments by stick growing. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>16(5):519-523, May 1994.

p-0027R. Nevatia and K. R Babu. Linear feature extraction and description. <i>Computer Graphics and Image Processing, </i>13:257-269, 1980.

p-0028P. Ramm, Y. Alexandrov, A. Cholewinski, Y. Cybuch, R. Nadon, and B. J. Soltys. Automated screening of neurite outgrow. <i>Journal of Bimolecular Screening, </i>8(1):7-18, February 2003.

p-0029P. Soille and H. Talbot. Directional morphological filtering. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>23(11):1313-1329, November 2001.

p-0030M. Sonka, M. D. Winniford, and S. M. Collins. Robust simultaneous detection of coronary borders in complex images. <i>IEEE Transactions on Medical Imaging, </i>14(1):151-161, March 1995.

p-0031J. Staal, M. D. Abramoff, M. Niemeijer, M. A. Viergever, and B. van Ginneken. Ridge-based vessel segmentation in color images of the retina. <i>IEEE Transactions on Medical Imaging, </i>23(4):501-509, April 2004.

p-0032C. Sun and B. Appleton. Multiple paths extraction in images using a constrained expanded trellis. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>27(12):1923-1933, December 2005.

p-0033Y. A. Tolias and S. M. Panas. A fuzzy vessel tracking algorithm for retinal images based on fuzzy clustering. <i>IEEE Transactions on Medical Imaging, </i>17(2):263-273, April 1998.

p-0034G. Van de Wouwer, R. Nuydens, T. Meert, and B. Weyn. Automatic quantification of neurite outgrowth by means of image analysis. In J.-A. Conchello, C. J. Cogswell, and T. Wilson, editors, <i>Proceedings of SPIE, Three</i>-<i>Dimensional and Multidimensional Microscopy: Image Acquisition and Processing XI, </i>volume 5324, pages 1-7, July 2004.

p-0035F. Zana and J.-C. Klein. Segmentation of vessel-like patterns using mathematical morphology and curvature evaluation. <i>IEEE Transactions on Image Processing, </i>10(7):1010-1019, July 2001.

p-0036Y. T. Zhou, V. Venkateswar, and R Chellappa. Edge detection and linear feature extraction using a 2-D random field model. <i>IEEE Transactions on Pattern Analysis and Machine Intelligence, </i>11(1):8495, January 1989.

DISCLOSURE OF THE INVENTION

p-0037It is an object of the preferred embodiments of the invention to provide a new linear feature detection algorithm for image analysis and pattern recognition that requires low computational time to achieve a result.

p-0038In accordance with a first aspect of the present invention, there is provided a method of extracting linear features from an image, the method including the steps of: (a) applying a linear non maximum suppression (NMS) filter to the image at different angles to produce a series of filtered image responses; (b) combining the filtered image responses into a combined image having extracted linear features.

p-0039Preferably, the method also includes the steps of: preprocessing the images by low pass filtering the image and of preprocessing the image to remove areas of low local variances. The method can also include the step of removing small objects from the combined image. Further, the maximum suppression filter can comprise a linear filter or line filter. The linear filter (or line filter) can be a linear window directed at angle 0, 45, 90 and 135 degrees. The combining step can comprise forming the union of the filtered image responses. Preferably, the number of filtered image responses can be 4 although other numbers are also possible. The method can also include the step of: joining closely adjacent extracted linear features.

p-0040In embodiments, one or more post-processing steps are preferably performed either after or in conjunction with the application the non maximum suppression filter. The post processing steps preferably can include one or more steps selected from: a. removing the linear features that do not have a symmetric profile; b. removing the linear features that are preferably small; c. identifying the end points of the linear features; d. identifying and link the end points of the linear features; e. joining closely adjacent extracted linear features; f. thinning the linear features; and g. identifying a feature boundary for the linear features.

p-0041The preprocessing steps preferably can also include one or more steps selected from: a. inverting the image; b. converting colour pixels of the image to an intensity value; c. converting colour pixels of the image to grey scale; d. removing of areas of the image having low local variances; e. masking out areas the image having low local variances; f. smoothing the image; and g. low pass filtering the image. The maximum suppression filter can comprise a line filter. The line filter can be a linear window directed at angle 0, 45, 90 and 135 degrees, or at other angles. The combining can comprise forming the unification of the series of filtered images. In some embodiments, the image can be a three-dimensional structure.

p-0042In accordance with a further aspect of the present invention, there is provided an apparatus for extracting linear features from an image, including: an input means; a processing system; and an output means; wherein the processing system is coupled to the input means; wherein the processing system is coupled to the output means; wherein the processing system is configured to accept an input image from the input device; wherein the processing system is further configured applying a non maximum suppression filter to the input image to produce a series of filtered images; wherein the processing system is further configured to combine the series of filtered images into an output image having extracted linear features; the processing system is further configured to present the image having extracted linear features to the output means.

p-0043In some embodiments, the processing system can be further configured to perform one or more post-processing steps applied to the series of filtered images, wherein the one or more post-processing steps are preferably performed either after or in conjunction with the application of the non maximum suppression filter.

p-0044In some embodiments, the processing system can be further configured to perform one or more post-processing steps applied to the output image having extracted linear features, wherein the one or more post-processing steps are preferably performed either after or in conjunction with the combining the filtered images into an output image having extracted linear features.

BRIEF DESCRIPTION OF THE DRAWINGS

p-0045A preferred embodiment of the invention will now be described, by way of example only, with reference to the accompanying drawings in which:

p-0046<figref idrefs="DRAWINGS">FIG. 1</figref> is a graphical depiction of the linear feature detection filter windows used in the multiple NMS technique according to the invention;

p-0047<figref idrefs="DRAWINGS">FIG. 2(</figref><i>a</i>) is a simple sample image for linear feature detection;

p-0048<figref idrefs="DRAWINGS">FIGS. 2(</figref><i>b</i>) and (<i>c</i>) are respectively the linear detection results of the NMS feature detection taken along a 45 degree and 0 degree directional orientations;

p-0049<figref idrefs="DRAWINGS">FIG. 3(</figref><i>a</i>) is a simple sample image for linear feature detection;

p-0050<figref idrefs="DRAWINGS">FIG. 3(</figref><i>b</i>) is the combined response resulting from the union of the of the NMS images of <figref idrefs="DRAWINGS">FIGS. 2(</figref><i>b</i>) and (<i>c</i>);

p-0051<figref idrefs="DRAWINGS">FIG. 4(</figref><i>a</i>) is a sample image for linear feature detection;

p-0052<figref idrefs="DRAWINGS">FIGS. 4(</figref><i>b</i>) to (<i>e</i>) are respectively the linear detection results of the NMS feature detection taken along 0, 45, 90 and 135 degree directional orientation;

p-0053<figref idrefs="DRAWINGS">FIG. 4(</figref><i>f</i>) is the combined response resulting from the union of the of the NMS images of <figref idrefs="DRAWINGS">FIGS. 4(</figref><i>b</i>) to (<i>e</i>);

p-0054<figref idrefs="DRAWINGS">FIG. 5</figref> represents an intensity profile along a linear window when examined across a linear feature;

p-0055<figref idrefs="DRAWINGS">FIG. 6</figref> represents an intensity profile along a linear widow when examined across a step feature;

p-0056<figref idrefs="DRAWINGS">FIG. 7(</figref><i>a</i>) is a sample image for linear feature detection;

p-0057<figref idrefs="DRAWINGS">FIGS. 7(</figref><i>b</i>) and (<i>c</i>) are the results of NMS linear feature detection on the image of <figref idrefs="DRAWINGS">FIG. 7(</figref><i>a</i>) with and without a symmetry check respectively;

p-0058<figref idrefs="DRAWINGS">FIG. 8</figref> illustrates a second intensity profine maximum along a linear window;

p-0059<figref idrefs="DRAWINGS">FIG. 9</figref> shows the detection of gap linkage between the end points of the detected features;

p-0060<figref idrefs="DRAWINGS">FIG. 10</figref> shows the results of gap linkage between the end points of the detected features;

p-0061<figref idrefs="DRAWINGS">FIG. 11</figref> and <figref idrefs="DRAWINGS">FIG. 12</figref> illustrate an extension to 3-D linear windowing;

p-0062<figref idrefs="DRAWINGS">FIG. 13</figref> represents a flow chart for one embodiment of the present invention;

p-0063<figref idrefs="DRAWINGS">FIG. 14</figref> is a Table showing the computational running times of the NMS feature detection algorithm in seconds;

p-0064<figref idrefs="DRAWINGS">FIG. 15</figref> is a Table showing the running times (in seconds) for different number of orientations of the NMS feature detection algorithm;

p-0065<figref idrefs="DRAWINGS">FIG. 16</figref> (<i>a</i>) is a sample image for linear feature detection;

p-0066<figref idrefs="DRAWINGS">FIGS. 16</figref> (<i>b</i>), (<i>c</i>) and (<i>d</i>) are the results form NMS linear feature detection performed using 2, 4, and 8 different orientations respectively;

p-0067<figref idrefs="DRAWINGS">FIG. 17</figref> illustrate an approach to obtain linear feature boundary;

p-0068<figref idrefs="DRAWINGS">FIG. 18</figref> shows and embodiment of a processing system configured to implement a method for extracting linear features from an image; and

p-0069<figref idrefs="DRAWINGS">FIGS. 19 to 22</figref> re further examples of linear feature detection on a variety of sample images.

PREFERRED EMBODIMENT OF THE INVENTION

p-0070The preferred embodiments of the invention provide a new algorithm for linear feature detection using multiple directional non-maximum suppression (NMS). This primarily involves NMS being applied to an image along multiple directions as defined by linear windows. Since there is no complicated filtering of the image involved, the algorithm is generally less computationally intensive than other methods.

Multiple Directional Non-Maximum Suppression

p-0071As described above, Hessian-based linear detection techniques or matched filters require a large amount of computational time and direct use of the image gradient is typically not very reliable. Instead of detecting the local direction of the linear features directly, the preferred embodiment detects the linear features by utilising the responses of multiple directional NMS.

p-0072Linear features are defined as a sequence of points where the image has an intensity maximum in the direction of the largest variance, gradient, or surface curvature. A directional local maximum is a pixel that is not surrounded by pixels of higher grey values (where the lowest grey value corresponds to black, and the highest grey value being white) in a linear window overlaid on the image. Linear features in the image may either correspond to local maximum or minimum in the grey value of the pixels and accordingly the linear feature may be bright (high grey value) or dark (low grey value) depending on the object property and the imaging techniques used. The following description of the preferred embodiments will describe the detection of bright linear features only for succinctness. It will be appreciated by those skilled in the art that dark features may be similarly detected by simply inverting the intensity of the input image, or using non-minimum suppression rather than non-maximum suppression. It will also be appreciated by those skilled in the art that colour images may be modified such that colours of interest are assigned a higher “grey value”.

p-0073NMS is a process for marking all pixels whose intensity is not maximal within a certain local neighbourhood as zero and all pixels whose intensity is maximal within the same local neighbourhood as one. In this process, the pixels belonging to a local neighbourhood are selected using a series of linear windows, whereby the pixels are specified geometrically with respect to the pixel under test. L_D_{<sub2>i </sub2>}is the result for non-maximum suppression in the direction D.

p-0074The outputs of directional non-maximum suppression are not greyscale images such as those obtained using most of the traditional directional filter methods, but rather, are binary images. The pixels of the binary image are each marked as 1 for the local directional maximum, and 0 (zero) otherwise.

Combining Multiple Directional Non-Maximum Suppression

p-0075Multiple NMS windows are applied and the union of the results from each window is constructed. To combine the outputs of the NMS images obtained at different angles, the following formula can be used:

p-0076<maths id="MATH-US-00001" num="00001"><math overflow="scroll"><mrow><mi>L</mi><mo>=</mo><mrow><munder><mover><mo>⋃</mo><msub><mi>N</mi><mi>D</mi></msub></mover><mrow><mi>i</mi><mo>=</mo><mn>1</mn></mrow></munder><mo></mo><msub><mi>L</mi><msub><mi>D</mi><mi>l</mi></msub></msub></mrow></mrow></math></maths>

p-0077where L_D_{<sub2>i </sub2>}is the result for non-maximum suppression in the direction D_i, N_Dis the number of directions used, and L is the combined result. The number of directions used, N_D, are typically 4, 8 or another multiple of 2.

p-0078<figref idrefs="DRAWINGS">FIG. 1</figref> displays a series of four linear windows <b>100</b> that may be used in one embodiment of the present invention, wherein the pixel under test <b>120</b> is marked. These four exemplary linear windows <b>110</b> though <b>113</b> are constructed at angles of 0°, 45°, 90°, and 135° respectively, and where the position X corresponds to the pixel under test <b>120</b> of the linear window. Additional directions, such as those at 22.5°, 67.5°, 112.5° and 157.5° may also be used (not shown). All pixels in an image are visited by the NMS process. Also, a maximum along a linear window is kept only if it lies in the centre of the linear window <b>120</b>.

p-0079<figref idrefs="DRAWINGS">FIG. 2</figref> displays a simple example of the application of only two directional NMS windows, <b>221</b> and <b>231</b>, as applied to a tri-level image <b>210</b>. The application of the first directional NMS windows <b>221</b> to the original image <b>210</b> produces Result A <b>220</b>. The application of the second directional NMS windows <b>231</b> to the original image <b>210</b> produces Result B <b>230</b>.

p-0080<figref idrefs="DRAWINGS">FIG. 3</figref> displays the union of results <b>320</b> of the previous example, wherein only two directional NMS windows, <b>221</b> and <b>231</b> are applied to a tri-level image <b>210</b>. The result for this example <b>320</b> is constructed from the union of Result A <b>220</b> and Result B <b>230</b>.

p-0081<figref idrefs="DRAWINGS">FIG. 4</figref> shows the NMS responses and the combined output as generated by the present invention. <figref idrefs="DRAWINGS">FIG. 4(</figref><i>a</i>) shows the input image displaying nitrite outgrowths, to which NMS is applied. <figref idrefs="DRAWINGS">FIGS. 4(</figref><i>b</i>) to (<i>e</i>) are the directional non-maximum suppression responses at angles 0°, 45°, 90°, and 135° respectively and using a linear window length of 11. <figref idrefs="DRAWINGS">FIG. 4(</figref><i>f</i>) shows the union of the multiple NMS responses given in <figref idrefs="DRAWINGS">FIGS. 4(</figref><i>b</i>) to (<i>e</i>).

p-0082<figref idrefs="DRAWINGS">FIG. 4(</figref><i>f</i>) represents an output, where it can be seen that the linear features present in <figref idrefs="DRAWINGS">FIG. 4(</figref><i>a</i>) are mostly shown in <figref idrefs="DRAWINGS">FIG. 4(</figref><i>f</i>). The outputs of directional non-maximum suppression are not greyscale images such as those obtained using most of the traditional directional filter methods, but rather, are binary images. The pixels of the binary image are each marked as 1 for the local directional maximum, and 0 (zero) otherwise.

Fast Linear Feature Detection

Background Processing

p-0083A small amount of smoothing of the image may reduce the possibility of multiple pixels being detected as a local maximum during processing. The smoothing window is not to be chosen excessively large—as this would degrade the fine-scale features of the image to the extent that they are not detected. A small window, of size e.g. 3×3 pixels is generally sufficient. Smoothing of the image can be achieved though low pass filtering.

p-0084To increase the processing speed, it is also possible to estimate the local variance of a point in the image and mask out those points with low local variance since a point with a low variance within a local window is unlikely to be on a linear feature.

Symmetry Check

p-0085It can generally be expected that the cross feature profile of a linear feature at a particular position in the image is near symmetrical around the maximum point. <figref idrefs="DRAWINGS">FIG. 5</figref> illustrates the parameters that can be used for checking whether a local maximum lies on a linear feature in the image. The amplitudes of the pixels selected across a linear feature <b>510</b> are shown. I_max<b>511</b> is the value of the local directional maximum. I_average1<b>520</b> and I_average2<b>530</b> are the average values of half windows on the two sides of the local maximum; and I_diff1<b>521</b> and I_diff2<b>531</b> are the differences between the maximum value and the two average values. For a linear feature in the image, both of the I_diff1and I_diff2values should be larger than a particular threshold value. If only one of the difference values is larger, the position of the local maximum resembles more like a step edge rather than a linear feature.

p-0086<figref idrefs="DRAWINGS">FIG. 6</figref> further illustrates the parameters that can be used for checking whether a local maximum lies on a linear feature in the image. The amplitudes of the pixels selected across a ‘step edge’ feature <b>610</b> are shown. I_max<b>611</b> is the value of the local directional maximum. I_{average 1}<b>620</b> and I_average2<b>630</b> are the average values of half windows on the two sides of the local maximum; and I_diff1<b>621</b> and I_diff2<b>631</b> are the differences between the maximum value and the two average values. For a ‘step edge’ feature in the image, one I_diffvalues will be smaller than a particular threshold value. If only one of the difference values is larger, the position of the local maximum resembles more like a step edge rather than a linear feature.

p-0087<figref idrefs="DRAWINGS">FIG. 7</figref> shows the effect of using a symmetry check especially to distinguish linear features from areas where step edges exist. <figref idrefs="DRAWINGS">FIG. 7(</figref><i>a</i>) is the input image; <figref idrefs="DRAWINGS">FIG. 7(</figref><i>b</i>) is the result without the symmetry check, and <figref idrefs="DRAWINGS">FIG. 7(</figref><i>c</i>) is the result with symmetry check. The effect of not performing the symmetry check is that the resulting image contains increased edge output and spurious noise artefacts as can be seen by comparing <figref idrefs="DRAWINGS">FIG. 7(</figref><i>b</i>) and <figref idrefs="DRAWINGS">FIG. 7(</figref><i>c</i>).

Extending to Multiple Local Maximums

p-0088In some images, certain linear features may be very close to each other. The NMS process described above may detect only one of the two or several linear features that are next to each other. <figref idrefs="DRAWINGS">FIG. 8</figref> shows an example where there are two local maxima <b>810</b>, <b>820</b> along this linear window. Only the indicated centre peak <b>810</b> may be detected. The NMS process can be extended so that multiple local maxima can be detected for any particular linear window. Once a local maximum has been found at the centre of a linear window, a second or more local maxima within the same window can be searched so that for a particular linear window direction, multiple local maxima, hence multiple linear features, can be detected. <figref idrefs="DRAWINGS">FIG. 8</figref> shows that a second local maximum <b>820</b> with intensity I_max2can be found. One may also check that the intensity I_max2should not be very different from I_maxfor the secondary maximum to be retained.

Linking Broken Linear Features

p-0089The resulting image from the union of the NMS for linear features may occasionally present gaps in the detected linear features. To correct for these gaps, a linking process may be employed to connect the end points of neighbouring linear features with short line segments if the distance between end points is small.

p-0090The end points of linear features can be detected by first pruning the end points within a binary image. The difference between the pruned and the input binary images then provides the end points. A thinning operation may be carried out before the pruning step so that better end points on linear features may be obtained.

p-0091<figref idrefs="DRAWINGS">FIG. 9</figref> shows two sub-regions of a binary image <b>910</b> and <b>920</b> with some end points on linear features marked, for example <b>911</b> and <b>921</b>. <figref idrefs="DRAWINGS">FIG. 9</figref> shows the sub-images contain gaps for a given linear feature, for example <b>912</b>, <b>913</b>, <b>914</b>, <b>915</b>, <b>922</b>, <b>923</b>, and <b>925</b>.

p-0092The distance between the identified end points is then calculated. If the distance of end points is within a threshold and the end points belong to different linear features, a short line segment can be drawn to connect the two end points.

p-0093<figref idrefs="DRAWINGS">FIG. 10</figref> shows two sub-regions of a binary image <b>1010</b> and <b>1020</b> with end points on linear features marked and that the sub-images contain gaps for a given linear feature. <figref idrefs="DRAWINGS">FIG. 10</figref> shows two sub-images where some of the gaps have been connected with short line segments <b>1012</b>, <b>1013</b>, <b>1014</b>, <b>1022</b> and <b>1023</b>. Some gaps remain non- connected because each end point belonged to the same feature <b>1025</b> (the sub-image does not show the connection of the same feature), or the distance between the end points were larger than the threshold <b>1015</b>.

p-0094For a pair of linear feature segments where each segment may have several end points detected, it is desirable to choose a best connection for the two segments using some of the end points that are close to each other. There may be several possible pairs of end points from the two segments to choose from. The connection can be chosen where the average image intensity on the connecting line segment is the maximum. Other information that can be used when selecting the best connection include the intensity values at the end points, the distance between the two connecting end points, the orientation information of each segment at each of the end points. While linking endpoints, the order with which this operation is carried out influences the final result: endpoints that have already been attributed in pairs are not free to be linked to another endpoint, even if the new link would result in a better overall outcome. By using a global optimization technique whereby all end-points are considered simultaneously using for example a graph-based algorithm, better results can be obtained.

p-0095Shortest path techniques may also be utilised to link linear features that have not been paired by the threshold process based on minimum distance. This is employed by ensuring that links between end points follow through the ridges of the detected linear features. The use of shortest path technique may add some cost to the feature detection process, but since the majority of the end points to be connected have already been connected by the distance threshold linking process, the additional processing time needed is typically minimal. Sometimes, it may be necessary to connect the end point of a linear feature to the middle of another linear feature, so that a “T” or “Y” shaped connection of two linear features can be made. In cases where there is a single pixel gap between neighbouring linear features but there is no end point on one of the features, this single pixel can be filled or linked.

Extension to 3D Images

p-0096The algorithms developed above can be easily extended for 3D images. For 3D images, the linear windows used need to be oriented in 3D space. <figref idrefs="DRAWINGS">FIG. 11</figref> shows an example of three 3D linear windows, <b>1110</b>, <b>111</b>, and <b>1112</b> of length <b>9</b> with their centre marked <b>1120</b>. <figref idrefs="DRAWINGS">FIG. 11(</figref><i>d</i>) shows the union of the three 3D linear windows <b>1130</b>. <figref idrefs="DRAWINGS">FIG. 12</figref> shows examples of the union of nine 3D linear windows <b>1200</b>. It will be appreciated by those skilled in the art that additional directions in 3D space could also be used.

p-0097An input image may be stored in a three-dimensional representation in a three dimensional array or as a plurality of 2 dimensional slices. It will be appreciated by those skilled in the art that other configurations are suitable for storing an image in a three-dimensional representation.

Algorithm Steps

p-0098The algorithm steps that may be included in the linear feature detection are represented in <figref idrefs="DRAWINGS">FIG. 13</figref>.

p-0099The steps of in the linear feature detection method <b>1300</b> can include: <ul><li id="ul0001-0001" num="0000"><ul><li id="ul0002-0001" num="0097">1. Inversion of the image intensity, if necessary, for example when detecting dark linear features, <b>1310</b>.</li><li id="ul0002-0002" num="0098">2. Pre-process the image including background extraction, creation of a low variance mask, and small window smoothing, <b>1320</b>.</li><li id="ul0002-0003" num="0099"><b>3</b>. Carry out 4 or 8 (3 or 9 for 3D case) directional non-maximum suppression on the images obtained in Step 2, and only carry out non-maximum suppression on foregrounds <b>1330</b>, and combine the multiple NMS outputs into one image <b>1340</b>. The symmetry check is performed during the NMS process <b>1350</b>. Multiple local maxima can also be found within a linear window.</li><li id="ul0002-0004" num="0100">4. Remove small features <b>1360</b>.</li><li id="ul0002-0005" num="0101">5. Obtain and link the end points of different linear objects when the end points are close <b>1370</b>.</li><li id="ul0002-0006" num="0102">6. Thin the obtained image if necessary <b>1380</b>.</li></ul></li></ul>

p-0100It will be appreciated by those skilled in the art that some of the above algorithm steps may be transposed or excluded.

Timings

p-0101The running time of the algorithm may depend on a number of factors. These factors include:

p-0102the threshold for background extraction;

p-0103the linear window size;

p-0104the number of linear windows used; and

p-0105the image contents.

p-0106If most regions of the image to be processed consists of relatively flat regions (in terms of the grey levels), then after the background extraction step of the algorithm, the number of remaining pixels surrounding the linear features of interest will be significantly reduced and the computational load will be correspondingly reduced.

p-0107The algorithm was tested on a variety of 8 bit images using a 2.40 GHz CPU Intel Pentium 4 computer running Linux. Colour input images were also tested, with only the intensity information being used for feature detection (individual channels of the colour image can also be used). The Table of <figref idrefs="DRAWINGS">FIG. 14</figref> shows the running time in seconds, including image IO, of the NMS method on different sized images and compared with the running time for linear feature detection using a traditional ranked filter approach. It can be seen that the NMS algorithm is approximately 5 times faster than the ranked filter approach.

p-0108The Table of <figref idrefs="DRAWINGS">FIG. 15</figref> gives the running times for the NMS algorithm with different number of orientations used. <figref idrefs="DRAWINGS">FIG. 15</figref> shows that for increases in the number of orientations, there is an increased computational cost. <figref idrefs="DRAWINGS">FIG. 16</figref> shows an example of the linear feature detection results achieved when using different number of orientation windows. <figref idrefs="DRAWINGS">FIG. 16(</figref><i>a</i>) is the input image and <figref idrefs="DRAWINGS">FIG. 16(</figref><i>b</i>), <b>16</b>(<i>c</i>), and <b>16</b>(<i>d</i>) are the results obtained with 2, 4, and 8 orientations respectively. <figref idrefs="DRAWINGS">FIG. 16</figref> shows an increased number of orientations, associated with an increased computational cost, gives a superior feature connection thereby reducing the need for post-process linking.

p-0109The NMS linear detection algorithm is very fast compared with the fastest algorithms reported in the literature. The main step is applying the non-maximum suppression at a number of different directions—typically 4—along linear windows. There is no computationally intensive operation such as matched filters or curvature calculation and the algorithm only requires a minimum number of tuning parameters, such as the linear window size, the number of directions/orientations, the threshold for background extraction, the threshold for symmetry check the maximum gap size; and the size of small objects/artefacts in the image to be removed. Furthermore, the NMS algorithm is not sensitive to background variation since no global threshold is used and all operations are local, keeping computations to a minimum.

p-0110The speed of the algorithm depends primarily on the content of the input images. If the intensity landscape is mostly flat, then the background estimation step will classify most regions of the image as background, resulting in a reduced amount of foreground pixels, that is, pixels close to linear features of interest that need to be processed. If a feature boundary is needed, extra steps are necessary to obtain the boundaries that can be used for width measurements. <figref idrefs="DRAWINGS">FIG. 17</figref> shows the linear feature skeleton <b>1710</b> bounded by two constructed boundaries <b>1720</b> and <b>1721</b> separated by distance <b>1730</b>. Shortest path or multiple paths techniques can be used for such purposes by first obtaining the linear feature skeleton <b>1710</b>. Then, the gradient perpendicular to the linear feature can be calculated. The gradient values on the left and the right sides of the skeleton can be added together and a shortest path next to the skeleton can be found. Multiple paths methods may also be used to extract the linear feature boundaries as discussed in Sun and Appleton.

p-0111<figref idrefs="DRAWINGS">FIG. 18</figref> shows an embodiment of a computer system <b>1800</b> used to implement a system for extracting linear features from an image, comprising a processing system <b>1810</b>, keyboard as an input device <b>1820</b> and a monitor as an output device <b>1830</b>. Further, the input image may be provided on a storage medium or via a computer network, input parameters may be pre-configured or entered at run time, the output image of linear features may be displayed on a monitor, saved to a storage medium or sent over a computer network. It will be appreciated by those skilled in the art that alternative or combinations of input devices or output devices are suitable for implementing alternative embodiments.

p-0112Further examples of the results obtained by the NMS linear detection feature algorithm can be seen in <figref idrefs="DRAWINGS">FIGS. 19 to 22</figref>. <figref idrefs="DRAWINGS">FIG. 19</figref> shows the NMS linear detection feature algorithm when applied to aerial or satellite images is capable of extracting geographical features. <figref idrefs="DRAWINGS">FIG. 20</figref> shows the NMS linear detection feature algorithm when applied to images of natural or constructed ground covering the algorithm is capable of extracting features. <figref idrefs="DRAWINGS">FIG. 21</figref> shows the NMS linear detection feature algorithm applied to biological images in which it is capable of extracting information including biometric information. <figref idrefs="DRAWINGS">FIG. 22</figref> shows the NMS linear detection feature algorithm applied to botanical images.

p-0113The linear feature detection method described herein, and/or shown in the drawings are presented by way of example only and are not limiting as to the scope of the invention. Unless otherwise specifically stated, individual aspects and components of the algorithm may be modified, or may have been substituted therefore known equivalents, or as yet unknown substitutes such as may be developed in the future or such as may be found to be acceptable substitutes in the future. The algorithm may also be modified for a variety of applications while remaining within the scope and spirit of the claimed invention, since the range of potential applications is great, and since it is intended that the present linear feature detection method be adaptable to many such variations.