Abstract:
A method and a device for computer-aided motion estimation in at least two temporally successive digital images with pixels to which coding information is assigned are provided, the motion being estimated on the basis of the spatial distribution of feature points.

Description:
BACKGROUND  
       [0001]     One embodiment of the invention relates to a method and a device for computer-aided motion estimation in at least two temporally successive digital images, a computer-readable storage medium and a computer program element.  
         [0002]     Development in the field of mobile radio telephones and digital cameras, together with the widespread use of mobile radio telephones and the high popularity of digital cameras, has led to modern mobile radio telephones often having built-in digital cameras. In addition, services such as, for example, the multimedia message service (MMS) are provided which enable digital image communications to be transmitted and received using mobile radio telephones suitable for this.  
         [0003]     Typically, the components of mobile radio telephones which enable digital images to be recorded do not afford high performance compared with commercially available digital cameras.  
         [0004]     The reasons for this are for example that mobile radio telephones are intended to be cost-effective and small in size.  
         [0005]     The resolution of digital images that can be recorded by means of mobile radio telephones with a built-in digital camera is too low for some purposes.  
         [0006]     By way of example, it is possible, in principle, to use a mobile radio telephone with a built-in digital camera to photograph printed text and to send it to another mobile radio telephone user in the form of an image communication by means of a suitable service, for example the multimedia message service (MMS), but the resolution of the built-in digital camera is insufficient for this in the case of a present-day commercially available device in a medium price bracket.  
         [0007]     However, it is possible to generate, from a suitable sequence of digital images which in each case represent a scene from a respective recording position, a digital image of the scene which has a higher resolution than that of the digital images of the sequence of digital images.  
         [0008]     This possibility exists for example when the positions from which digital images of a sequence of digital images of the scene have been recorded differ in a suitable manner.  
         [0009]     The recording positions, that is, the positions from which the digital images of the sequence of digital images of the scene have been recorded, may differ in a suitable manner for example when the plurality of digital images has been generated by recording a plurality of digital images by means of a digital camera held manually over a printed text.  
         [0010]     In this case, the differences in the recording positions that are generated as a result of the slight movement of the digital camera that arises as a result of shaking of the hand typically suffice to enable the generation of a digital image of the scene with high resolution.  
         [0011]     However, this necessitates calculation of the differences in the recording positions.  
         [0012]     If a first digital image is recorded from a first recording position and a second digital image is recorded from a second recording position, an image content constituent, for example an object of the scene, is represented in the first digital image at a first image position and in a first form, which is taken to mean the geometrical form hereinafter, and is represented in the second digital image at a second image position and in a second form.  
         [0013]     The change in the recording position from the first recording position to the second recording position is reflected in the change in the first image position to the second image position and the first form to the second form.  
         [0014]     Therefore, a calculation of a recording position change which is necessary for generating a digital image having a higher resolution than that of the digital images of the sequence of digital images can be effected by calculating the change in the image position at which image content constituents are represented and the form in which image content constituents are represented.  
         [0015]     If an image content constituent is represented in a first image at a first (image) position and in a first form and is represented in a second image at a second position and in a second form, then a motion of the image content constituent or an image motion will be how this is referred to hereinafter.  
         [0016]     Not only is it possible for the position of the representation of an image content constituent to vary in successive images, but the representation may also be distorted or its size may change.  
         [0017]     Moreover, the representation of an image content constituent may change from one digital image of the sequence of digital images to another digital image of the sequence of digital images, for example the brightness of the representation may change.  
         [0018]     Only the temporal change in the image data can be utilized for determining the image motion. However, this temporal change is caused not just by the motion of objects in the vicinity observed and by the observer&#39;s own motion, but also by the possible deformation of objects and by changing illumination conditions in natural scenes.  
         [0019]     In addition, disturbances have to be taken into account, for example, vibration of the camera or noise in the processing hardware.  
         [0020]     Therefore, the pure image motion can only be obtained with knowledge of the additional influences or be estimated from assumptions about the latter.  
         [0021]     For the generation of a digital image having a higher resolution than that of the digital images of the sequence of digital images, it is very advantageous for the calculation of the motion of the image contents from one digital image of the sequence of digital images to another digital image of the sequence of digital images to be effected with subpixel accuracy.  
         [0022]     Subpixel accuracy is to be understood to mean that the motion is accurately calculated over a length shorter than the distance between two locally adjacent pixels of the digital images of the sequence of digital images.  
         [0023]     Hereinafter an image is always to be understood to mean a digital image.  
         [0024]     One conventional method for carrying out motion estimation with subpixel accuracy is the determination of the optical flow J. J. Gibson,  The Perception of the Visual World,  Boston, 1950.  
         [0025]     The optical flow relates to the image changes, that is, to the changes in the representation of image contents of an image of the sequence of digital images with respect to the temporally succeeding or preceding image of the sequence of digital images which arise from the motion of the objects and the observer&#39;s own motion. The image motions generated can be interpreted as velocity vectors which are attached to the pixels. The optical flow is understood to mean the vector field of these vectors. In order to determine the motion components, assumptions about the temporal change in the image values are usually made.  
         [0026]     I(x, y, t) designates the time-dependent, two-dimensional image. I(x, y, t) is an item of coding information which is assigned to the pixel at the location (x, y) of the image at the instant t.  
         [0027]     Coding information is to be understood hereinafter to mean an item of brightness information (luminance information) and/or an item of color information (chrominance information) which is assigned in each case to one pixel or a plurality of pixels.  
         [0028]     A sequence of digital images is expressed as a single, time-dependent image, that is, that the first image of the sequence of digital images corresponds to a first instant t 1 , the second image of the sequence of digital images corresponds to a second instant t 2  and so on.  
         [0029]     I(x, y, t 1 ) is thus for example the gray-scale value of an image at the location (x, y) of the image of the sequence of digital images which corresponds to the first instant t 1 ; for example, it was recorded at the first instant t 1 .  
         [0030]     The change for a pixel which the latter experiences in the time dt with rate (dx, dy) can be expressed by means of development into a Taylor series  
               I   ⁡     (       x   +     ⅆ   x       ,     y   +     ⅆ   y       ,     t   +     ⅆ   t         )       =       I   ⁡     (     x   ,   y   ,   t     )       +         ∂   I       ∂   x       ⁢     ⅆ   x       +         ∂   I       ∂   y       ⁢     ⅆ   y       +         ∂   I       ∂   t       ⁢     ⅆ   t       +   …             (   1   )             
 
 For the determination of the optical flow, the assumption is made that the image values remain constant along the motion direction. This is formulated by the equation 
 
 I ( x+dx, y+dy, t+dt )= I ( x, y, t )   (2) 
 
 from which the equation  
                     ∂   I       ∂   x       ⁢     ⅆ   x       +         ∂   I       ∂   y       ⁢     ⅆ   y       +         ∂   I       ∂   t       ⁢     ⅆ   t       +   …     =   0           (   3   )             
 
 follows, wherein as in equation (1), the three dots symbolize the terms having higher derivatives than the first partial derivatives of the function I. 
 
         [0031]     If equation (3) is divided by the expression dt and the terms having higher derivatives than the first partial derivatives of I are disregarded, this results in the equation  
                 ∂   I       ∂   t       =         -       ∂   I       ∂   x         ⁢       ⅆ   x       ⅆ   t         -         ∂   I       ∂   y       ⁢         ⅆ   y       ⅆ   t       .                 (   4   )             
 
         [0032]     Disregarding the higher derivatives leads to errors if the image motion is large in relation to the pixel grid.  
         [0033]     The vector  
       [         ⅆ   x       ⅆ   t       ,       ⅆ   y       ⅆ   t         ]       
 
 specifies the components of the optical vector field, and is usually designated by [u, v]. 
 
         [0034]     The following thus holds true for equation (4)  
                 ∂   I       ∂   t       =         -       ∂   I       ∂   x         ⁢   u     -         ∂   I       ∂   y       ⁢     v   .                 (   5   )             
 
         [0035]     This equation is deemed to be a fundamental equation of the optical flow.  
         [0036]     In order that u and v can be determined unambiguously, it is known to make further assumptions about the temporal change in the image data.  
         [0037]     In accordance with B. K. P Horn and B. G. Schunck,  Determining Optical Flow,  Artificial Intelligence, 1981, an additional assumption taken in this respect is that the optical flow is smooth.  
         [0038]     Both assumptions together lead to a minimization problem, which is formulated as follows:  
               (       u   ^     ,     v   ^       )     =         arg   ⁢           ⁢   min       u   ,   v       ⁢         ∫   ∫       x   ,   y       ⁡     [       (         ∂   I       ∂   t       +         ∂   I       ∂   x       ⁢   u     +         ∂   I       ∂   y       ⁢   v       )     +     λ   ⁡     (         (       ∂   u       ∂   x       )     2     +       (       ∂   u       ∂   y       )     2     +       (       ∂   v       ∂   x       )     2     +       (       ∂   v       ∂   y       )     2       )         ]       ⁢           ⁢     ⅆ   x     ⁢           ⁢     ⅆ   y               (   6   )             
 
 The first term of the integral corresponds to the fundamental equation of the optical flow (5) and the second term represents the smoothness condition in accordance with B. K. P Horn and B. G. Schunck,  Determining Optical Flow,  Artificial Intelligence, 1981. 
 
         [0039]     This clearly means that the first term has the effect that the vector field which solves the minimization problem given by equation (6) fulfils equation (5) as well as possible. The smoothness condition has the effect that the partial derivatives of the vector field which solves the minimization problem given by equation (6) with respect to the position variables x and y are as small as possible.  
         [0040]     The minimization problem given by equation (6) can be solved by means of a variation calculation approach.  
         [0041]     A system of linear equations is solved in this case, the number of unknowns in the system of linear equations being twice the number of pixels.  
         [0042]     In B. K. P Horn and B. G. Schunck,  Determining Optical Flow,  Artificial Intelligence, 1981, an iterative procedure in accordance with the so-called Gauss-Seidel method is proposed for solving the system of linear equations.  
         [0043]     In accordance with B. Lucas, T. Kanade,  An Iterative Image Registration Technique with an Application to Stero Vision,  7th International Joint Conf. on Artificial Intelligence (IJCAI), pp. 674-679, 1981, the condition that adjacent pixels must have the same motion vector is used as a second assumption for the determination of the optical flow.  
         [0044]     It can be deduced from equation (5) that this assumption must be fulfilled for at least two points.  
         [0045]     However, a small local neighborhood of a pixel is usually used.  
         [0046]     The determination of u, v can be formulated under this assumption as a least squares problem:  
                 ∑     x   ,   y       ⁢           ⁢       (         ∂   I       ∂   t       +         ∂   I       ∂   x       ⁢   u     +         ∂   I       ∂   y       ⁢   v       )     2       →     min   .             (   7   )             
 
         [0047]     This leads to the system of equations:  
                   ∑     x   ,   y       ⁢           ⁢         (       ∂   I       ∂   x       )     2     ⁢   u       +       ∑     x   ,   y       ⁢           ⁢       (         ∂   I       ∂   x       ⁢       ∂   I       ∂   y         )     ⁢   v         =     -       ∑     x   ,   y       ⁢           ⁢     (         ∂   I       ∂   t       ⁢       ∂   I       ∂   x         )                 (   8   )                     ∑     x   ,   y       ⁢           ⁢       (         ∂   I       ∂   x       ⁢       ∂   I       ∂   y         )     ⁢   u       +       ∑     x   ,   y       ⁢           ⁢         (       ∂   I       ∂   y       )     2     ⁢   v         =     -       ∑     x   ,   y       ⁢           ⁢     (         ∂   I       ∂   t       ⁢       ∂   I       ∂   y         )                 (   9   )             
 
 The sums in the equations (7), (8) and (9) proceed over all x, y from the local neighborhood of the pixel that is used. 
 
         [0048]     Through the evaluation of a local neighborhood, an optical flow vector is determined with subpixel accuracy in both of the methods explained above.  
         [0049]     However, the following problems occur in both methods: 
        no motion can be determined in homogeneous regions since the required smoothness or the group formation does not yield additional information.     in both methods, the local and temporal derivatives are approximated by discrete differences, which may lead to a low accuracy.     problems arise if the motion is large in relation to the sampling time of the images. It is then no longer possible to straightforwardly disregard the higher derivatives in the Taylor series development. In this case, so-called block matching methods that are based on a correlation analysis often even lead to better results. These methods are comparable, in principle, with the approaches in accordance with B. Lucas, T. Kanade,  An Iterative Image Registration Technique with an Application to Stero Vision,  7th International Joint Conf. on Artificial Intelligence (IJCAI), pp. 674-679, 1981.     the assessment of small local neighborhoods leads to a further problem for example for photographing text documents. Even if the intensity pattern is high in contrast within the local neighborhood, ambiguities can arise because the pattern is repeated in the vicinity of the local neighborhood. This occurs, for example, in the case of texts since there are no intensity differences between the letters and letters are formed from the same geometrical forms. The correlation analysis especially leads to errors here.     contrary to the assumption in accordance with B. K. P Horn and B. G. Schunck,  Determining Optical Flow,  Artificial Intelligence, 1981; B. Lucas, T. Kanade,  An Iterative Image Registration Technique with an Application to Stero Vision,  7th International Joint Conf. on Artificial Intelligence (IJCAI), pp. 674-679, 1981, it cannot be expected in general (if, for example, moving objects are present in the image) that the optical flow will proceed in locally constant or smooth fashion. Rather, it must be regarded as smooth in portions since discontinuities occur, for example, at object boundaries. These discontinuities have to be taken into account in the determination of the optical flow.        
 
         [0055]     Numerous studies are concerned with the problem of the optical flow taking account of discontinuities, obscurations, etc.  
         [0056]     For the above-described application in which a high resolution image is to be generated from a low resolution sequence of digital images generated by means of a digital camera, these approaches are not necessary, however, since, in the case of the above application, only a motion of the camera is present and the above assumptions are thereforee fulfilled to an approximation.  
         [0057]     Consequently, the use of these approaches for a method would lead to unnecessarily high complexity of the method and thus to a low efficiency of the method.  
         [0058]     The problem of unambiguously determining the optical flow in homogeneous image sections or along extended horizontal or vertical edges can be avoided by the optical flow not being implemented at all pixels but rather at points with significant image values (see, for example, J. Shi, C. Tomasi,  Good Features to Track,  IEEE Conf. on Computer Vision and Pattern Recognition (CVPR94), 1994).  
         [0059]     This has the effect that only a thinned optical flow is present. The problem of the approximation of the derivatives by discrete differences in the case of fast motions can be reduced by using image pyramids (see, for example, W. Enkelmann,  Investigations of multigrid algorithms for the estimation of optical flow fields in image sequences,  Computer Vision, Graphics and Image Processing, 150-177, 1988).  
         [0060]     In the case of the above methods, motion vectors are determined as individual pixels on the basis of the fundamental equation of the optical flow (5). A local vicinity is taken into account in the determination. The calculations of the motion vectors at the individual pixels are effected independently of one another.  
         [0061]     This enables the determination of different motions generated by different objects.  
         [0062]     It is furthermore known, under the presumption that the scene represented by the images of the sequence of digital images is static and the image motion is caused only by the observer, to determine a motion model for all pixels on the basis of the fundamental equation of the optical flow (5) from the coding information of the images.  
         [0063]     This is explained below.  
         [0064]     If u(x, y, t) and v(x, y, t) designate the motion at a pixel (x, y) at an instant t then the following holds true: 
 
 I   x ( x, y, t ) u ( x, y, t )+ I   y ( x, y, t ) v ( x, y, t )+ I   t ( x, y, t )=0   (10) 
 
 (see Equation (5)). 
 
         [0065]     I x (x, y, t), I y (x, y, t), I t (x, y, t) designate the partial derivatives of the function I(x, y, t) with respect to the variable x and the variable y and the variable z at the location (x, y) at the instant t.  
         [0066]     Various motion models can be used for u(x, y, t) and v(x, y, t) in order to model the motion sought in the image as well as possible.  
         [0067]     For an affine motion model, the following holds true for example: 
 
 u ( x, y, t )= a   0   x+a   1   y+a   2    (11) 
 
 v ( x, y, t )= a   3   x+a   4   y+a   5    (12) 
 
 The determination of u(x, y, t) and v(x, y, t) with the aid of equation (10) may be formulated, for example, as the minimization of a square error:  
                 a   ^     _     =         arg   ⁢           ⁢   min       a   _       ⁢           ⁢       ∑   y     ⁢           ⁢       ∑   x     ⁢           ⁢       (         I   x     ⁡     (         a   0     ⁢   x     +       a   1     ⁢   y     +     a   2       )       +       I   y     ⁡     (         a   3     ⁢   x     +       a   4     ⁢   y     +     a   5       )       +     I   t       )     2                   (   13   )             
 
         [0068]     In the solution of the minimization problem given by (13), the parameters a 0 , a 1 , a 2 , a 3 , a 4  and a 5  from equations (11) and (12) are determined in the form of a optimum parameter vector  
           a   ^     _     =       (           a   0             ⋮             a   5           )     .         
 
         [0069]     This method also leads to poor results in the case of large motions because ignoring the higher derivatives in the Taylor development leads to errors.  
         [0070]     Therefore, in accordance with the prior art, a hierarchical procedure is employed in the case of this method, too.  
         [0071]     Firstly, the motion is determined at a low resolution level since the size of the motion is also reduced by reducing the resolution. The resolution is then increased progressively up to the original resolution.  
         [0072]     Moreover, the quality of the motion estimation is improved by means of an iterative procedure.  
         [0073]      FIG. 11  illustrates a flow diagram  1100  of a method for parametric motion determination that is disclosed in Y. Altunbasak, R. M. Mersereau, A. J. Patti,  A Fast Parametric Motion Estimation Algorithm with Illumination and Lens Distortion Correction,  IEEE Transactions on Image Processing, 12(4), pp. 395-408, 2003, for example.  
         [0074]     For a plurality of images  1101 , a first loop  1102  is implemented over all resolution levels, that is, over all resolution stages.  
         [0075]     Within each iteration of the first loop  1102 , the image in the current resolution level is subjected to low-pass filtering in step  1103  and is subsequently subsampled in step  1104 .  
         [0076]     The local gradients, that is, clearly the image directions with the greatest rise in brightness, are determined in step  1105 .  
         [0077]     Afterward, a second loop  1106  is implemented within each iteration of the first loop  1102 .  
         [0078]     Within each pass through the second loop, firstly step  1107  involves calculating the temporal gradient, that is, clearly the change in brightness at a pixel from the image that was recorded at the instant t to the image that was recorded at the instant t+1.  
         [0079]     In step  1108 , within the first iteration of the second loop  1106  for the first resolution level, a first parameter vector  a   0  is calculated from the coding information I(x, y, t) of the image that was recorded at the instant t and the coding information I(x, y, t+1) of the image that was recorded at the instant t+1, for example by means of a least squares estimation as in the case of the method described above, which determines the parametric motion model.  
         [0080]     The quality of the current motion model which is determined by the currently calculated parameter vector is measured in step  1109 .  
         [0081]     If the quality has not improved, the current iteration of the second loop  1106  is ended.  
         [0082]     If the quality has improved, in step  1111 , from the coding information I(x, y, t+1) of the image that was recorded at the instant t+1, by means of the motion model which is determined by the currently calculated parameter vector, a compensated coding information item I 1 (x, y, t+1) of the image that was recorded at the instant t+1 is determined by compensation and the currently calculated parameter vector is accepted in step  1112 .  
         [0083]     The current iteration of the second loop  1106  is subsequently ended.  
         [0084]     In all subsequent iterations of the second loop  1106 , the procedure is analogous to the first iteration of the second loop  1106  for the first resolution level except that in each case instead of the coding information I 1 (x, y, t+1) the compensated coding information from the last iteration I 1 (x, y, t+1), I 2 (x, y, t+1), . . . is used in order to determine parameter vectors a 1 , a 2 , . . . .  
         [0085]     The loop  1106  is implemented until a predetermined termination criterion is satisfied, for example the least squares error lies below a predetermined limit.  
         [0086]     If the determination criterion is satisfied, the current iteration of the first loop  1102  is ended.  
         [0087]     If an iteration of the first loop  1102  has been implemented for each desired resolution level, then a parameter vector a is calculated from the calculated parameter vectors  a   0 ,  a   1 ,  a   2 , . . . , and the motion  1113  is deemed to have been determined.  
         [0088]     This method is known by the designation “parametric motion determination”.  
         [0089]     It is furthermore known to determine image motion by means of temporal tracking of objects.  
         [0090]     Numerous methods exist which presuppose explicit model knowledge about the object and for which a preceding step of object detection is necessary (see, for example, D. Noll, M. Werner, and W. von Seelen.,  Real - Time Vehicle Tracking and Classification,  Proceedings of the Intelligent Vehicles &#39;95, pp. 101-106, 1995).  
         [0091]     However, these methods are not suitable for applications such as the one described above, in which a high resolution image is to be generated from a low resolution sequence of digital images generated by means of a digital camera, since these methods have a severe limitation of the variation possibilities, that is, of the image alterations that can be ascertained.  
         [0092]     Another group of methods uses a contour of an object to be tracked. These methods are known by the key words “active contours” or “snakes”.  
         [0093]     These approaches are not suitable either for applications such as the one described above, in which a high resolution image is to be generated from a low resolution sequence of digital images generated by means of a digital camera, since in general no object contour is present.  
         [0094]     A further group of customary methods for object tracking uses a representation of objects by feature points and tracks these points over time, that is, over the sequence of digital images.  
         [0095]     The points are firstly tracked independently of one another.  
         [0096]     A motion model that enables the displacements of the individual points is subsequently ascertained.  
         [0097]     For the motion of the individual object points, it is possible to use methods for determining the optical flow. The disadvantages of the optical flow that have already been discussed thus occur here with the addition that the evaluation of homogeneous regions is avoided by the selection of feature points.  
         [0098]     It is also possible to determine a uniform motion for all object points.  
         [0099]     In contrast to the methods based on the optical flow, there is the problem here that the parameters of the motion model can no longer be determined directly by means of a system of linear equations, rather an optimization over the entire parameter range is necessary.  
         [0100]     In the case of the method by Werner et al. disclosed in Martin Werner,  Objektverfolgung und Objekterkennung mittels der partiellen Hausdorffdistanz,  [Object tracking and object identification by means of the partial Hausdorff distance], Faculty for Electrical Engineering, Bochum, Ruhr University, 1998, a motion model is determined by means of a minimization of the Hausdorff distance. This requires carrying out a minimization over all motion parameters, which leads to a considerable computational complexity.  
         [0101]     An alternative approach, described by Capel et al. in D. Capel, A. Zisserman,  Computer vision applied to super resolution,  Signal Processing Magazine, IEEE, May 2003, pages 75-86, Vol. 20, Issue: 3, ISSN: 1053-5888, consists in dividing the object features into small subsets. For each of these subsets, firstly a dedicated motion model is determined in which corresponding object features are sought at the instant t 1  and t 2 . Corresponding object features are determined by comparing the intensity patterns. By means of these corresponding points, a motion model can be determined directly by means of least squares approach. The model which permits the best assignment for all object features is ultimately selected from the motion models of the subsets. An assessment for the best assignment is for example the minimization of the sum of absolute image differences.  
         [0102]     In order to reduce the complexity for determining corresponding points, it is necessary in the case of this method, however, to determine a minimum number of subsets with a minimum number of feature points. Therefore, inaccuracies and ambiguities such as have already been described above with reference to methods based on the optical flow occur in this method.  
         [0103]     A further possibility for determining a motion model for an object representation by feature points is described in A. Techmer,  Contour - Based Motion Estimation and Object Tracking for Real - Time Applications,  IEEE International Conference on Image Processing (ICIP 2001), pp. 648-651, 2001. A contour-based determination of the image motion is presented here. The motion is calculated by a comparison of contour point positions and contour forms. The approach may additionally be extended to object tracking. The method is based solely on the evaluation of distances and thus on the evaluation of the geometrical form of objects. This makes the method less sensitive toward illumination or exposure changes in comparison with approaches that assess the intensity pattern. The determination of the motion model only requires a variation over the two translation components of the motion. The remaining parameters can be determined directly by means of a least squares estimation. This achieves a substantial reduction of the computational complexity in comparison with methods requiring a variation over all parameters of the motion model (see, for example, Martin Werner,  Objektverfolgung und Objekterkennung mittels der partiellen Hausdorffdistanz,  [Object tracking and object identification by means of the partial Hausdorff distance], Faculty for Electrical Engineering, Bochum, Ruhr University, 1998).  
         [0104]     These approaches have the disadvantage that motion models such as affine transformation models, for example, have to be determined which have a high number of degrees of freedom.  
         [0105]     William H. Press, et al.,  Numerical Recipes  in C, ISBN: 0-521-41508-5, Cambridge University Press discloses interpolating a function bicubically.  
       SUMMARY  
       [0106]     One embodiment of the invention ascertains the image motion in at least two temporally successive digital images efficiently and with high accuracy. One embodiment includes a method and a system for ascertaining the image motion of at least two temporally successive digital images, a computer-readable storage medium and a computer program element.  
         [0107]     Image motion in a first image and a second image that temporally follows the first image is to be understood to mean that an image content constituent is represented in the first image at a first (image) position and in a first form and is represented in the second, following image at a second position and in a second form, the first position and the second position or the-first form and the second form being different.  
         [0108]     Efficiently means, for example, that the calculation can be implemented by means of simple and cost-effective hardware in a short time.  
         [0109]     By way of example, the intention is that the hardware required for the calculation can be provided in a cost-effective mobile radio telephone.  
         [0110]     As mentioned above, coding information is to be understood to mean brightness information (luminescence information) and/or an item of color information (chrominance information) which is in each assigned to a pixel.  
         [0111]     In one embodiment, a method for computer-aided motion estimation in at least two temporally successive digital images with pixels to which coding information is assigned, is provided 
        a set of feature points of the first image being determined using a first selection criterion, a feature point of the first image being a pixel of the first image, in the case of which the coding information which is assigned to the pixel and the coding information which is assigned in each case to the pixels in a vicinity of the pixel satisfy the first selection criterion;     a set of feature points of the second image being determined using a second selection criterion, a feature point of the second image being a pixel of the second image, in the case of which the coding information which is assigned to the pixel and the coding information which is assigned in each case to the pixels in a vicinity of the pixel satisfy the second selection criterion;     an assignment of each feature point of the first image to a respective feature point of the second image is determined on the basis of the spatial distribution of the set of feature points of the first image and on the basis of the spatial distribution of the set of feature points of the second image; and     the motion is estimated with subpixel accuracy on the basis of the assignment.        
 
         [0116]     Furthermore, in one embodiment, a computer program element is provided, which, after it has been loaded into a memory of a computer, has the effect that the computer carries out the above method.  
         [0117]     Furthermore, in one embodiment, a computer-readable storage medium is provided, on which a program is stored which enables a computer, after it has been loaded into a memory of the computer, to carry out the above method.  
         [0118]     Furthermore, in one embodiment, a device is provided, which is set up such that the above method is carried out.  
         [0119]     In one case, the motion determination is effected by means of a comparison of feature positions.  
         [0120]     In one embodiment, features are determined in two successive images and assignment is determined by attempting to determine those features in the second image to which the features in the first image respectively correspond. If that feature in the second image to which a feature in the first image corresponds has been determined, then this is interpreted such that the feature in the first image has migrated to the position of the feature in the second image and this position change, which corresponds to an image motion of the feature, is calculated. Furthermore, a uniform motion model which models the position changes as well as possible is calculated on the basis of the position changes of the individual features.  
         [0121]     Therefore, an assignment is fixedly chosen and a motion model is determined which best maps all feature points of the first image onto the feature points—respectively assigned to them—of the second image in a certain sense, for example in a least squares sense as described below.  
         [0122]     In one case, a distance between the set of feature points of the first image that is mapped by means of the motion model and the set of the feature points of the second image is not calculated for all values of the parameters of the motion model. Consequently, a low computational complexity is achieved in the case of the method provided.  
         [0123]     Features are points of the image which are significant in a certain predetermined sense, for example edge points.  
         [0124]     An edge point is a point of the image at which a great local change in brightness occurs; for example, a point whose neighbor on the left is black and whose neighbor on the right is white is an edge point.  
         [0125]     Formally, an edge point is determined as a local maximum of the image gradient in the gradient direction or is determined as a zero crossing of the second derivative of the image information.  
         [0126]     Further image points that can be used as feature points in the method provided are, for example: 
        gray-scale value corners, that is, pixels which have a local maximum of the image gradient in the x and y direction.     corners in contour profiles, that is, pixels at which a significant high curvature of a contour occurs.     pixels with a local maximum filter response in the case of filtering with local filter masks (for example, sobel operator, gabor functions, etc.).     pixels which characterize the boundaries of different image regions. These image regions are generated, for example, by image segmentations such as “region growing” or “watershed segmentation”.     pixels which describe centroids of image regions, as are generated for example by the image segmentations mentioned above.        
 
         [0132]     The positions of a set of features are determined by a two-dimensional spatial feature distribution of an image.  
         [0133]     In the determination of the motion of a first image and a second image in accordance with one method provided, the spatial feature distribution of the first image is compared with the spatial feature distribution of the second image.  
         [0134]     In contrast to a method based on the optical flow, in the case of one method provided the motion is not calculated on the basis of the brightness distribution of the images, but rather on the basis of the spatial distribution of significant points.  
         [0135]     In addition to the above-described “super-resolution”, that is, the generation of high resolution images from a sequence of low resolution images, the motion estimation method provided may furthermore be used 
        for structure-from-motion methods that serve to infer the 3D geometry of the vicinity from a sequence of images recorded by a moving camera;     for methods for generating mosaic images in which a large high resolution image is assembled from individual high resolution smaller images; and     for video compression methods in which an improved compression rate can be achieved by means of a motion estimation.        
 
         [0139]     One embodiment of the method provided is distinguished by its high achievable accuracy and by its simplicity.  
         [0140]     It is not necessary for any spatial and temporal derivatives to be approximated, which is computationally intensive and typically leads to inaccuracies.  
         [0141]     On account of the simplicity of one method provided, it is possible to implement the method in a future mobile radio telephone, for example, without the latter having to have a powerful and cost-intensive data processing unit.  
         [0142]     In one case of the above method, the motion estimation based on the spatial distribution of the set of feature points of the first image and on the spatial distribution of the set of feature points of the second image is carried out by a feature point from the set of features of the second image being assigned to each feature point from the set of feature points of the first image.  
         [0143]     In one embodiment, a feature point from the set of feature points of the first image is assigned to a feature point from the set of feature points of the second image with respect to which the feature point from the set of feature points of the first image has a minimum spatial distance which is determined from the coordinates of the feature point from the set of feature points of the first image and the coordinates of the feature point from the set of feature points of the second image.  
         [0144]     In one case, the motion estimation can be carried out with low computational complexity. By way of example, the abovementioned assignment can be carried out with the aid of a distance transformation, for which efficient methods are known.  
         [0145]     In one case of the above method, the determination of the set of feature points of the first image, the determination of the set of feature points of the second image and the motion estimation are effected with subpixel accuracy.  
         [0146]     In one case of the above method, the motion estimation is effected by determination of a motion model.  
         [0147]     In one case, a translation is determined prior to the determination of the motion model.  
         [0148]     In one embodiment, which is described below, a translation is determined before the above-described assignment of the feature points of the first image to feature points of the second image is determined.  
         [0149]     By determining a translation prior to the determination of the motion model, it is possible to increase the accuracy of the motion estimation with low computational complexity.  
         [0150]     The translation can be determined with low computational complexity since a translation can be determined by means of a small number of motion parameters.  
         [0151]     In one case, an affine motion model or a perspective motion model is determined.  
         [0152]     The motion model is in one case determined iteratively.  
         [0153]     In this case, in each iteration the assignment of each feature point of the first image to a respective feature point of the second image is fixedly chosen, but the assignments that are used in different iterations may be different.  
         [0154]     As a result, it is possible to obtain a high accuracy.  
         [0155]     In one case of the above method, the first selection criterion and the second selection criterion being chosen such that the feature points from the set of feature points of the first image are edge points of the first image and the feature points from the set of feature points of the second image are edge points of the second image.  
         [0156]     In one case, the above method is used, in the case of a structure-from-motion method, in the case of a method for generating mosaic images, in the case of a video compression method or a super-resolution method. 
     
    
     BRIEF DESCRIPTION OF THE DRAWINGS  
       [0157]     The accompanying drawings are included to provide a further understanding of the present invention and are incorporated in and constitute a part of this specification. The drawings illustrate the embodiments of the present invention and together with the description serve to explain the principles of the invention. Other embodiments of the present invention and many of the intended advantages of the present invention will be readily appreciated as they become better understood by reference to the following detailed description. The elements of the drawings are not necessarily to scale relative to each other. Like reference numerals designate corresponding similar parts.  
         [0158]      FIG. 1  illustrates an arrangement in accordance with one exemplary embodiment of the invention.  
         [0159]      FIG. 2  illustrates a flow diagram of a method in accordance with one exemplary embodiment of the invention.  
         [0160]      FIG. 3  illustrates a flow diagram of a determination of a translation in accordance with one exemplary embodiment of the invention.  
         [0161]      FIG. 4  illustrates a flow diagram of a determination of an affine motion in accordance with one exemplary embodiment of the invention.  
         [0162]      FIG. 5  illustrates a flow diagram of a method in accordance with a further exemplary embodiment of the invention.  
         [0163]      FIG. 6  illustrates a flow diagram of an edge detection in accordance with one exemplary embodiment of the invention.  
         [0164]      FIG. 7  illustrates a flow diagram of an edge detection with subpixel accuracy in accordance with one exemplary embodiment of the invention;  
         [0165]      FIG. 8A  and  FIG. 8B  illustrate the results of a performance comparison of one embodiment of the invention with known methods.  
         [0166]      FIG. 9  illustrates a flow diagram of a method in accordance with a further exemplary embodiment of the invention.  
         [0167]      FIG. 10  illustrates a flow diagram of a determination of a perspective motion in accordance with one exemplary embodiment of the invention.  
         [0168]      FIG. 11  illustrates a flow diagram of a known method for parametric motion determination. 
     
    
     DETAILED DESCRIPTION  
       [0169]     In the following Detailed Description, reference is made to the accompanying drawings, which form a part hereof, and in which is shown by way of illustration specific embodiments in which the invention may be practiced. In this regard, directional terminology, such as “top,” “bottom,” “front,” “back,” “leading,” “trailing,” etc., is used with reference to the orientation of the Figure(s) being described. Because components of embodiments of the present invention can be positioned in a number of different orientations, the directional terminology is used for purposes of illustration and is in no way limiting. It is to be understood that other embodiments may be utilized and structural or logical changes may be made without departing from the scope of the present invention. The following detailed description, therefore, is not to be taken in a limiting sense, and the scope of the present invention is defined by the appended claims.  
         [0170]      FIG. 1  illustrates an arrangement  100  in accordance with one exemplary embodiment of the invention.  
         [0171]     A low resolution digital camera  101  is held over a printed text  102  by a user (not shown).  
         [0172]     Low resolution is to be understood to mean a resolution which does not suffice to enable a digital image with this resolution of the printed text  102  which has been recorded by means of the digital camera  101  and is displayed on a screen to represent the text with sufficiently high resolution such that it can be read in a simple manner by a user or can be automatically processed further in a simple manner, for example in the case of optical pattern recognition, optical character recognition.  
         [0173]     The printed text  102  may be for example a text printed on paper which the user wishes to send to another person.  
         [0174]     The digital camera  101  is coupled to a (micro)processor  107 .  
         [0175]     The digital camera  101  generates a sequence of low resolution digital images  105  of the printed text  102 . The recording positions of the digital images from the sequence of low resolution digital images  105  of the printed text  102  are different since the user&#39;s hand is not completely motionless.  
         [0176]     The sequence of low resolution digital images  105  is fed to the processor  107 , which calculates a high resolution digital image  106  from the sequence of low resolution digital images  105 .  
         [0177]     For this purpose, the processor  107  uses a method for ascertaining the image motion such as is described further below in exemplary embodiments.  
         [0178]     The high resolution digital image  106  is displayed on a screen  103  and can be transmitted to another person by the user by means of a transmitter  104 .  
         [0179]     In an exemplary embodiment, the digital camera  101 , the processor  107 , the screen  103  and the transmitter  104  are contained in a mobile radio telephone.  
         [0180]      FIG. 2  illustrates a flow diagram  200  of a method in accordance with one exemplary embodiment of the invention.  
         [0181]     The method explained below serves for calculating the motion in the sequence of low resolution images  105  that have been recorded by means of the digital camera  101 . Each image of the sequence of low resolution images  105  is expressed by a function I(x, y, t), where t is the instant at which the image was recorded and I(x, y, t) specifies the coding information of the image at the location (x, y) which was recorded at the instant t.  
         [0182]     It is assumed in this exemplary embodiment that no illumination fluctuations or disturbances in the process hardware occurred during the recording of the digital images.  
         [0183]     Under this assumption, the following equation holds true for two successive digital images in the sequence of low resolution images  105  with the coding information I(x, y, t) and I(x, y, t+dt) occur respectively: 
 
 I ( x+dx, y+dy, t+dt )= I ( x, y, t )   (14) 
 
 In this case, dt is the difference between the recording instants of the two successive digital images in the sequence of low resolution images  105 . 
 
         [0184]     Under the assumption that only one cause of motion exists equation (14) can also be formulated by 
 
 I ( x, y, t+dt )= I (Motion( x, y, t ),  t )   (15) 
 
 where motion (x, y, t) describes the motion of the pixels. 
 
         [0185]     The image motion can be modeled for example by means of an affine transformation  
                 [           x   ⁡     (     t   +   dt     )                 y   ⁡     (     t   +   dt     )             ]     =         [           m     x   ⁢           ⁢   0             m     x   ⁢           ⁢   1                 m     y   ⁢           ⁢   0             m     y   ⁢           ⁢   1             ]     ⁡     [           x   ⁡     (   t   )                 y   ⁡     (   t   )             ]       +     [           t   x               t   y           ]         ,           (   16   )             
 
 An image of the sequence of low resolution digital images  105  is provided in step  201  of the flow diagram  200 . 
 
         [0186]     It is assumed that the digital image was recorded by means of the digital camera  101  at an instant t+1.  
         [0187]     An image that was recorded at an instant τ is designated hereinafter as image τ for short.  
         [0188]     Consequently, by way of example, the image that was recorded by means of the digital camera  101  at an instant t+1 is designated as image t+1.  
         [0189]     It is furthermore assumed that a digital image that was recorded at an instant t is present, and that the image motion from the image t to the image t+1 is to be determined.  
         [0190]     The feature detection, that is, the determination of feature points and feature positions, is prepared in step  202 .  
         [0191]     By way of example, the digital image is preprocessed by means of a filter for this purpose.  
         [0192]     A feature detection with a low threshold is carried out in step  202 .  
         [0193]     This means that during the feature detection, a value is assigned to each pixel, and a pixel belongs to the set of feature points only when the value assigned to it lies above a certain threshold value.  
         [0194]     In the case of the feature detection carried out in step  202 , said threshold value is low, where “low” is to be understood to mean that the value is less than the threshold value of the feature detection carried out in step  205 .  
         [0195]     A feature detection in accordance with one embodiment of the invention is described further below.  
         [0196]     The set of feature points that is determined during the feature detection carried out in step  202  is designated by P t+1   K :  
               P     t   +   1     K     =     {         [         P       t   +   1     ,   x       ⁡     (   k   )       ,       P       t   +   1     ,   y       ⁡     (   k   )         ]     T     ,     0   ≤   k   ≤     K   -   1         }             (   17   )             
 
         [0197]     In this case,  P   t+1 =[P t+1,x (k) P t+ 1,y(k)] T  designates a feature point with the index k from the set of feature points P t+1   K  in vector notation.  
         [0198]     The image information of the image t is written as function I(x, y, t) analogously to above.  
         [0199]     A global translation is determined in step  203 .  
         [0200]     This step is described below with reference to  FIG. 3 .  
         [0201]     Affine motion parameters are determined in step  204 .  
         [0202]     This step is described below with reference to  FIG. 4 .  
         [0203]     A feature detection with a high threshold is carried out in step  205 .  
         [0204]     In other words, the threshold value is high during the feature detection carried out in step  205 , where high is to be understood to mean that the value is greater than the threshold value of the feature detection with a low threshold value that is carried out in step  202 .  
         [0205]     As mentioned, a feature detection in accordance with one embodiment of the invention is described further below.  
         [0206]     The set of feature points determined during the feature detection carried out in step  205  is designated by O t+1   N :  
               O     t   +   1     N     =     {         [         O       t   +   1     ,   x       ⁡     (   n   )       ,       O       t   +   1     ,   y       ⁡     (   n   )         ]     T     ,     0   ≤   n   ≤     N   -   1         }             (   18   )             
 
         [0207]     In this case,  O   t+1 (n)=[O t+1,x (n), O t+1,y (n)] T  designates the n-th feature point of the set O t+1   N  in vector notation. The feature detection with a high threshold that is carried out in step  205  does not serve for determining the motion from image t to image t+1, but rather serves for preparing for the determination of motion from image t+1 to image t+2.  
         [0208]     Accordingly, it is assumed hereinafter that a feature detection with a high threshold for the image t to step  205  was carried out in which a set of feature points  
               O   t   N     =     {         [         O     t   ,   x       ⁡     (   n   )       ,       O     t   ,   y       ⁡     (   n   )         ]     T     ,     0   ≤   n   ≤     N   -   1         }             (   19   )             
 
 was determined. 
 
         [0209]     Step  203  and step  204  are carried out using the set of feature points O t   N .  
         [0210]     In step  203  and step  204 , a suitable affine motion determined by a matrix  {circumflex over (M)}   t  and a translation vector  {circumflex over (T)}   t  is calculated, so that for 
 
   Ô     t+1   N   = {circumflex over (M)}     t     O     t   N   + {circumflex over (T)}     t    (20) 
 
 the relationship 
 
Ô t+1   N ⊂P t+1   N    (21) 
 
 holds true, where Ô t+1   N  is the set of column vectors of the matrix  Ô   t+1   N . 
 
         [0211]     In this case,  O   t   N  designates the matrix whose column vectors are the vectors of the set O t   N .  
         [0212]     This can be interpreted such that a motion is sought which maps the feature points of the image t onto feature points of the image t+1.  
         [0213]     The determination of the affine motion is made possible by the fact that a higher threshold is used for the detection of the feature points from the set O t   N  than for the detection of the feature points from the set P t+1   K .  
         [0214]     If the same threshold is used for both detections, there is the possibility that some of the pixels corresponding to the feature points from O t   N  will not be detected as feature points at the instant t+1.  
         [0215]     The pixel in image t+1 that corresponds to a feature point in image t is to be understood as the pixel at which the image content constituent represented by the feature point in image t is represented in image t+1 on account of the image motion.  
         [0216]     In general,  {circumflex over (M)}   t  and  {circumflex over (T)}   t  cannot be determined such that (21) holds true, thereforee  {circumflex over (M)}   t  and  {circumflex over (T)}   t  are determined such that O t   N  is mapped onto P t+1   K  as well as possible by means of the affine motion in a certain sense as is defined below.  
         [0217]     In this embodiment, the minimum distances of the points from Ô t   N  to the set P t+1   N , are used for a measure of the quality of the mapping of O t   N  onto P t+1   K .  
         [0218]     The minimum distance  
              D     min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )              
 
 of a point (x, y) from the set P t+1   K  is defined by  
                      D     min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )            =       min   k     ⁢           ⁢              [     x   ,   y     ]     T     -       P     t   +   1       ⁡     (   k   )                        (   22   )             
 
 The minimum distances of the points from O t   N  from the set P t+1   K  can be determined efficiently for example with the aid of a distance transformation, which is a morphological operation (see G. Borgefors,  Distance Transformation in Digital Images,  Computer Vision, Graphics and Image Processing, 34, pp. 344-371, 1986). 
 
         [0219]     In the case of a distance transformation such as is described in G. Borgefors,  Distance Transformation in Digital Images,  Computer Vision, Graphics and Image Processing, 34, pp. 344-371, 1986, a distance image is generated from an image in which feature points are identified, in which distance image the image value at a point specifies the minimum distance to a feature point.  
         [0220]     Clearly,  
              D     min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )              
 
 specifies for a point the distance to the point from P t+1   K  with respect to which the point (x, y) has the smallest distance. 
 
         [0221]     The affine motion is determined in the two steps  203  and  204 .  
         [0222]     For this purpose, the affine motion formulated in (20) is decomposed into a global translation and a subsequent affine motion: 
 
   Ô     t+1   N   = {circumflex over (M)}     t (   O     t   N   + {circumflex over (T)}     t   0 )+   {circumflex over (T)}     t   1    (23) 
 
         [0223]     The translation vector  {circumflex over (T)}   t   0  determines the global translation and the matrix  {circumflex over (M)}   t  and the translation vector  {circumflex over (T)}   t   1  determine the subsequent affine motion.  
         [0224]     Step  203  is explained below with reference to  FIG. 3 .  
         [0225]      FIG. 3  illustrates a flow diagram  300  of a determination of a translation in accordance with an exemplary embodiment of the invention.  
         [0226]     In step  203 , which is represented by step  301  of the flow diagram  300 , the translation vector is determined using P t+1   K  and O t   N  such that  
                 T   ^     t   0     =         arg   ⁢           ⁢   min       T   t   0       ⁢           ⁢       ∑   n     ⁢           ⁢            D     min   ,     P     t   +   1     K         ⁡     (           O   tx     ⁡     (   n   )       +     T   tx   0       ,         O   ty     ⁡     (   n   )       +     T   ty   0         )                        (   24   )             
 
         [0227]     Step  301  has steps  302 ,  303 ,  304  and  305 .  
         [0228]     For the determination of {circumflex over (T)} t   0 , such that equation (24) holds true, step  302  involves choosing a value T y   0  in an interval [{circumflex over (T)} y0   0 , {circumflex over (T)} y1   0 ].  
         [0229]     Step  303  involves choosing a value T x   0  in an interval [{circumflex over (T)} x0   0 , {circumflex over (T)} x1   0 ].  
         [0230]     Step  304  involves determining the value sum (T x   O , T y   0 ) in accordance with the formula  
               sum   ⁡     (         T   _     x   0     ,       T   _     y   0       )       =       ∑   n     ⁢           ⁢              D   _       min   ,     P     t   +   1     K         ⁡     (             O   _     tx     ⁡     (   n   )       +       T   _     tx   0       ,           O   _     ty     ⁡     (   n   )       +       T   _     ty   0         )                      (   25   )             
 
 for the chosen values  T   x   0  and  T   y   0 . 
 
         [0231]     Steps  302  to  304  are carried out for all chosen pairs of values  
         T   y   0     ∈       [         T   ^       y   ⁢           ⁢   0     0     ,       T   ^       y   ⁢           ⁢   1     0       ]     ⁢           ⁢   and   ⁢           ⁢     T   x   0       ∈       [         T   ^       x   ⁢           ⁢   0     0     ,       T   ^       x   ⁢           ⁢   1     0       ]     .         
 
         [0232]     In step  305 ,  {circumflex over (T)}   y   0  and  {circumflex over (T)}   x   0  are determined such that sum ( {circumflex over (T)}   x   0 ,  {circumflex over (T)}   y   0 ) is equal to the minimum of all sums calculated in step  304 .  
         [0233]     The translation vector  {circumflex over (T)}   t   0  is given by  
                   T   ^     _     t   0     =     [           T   ^     _     x   0     ,         T   ^     _     y   0       ]             (   26   )             
 
 Step  204  is explained below with reference to  FIG. 4 . 
 
         [0234]      FIG. 4  illustrates a flow diagram  400  of a determination of an affine motion in accordance with an exemplary embodiment of the invention.  
         [0235]     Step  204 , which is represented by step  401  of the flow diagram  400 , has steps  402  to  408 .  
         [0236]     Step  402  involves calculating the matrix 
 
   O ′   t   N   = O     t   N   + {circumflex over (T)}     t   0    (27) 
 
 whose column vectors form a set of points O′ t   N . 
 
         [0237]     A distance vector  
           D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )         
 
 is determined for each point (x, y) from the set O′ t   N . 
 
         [0238]     The distance vector is determined such that it points from the point (x, y) to the point from P t+1   K  with respect to which the distance of the point (x, y) is minimal.  
         [0239]     The determination is thus effected in accordance with the equations  
               k   min     =         arg   ⁢           ⁢   min     k     ⁢           ⁢              [     x   ,   y     ]     T     -       P     t   +   1       ⁡     (   k   )                        (   28   )             
 
                   D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )       =         [     x   .   y     ]     T     -       P     t   +   1       ⁡     (     k   min     )                 (   29   )             
 
 The distance vectors can also be calculated from the minimum distances which are present in the form of a distance image, for example, in accordance with the following formula:  
                   D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )       =                D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )            ⁡     [             ∂              D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )                ∂   x                   ∂              D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )                ∂   y             ]               (   30   )             
 
 In steps  403  to  408 , assuming that the approximation  
                   O   _       t   +   1     N     ≈         O   ~     _       t   +   1     N       =         O   _     t     ′   N       +         D   _       min   ,     P     t   +   1     K         ⁡     (     O   t     ′   N       )                 (   31   )             
 
 holds true for the feature point set P t+1   N , the affine motion is determined by means of a least squares estimation, that is, that the matrix  {circumflex over (M)}   t  and the translation vector  {circumflex over (T)}   t   1  are determined such that the term  
               ∑   n     ⁢           ⁢       (             O   ~     _       t   +   1       ⁡     (   n   )       -     (             M   ^     _     t   1     ⁢         O   _     t   ′     ⁡     (   n   )         +         T   ^     _     t   1       )       )     2             (   32   )             
 
 is minimal, which is the case precisely when the term  
               ∑   n     ⁢           ⁢         (       (           O   _     t   ′     ⁡     (   n   )       +         D   _       min   ,     P     t   +   1     K         ⁡     (         O   _     t   ′     ⁡     (   n   )       )         )     -     (             M   ^     _     t   1     ⁢         O   _     t   ′     ⁡     (   n   )         +         T   ^     _     t   1       )       )     2     .             (   33   )             
 
 is minimal. 
 
         [0240]     In this case, the n-th column of the respective matrix is designated by  O ′ t (n) and  Õ   t+1 (n).  
         [0241]     The use of the minimum distances in equation (33) can be interpreted such that it is assumed that a feature point in image t corresponds to the feature point in image t+1 which lies nearest to it, that is, that the feature point in image t has moved to the nearest feature point in image t+1.  
         [0242]     The least squares estimation is iterated in this embodiment.  
         [0243]     This is effected in accordance with the following decomposition of the affine motion: 
 
   {circumflex over (M)} O+ {circumflex over (T)} = {circumflex over (M)}     L (   {circumflex over (M)}     L-1 ( . . . ( {circumflex over (M)}   1 ( O+ {circumflex over (T)}     0 )+   {circumflex over (T)}     1 ) . . . )+   {circumflex over (T)}     L-1 )+   {circumflex over (T)}     L    (34) 
 
 The temporal dependence has been omitted in equation (34) for the sake of simplified notation. 
 
         [0244]     That is, that L affine motions are determined, the L-th affine motion being determined in such a way that it maps the feature point set which arises as a result of progressive application of the 1 st , 2 nd , . . . and the (l-1)th affine motion to the feature point set O′ t   N  onto the set P t+1   K , as well as possible, in the above-described sense of the least squares estimation.  
         [0245]     The l-th affine motion is determined by the matrix  {circumflex over (M)}   t   l  and the translation vector  {circumflex over (T)}   t   1 .  
         [0246]     At the end of step  402 , the iteration index 1 is set to zero and the procedure continues with step  403 .  
         [0247]     In step  403 , the value of 1 is increased by 1 and a check is made to ascertain whether the iteration index 1 lies between 1 and L.  
         [0248]     If this is the case, the procedure continues with step  404 .  
         [0249]     Step  404  involves determining the feature point set O′ 1  that arises as a result of the progressive application of the 1 st , 2 nd , and the (l-1)-th affine motion to the feature point set O′ t   N .  
         [0250]     Step  405  involves determining distance vectors analogously to equations (28) and (29) and a feature point set analogously to (31).  
         [0251]     Step  406  involves calculating a matrix  {circumflex over (M)}   t   l  and a translation vector  {circumflex over (T)}   t   l , which determine the l-th affine motion.  
         [0252]     Moreover, a square error is calculated analogously to (32).  
         [0253]     Step  407  involves checking whether the square error calculated is greater than the square error calculated in the last iteration.  
         [0254]     If this is the case, in step  408  the iteration index 1 is set to the value L and the procedure subsequently continues with step  403 .  
         [0255]     If this is not the case, the procedure continues with step  403 .  
         [0256]     If the iteration index is set to the value L in step  408 , then in step  403  the value of 1 is increased to the value L+1 and the iteration is ended.  
         [0257]     In a one embodiment, steps  202  to  205  of the flow diagram  200  illustrated in  FIG. 2  are carried out with subpixel accuracy.  
         [0258]      FIG. 5  illustrates a flow diagram  500  of a method in accordance with a further exemplary embodiment of the invention.  
         [0259]     In this embodiment, a digital image that was recorded at the instant  0  is used as a reference image, which is designated hereinafter as reference window.  
         [0260]     The coding information  502  of the reference window  501  is written hereinafter as function I(x, y,  0 ) analogously to the above.  
         [0261]     Step  503  involves carrying out an edge detection with subpixel resolution in the reference window  501 .  
         [0262]     A method for edge detection with subpixel resolution in accordance with one embodiment is described below with reference to  FIG. 7 .  
         [0263]     In step  504 , a set of feature points ON of the reference window is determined from the result of the edge detection.  
         [0264]     For example, the significant edge points are determined as feature points.  
         [0265]     The time index t is subsequently set to the value zero.  
         [0266]     In step  505 , the time index t is increased by one and a check is subsequently made to ascertain whether the value of t lies between 1 and T.  
         [0267]     If this is the case, the procedure continues with step  506 .  
         [0268]     If this is not the case, the method is ended with step  510 .  
         [0269]     In step  506 , an edge detection with subpixel resolution is carried out using the coding information  511  of the t-th image, which is designated as image t analogously to the above.  
         [0270]     This yields, as is described in greater detail below, a t-th edge image, which is designated hereinafter as edge image t, with the coding information e h (x, y, t) with respect to the image t.  
         [0271]     The coding information e h (x, y, t) of the edge image t is explained in more detail below with reference to  FIG. 6  and  FIG. 7 .  
         [0272]     Step  507  involves carrying out a distance transformation with subpixel resolution of the edge image t.  
         [0273]     That is, that a distance image is generated from the edge image t, in the case of which distance image the image value at a point specifies the minimum distance to an edge point.  
         [0274]     The edge points of the image t are the points of the edge image t in the case of which the coding information e h (x, y, t) has a specific value.  
         [0275]     This is explained in more detail below.  
         [0276]     The distance transformation is effected analogously to the embodiment described with reference to  FIG. 2 ,  FIG. 3  and  FIG. 4 .  
         [0277]     In this case, use is made of the fact that the positions of the edge points of the image t were determined with subpixel accuracy in step  506 .  
         [0278]     The distance vectors are calculated with subpixel accuracy.  
         [0279]     In step  508 , a global translation is determined analogously to step  203  of the exemplary embodiment described with reference to  FIG. 2 ,  FIG. 3  and  FIG. 4 .  
         [0280]     The global translation is determined with subpixel accuracy.  
         [0281]     Parameters of an affine motion model are calculated in the processing block  509 .  
         [0282]     The calculation is effected analogously to the flow diagram illustrated in  FIG. 4  as explained above.  
         [0283]     The parameters of an affine motion model are calculated with subpixel accuracy.  
         [0284]     After the end of the processing block  509 , the procedure continues with step  505 .  
         [0285]     The method is ended if t=T, that is, if the motion of the image content between the reference window and the t-th image has been determined.  
         [0286]      FIG. 6  illustrates a flow diagram  600  of an edge detection in accordance with an exemplary embodiment of the invention.  
         [0287]     The determination of edges represents an expedient compromise for the motion estimation with regard to concentration on significant pixels during the motion determination and obtaining as many items of information as possible.  
         [0288]     Edges are usually determined as local maxima in the local derivative of the image intensity. The method used here is based on the papers by Canny J. Canny,  A Computational Approach to Edge Detection,  IEEE Transactions on Pattern Analysis and Machine Intelligence, 6, 1986.  
         [0289]     In step  602 , a digital image in the case of which edges are intended to be detected is fitted by means of a Gaussian filter.  
         [0290]     This is effected by convolution of the coding information  601  of the image, which is given by the function I(x, y), using a Gaussian mask designated by gmask.  
         [0291]     Step  603  involves determining the partial derivative with respect to the variable x of the function I g (x, y).  
         [0292]     Step  604  involves determining the partial derivative with respect to the variable y of the function I g (x, y).  
         [0293]     In step  605 , a decision is made as to whether an edge point is present at a point (x, y).  
         [0294]     For this purpose, two conditions have to be met at the point (x, y).  
         [0295]     The first condition is that the sum of the squares of the two partial derivatives determined in step  603  and step  604  at the point (x, y), designated by I g, x, y (x, y) lies above a threshold value.  
         [0296]     The second condition is that I g, x, y (x, y) has a local maximum at the point (x, y).  
         [0297]     The result of the edge detection is combined in an edge image whose coding information  606  is written as a function and designated by e(x, y).  
         [0298]     The function e(x, y) has the value I g, x, y (x, y) at a location (x, y) if it was decided with regard to (x, y) in step  605  that (x, y) is an edge point, and has the value zero at all other locations.  
         [0299]     The approach for detecting gray-scale value corners as illustrated in  FIG. 6  affords the possibility of controlling the number and the significance of the edges by means of a threshold.  
         [0300]     It can thus be ensured that O t+1   N  is contained in P t+1   K .  
         [0301]     The point sets O t+1   N  and P t+1   K  can be read from the edge image having the coding information e(x, y).  
         [0302]     If the method illustrated in  FIG. 6  is used in the exemplary embodiment illustrated in  FIG. 2 , then for generating P t+1   K  from e(x, y) the threshold used in step  605  corresponds to the “low threshold” used in step  205 .  
         [0303]     For determining O t+1   N  using the “high threshold” used in  205 , a selection is made from the edge points given by e(x, y).  
         [0304]     This is effected for example analogously to the checking of the first condition from step  605  as explained above.  
         [0305]      FIG. 7  illustrates a flow diagram  700  of an edge detection with subpixel accuracy in accordance with an exemplary embodiment of the invention.  
         [0306]     Steps  702 ,  703  and  704  do not differ from steps  602 ,  603  and  604  of the edge detection method illustrated in  FIG. 6 .  
         [0307]     In order to achieve a detection with subpixel accuracy, the flow diagram  700  has a step  705 .  
         [0308]     Step  705  involves extrapolating the partial derivatives in the x direction and y direction determined in step  703  and step  704 , which are designated as local gradient images with coding information I gx (x, y) and I gy (x, y), to a higher image resolution.  
         [0309]     The missing image values are determined by means of a bicubic interpolation. The method of bicubic interpolation is explained, for example, in William H. Press, et al.,  Numerical Recipes in C,  ISBN: 0-521-41508-5, Cambridge University Press.  
         [0310]     The coding information of the resulting high resolution gradient images is designated by I hgx (x, y) and I hgy (X, y).  
         [0311]     Step  706  is effected analogously to step  605  using the high resolution edge images.  
         [0312]     The coding information  707  of the edge image generated in step  706  is designated by e h (x, y), where the index h is intended to indicate that the edge image likewise has a high resolution.  
         [0313]     The function e h (x, y) generated in step  707 , in contrast to that in step  706 , in this exemplary embodiment does not have the value I g, x, y (x, y) if it was decided that an edge point is present at the location (x, y), but rather the value 1.  
         [0314]     The results of a performance comparison between the method provided and known methods are explained below.  
         [0315]      FIG. 8A  and  FIG. 8B  illustrate the results of a performance comparison between an embodiment of the invention and known methods.  
         [0316]     In order to generate reference data for the evaluation of the motion estimation, “camera shake” was simulated.  
         [0317]     For this purpose, different views, that is, recordings from different camera positions, were generated by means of simulation from a high resolution image using affine transformations.  
         [0318]     These views were subsequently filtered by means of a low-pass filter and subsampled. The sequence of digital images thus generated, which was used as an example of a sequence of digital images recorded by a moving camera, was processed by means of various methods for motion estimation.  
         [0319]     The following reference methods were used:  
         [0320]     1. An optical flow method based on the papers by Lucas and Kanade (see B. Lucas, T. Kanade,  An Iterative Image Registration Technique with an Application to Stero Vision,  7th International Joint Conf. on Artificial Intelligence (IJCAI), pp. 674-679, 1981), using gray-scale value comers with a resolution with subpixel accuracy. The method additionally uses a resolution pyramid in order to avoid problems in the case of fast motions.  
         [0321]     This method corresponds to the dotted line in  FIG. 8A  and  FIG. 8B .  
         [0322]     2. A parametric motion estimation method based on the optical flow.  
         [0323]     This method corresponds to the dash-dotted line in  FIG. 8A  and  FIG. 8B .  
         [0324]     3. A method for distance-based motion estimation without improvement of the subpixel accuracy.  
         [0325]     This method corresponds to the dashed line in  FIG. 8A  and  FIG. 8B .  
         [0326]      FIG. 8A  illustrates the profiles of the average error of the motion estimation in an embodiment of the method provided with subpixel accuracy and the three reference methods.  
         [0327]     The deviation between the simulated displacement and the measured displacement vectors was averaged over all pixels.  
         [0328]     The motion of the camera was firstly simulated as a pure translation assuming ideal conditions.  
         [0329]      FIG. 8B  illustrates the profiles of the average error of the motion estimation in an embodiment of the method provided with subpixel accuracy and the three reference methods for the simulation of an affine transformation as camera motion.  
         [0330]     The error profiles illustrated in  FIG. 8A  and  FIG. 8B  illustrate that the greatest accuracy is obtained with an embodiment of the method provided.  
         [0331]     An overview of the required number of additions and multiplications in an embodiment of the method provided which generated the results illustrated in  FIG. 8A  and  FIG. 8B  is given below.  
         [0332]     In addition, typical values for the number of additions and multiplications are specified for the example of a QVGA resolution.  
                                                                                           Processing Steps                Addi-   Multipli-   Add. for   Mult. for           tions   cations   QVGA   QVGA                        Gaussian   (s g  − 1)2 r c   s g 2 r c   1 075 200     153 600       filter       Gradient   2 r c       153 600       0       filter       Edge   rc + 4n e     2rc +   107 520   192 000       detection       5n e         Edges with   9n e   ·   9n e  ·   29 859 840     13 271 040         subpixel   (20 +   (12 +       accuracy   u x u y 103)   u x u y 45)       Distance   8 r u y   c u x         2 457 600         0       transfor-       mation       optimal   s x s y N       232 320       0       trans-       lation       Affine   L (30 +   L (56 +   460 950   269 080       transfor-   48N)   28N)       mation           Total           34 347 030     13 885 720                    
 
         [0333]     The definitions of the variables for the assessment of the computation time are specified in the table below.  
                                       Variable   Meaning   Typical value                   s g     Magnitude of the Gaussian mask   s g  = 7       r, c   Number of rows (r), number of   r = 240, c = 320           columns (c)       n e     Number of edge points   n e   = 0.1 r c       u x , u y     Scaling factors in x and y direction   u x  = 2, u y  = 2       s x , s y     Search range for optimal translation   s x  = 11, s y  = 11           in x and y direction       L   Number of iterations for the deter-   L = 5           mination of the affine transformation       N   Number of object points (N &lt; n e )   N = 0.25               n e  = 0.025 r c                  
 
         [0334]     It is evident that the complexity for the actual motion determination is low in relation to the feature extraction with subpixel accuracy.  
         [0335]     A feature extraction with subpixel accuracy is also required for example for the reference method specified under 3.  
         [0336]     For comparison of the number of operations, the assessment was likewise carried out for the method described with reference to  FIG. 11 .  
         [0337]     It was assumed in this case that 3 pyramid levels were used. On average 5 iterations were performed for each pyramid level.  
         [0338]     It was additionally taken into account that the optical flow is only carried out at points with high significance (for example, gray-scale value edges).  
         [0339]     The complexity for determining the significant pixels was not taken into account.  
         [0340]     The table below shows the results of the assessment of the required number of operations.  
                                                       Processing                   Add.   Mul. for       Steps   Iteration   Pyr.   Additions   Multiplications   QVGA   QVGA                   Low-pass       X   4 r c   2 r c   403 200   201 600       filter       Sampling       X               0       0       Local       X   2 r c       201 600       0       gradients       Temporal   x   X   rc       504 000       0       gradient               Motion determination   x   X             42   ⁢           ⁢     n   e       +       n   a   3     6                       46   ⁢           ⁢     n   e       +       n   a   3     6             2 117 880     2 319 480                 Quality   x   X   8n e     7n e     403 200   352 800       measurement       Motion   x   X   103rc   45rc   51 912 000     22 680 000         compensation       Update motion   x   X   8   12      120      180       parameter           Total                   55 542 000     25 554 060                    
 
         [0341]     It is striking that both methods are dominated by the complexity for the interpolation of image data.  
         [0342]     In the approach presented here, the interpolation is necessary for the edge detection with subpixel accuracy; in the reference method mentioned under 3, an interpolation is necessary for the motion compensation. A bicubic interpolation was used in both implementations.  
         [0343]     It is evident from the assessment of the computation times that the method provided, in one embodiment, is not more complex than previously known methods even though a higher accuracy can be achieved.  
         [0344]     The computation time for the novel method may additionally be significantly reduced if the detection of the features with subpixel accuracy is reworked.  
         [0345]     In one embodiment, the gradient images in x and y are converted into a higher image resolution by means of a linear interpolation. In contrast to the reference method by means of optical flow, this is appropriate here since the gradient images are locally smooth on account of the low-pass filter character of the gradient filter.  
         [0346]     In another embodiment, the interpolation is only performed at feature positions to be expected.  
         [0347]      FIG. 9  illustrates a flow diagram  900  of a method in accordance with a further exemplary embodiment of the invention.  
         [0348]     This exemplary embodiment differs from that explained with reference to  FIG. 2  in that a perspective motion model is used instead of an affine motion model such as is given by equation (16), for example.  
         [0349]     Since a camera generates a perspective mapping of the three-dimensional environment onto a two-dimensional image plane, an affine model yields only an approximation of the actual image motion which is generated by a moving camera.  
         [0350]     If an ideal camera, that is, without lens distortions, is assumed, the motion can be described by a perspective motion model such as is given by the equation below, for example.  
                     [           x   ⁡     (     t   +   dt     )                 y   ⁡     (     t   +   dt     )             ]     =       Motion   pers     ⁡     (       M   _     ,     x   ⁡     (   t   )       ,     y   ⁡     (   t   )         )                   =     [                 a   1     ⁢     x   ⁡     (   t   )         +       a   2     ⁢     y   ⁡     (   t   )         +     a   3             n   1     ⁢     x   ⁡     (   t   )         +       n   2     ⁢     y   ⁡     (   t   )         +     n   3                         b   1     ⁢     x   ⁡     (   t   )         +       b   2     ⁢     y   ⁡     (   t   )         +     b   3             n   1     ⁢     x   ⁡     (   t   )         +       n   2     ⁢     y   ⁡     (   t   )         +     n   3               ]                   (   35   )             
 
  M  designates the parameter vector for the perspective motion model. 
 
 M =[a 1 , a 2 , a 3 , b 1 , b 2 , b 3 , n 1 , n 2 , n 3 ]  (36) 
 
         [0351]     The method steps of the flow diagram  900  are analogous to those of the flow diagram  200 ; therefore, only the differences are discussed below.  
         [0352]     As in the case of the method described with reference to  FIG. 2 , a feature point set  
               O   t   N     =     {         [         O   tx     ⁡     (   n   )       ,       O   ty     ⁡     (   n   )         ]     T     ,     0   ≤   n   ≤     N   -   1         }             (   37   )             
 
 is present. 
 
         [0353]     This feature point set represents an image excerpt or an object of the image which was recorded at the instant t.  
         [0354]     The motion that maps O 1   N  onto the corresponding points of the image that was recorded at the instant t+1 is now sought.  
         [0355]     In contrast to the method described with reference to  FIG. 2 , the parameters of a perspective motion model are determined in step  904 .  
         [0356]     The motion model according to equation (36) has nine parameters but only eight degrees of freedom, as can be seen from the equation below.  
                     [           x   ⁡     (     t   +   dt     )                 y   ⁡     (     t   +   dt     )             ]     =     [                 a   1     ⁢     x   ⁡     (   t   )         +       a   2     ⁢     y   ⁡     (   t   )         +     a   3             n   1     ⁢     x   ⁡     (   t   )         +       n   2     ⁢     y   ⁡     (   t   )         +     n   3                         b   1     ⁢     x   ⁡     (   t   )         +       b   2     ⁢     y   ⁡     (   t   )         +     b   3             n   1     ⁢     x   ⁡     (   t   )         +       n   2     ⁢     y   ⁡     (   t   )         +     n   3               ]                 =       [               n   3       n   3       ⁢             a   1       n   3       ⁢     x   ⁡     (   t   )         +         a   2       n   3       ⁢     y   ⁡     (   t   )         +       a   3       n   3                 n   1       n   3       ⁢     x   ⁡     (   t   )         +         n   2       n   3       ⁢     y   ⁡     (   t   )         +   1                       n   3       n   3       ⁢             b   1       n   3       ⁢     x   ⁡     (   t   )         +         b   2       n   3       ⁢     y   ⁡     (   t   )         +       b   3       n   3                 n   1       n   3       ⁢     x   ⁡     (   t   )         +         n   2       n   3       ⁢     y   ⁡     (   t   )         +   1               ]     =     [                 a   1   ′     ⁢     x   ⁡     (   t   )         +       a   2   ′     ⁢     y   ⁡     (   t   )         +     a   3   ′             n   1   ′     ⁢     x   ⁡     (   t   )         +       n   2   ′     ⁢     y   ⁡     (   t   )         +   1                       b   1   ′     ⁢     x   ⁡     (   t   )         +       b   2   ′     ⁢     y   ⁡     (   t   )         +     b   3   ′             n   1   ′     ⁢     x   ⁡     (   t   )         +       n   2   ′     ⁢     y   ⁡     (   t   )         +   1             ]                     (   38   )             
 
         [0357]     The parameters of the perspective model can be determined like the parameters of the affine model by means of a least squares estimation by minimizing the term  
                 E   pers     ⁡     (       a   1   ′     ,     a   2   ′     ,     a   3   ′     ,     b   1   ′     ,     b   2   ′     ,     b   3   ′     ,     n   1   ′     ,     n   2   ′       )       =     ∑           ⁢     (         (         n   1   ′     ⁢         O   _     x   ′     ⁡     (   n   )         +       n   2   ′     ⁢         O   _     y   ′     ⁡     (   n   )         +   1     )     ⁢     (           O   _     x   ′     ⁡     (   n   )       +     d     n   ,   x         )     ⁢             -     (         a   1   ′     ⁢         O   _     y   ′     ⁡     (   n   )         +       a   2   ′     ⁢         O   _     y   ′     ⁡     (   n   )         +     a   3   ′       )       )     2       +       (         (         n   1   ′     ⁢         O   _     x   ′     ⁡     (   n   )         +       n   2   ′     ⁢         O   _     y   ′     ⁡     (   n   )         +   1     )     ⁢     (           O   _     y   ′     ⁡     (   n   )       +     d     n   ,   y         )       -     (         b   1   ′     ⁢         O   _     x   ′     ⁡     (   n   )         +       b   2   ′     ⁢         O   _     y   ′     ⁡     (   n   )         +     b   3   ′       )       )     2                   (   39   )             
 
         [0358]     In this case,  O ′ is defined in accordance with equation (27) analogously to the embodiment described with reference to  FIG. 2 .  
         [0359]      O ′ x (n) designates the first component of the n-th column of the matrix  O ′ and  O ′ y (n) designates the second component of the n-th column of the matrix  O ′.  
         [0360]     The minimum distance vector  
           D   _       min   ,     P     t   +   1     K         ⁡     (     x   ,   y     )         
 
 calculated in accordance with equation (29) is designated in abbreviated fashion as [d n,x d n,y ] T . 
 
         [0361]     The time index t has been omitted in formula (39) for the sake of simpler representation.  
         [0362]     Analogously to the method described with reference to  FIG. 2 , in which an affine motion model is used, the accuracy can be improved for the perspective model, too, by means of an iterative procedure.  
         [0363]      FIG. 10  illustrates a flow diagram  1000  of a determination of a perspective motion in accordance with an exemplary embodiment of the invention.  
         [0364]     Step  1001  corresponds to step  904  of the flow diagram  900  illustrated in  FIG. 9 .  
         [0365]     Steps  1002  to  1008  are analogous to steps  402  to  408  of the flow diagram  400  illustrated in  FIG. 4 .  
         [0366]     The difference lies in the calculation of the error E pers , which is calculated in accordance with equation (39) in step  1006 .  
         [0367]     Although specific embodiments have been illustrated and described herein, it will be appreciated by those of ordinary skill in the art that a variety of alternate and/or equivalent implementations may be substituted for the specific embodiments shown and described without departing from the scope of the present invention. This application is intended to cover any adaptations or variations of the specific embodiments discussed herein. Therefore, it is intended that this invention be limited only by the claims and the equivalents thereof.