(19)
(11)EP 2 856 429 B1

(12)EUROPEAN PATENT SPECIFICATION

(45)Mention of the grant of the patent:
29.07.2020 Bulletin 2020/31

(21)Application number: 13735059.1

(22)Date of filing:  28.05.2013
(51)International Patent Classification (IPC): 
G06T 7/42(2017.01)
(86)International application number:
PCT/GB2013/051412
(87)International publication number:
WO 2013/179021 (05.12.2013 Gazette  2013/49)

(54)

METHODS AND APPARATUS FOR IMAGE PROCESSING, AND LASER SCANNING OPHTHALMOSCOPE HAVING AN IMAGE PROCESSING APPARATUS

VERFAHREN UND VORRICHTUNG ZUR BILDVERARBEITUNG SOWIE LASERABTASTUNGS-OPHTHALMOSKOP MIT EINER BILDVERARBEITUNGSVORRICHTUNG

PROCÉDÉS ET APPAREIL POUR LE TRAITEMENT D'IMAGES, ET OPHTALMOSCOPE À BALAYAGE LASER AYANT UN APPAREIL DE TRAITEMENT D'IMAGES


(84)Designated Contracting States:
AL AT BE BG CH CY CZ DE DK EE ES FI FR GB GR HR HU IE IS IT LI LT LU LV MC MK MT NL NO PL PT RO RS SE SI SK SM TR

(30)Priority: 28.05.2012 GB 201209390

(43)Date of publication of application:
08.04.2015 Bulletin 2015/15

(73)Proprietor: Optos Plc
Dunfermline, Scotland, KY11 8GR (GB)

(72)Inventor:
  • CLIFTON, David
    Shandon Edinburgh EH11 1SR (GB)

(74)Representative: Hoffmann Eitle 
Patent- und Rechtsanwälte PartmbB Arabellastraße 30
81925 München
81925 München (DE)


(56)References cited: : 
  
  • SURYA PRAKASH ET AL: "Segmenting Multiple Textured Objects Using Geodesic Active Contour and DWT", 18 December 2007 (2007-12-18), PATTERN RECOGNITION AND MACHINE INTELLIGENCE; [LECTURE NOTES IN COMPUTER SCIENCE], SPRINGER BERLIN HEIDELBERG, BERLIN, HEIDELBERG, PAGE(S) 111 - 118, XP019084438, ISBN: 978-3-540-77045-9 abstract figure 2 Section 3
  • YINGEN XIONG ET AL: "Hand Motion Gesture Frequency Properties and Multimodal Discourse Analysis", INTERNATIONAL JOURNAL OF COMPUTER VISION, KLUWER ACADEMIC PUBLISHERS, BO, vol. 69, no. 3, 1 April 2006 (2006-04-01), pages 353-371, XP019410160, ISSN: 1573-1405, DOI: 10.1007/S11263-006-8112-5
  
Note: Within nine months from the publication of the mention of the grant of the European patent, any person may give notice to the European Patent Office of opposition to the European patent granted. Notice of opposition shall be filed in a written reasoned statement. It shall not be deemed to have been filed until the opposition fee has been paid. (Art. 99(1) European Patent Convention).


Description

FIELD



[0001] The present invention relates to image processing, more particularly to methods and apparatus for image processing. In one example, the present invention relates to but not exclusively to methods of textural analysis applied in the field of ocular imaging.

BACKGROUND



[0002] It is well known to capture image data using digital image sensors. An array of light sensitive picture elements (pixels) is provided, which can be manufactured as charge coupled devices or using complementary metal oxide semiconductor (CMOS) techniques. Incident radiation causes a charge to be generated at each pixel. This charge is converted to a voltage whose value is then digitised. The value of the voltage depends upon the intensity of the illumination incident on the pixel during an integration time, when the pixels are set in a light sensitive mode. A pixel array can be formed in one, two or three dimensions. The most common form is a two dimensional pixel array as found in everyday cameras for example. A one dimensional pixel array is typically referred to as a "linear array", and an image sensor with such an array can be termed a linear sensor or a linear scanner.

[0003] The set of intensity values derived from the pixel array is known as image data. The "raw" image data output by the pixel array may be subjected to various post- processing techniques in order to reproduce an image either for viewing by a person or for processing by a machine. These post-processing techniques can include various statistical methods for image analysis and for performing various camera and image processing functions.

[0004] One example of such techniques is the recognition and/or classification of textures within an image. The texture can represent surface characteristics of an object or region, and can be used as a basis for identifying different objects or different regions of an object within an image. Modelling texture is usually characterised by variations in signal intensity, and sometimes the spatial relationship (local neighbourhood properties) of these variations in an image.

[0005] It is known to use two dimensional wavelet transforms in a method for modelling texture. Wavelet transforms are useful because they give the ability to construct a time-frequency representation of a signal that offers very good time and frequency localisation.

[0006] An introduction to wavelets can be found from US 2010/0014761 and is provided below in the detailed description section.

[0007] Further background is provided in the article titled "Segmenting Multiple Textured Objects Using Geodesic Active Contour and DWT" by S. Prakash et al., published in the International Conference on Pattern Recognition and Machine Intelligence (PReMI 2007), at pages 111-118 on 18 December 2007. The article addresses the issue of segmenting multiple textured objects in presence of a background texture. The proposed technique is based on Geodesic Active Contour (GAC) and can segment multiple textured objects from the textured background. For an input texture image, a texture feature space is created using scalogram obtained from discrete wavelet transform (DWT). Then, a 2-D Riemannian manifold of local features is extracted via the Beltrami framework. The metric of this surface provide a good indicator of texture changes, and therefore, is used in GAC algorithm for texture segmentation. The main contribution in this work lie in the development of new DWT and scalogram based texture features which have a strong discriminating power to define a good texture edge metric which is used in GAC technique. The described technique is validated using a set of synthetic and natural texture images.

[0008] However these methods are not tolerant to fragmentation of the textural features being measured, that is, when the textural features comprise diffuse, irregular, broken or spaced patterns in the image. There are also limits to the resolution and resolving power of existing techniques.

[0009] There is therefore a need for a method that is more robust to fragmentation in the textural images, and/or that has improved resolution, and/or that has an improved resolving power with respect to existing wavelet based techniques.

SUMMARY



[0010] The present invention provides a method of image processing according to claim 1, an image processing apparatus according to claim 11, and a computer program product according to claim 12. Optional features are set out in the remaining claims.

BRIEF DESCRIPTION OF THE DRAWINGS



[0011] Embodiments of the present invention are described, by way of example only, with reference to the accompanying drawings in which:

Figure 1 illustrates steps of an exemplary algorithm according to a first embodiment of the present invention;

Figure 2 illustrates aspects of the step of mapping an image along slice coordinates illustrated in Figure 1;

Figure 3 illustrates aspects of the step of computing a wavelet scalogram of a slice map as illustrated in Figure 1;

Figure 4 shows example scalogram ridges, for illustration purposes;

Figure 5 illustrates aspects of the step of mapping ridge features on a scalogram surface as illustrated in Figure 1;

Figure 6 illustrates aspects of the steps of superimposing ridge features for each slice into a ridge composite mapping, generating a features position histogram and thresholding the histogram, to determine a frequency of upper components as illustrated in Figure 1;

Figure 7 illustrates aspects of the step of outputting a feature count over an image's selected region as illustrated in Figure 1;

Figure 8 illustrates an exemplary embodiment for wavelet based surface texture mapping for a retinal image region showing a step where an image is mapped in radial slices;

Figure 9 illustrates the superimposition of slices shown in Figure 8 in an orthogonal texture image with equal sample intervals;

Figure 10 illustrates the condensing of the information from Figure 9 into the frequency domain in the form of overlaid scalograms for each slice of the orthogonal texture image;

Figure 11 illustrates a peak in wavelet scale entropy from the diagram of Figure 10, used to parameterise the retinal texture;

Figure 12 illustrates a further aspect of wavelet based surface texture mapping for the retinal region showing a step where a radial scan is imaged to map the extent of retinal texture characteristics;

Figure 13 illustrates the scalogram entropy peak value plotted against scalogram entropy peak position, showing how retinal texture can be partitioned between a set of lid textures and a set of lash textures; and

Figure 14 illustrates an exemplary application of a textural classification method according to an embodiment of the present invention.


DETAILED DESCRIPTION



[0012] There is described in the following a method of image processing comprising:
  1. (i) mapping an image along a one dimensional slice;
  2. (ii) computing a wavelet scalogram of the slice;
  3. (iii) mapping ridge features within the wavelet scalogram;
repeating steps (i), (ii) and (iii) for one or more mapped image slices;
superimposing the mapped ridge features from the slices; and deriving textural information from said superimposed mapped ridge features.

[0013] Because the textural analysis is performed based on the extracted ridge features only, the method provides robust, reliable performance in cases where the textural features being measured are fragmented.

[0014] The step of deriving textural information from said superimposed mapped ridge features comprises thresholding or deriving statistics from a histogram representation of the 2D frequency and scale space spanning the mapped ridge features.

[0015] Optionally, said step of mapping an image along a one dimensional slice comprises selecting a row or column of image data from the image.

[0016] Alternatively, said step of mapping an image along a one dimensional slice comprises mapping a path along a set of image data pixels to a straight line representation.

[0017] Optionally, said path along a set of image data pixels comprises a straight line which is angled away from the horizontal or vertical axes defined by the rows and columns of the image data pixels.

[0018] Optionally, said path along a set of image data pixels comprises a circular path, or part thereof.

[0019] Optionally, said selected row or column or said path is chosen to have a directionality which corresponds to a directionality of a characteristic of the image which is expected or is being searched for.

[0020] Optionally, said characteristic of the image is a feature associated with a known pathology or medical condition.

[0021] Examples of these pathologies or conditions include retinal ischemic areas, macular degeneration, or other regions of retinal pathologies.

[0022] Optionally, said derived textural information is used to classify regions of a retinal image as comprising either retinal texture, or texture which comprises lid or lashes.

[0023] Optionally, the regions comprising lid or lash textures are excluded from a subsequent further image analysis and/or examination procedure.

[0024] Optionally, said derived textural information is used to identify scanner-specific image anomalies.

[0025] Examples include subtle, periodic variations in image intensity that occur due to inaccuracies in optical scan components.

[0026] Optionally, the step of computing a wavelet scalogram of the slice comprises application of a continuous wavelet transform.

[0027] Optionally, the step of computing a wavelet scalogram comprises selecting a characteristic frequency of the applied wavelet transform to match a chosen shape.

[0028] Optionally, the step of computing a wavelet scalogram comprises selecting a scale of the applied wavelet transform to match a chosen feature size.

[0029] The selection of one of both of a characteristic frequency and a scale of the applied wavelet transform allows it to be tuned to latch on to features of interest, which are expected to be present in an image or which are being searched for.

[0030] Optionally, a characteristic frequency and a scale are chosen to match the size and shape of rods and/or cones in a retina.

[0031] There is also described in the following an image processing apparatus comprising: means for receiving image data from an image sensor; and
a processor arranged to:
  1. (i) map an image along a one dimensional slice;
  2. (ii) compute a wavelet scalogram of the slice;
  3. (iii) map ridge features within the wavelet scalogram;
repeat steps (i), (ii) and (iii) for one or more mapped image slices;
superimpose the mapped ridge features from the slices; and
derive textural information from said superimposed mapped ridge features.

[0032] Optionally, the apparatus is arranged to perform the methods of the features defined above.

[0033] There is also described a laser scanning ophthalmoscope having an image processing apparatus comprising:

means for receiving image data from an image sensor; and

a processor arranged to:

  1. (i) map an image along a one dimensional slice;
  2. (ii) compute a wavelet scalogram of the slice;
  3. (iii) map ridge features within the wavelet scalogram;

repeat steps (i), (ii) and (iii) for one or more mapped image slices;

superimpose the mapped ridge features from the slices; and

derive textural information from said superimposed mapped ridge features.



[0034] There is also described a computer program product encoded with instructions that when run on a computer, cause the computer to receive image data and perform a method of image processing comprising:
  1. (i) mapping an image along a one dimensional slice;
  2. (ii) computing a wavelet scalogram of the slice;
  3. (iii) mapping ridge features within the wavelet scalogram;
repeating steps (i), (ii) and (iii) for one or more mapped image slices;
superimposing the mapped ridge features from the slices; and
deriving textural information from said superimposed mapped ridge features.

[0035] The computer program product may be stored on or transmitted over as one or more instructions or code on a computer-readable medium. Computer-readable media includes both computer storage media and communication media including any medium that facilitates transfer of a computer program from one place to another. A storage media may be any available media that can be accessed by a computer. By way of example such computer-readable media can comprise RAM, ROM, EEPROM, CD-ROM or other optical disk storage, magnetic disk storage or other magnetic storage devices, or any other medium that can be used to carry or store desired program code in the form of instructions or data structures and that can be accessed by a computer. Also, any connection is properly termed a computer-readable medium. For example, if the software is transmitted from a website, server, or other remote source using a coaxial cable, fiber optic cable, twisted pair, digital subscriber line (DSL), or wireless technologies such as infrared, radio, and microwave, then the coaxial cable, fiber optic cable, twisted pair, DSL, or wireless technologies such as infrared, radio, and microwave are included in the definition of medium. Disk and disc, as used herein, includes compact disc (CD), laser disc, optical disc, digital versatile disc (DVD), floppy disk and Blu-ray disc where disks usually reproduce data magnetically, while discs reproduce data optically with lasers. Combinations of the above should also be included within the scope of computer-readable media. The instructions or code associated with a computer-readable medium of the computer program product may be executed by a computer, e.g., by one or more processors, such as one or more digital signal processors (DSPs), general purpose microprocessors, ASICs, FPGAs, or other equivalent integrated or discrete logic circuitry.

[0036] The following introductory description of wavelets is taken from US 2010/0014761, with some minor modifications:

Introduction to Wavelet Transform



[0037] The continuous wavelet transform of a signal x(t) in accordance with the present disclosure may be defined as

where ψ*(t) is the complex conjugate of the wavelet function ψ(t), a is the dilation parameter of the wavelet and b is the location parameter of the wavelet. The transform given by equation (1) may be used to construct a representation of a signal on a transform surface. The transform may be regarded as a time-scale representation. Wavelets are composed of a range of frequencies, one of which may be denoted as the characteristic frequency of the wavelet, where the characteristic frequency associated with the wavelet is inversely proportional to the scale a. One example of a characteristic frequency is the dominant frequency. Each scale of a particular wavelet may have a different characteristic frequency.

[0038] The continuous wavelet transform decomposes a signal using wavelets, which are generally highly localized in time. The continuous wavelet transform may provide a higher resolution relative to discrete transforms, thus providing the ability to garner more information from signals than typical frequency transforms such as Fourier transforms (or any other spectral techniques) or discrete wavelet transforms. Continuous wavelet transforms allow for the use of a range of wavelets with scales spanning the scales of interest of a signal such that small scale signal components correlate well with the smaller scale wavelets and thus manifest at high energies at smaller scales in the transform. Likewise, large scale signal components correlate well with the larger scale wavelets and thus manifest at high energies at larger scales in the transform. Thus, components at different scales may be separated and extracted in the wavelet transform domain. Moreover, the use of a continuous range of wavelets in scale and time position allows for a higher resolution transform than is possible relative to discrete techniques.

[0039] The energy density function of the wavelet transform, the scalogram, is defined as

where '∥' is the modulus operator. The scalogram may be rescaled for useful purposes. One common resealing is defined as

and is useful for defining ridges in wavelet space when, for example, the Morlet wavelet is used. Ridges are defined as the locus of points of local maxima in the plane.

[0040] For implementations requiring fast numerical computation, the wavelet transform may be expressed as an approximation using Fourier transforms. Pursuant to the convolution theorem, because the wavelet transform is the cross-correlation of the signal with the wavelet function, the wavelet transform may be approximated in terms of an inverse FFT of the product of the Fourier transform of the signal and the Fourier transform of the wavelet for each required a scale and then multiplying the result by



[0041] In the discussion of the technology which follows herein, the "scalogram" may be taken to include all suitable forms of resealing including, but not limited to, the original unsealed wavelet representation, linear rescaling, any power of the modulus of the wavelet transform, or any other suitable resealing. In addition, for purposes of clarity and conciseness, the term "scalogram" shall be taken to mean the wavelet transform, T(a,b) itself, or any part thereof. For example, the real part of the wavelet transform, the imaginary part of the wavelet transform, the phase of the wavelet transform, any other suitable part of the wavelet transform, or any combination thereof is intended to be conveyed by the term "scalogram".

[0042] A scale, which may be interpreted as a representative temporal period, may be converted to a characteristic frequency of the wavelet function. The characteristic frequency associated with a wavelet of arbitrary a scale is given by:

where fc, the characteristic frequency of the mother wavelet (i.e., at a=1), becomes a scaling constant and f is the representative or characteristic frequency for the wavelet at arbitrary scale a.

[0043] Any suitable wavelet function may be used in connection with the present disclosure. One of the most commonly used complex wavelets, the Morlet wavelet, is defined as:

where f0 is the central frequency of the mother wavelet. The second term in the parenthesis is known as the correction term, as it corrects for the non-zero mean of the complex sinusoid within the Gaussian window. In practice, it becomes negligible for values of f0>>0 and can be ignored, in which case, the Morlet wavelet can be written in a simpler form as



[0044] This wavelet is a complex wave within a scaled Gaussian envelope. While both definitions of the Morlet wavelet are included herein, the function of equation (6) is not strictly a wavelet as it has a non-zero mean (i.e., the zero frequency term of its corresponding energy spectrum is non-zero). However, it will be recognized by those skilled in the art that equation (6) may be used in practice with f0>>0 with minimal error and is included (as well as other similar near wavelet functions) in the definition of a wavelet herein.

An Introduction to Image Forming



[0045] As described above, the present disclosure provides for the use of a wavelet transform for the parameterisation of the rate of occurrence, or density, of a certain feature (or features of similar characteristics) in an image. In one embodiment, an overview of an algorithm for performing an exemplary method of the disclosure is shown in Figure 1. Where the conventional use of wavelet analysis refers to time-varying waveforms and temporal frequencies, the skilled reader will understand that the same technique can be adapted to process images having spatial frequency characteristics.

[0046] The method illustrated in Figure 1 comprises the steps of loading an image 100, selecting a region for analysis 102, beginning a slice loop 104, mapping the image along slice coordinates 106, computing a wavelet scalogram of a slice map 108, mapping ridge features on the scalogram surface 110, and superimposing ridge features for each slice into a ridge composite mapping 112. Steps 104 - 112 are then looped N times, for N slices, following the end of the slice loop, i.e. when a counter has reached N. The method then proceeds with a step of generating a featured position (in spatial frequency) histogram from the composite mapping 114, thresholding the histogram to determine frequency of upper components 116 and outputting a feature count over the image selected region 118. It is understood by the skilled person that it is not necessary to arrange the superimposing step within the slice loop. For example, at the end of the slice loop, the ridge features are mapped on the scalogram surface of each image slice. With the mapped ridge features from each image slice, the superimposing step can be performed after each slice within the slice loop or in a separate loop to obtain a ridge composite mapping.

[0047] Figures 2 - 6 illustrate various aspects of these steps, applied for the specific example of computing a density estimation for retinal rod and cone features. This is achieved by tuning the wavelet characteristics to latch to the required features characteristic of rods and cones.

[0048] The tuning of a wavelet is achieved through the adjustment of its characteristic frequency, fc (which defines the shape that can be latched), and a selection of the scale, a (which defines, in effect, the feature size that can be latched). These parameters are explained above with reference to equation (4). Particular values of fc and a can be chosen based on known typical size and shape of rods and cones in a typical eye. It is to be appreciated that in general the wavelet can be selected to be tuned for other textural features such as for example, image artefacts arising from the physical structure of the image scanning arrays or a lid / lash classification or retinal pathology classification.

[0049] Figure 2 illustrates aspects of step 106; that of mapping an image along slice coordinates. An image 200 is shown from which a region 202 is selected for analysis. Here the region 202 is represented by a two dimensional pixel array. In one embodiment, each of the slices may comprise successive horizontal rows of the selected region 202. Figure 2 also shows a graph plotting the signal intensity on the y axis versus the horizontal position along the slice, on the x axis.

[0050] The mapping of the slices can be at any orientation so long as the slices are parallel to each other. Sliced positions are incremented by a pixel spacing of at least half the pixel width of the cone/rod features. This minimum spacing corresponds to the Nyquist frequency of the patterns of interest.

[0051] During or following on from the looping of the slices (repetition of steps 104 through 112 of Figure 1), the data from successive slices can be superimposed and can be plotted on the same graph.

[0052] Figure 3 illustrates aspects of step 108 illustrated in Figure 1, that of the computation of the wavelet scalogram of the slice map. Plot 300 is a representation of the mapping of the slices illustrated at the top of the region 202 of image 200. The slices (each showing intensity against x position) are shown plotted on an N axis, where N is the integer count number of the slice. Exemplary wavelet scalograms 302, 304 and 306 are shown for first, second and third slices illustrated in representative plot 300.Similar plots are made for each of the N slices present.

[0053] In these plots, the x-axis is pixel position (or location along the line of the mapping, commonly given the symbol b); the y-axis is scale (or frequency, commonly given the symbol a). The z-axis is a colour scale where different colours indicate the values of the wavelet transform coefficients (commonly given by notation: T(a,b)), and indicated on a scale running from blue (low values) to red (high values) colour mapping. T(a,b) is given by equation (1) above.

[0054] Figure 4 illustrates the ridges of the scalogram surfaces. In this diagram two ridges are plotted as black lines, representing the local maxima, at which points δT/δb = 0.

[0055] Figure 5 illustrates aspects of step 110 of the algorithm shown in Figure 1, that of mapping ridge features on the scalogram surfaces. An example scalogram 500 is shown, which is one from the N scalograms illustrated in Figure 3. In this step, the ridges (illustrated at 502 for example) of the scalogram are mapped and identified, as shown in the bottom half of the figure.

[0056] Figure 6 illustrates aspects of steps 112, 114, 116 illustrated in Figure 1, those of superimposition of ridge features for each slice into a ridge composite mapping, generation of a featured position histogram from the composite mapping and thresholding from the histogram to determine frequency of upper components. Figure 6 illustrates wavelet scalograms 600, 602, 604 with mapped wavelet ridges (represented by 606, 608, 610 for example). The ridges can be tagged and classified according to a band index (a ridge ID), and their mean positions in time and frequency.

[0057] The ridges are then extracted from the scalograms 600, 602, 604 forming ridge data 612, 614, 616 for each slice. The ridge data are then superimposed (it will be appreciated that the ridge data can be superimposed directly from the scalograms if preferred, there does not have to be a separate extraction step before the superimposition. In either case, the superimposition is performed with the ridge data only, rather than the entire scalogram data). A representation of the superimposed ridges is shown by plot 618 which plots the frequency of the features on ay axis against spatial position in an x axis. It can be seen from the diagram 618 that there is a features boundary at approximately 1.4 on the spatial frequency scale (vertical axis). The superimposed ridges can also be plotted as a histogram 620 which can be thresholded to determine the spatial frequency of the upper components. The frequency of occurrence of the uppermost ridges cluster is used to estimate the spatial density of coherent features that are resolvable just below a base noise level.

[0058] These features will correspond to the features for which the wavelet characteristics have been tuned. In one example, the wavelet characteristics can be tuned for the correlation of rod and cone textural features. The characteristic frequency, fc, of the wavelet is tuned to give the maximum resonance (i.e. matched in shape) and the scale range, a, is set so that the size of the features (or size range) is correctly matched.

[0059] Figure 7 illustrates further aspects of the method shown in Figure 1, the outputting of a feature count over the image selected region. Figure 7A shows the image 200 with selected region 202. In this instance, all available slices 700 from the region 202 have now been processed. Figure 7B shows a selected image segment 702. Figure 7C illustrates a multi slice view, showing the intensity versus spatial position for each of the N slices. Figure 7D illustrates the composite scalogram 706 and scalogram profile 708. Figure 7E illustrates a specific slice 710 and its scalogram 712. Figure 7F illustrates the superimposed ridge sets 714 and the histogram 716 of the superposition. The boundary of feature counts illustrated at 718 is used to deduce the feature density and therefore to provide a count of rod/cone shape occurrences in the image.

[0060] Figures 8 through 13 illustrate an alternative embodiment, wherein an image is mapped in radial slices rather than horizontal slices. As shown in Figure 8, an image 800 can comprise an annular region of interest 802 which is mapped in a number of different radial slices 804.

[0061] The slices 804 are superimposed to form an orthogonal texture image 900 as shown in Figure 9 with equal sample intervals. The diagram represents the concentric circle mappings of Figure 8 represented as a surface plot, with the y axis representing the radius from the centre of each concentric mapping circle and x axis representing the angular position of mapped points on each of the circles. Oversampling can be used to make sure the 'length' of each map is the same and does not increase with the increasing circumference as the mapping moves away from the centre.

[0062] The information is then condensed into the frequency domain in the form of overlaid scalograms for each slice of the orthogonal texture image. A plot of the scalogram composite 1000 is shown in Figure 10. Entropy is then mapped across the scale, as illustrated in plot 1100 in Figure 11. The peak 1102 in wavelet scale entropy is then used to parameterise the texture of the image, i.e. the retinal texture in this example.

[0063] Figure 12 illustrates a plot 1200 illustrating a further embodiment of a wavelet based surface texture mapping for the retinal region. The regions mapped are shown by the white lines and represent a subset of the radial mappings shown in Figure 8, defining a limited angular region. From this plot the extent of retinal texture partitioned from lids and lashes textures can be derived, as shown in the plot 1300 of Figure 13 which illustrates the scalogram entropy peak energy on the y axis versus the scalogram entropy peak position (on a frequency scale and so measured in Hz) on the x axis. If a composite scalogram for each mapped region is taken and the scalogram entropy peak (along scale bands) is found, and then the peak value is plotted against peak position, a set of points is obtained which can discriminate between lashes regions (points 1302 above the dotted line) and retinal texture (points 1304 below the dotted line). This could be used as an automated method of finding the extent of useful retinal texture in an image.

[0064] Figure 14 shows an alternative example where the general methods of the algorithm shown in Figure 1 are used to classify a scanner induced imaging artefact. One of the example instruments in which the present disclosure may be used is a scanning laser ophthalmoscope of the type disclosed in EP0730428. Another example instrument in which the present disclosure may be used is a scanning laser ophthalmoscope of the type disclosed in type is WO 2008/009877. In one example implementation of this type of arrangement, an image scan is carried out by a linear scanner having sixteen individual mirrors (facets) in a polygon arrangement. A reflection along each of these mirrors represents one vertical scan line in the image. A vertical block of sixteen lines represents a single turn of the polygon mirror. If a reflection timing error occurs due to the fact that the inter-mirror angles are not sufficiently accurate, then the vertical position of the line associated with that facet will be displaced in a vertical direction (relative to lines on either side). Timing anomalies from more than one mirror can also occur. This will result in a pattern of timing anomalies in the image (i.e. a displacement pattern of vertical positioning errors) that repeats every sixteen lines across the image. The phenomena is often referred to as 'scanning jitter' (a sixteen-pixel jitter in this particular case). The method of the present disclosure can be used to measure the occurrence of and the intensity of these classes of anomalies.

[0065] Figure 14 shows how these effects can be identified, so that they can be subsequently measured. The diagram shows an image 1400, a selected region 1402 and a composite scalogram 1406 of that region. Plot 1408 then illustrates the energy, summed along each scale band (i.e. summed along rows in the scalogram plot 1406) plotted on they axis versus pixel position on the x axis. The measurement 1410 of the peak 1412 with respect to the linear gradient 1414 gives the pixel jitter metric. The linear ramp could be subtracted from the scalogram summed energy signal give a result that looks more like a conventional spectrum plot. The peak in textural response at the pixel spacing of sixteen in this example is indicative of the degree of sixteen-pixel jitter or timing anomalies information contained within an image.

[0066] The various features and methods discussed herein enable automated characterisation of "textural" image features. In general, the algorithm will condense spatial information to provide a compact quantitative description of textural information. The algorithm may have specific application in for example lids and lashes classification and the counting of features in a common class. An example of this could be quantifying (sixteen) pixel jitter in images and the "counting" of rods and cone features in a retinal image. In this latter example the measurement could be relevant as an early indicator to conditions such as Alzheimer's dementia.

[0067] The measurement of pixel distortion components or scanner artefacts in a retinal image could also be used for example as a fail/pass criterion in polygon manufacture, where the polygon is of the type disclosed in the European Patent mentioned above.

[0068] The use of ridge and scalogram supposition ameliorates noise performance of the measurements made.

[0069] Various improvements and modifications may be made to the above without departing from the scope of the disclosure.

[0070] For example, "rods and cones" are generally mentioned together in the present disclosure, and are generally of similar size so can be detected as taught above. However it may be possible to further "fine tune" the methods and apparatus of the disclosure to treat rods and cones independently, as the sizing of rods and cones may in some cases be different from each other.

[0071] Furthermore, the use of the continuous wavelet transform will in general be preferred, and is described above. However it is possible to use non-continuous wavelet transforms. The above and other variations can be made by the skilled person without departing from the scope of the invention as defined in the appended claims.


Claims

1. A computer-implemented method of image processing comprising:

(i) selecting (106) a one-dimensional slice (204, 804) from an image;

(ii) computing (108) a wavelet scalogram of the slice;

(iii) identifying (110) ridge features within the wavelet scalogram;
repeating steps (i), (ii) and (iii) for one or more further slices;
concatenating (112) the identified ridge features from the slices such that they are representable in a two-dimensional ridge composite mapping (618); and
determining (116) a spatial frequency of the uppermost components of said concatenated identified ridge features by thresholding a histogram representation (620) of the ridge composite mapping (618).


 
2. The method of claim 1, wherein said step of selecting (106) a one-dimensional slice from the image comprises selecting a row or column of image data from the image.
 
3. The method of claim 1, wherein said step of selecting (106) a one-dimensional slice from the image comprises mapping a path along a set of image data pixels to a straight line representation.
 
4. The method of claim 3, wherein said path along a set of image data pixels comprises:

a straight line which is angled away from the horizontal or vertical axes defined by the rows and columns of the image data pixels; or

at least a part of a circular path.


 
5. The method of any of claims 2 to 4, wherein said selected row or column or said path is chosen to have a directionality which corresponds to a directionality of a characteristic of the image which is expected or is being searched for.
 
6. The method of claim 5, wherein said characteristic of the image is a feature associated with a known pathology or medical condition.
 
7. The method of any preceding claim, wherein the step of computing (108) a wavelet scalogram of the slice comprises application of a continuous wavelet transform.
 
8. The method of claim 7, wherein the step of computing (108) a wavelet scalogram comprises selecting a characteristic frequency of the applied wavelet transform to match a chosen shape.
 
9. The method of any preceding claim, wherein the step of computing (108) a wavelet scalogram comprises:

selecting a scale of an applied wavelet transform to match a chosen feature size; or

selecting a characteristic frequency of an applied wavelet transform to match a chosen shape.


 
10. The method of any of claims 1 to 8, wherein the step of computing (108) a wavelet scalogram comprises selecting a scale and a characteristic frequency a of an applied wavelet transform to match the size and shape of rods and/or cones in a retina.
 
11. An image processing apparatus configured to perform the methods of any of claims 1 to 10.
 
12. A computer program product encoded with instructions that when run on a computer, cause the computer to perform the method of any of claims 1 to 10.
 


Ansprüche

1. Ein computerimplementiertes Verfahren zur Bildverarbeitung, umfassend:

(i) Auswählen (106) einer eindimensionalen Schicht (204, 804) aus einem Bild;

(ii) Berechnen (108) eines Wavelet-Skalogramms der Schicht;

(iii) Identifizieren (110) von Kammmerkmalen innerhalb des Wavelet-Skalogramms;
Wiederholen der Schritte (i), (ii) und (iii) für eine oder mehrere weitere Schichten;
Verketten (112) der identifizierten Kammmerkmale aus den Schichten, so dass sie in einer zweidimensionalen zusammengesetzten Kammkartierung (618) darstellbar sind; und
Bestimmen (116) einer räumlichen Häufigkeit der obersten Komponenten der genannten verketteten identifizierten Kammmerkmale durch Schwellenwertbestimmung einer Histogrammdarstellung (620) der zusammengesetzten Kammkartierung (618).


 
2. Das Verfahren nach Anspruch 1, wobei der Schritt des Auswählens (106) einer eindimensionalen Schicht aus dem Bild das Auswählen einer Zeile oder Spalte von Bilddaten aus dem Bild umfasst.
 
3. Das Verfahren nach Anspruch 1, wobei der Schritt des Auswählens (106) einer eindimensionalen Schicht aus dem Bild das Abbilden eines Pfades entlang eines Satzes von Bilddatenpixeln auf eine geradlinige Darstellung umfasst.
 
4. Das Verfahren nach Anspruch 3, wobei der Pfad entlang eines Satzes von Bilddatenpixeln umfasst:

eine gerade Linie, die von den durch die Zeilen und Spalten der Bilddatenpixel definierten horizontalen oder vertikalen Achsen weggewinkelt ist; oder

zumindest einen Teil einer Kreisbahn.


 
5. Das Verfahren nach einem der Ansprüche 2 bis 4, wobei die ausgewählte Zeile oder Spalte oder der Pfad so gewählt wird, dass er eine Richtungsabhängigkeit aufweist, die einer Richtungsabhängigkeit einer Eigenschaft des Bildes entspricht, die erwartet wird oder nach der gesucht wird.
 
6. Das Verfahren nach Anspruch 5, wobei das Merkmal des Bildes ein Merkmal ist, das mit einer bekannten Pathologie oder einem bekannten medizinischen Zustand verbunden ist.
 
7. Das Verfahren nach einem der vorhergehenden Ansprüche, wobei der Schritt des Berechnens (108) eines Wavelet-Skalogramms der Schicht die Anwendung einer kontinuierlichen Wavelet-Transformation umfasst.
 
8. Das Verfahren nach Anspruch 7, wobei der Schritt des Berechnens (108) eines Wavelet-Skalogramms die Auswahl einer charakteristischen Frequenz der angewandten Wavelet-Transformation zur Anpassung an eine gewählte Form umfasst.
 
9. Das Verfahren nach einem der vorhergehenden Ansprüche, wobei der Schritt des Berechnens (108) eines Wavelet-Skalogramms umfasst:

Auswählen einer Skala einer angewandten Wavelet-Transformation, um einer gewählten Merkmalsgröße zu entsprechen; oder

Auswählen einer charakteristischen Frequenz einer angewandten Wavelet-Transformation zur Anpassung an eine gewählte Form.


 
10. Das Verfahren nach einem der Ansprüche 1 bis 8, wobei der Schritt des Berechnens (108) eines Wavelet-Skalogramms das Auswählen einer Skala und einer charakteristischen Frequenz a einer angewandten Wavelet-Transformation umfasst, um die Größe und Form von Stäbchen und/oder Zapfen in einer Netzhaut anzupassen.
 
11. Eine Bildverarbeitungsvorrichtung konfiguriert, die Verfahren nach einem der Ansprüche 1 bis 10 durchzuführen.
 
12. Ein Computerprogrammprodukt, das mit Anweisungen kodiert ist, die, wenn sie auf einem Computer ausgeführt werden, den Computer veranlassen, das Verfahren nach einem der Ansprüche 1 bis 10 auszuführen.
 


Revendications

1. Procédé mis en œuvre par ordinateur de traitement d'images comprenant :

(i) la sélection (106) d'une tranche unidimensionnelle (204, 804) à partir d'une image ;

(ii) le calcul (108) d'un scalogramme d'ondelettes de la tranche ;

(iii) l'identification (110) de caractéristiques de crête dans le scalogramme d'ondelettes ;
la répétition des étapes (i), (ii) et (iii) pour une ou plusieurs autres tranches ;
la concaténation (112) des caractéristiques de crête identifiées à partir des tranches de manière qu'elles soient représentables dans un mappage composite de crêtes en deux dimensions (618) ; et
la détermination (116) d'une fréquence spatiale des composantes les plus hautes desdites caractéristiques de crête identifiées concaténées par seuillage d'une représentation d'histogramme (620) du mappage composite de crêtes (618).


 
2. Procédé selon la revendication 1, dans lequel ladite étape de sélection (106) d'une tranche unidimensionnelle à partir de l'image comprend la sélection d'une ligne ou colonne de données d'image à partir de l'image.
 
3. Procédé selon la revendication 1, dans lequel ladite étape de sélection (106) d'une tranche unidimensionnelle à partir de l'image comprend le mappage d'une trajectoire suivant un ensemble de pixels de données d'image à une représentation par ligne droite.
 
4. Procédé selon la revendication 3, dans lequel ladite trajectoire suivant un ensemble de pixels de données d'image comprend :

une ligne droite qui est inclinée en s'éloignant des axes horizontaux ou verticaux définis par les lignes et les colonnes des pixels de données d'image ; ou

au moins une partie d'une trajectoire circulaire.


 
5. Procédé selon l'une quelconque des revendications 2 à 4, dans lequel ladite ligne ou colonne sélectionnée ou ladite trajectoire est choisie pour avoir une directionnalité qui correspond à une directionnalité d'une caractéristique de l'image qui est prévue ou est actuellement recherchée.
 
6. Procédé selon la revendication 5, dans lequel ladite caractéristique de l'image est une caractéristique associée à une pathologie ou condition médicale connue.
 
7. Procédé selon une quelconque revendication précédente, dans lequel l'étape de calcul (108) d'un scalogramme d'ondelettes de la tranche comprend l'application d'une transformée en ondelettes continue.
 
8. Procédé selon la revendication 7, dans lequel l'étape de calcul (108) d'un scalogramme d'ondelettes comprend la sélection d'une fréquence caractéristique de la transformée en ondelettes appliquée pour correspondre à une forme choisie.
 
9. Procédé selon une quelconque revendication précédente, dans lequel l'étape de calcul (108) d'un scalogramme d'ondelettes comprend :

la sélection d'une échelle d'une transformée en ondelettes appliquée pour correspondre à une taille de caractéristique choisie ; ou

la sélection d'une fréquence caractéristique d'une transformée en ondelettes appliquée pour correspondre à une forme choisie.


 
10. Procédé selon l'une quelconque des revendications 1 à 8, dans lequel l'étape de calcul (108) d'un scalogramme d'ondelettes comprend la sélection d'une échelle et d'une fréquence caractéristique d'une transformée en ondelettes appliquée pour correspondre à la taille et la forme de tiges et/ou de cônes dans une rétine.
 
11. Appareil de traitement d'images configuré pour effectuer les procédés selon l'une quelconque des revendications 1 à 10.
 
12. Produit de programme informatique codé avec des instructions qui, quand elles sont exécutées sur un ordinateur, amènent l'ordinateur à effectuer le procédé selon l'une quelconque des revendications 1 à 10.
 




Drawing





















































Cited references

REFERENCES CITED IN THE DESCRIPTION



This list of references cited by the applicant is for the reader's convenience only. It does not form part of the European patent document. Even though great care has been taken in compiling the references, errors or omissions cannot be excluded and the EPO disclaims all liability in this regard.

Patent documents cited in the description




Non-patent literature cited in the description