Laser scanning system and method

09788717 · 2017-10-17

Assignee

Inventors

Cpc classification

International classification

Abstract

A method of reducing jitter error in a laser scanning system adapted to produce a scanned image including a number of lines of an object. The method includes the steps of providing a reference object arranged such that the scanned image produced by the laser scanning system includes a reference image of the reference object, processing the reference image to calculate an error arising from non-repeatable displacement of the lines of the reference image, and adjusting at least one operating parameter of the laser scanning system in response to the calculated error.

Claims

1. A method of reducing jitter error in a laser scanning system adapted to produce a scanned image of an object, the scan comprising a number of lines, the method comprising the steps of: providing a reference object such that the scanned image produced by the laser scanning system includes a reference image of the reference object; processing the reference image to calculate an error arising in the laser scanning system; and adjusting at least one operating parameter of the laser scanning system in response to the calculated error; wherein the reference object is periodic in the direction of the line of scan of the laser scanning system.

2. The method according to claim 1 wherein the reference image comprises a number of pixels.

3. The method according to claim 1, wherein the reference object includes a series of alternate black and white bands.

4. The method according to claim 1, wherein the reference object includes a series of alternate black and white bands, and wherein the alternate black and white bands are equally spaced.

5. The method according to claim 4, wherein the width of each band is between 0.25 mm and 0.75 mm.

6. The method according to claim 1, wherein the laser scanning system forms part of a scanning laser ophthalmoscope.

7. The method according to claim 6, wherein the method is carried out as a calibration before imaging a patient's eye.

8. The method according to claim 7, wherein the reference object appears in an image of the patient's eye.

Description

(1) Embodiments of the present invention will now be described, by way of example only, with reference to the accompanying drawings, in which:—

(2) FIG. 1 is a schematic view of the scanning and reflecting elements of a laser scanning ophthalmoscope;

(3) FIG. 2 is a first partial schematic view of the elements of FIG. 1;

(4) FIG. 3 is a second partial schematic view of the elements of FIG. 1;

(5) FIG. 4 is a flow chart detailing the method of reducing jitter error in the laser scanning system in accordance with the present invention;

(6) FIG. 5 is an example of a scanned image produced by the laser scanning system of FIG. 1;

(7) FIG. 6a is a graph illustrating two column data signals from the reference image of FIG. 5;

(8) FIG. 6b is a graph illustrating the two column signals of FIG. 6a after an approximate shift has been applied to the second signal;

(9) FIG. 7 is a pair-wise distance matrix between a reference column and a second column of the reference image;

(10) FIG. 8 is a schematic representation of the approximated column shift error in the reference image;

(11) FIG. 9 is an example of a corrected scanned image produced by the laser scanning system of FIG. 1; and

(12) FIG. 10 is a schematic view of an alternative embodiment of the reflecting and scanning elements of the laser scanning ophthalmoscope of FIG. 1.

(13) FIG. 1 is a schematic view of the scanning and reflecting elements of the laser scanning system 1, as applied to a scanning laser ophthalmoscope (SLO). As illustrated, the main components are a high speed rotating polygon mirror 10 (an example of a first scanning element), a slow speed mirror 12 (an example of a second scanning element), a first ellipsoidal mirror (slit mirror) 14 (an example of scan relay means) and a second ellipsoidal mirror (main mirror) 16 (an example of a scan transfer means).

(14) The slow speed mirror 12 may be an oscillating plane mirror, such as a galvanometer mirror.

(15) The rotating polygon mirror 10 is located at the first focal point of the slit mirror 14 and the slow speed mirror 12 is located at the second focal point of the slit mirror 14. The slow speed mirror 12 is also located at the first focal point of the main mirror 16. A patient's eye 18 is located at the second focal point of the main mirror 16.

(16) A source of collimated light 20a produces a laser light beam 20. The laser light beam 20 is then reflected from the polygon 10 to the patients eye 18 via the slit mirror 14, the slow speed mirror 12 and the main mirror 16.

(17) The slit mirror 14 acts as an optical amplifier, as it amplifies the scanning aperture of the beam 20. For example, with reference to FIG. 2, if the polygon 10 has 16 facets, it will rotate 22.5 degrees between facets and the resulting scanning aperture cc of the beam 20 on the slit mirror 14 will be 45 degrees. After the beam 20 is scanned across the slit mirror 14, the scanning aperture β is 120 degrees.

(18) With reference to FIG. 3, the beam 20 is then scanned across the main mirror 16 by the slow speed mirror 12 and reflected into the patient's eye 18. The source of collimated light 20a, the polygon 10 and the slow speed mirror 12 combine to provide a two-dimensional collimated light scan. The laser scanning system 1 therefore produces a scanned image comprising a number of lines of an object, e.g. a human retina.

(19) The source of collimated light 20a, the polygon 10 and the slow speed mirror 12 may also combine to provide a two-dimensional collimated light scan from an apparent point source, and the apparent point source may be provided at the first focus point of the main mirror 16, such that the main mirror 16 transfers the two-dimensional collimated light scan from the apparent point source into the patient's eye 18.

(20) The main mirror 16 may comprise an elliptical mirror, an aspherical mirror, an ellipsoidal mirror, a pair of parabola mirrors or a pair of paraboloidal mirrors.

(21) The slit mirror 14 may comprise an elliptical mirror, an aspherical mirror, an ellipsoidal mirror, a pair of parabola mirrors or a pair of paraboloidal mirrors.

(22) This type of laser scanning ophthalmoscope is described in the Applicant's European Patent Nos. 0730428 and 07733214.6.

(23) FIG. 1 also illustrates a reference object 22. In the embodiment described and illustrated here the reference object 22 is a striped target comprising a series of alternate black and white bands. The bands are equally spaced with the width of each band being between approximately 0.25 mm and 0.75 mm. Ideally the width of each band is approximately 0.5 mm.

(24) The reference object 22 is placed in the laser scanning system 1 before the system is used to obtain an image of the patient's retina.

(25) The reference image is positioned between the polygon 10 and the slit mirror 14 such that is in the line of scan of the laser scanning system 1, i.e. the laser beam 20 moves across the reference object 22 during operation of the laser scanning system 1. The reference object 22 is therefore periodic in the direction of the line of scan of the laser scanning system 1.

(26) With the reference object 22 positioned in the optical path after the polygon 10, i.e. between the polygon 10 and the slow speed mirror 12, the laser beam 20 repeatedly moves across the same portion of the reference object 22 to produce the reference image 24 (see FIG. 5 below). That is, the reference image 24 is constructed of a plurality of repeat image scans of the same portion of the reference object 22.

(27) As illustrated in FIG. 5, the reference image 24 comprises a plurality of columns 25 of alternate black and white bands. The reference image 24 may comprise a number of pixels. Each column (scan line) 25 contains a jitter error that manifests itself as a vertical shift. That is, each column 25 of the reference image 24 contains a displacement along the scan direction of the polygon 10. This displacement is quasi-identical in columns occurring from the same polygon facet. That is, the displacement values are quasi-periodic. Small differences of less than % pixel may result from noise in the system which prevents the displacements occurring from the same polygon facet from being exactly identical.

(28) Typically, the reference image 24 comprises 3900 columns (and 3072 lines). However, it should be appreciated that the scanned reference image 24 may comprise any suitable number of columns 25.

(29) The method for reducing the jitter error in the laser scanning system 1 will now be described with reference to FIGS. 4 to 8.

(30) With reference to FIG. 4, in step 1 the reference image 24 is captured by the laser scanning system 1. As described above, the laser beam 20 repeatedly moves across the striped target of the reference object 22 to produce the reference image 24 of FIG. 5. In obtaining the reference image 24 the laser scanning system 1 uses the green channel (i.e. wavelength of approximately 510 nm).

(31) The reference image 24 may be a pixelated image with m rows of n columns (i.e. m×n). (For example, m=3072, n=3900).

(32) Steps 2 to 6 describe how the reference image 24 is processed to calculate the error arising from non-repeatable displacement of the columns 25 of the reference image 24.

(33) In step 2 a reference column 26 from the columns 25 of the reference image 24 is determined. The reference column 26 is determined using the formula n/2, where n is the total number of columns 25.

(34) In step 3 the similarity between the reference column 26 and each of the other columns 25 is determined. The similarity between the reference column 26 and each of the other columns 25 is determined by using a dynamic time warping (DTW) algorithm or a derivative dynamic time warping (DDTW) algorithm. DTW and DDTW are known techniques for efficiently finding an alignment between two signals (or sequences). DDTW and DTW compare vectors locally, looking at local trends and/or variations between vectors.

(35) In the embodiment described and illustrated here a DDTW algorithm is used to find alignment between the reference column signal and every other column signal. However, it should be appreciated that DTW may alternatively be used.

(36) The warping between the reference column 26 and each of the other columns 25 is determined by firstly generating data signals (vector sequences) for each column 25 (including the reference column 26) which are representative of the image information of the column. The image information of the column may include the light intensity, or brightness, of the image and the data signals may include values which are representative of this intensity, or brightness. Where the reference image 24 is a pixelated image, the light intensity, or brightness, may include the light intensity or brightness of each pixel in the reference image.

(37) The data signals may be represented as vectors:
R=r.sub.1,r.sub.2, . . . ,r.sub.i, . . . ,r.sub.n,
C=c.sub.1,c.sub.2, . . . ,c.sub.j, . . . ,c.sub.m, where R is the reference column signal C is one of the other column signals n and m are the lengths of the signals. (In the present embodiment m=n, as all columns have the same length=3072.)

(38) The data signals are illustrated in FIG. 6a. As illustrated in FIG. 6a, the two signals R and C are out of alignment with one another. As explained above, this is a result of the jitter error.

(39) To find the alignment between the two signals using DDTW an n by m matrix is constructed where the (i.sup.th,j.sup.th) element of the matrix contains the pair-wise distance d(r.sub.i,c.sub.j) between the two points r.sub.i and c.sub.j. Using DDTW the distance d(r.sub.i,c.sub.j) is not Euclidean but rather the square of the difference of the estimated derivatives of r.sub.i and c.sub.j. The distance matrix is therefore:
D[i,j]=((r.sub.i)′−(c.sub.j)′)((r.sub.i)′−(c.sub.j)′) with 1≦i≦m and 1≦j≦m, where (r.sub.i)′, (r.sub.j)′, (c.sub.i)′ and (c.sub.j)′ are the derivatives of R and C at points i and j.

(40) Each matrix element (i,j) corresponds to the alignment between the points r.sub.i and c.sub.j. The distance matrix D[i,j] between the data signals R and C is illustrated in FIG. 7.

(41) Once the matrix is constructed, the DDTW algorithm finds the alignment path which runs through the “low cost” areas, or “valleys”, on the matrix. This alignment path is an example of a function representative of the minimal distance path of the matrix. The minimal distance path may be termed the “alignment path”, the “optimal path”, the “minimal cost path”, the “warping path”, or the “warping function”.

(42) The alignment path is a contiguous (in the sense stated below) set of matrix column and row indexes (or coordinates) that defines a mapping between R and C.

(43) The alignment path constructed by the DDTW algorithm may be expressed as:
P=(p.sub.1,p.sub.2, . . . ,p.sub.k, . . . p.sub.K) with max(m,n)≦K(m+n−1, where the k.sup.th element of P is defined as p.sub.k=(px,py).sub.k.

(44) The alignment path P is subject to certain constraints:

(45) (a) Boundary conditions: p.sub.1=(1,1) and p.sub.K=(m,n). This requires the alignment path to start and finish in diagonally opposite corner cells of the matrix, i.e. the starting and end points of the alignment path must be the first and last points of the sequences.
(b) Monotonicity: Given p.sub.k=(a,b), then p.sub.k-1=(a′,b′) where a−a′≧0 and b−b′≧0. This forces the points in the alignment path P to be monotonically spaced in time, i.e. the coordinates of the alignment path must be increasing.
(c) Continuity: Given p.sub.k=(a,b), then p.sub.k-1=(a′,b′) where a−a′≦1 and b−b′≦1. This restricts the allowable steps in the alignment path P to adjacent cells (including diagonally adjacent cells), i.e. the path can only go one step up, or left, or both.

(46) There may be many alignment paths which satisfy the above criteria. However, it is necessary to determine the path which minimizes the warping cost, i.e. the optimal path.

(47) The optimal alignment path should minimize the warping cost:

(48) D D T W ( R , C ) = min ( .Math. k = 1 K p k K ) , where K in the denominator is used to compensate for the fact that the alignment paths may have different lengths.

(49) The alignment path can be found by evaluating the following recurrence which defines the cumulative distance γ(i,j) as the distance d(i,j) found in the current cell and the minimum of the cumulative distances of the adjacent elements:
γ(i,j)=d(r.sub.i,c.sub.i)+min{γ(i−1,j−1)γ(i−1,j),γ(i,j−1)}.

(50) The alignment path is illustrated in FIG. 7, referenced as 28. FIG. 7 also illustrates the reference column and other column data signals R and C.

(51) The pseudo code for determining the matrix and the alignment path are detailed in Annexes 1 and 2.

(52) FIG. 6b illustrates the mapping of signal C onto signal R. As illustrated, once the alignment path is known it can be used to produce a close mapping of signals.

(53) In step 4 the warping between the reference column data signals and the other column signals is approximated with a shift (or offset). The shift between signals corresponds to a column shift in the reference image 24. The column shift corresponds to the error arising from the non-repeatable displacement of the columns 25 of the reference image 24.

(54) The alignment path has the form:
P=(p.sub.1,p.sub.2, . . . ,p.sub.K) with p.sub.k=(p.sub.kx,p.sub.ky) are the coordinates of each point in the path.

(55) This is therefore the parametric representation of the curve of the alignment path.

(56) The parametric representation is then transformed to an explicit form:
p.sub.ky=f(p.sub.kx) with 1≦k≦K.

(57) The transformation is obtained by removing duplicates of the abscissas of p, obtaining their corresponding coordinates and then interpolating the results over the frame [1,m]. In the embodiment described here, m=n is the length of the two original signals R and C.

(58) The pseudo code for removing the abscissas of p and interpolating the results over the frame [1, m] are detailed in Annexes 3 and 4.

(59) The shift is obtained by comparing the implicit form of the alignment path with the diagonal reference path of the matrix.

(60) The reference path of the matrix has a parametric representation X=t, Y=t, t in [1,N] where N is the signal length. The implicit form of the reference path is Y=X. The reference path would be the same as the alignment path if the two compared signals were identical. That is, the shift is obtained by comparing the discrete form of the alignment path with an “ideal path”, i.e. the path that would exist if there were no warping between the two data signals.

(61) The alignment path 28 and the reference path 30 are illustrated in FIG. 7. The pseudo code for comparing the discrete form of the alignment path with the diagonal distance path of the matrix and obtaining the shift between each signal is detailed in Annex 5. The shift signal may be referenced as s (see FIG. 4).

(62) It should be noted that the above-described shift is determined between the reference column 26 and every other column 25. The result of this is that n shift values for each column 25 (n=3900).

(63) In step 5 the column shift for each facet of the polygon 10 is determined. Given that most of the jitter present in the reference image 24 is polygon facet dependent, the average of all the shifts of the columns 25 originating from the same polygon facet is averaged. In the case where the polygon 10 has 16 facets, these columns (or pixels) are separated by 16 columns (or pixels).

(64) The column shift for each facet of the polygon 10 may be determined by:

(65) Offset ( j ) = .Math. k , n - 16 k + j l k , n - 16 k + j n s _ [ ( n - 16 k ) + j ] _ , where Offset(j) is the column shift for each facet and k k is an index representing image columns n number of columns in the image (n=3900) j variable corresponding to each facet. Here 1≦j≦16

(66) The approximated column shift error in the reference image 24 is illustrated in FIG. 8. Columns which originate from the same facet of the polygon 10 are referenced with the number 32. The reference column 26 is also illustrated.

(67) FIG. 8 thus illustrates the averaged calculated error arising from non-repeatable displacement of the columns of the reference image 24. From FIG. 8 it can clearly be seen that the same column shift value has been assigned to each column 25 originating from the same facet of the polygon 10.

(68) In step 6 the column shift (or offset) determined in step 5 is approximated to 0.25 pixel accuracy. This is due to the fact that the sampling of the scan has 0.25 pixel rate.

(69) The approximated offset can be calculated by:
Offset(j)=0.25×(round(Offset(j))).

(70) Thus, for example, if the offset value was calculated as having a value of 1.6, it would be approximated to 0.25 pixel accuracy to 1.5.

(71) In step 7 a scanning time delay (capture time) of the polygon 10 of the laser scanning system 1 is determined based on the approximated offset value determined in step 6. This scanning time delay of the polygon of the laser scanning system 1 is an example of an operating parameter of the laser scanning system 1.

(72) In the embodiment described and illustrated here, each pixel of the reference image 24 is scanned at a frequency of 132 MHz, i.e. 7.57 ns.

(73) This yields a linear relationship between the time delay and the offset of step 6 as:
Time Delay(ns)=7.57×Offset(j).

(74) This time delay is applied to each facet of the polygon 10 of the laser scanning system 1. That is, the time delay, or capture time, for each facet of the polygon is determined by advancing or delaying a reference capture time by a factor dependent upon the average approximate signal shift for that specific facet of the polygon.

(75) The system may start imaging with the longest facet of the polygon 10. In order to know which shift corresponds to which facet of the polygon 10, the system 1 is forced to always start scanning using the longest facet. The facet data are recorded during an initial stage which imaging is started with facet 1, for example. The DDTW algorithm is run during this stage to obtain 16 shift values. s1 corresponding to facet 1, s2 corresponding to facet 2 etc. If a second image is captured, there is no way of knowing if the imaging will start with facet 1, and therefore no guarantee that the offsets will be applied correctly. To overcome this, imaging is started with the longest facet, as determined by the initial stage. Since the polygon always rotates in the same direction, and knowing the longest facet (or reference facet) and its offset, it is possible to attribute the remaining offsets correctly.

(76) Once the time delay values (operating parameters) of the polygon 10 of the laser scanning system 1 have been adjusted in response to the calculated errors, the reference object 22 is scanned again. An example of the reduced jitter error reference object image 24′ is illustrated in FIG. 9. It can clearly be seen that the jitter error has been reduced markedly from the original reference image 24 of FIG. 5.

(77) Once the start of scan time delay/advance values of the polygon 10 of the laser scanning system 1 have been adjusted in response to the calculated errors, the reference object 22 is removed from the system and the system is used to obtain an image of the patient's retina in the normal manner. With this adjustment, the image of the patient's retina has a greatly reduced jitter error.

(78) Steps 1 to 7 are preferably carried out by a computer, processing means, or the like. The steps 1 to 7 may be incorporated as part of the existing computer control software of the laser scanning system 1. The steps 1 to 7 may be provided on a separate computer, or the like, from the existing computer of the laser scanning system 1.

(79) Step 7 is carried out by a computer, processing means, or the like, which is used to control the scanning elements 10, 12 of the laser scanning system 1. The software of the computer of the laser scanning system 1 controls the start of scan time 12. The software is reprogrammed to adjust the operating parameters of the laser scanning system 1 in response to the calculated error to correct for the non-repeatable displacement of the lines of the reference image 24.

(80) FIG. 10 is an alternative embodiment of the laser scanning system 1. The second embodiment is similar to the first embodiment, with the exception that the reference object 22′ is located on the slit mirror 14.

(81) The reference object 22′ is placed on the edge of the slit mirror 14. Therefore, the scan laser beam 20 moves across the reference object 22′ and the slit mirror 14 during operation of the laser scanning system 1′.

(82) This means that the scanned image produced by the laser scanning system 1′ includes an image of the patient's eye 18 and a reference image of the reference object 22′. That is, the scanned image is a two dimensional image which includes a portion of the patient's eye 18 and a portion of the reference object 22′.

(83) In the embodiment described and illustrated here the reference object 22′ covers approximately 2.5% of the scanning aperture and is imaged at the beginning of every scan line.

(84) The method of reducing the jitter error in the laser scanning system 1′ is the same as that described above in respect of the first embodiment of the laser scanning system 1.

(85) The method of the present invention therefore obviates or mitigates the disadvantages of previous proposals by correcting for the jitter error introduced to the scanned image from the imperfect timing of the polygon facet position, the random noise introduced by the triggering system of the facet detector, the variations in the cut-depth of the facets of the polygon, the variations in the flatness of the facet of the polygon and the variations in the rotational speed of the polygon. The benefits of the method are increased jitter measurement accuracy, improved visual appearance of the scanned image, reduced facet detector accuracy requirements, reduced polygon building requirements, reduced polygon and facet detector failures in the field and improved polygon performance monitoring.

(86) The method of the present invention exploits the linear dependency between imaging signal delay and pixel shifts and incorporates the shifts as timing delays when capturing the image. The method does not require any image processing, which reduces the computational cost. Since the method may be performed as a calibration before imaging of the patient's eye, the time taken to obtain the image of the patient's eye is unchanged, as the shift errors have already been corrected for beforehand.

(87) Modifications and improvements can be made to the above without departing from the scope of the present invention. For example, although the laser scanning systems 1, 1′ have been illustrated above as comprising slit mirrors 14, it is not essential that the system includes a slit mirror. In this case the reference object 22 is placed in the optical path after the polygon 10 or located on the main mirror 16.

(88) Furthermore, it is also possible that the reference object 22 is positioned in the optical path after the slow speed mirror 12 and the polygon 10, i.e. the reference object 22 may be positioned in the system such that the polygon 10 and the slow speed mirror 12 combine to produce a two-dimensional raster scan light pattern across the reference object 22 to produce the reference image.

(89) Also, the step of generating data signals for each column 25 may not necessarily include all columns in the reference image 24. It is possible to only generate data signals for a portion of these columns, i.e. a portion of the reference image.

(90) Furthermore, although the method of reducing the jitter error in a laser scanning system has been described above with respect to the rotating polygon scanning element, it should be appreciated that the method could be applied to correcting for non-repeatable line displacement errors in the reference image caused by other scanning elements, e.g. oscillating plane mirrors (slow speed mirrors), or any other rotating and/or oscillating reflective elements used in the scanning system.

(91) TABLE-US-00001 Annex 1 Accumulated Cost Matrix(X;Y;D) //X and Y are input signals, D is the distance matrix defined as piece- //wise squared distances between X and Y in the case of DTW, or //piece- wise distances between the first derivatives of X and Y in the //case of DDTW ---------------------------------------------------------------------------------------------------  1: m ←[X] //m is the size of X  2: m ←[Y] //m is the size of Y. In our case X and Y have the same length  3: ddtw[ ] ← new [m ,m] //allocate space for ddtw:the Accumulated Cost Matrix  4: ddtw(0,0) ← 0 //initialize temporary variable  5: for i = 1; i ≦ m; i ++, do  6: ddtw(i; 1) ← ddtw(i−1; 1) + D(i; 1) //fill first column of ddtw  7: end for  8: for j = 1; j ≦ m; j ++, do  9: ddtw(1; j) ← ddtw(1; j−1) + D(1; j) //fill first column of ddtw 10: end for 11: for i = 1; i ≦ m; i ++, do 12: for j = 1; j ≦ m; j ++, do 13: ddtw(i; j) ← D(i; j)+min {ddtw(i − 1; j); ddtw(i; j − 1); ddtw(i − 1; j − 1)}//fill rest of matrix 14: end for 15: end for 16: return ddtw //return accumulated cost matrix ddtw ---------------------------------------------------------------------------------------------------

(92) TABLE-US-00002 Annex 2 Optimal Warping Path(ddtw) ---------------------------------------------------------------------------------------------------  1: path[ ] ← new array //allocate space for the warping path  2: i = rows(ddtw) //i is an index for ddtw rows  3: j = columns(ddtw) //j is index for ddtw columns //The following steps fill the path by a backtracking algorithm from point (m,m) to point (1,1)  4: while (i > 1) & (j > 1) do  5: if i == 1 then  6: j = j − 1  7: else if j == 1 then  8: i = i − 1  9: else 10:  if ddtw(i−1; j) == min {ddtw(i − 1; j); ddtw(i; j − 1); ddtw(i − 1; j − 1)} 11:  then 12: i = i − 1 13: else if ddtw(i; j−1) == min {ddtw(i − 1; j); ddtw(i; j − 1); ddtw(i − 1; j − 1)} 14: then 15: j = j − 1 16:  else 17: i = i − 1; j = j − 1 18:  end if 19: path:add((i; j)) 20: end if 21: end while 22: return path ---------------------------------------------------------------------------------------------------

(93) TABLE-US-00003 Annex 3 Remove Abscissa Duplicates(path) //the path obtained from the Annex 2 is in the form (xk,yk) 1≦k≦K, //we need to write this path in the form Yk=f(Xk) 1≦k≦K. To be able to do so, we need remove duplicates from xk 1≦k≦K. ---------------------------------------------------------------------------------------------------  1: new_path[ ] ← new array  2: K = length(pathx) //K is the length of the x-coordinates (abscissas) of the path  3: new_pathx(1)=pathx(1)  4: j=1 //initialize index of elements of newpath  5: for i = 1; i ≦ K−1; i ++, do  6: if pathx(i+1)!=pathx(i)  7: new_pathx(j+1)=pathx(i+1)  8: new_pathy(j+1)=pathy(i+1)  9: j=j+1 10: end if 11: end for ---------------------------------------------------------------------------------------------------

(94) TABLE-US-00004 Annex 4 Interpolate New Path(x[1...m],xi[1..Ni], yi[1..Ni],)) //the previous step Annex3 gave a modified path without repetition. // The new path doesn't necessarily cover all the range [1,m] where m is the size of the input signals. We therefore interpolate the path to cover each point in the range [1,m]. In this case x[1...m]=[1,2,3...m] The output will be the y-coordinate corresponding to each element of [1,2,3...m] ---------------------------------------------------------------------------------------------------  1: For i = 1; i ≦ m; i ++, do  2: if (x[i] ≦xi[1]) do  3: y[i] = yi[1];  4: end if  5: if (x[i] ≧xi[Ni]) do  6: y[i] = yi[Ni];  7: end if  8: j = 1;  9: while (j ≦m) 10: if (xi[j] >= x[i]) break; 11: j = j + 1; 12: end while 13: newPathInt[i] = yi[j−1] + (yi[j] − yi[j−1])*(x[i] − xi[j−1])/(xi[j]−xi[j−1]); 14: end for ---------------------------------------------------------------------------------------------------

(95) TABLE-US-00005 Annex 5 Distance To Diagonal(new_pathInt) //this algorithm approximates the warp with a shift. //the shift is the distance between the path and the diagonal (defined as y=x, or x=i, y=i 1≦i≦m) --------------------------------------------------------------------------------------------------- 1: Shift=0 2: for i = 1; i ≦ m; i ++, do 3: shift=shift+(new_path(i)−i) 4: end for 5: shift=shift/m ---------------------------------------------------------------------------------------------------