System for Computing the Self-Motion of Moving Images Devices

The present invention pertains to imaging and particularly to three-dimensional interpretation of two-dimensional image sequences. More particularly, the invention pertains to determination of three-dimensional self-motion of moving imaging devices.

Visual information is an indispensable clue for the successful operation of an autonomous land vehicle. Even with the use of sophisticated inertial navigation systems, the accumulation of position error requires periodic corrections. Operation in unknown environments or mission tasks involving search, rescue, or manipulation critically depend upon visual feedback.

Assessment of scene dynamics becomes vital when moving objects may be encountered, e.g., when the autonomous land vehicle follows a convoy, approaches other vehicles, or has to detect moving threats. For the given case of a moving camera, such as one mounted on the autonomous land vehicle, image motion can supply important information about the spatial layout of the environment ("motion stereo") and the actual movements of the land vehicle.

Previous work in motion analysis has mainly concentrated on numerical approaches for the recovery of three-dimensional (3-D) motion and scene structure from two-dimensional (2-D) image sequences. The most common approach is to estimate 3-D structure and motion in one computational step by solving a system of linear or non-linear equations. This technique is characterized by several severe limitations. First, it is known for its notorious noise-sensitivity. To overcome this problem, some researchers have extended this technique to cover multiple frames. Secondly, it is designed to analyze the relative motion and 3-D structure of a single rigid object. To estimate the egomotion of an autonomous land vehicle (ALV), having the imaging device or camera, and the accompanying scene structure, the environment would have to be treated as a large rigid object. However, rigidness of the environment cannot be guaranteed due to the possible presence of moving objects in the scene. The consequence of accidentally including a moving 3-D point into the system of equations, representing the imaged environment, in the best case, would be a solution (in terms of motion and structure) exhibiting a large residual error, indicating some non-rigid behavior. The point in motion, however, could not be immediately identified from this solution alone. In the worst case (for some forms of motion), the system may converge towards a rigid solution (with small error) in spite of the actual movement in the point set. This again shows another (third) limitation: there is no suitable means of expressing the ambiguity and uncertainty inherent to dynamic scene analysis. The invention, that solves the aforementioned problems, is novel in two important aspects. The scene structure is not treated as a mere by-product of the motion computation but as a valuable means to overcome some of the ambiguities of dynamic scene analysis. The key idea is to use the description of the scene's 3-D structure as a link between motion analysis and other processes that deal with spatial perception, such as shape-from-occlusion, stereo, spatial reasoning, etc. A 3-D intepretation of a moving scene can only be correct if it is acceptable by all the processes involved.

Secondly, numeral techniques are largely replaced by a qualitative strategy of reasoning and modeling. Basically, instead of having a system of equations approaching a single rigid (but possibly incorrect) numerical solution, multiple qualitative interpretations of the scene are maintained. All the presently existing interpretations are kept consistent with the observations made in the past. The main advantage of this approach of the present invention is that a new interpretation can be supplied immediately when the currently favored interpretation turns out to be unplausible.

The problem of determining the motion parameters of a moving camera relative to its environment from a sequence of images is important for applications for computer vision in mobile robots. Short-term control, such as steering and braking, navigation, and obstacle detection/avoidance are all tasks that can effectively utilize this information.

The present invention deals with the computation of sensor platform motion from a set of displacement vectors obtained from consecutive pairs of images. It is directed for application to autonomous robots and land vehicles. The effects of camera rotation and translation upon the observed image are overcome. The new concept of "fuzzy" focus of expansion (FOE), which marks the direction of vehicle heading (and provides sensor rotation), is exploited. It is shown that a robust performance for FOE location can be achieved by computing a 2-D region of possible FOE-locations (termed "Fuzzy FOE") instead of looking for a single-point FOE. The shape of this FOE is an explicit indicator of the accuracy of the result. Given the fuzzy FOE, a number of very effective inferences about the 3-D scene structure and motion are possible and the fuzzy FOE can be employed as a practical tool in dynamic scene analysis. The results are realized in real motion sequences.

The problem of understanding scene dynamics is to find consistent and plausible 3-D interpretations for any change observed in the 2-D image sequence. Due to the motion of the autonomous land vehicle (ALV), containing the scene sensing device, stationary objects in the scene generally do not appear stationary in the image, whereas moving objects are not necessarily seen in motion. The three main tasks of the present approach for target motion detection and tracking are: (1) to estimate the vehicle's motion; (2) to derive the 3-D structure of the stationary environment; and (b 3) to detect and classify the motion of individual targets in the scene. These three tasks are interdependent. The direction of heading (i.e., translation) and rotation of the vehicle are estimated with respect to stationary locations in the scene. The focus of expansion (FOE) is not determined as a particular image location, but as a region of possible FOE-locations called the fuzzy FOE. We present a qualitative strategy of reasoning and modeling for the perception of 3-D space from motion information. Instead of refining a single quantitative description of the observed environment over time, multiple qualitative interpretations are maintained simultaneously. This offers superior robustness and flexibility over traditional numerical techniques which are often ill-conditioned and noise-sensitive. A rule-based implementation of this approach is discussed results on real ALV imagery are presented.

The system of the present invention tracks stationary parts of the visual environment in the image plane, using corner points, contour segments, region boundaries and other two-dimensional tokens as references. This results in a set of 2-D displacement vectors for the selected tokens for each consecutive pair of camera images. The self-motion of the camera is modeled as two separate rotations about horizontal and the vertical axes passing through the lens center and a translation in 3-D space. If the camera performs pure translation along a straight line in 3-D space, then (theoretically) all the displacement vectors extend through one particular location in the image plane, called the focus of expansion (FOE) under forward translation or focus of contraction (FOC) under backward translation. The 3-D vector passing through the lens center and the FOE (on the image plane) corresponds to the direction of camera translation in 3-D space.

The invention can provide the directions of instantaneous heading of the vehicle (FOE) within 1.degree. and self motions of moving imaging devices can be accurately obtained. This includes rotations of .+-.5.degree. or larger in horizontal and vertical directions. To cope with the problems of noise and errors in the displacement field, a region of possible FOE-locations (i.e., the fuzzy FOE) is determined instead of a single FOE.

In practice, however, imaging noise, spatial discretization errors, etc., make it impractical to determine the FOE as an infinitesimal image location. Consequently, a central strategy of our method is to compute a region of possible FOE-locations (instead of a single position) which produces more robust and reliable results than previous approaches.

FIG. 1 is a functional block diagram of the present invention.

FIG. 2 is a diagram showing an extended application of the invention to three-dimensional scene constructing.

FIG. 3 shows a camera model and corresponding coordinate system.

FIG. 4 illustrates a successive application of horizontal and vertical rotation of the camera.

FIG. 5 diagrams the effect of a focus of expansion location (FOE) for pure camera translation.

FIG. 6 illustrates the concept of the FOE for discrete time steps during the motion of a vehicle having the camera.

FIG. 7 shows the amount of expansion from the FOE for discrete time steps.

FIGS. 8a and 8b display a displacement field caused by horizontal and vertical rotation and translation of the camera and a derotated displacement field, respectively.

FIGS. 9a and 9b show an image plane and a rotation space, respectively.

FIGS. 10a and 10b show mappings of a polygon from rotation space to an image plane and vice versa, respectively.

FIGS. 11a-f illustrate a changing rotation polygon.

FIG. 12 illustrates an intersection of displacement vectors with a vertical line which, if moved, changes the variance of intersection.

FIG. 13 shows a displacement field used to evaluate various error functions.

FIGS. 14a-d reveal the standard deviation of intersection at a vertical cross section at position x for different amounts of vertical rotation.

FIGS. 15a-d reveal the standard deviation of intersection (square root) at a vertical cross section at position x for different amounts of vertical rotation with no horizontal rotation and with no pixel noise applied to the image locations.

FIGS. 16a-d reveal the standard deviation of intersection (square root) at a vertical cross section at position x for different amounts of vertical rotation with no horizontal rotation and with .+-.1 pixels of noise applied to the image locations.

FIGS. 17a-d reveal the standard deviation of intersection (square root) at a vertical cross section at position x for different amounts of vertical rotation with no horizontal rotation and with .+-.2 pixels of noise applied to the image locations.

FIGS. 18a-d show the location of minimum intersection standard deviation under varying horizontal rotation with the horizontal location of the FOE marked xf.

FIGS. 19a-d show the amount of minimum intersection standard deviation under varying horizontal rotation and with no noise added to the image locations.

FIGS. 20a-d show the amount of minimum intersection deviation under varying horizontal rotation and with .+-.2 pixels noise added to the image locations.

FIG. 21 illustrates intersecting displacement vectors with two vertical lines which lie on the same side of the FOE.

FIGS. 22a-d show the correlation coefficient for the intersection of displacement vectors at two vertical lines under varying horizontal rotations with no noise added.

FIGS. 23a-d show the correlation coefficient for the intersection of displacement vectors at two vertical lines under varying horizontal rotations and with .+-.2 pixels of noise added to the image locations.

FIG. 24 illustrates displacement vectors and measurement of their error.

FIG. 25 shows how to determine the optimum two-dimensional shift for a set of displacement vectors.

FIG. 26 reveals how FOE locations are dismissed if the displacement field resulting from the application of the optimal shift results in a vector not pointing away from the FOE.

FIGS. 27a-d illustrate the displacement field and minimum error at selected FOE locations.

FIGS. 28a-e display the effects of increasing the average length of displacement vectors upon the shape of the error function.

FIGS. 29a-e display the effects of increasing residual rotation in a horizontal direction upon the shape of the error function for relatively short vectors.

FIGS. 30a-e display the effects of increasing residual rotation in a vertical direction upon the shape of the error function for relatively short vectors.

FIGS. 31a-e show the effects of increasing residual rotation in horizontal and vertical directions upon the shape of the error function for relatively short vectors.

FIGS. 31f-j indicate the amount of optimal linear shift obtained under the same conditions in FIGS. 31a-e.

FIGS. 32a-e show the effects of the uniform noise applied to image point coordinates for a constant average vector length.

FIGS. 33a and 33b reveal the different effects of uniform noise applied to image point coordinates for shorter and longer average vector lengths, respectively.

FIG. 34 reveals a side view of a camera traveling parallel to a flat surface.

FIGS. 35a-i show an original image sequence taken from a moving vehicle after edge detection and point detection with the selected points located at the lower-left corners of their marks.

FIG. 35j-p show the original image sequence after edge detection and point selection.

FIG. 36a-p illustrate the displacement vectors and estimates of vehicles motion for the image sequence shown in FIGS. 35a-p.

FIG. 37 illustrates a specific hardware implementation of the embodiment of the invention.

FIGS. 1 and 2 reveal the main portions of the present invention. FIG. 2 also expands on 3-D interpretations of 2-D images in items 118 and 120. Item 116 incorporates portions of items 122 and 124 of FIG. 1. First, significant features (points, boundaries, corners, etc.) are extracted by feature extraction and tracking 114 from the image data 112 and the 2-D displacement vectors are computed for this set of features. For the examples shown here, points were selected and tracked between individual frames by item 114. Automatic techniques suitable for this task are in the related art. In the second step, the vehicle's direction of translation, i.e. the focus of expansion (FOE), and the amount of rotation in space are determined. The effects of vehicle motion on the FOE computation is described below. Nearly all the necessary numerical computation is perforated in the FOE computation stage, which also is described below. The third step (at 2-D change analysis 118) constructs an internal 3-D model of the scene. Also disclosed are the concepts and operation of the qualitative scene model 120. Experiments with the present invention on real imagery taken from a moving ALV are discussed below.

System 110 contains the following three main components--token tracking 114, FOE seeker 124 and optimal derotation 122--in FIG. 1. The 2-D displacement vectors for selected image tokens are determined in the first stage (token tracking 114). Since those original displacement vectors are caused by some arbitrary and (at this point) unknown camera motion, including camera rotations, they do not yet exhibit the characteristic radial pattern of pure camera translation.

The second component (FOE seeker 124) selects a set of candidate locations for the FOE and forms a connected image region of feasible FOE-locations plus the range of corresponding camera rotations, based on the results from the third component (optimal derotation 122). A particular FOE-location is feasible, if the corresponding error value (computed by the optimal derotation module 122) is below some dynamically adjusted threshold. The size of the final FOE-region reflects the amount of uncertainty contained in the visual information (large regions reflect high uncertainty).

The third component (optimal derotation 122) determines the optimal 3-D camera rotations for a particular (hypothesized) FOE-location. This is accomplished by simulating the effects of reverse camera rotations upon the given set of displacement vectors. The camera is virtually rotated until the modified displacement field is closest to a radial expansion pattern with respect to the selected FOE-location. Module 122 returns the necessary amount of reverse rotation and the deviation from a radial displacement field (i.e., an error value).

Component 123, comprising FOE seeker 124 and optimal derotation 122, represents the fuzzy FOE means which outputs the region of possible FOE locations.

The first step of the present invention is to estimate the vehicle's motion relative to the stationary environment using visual information. Arbitrary movement of an object in 3-D space and thus the movement of the vehicle itself can be described as a combination of translation and rotation. While knowledge about the composite vehicle motion is essential for control purposes, only translation can supply information about the spatial layout of the 3-D scene (motion stereo). This, however, requires the removal of all image effects results from vehicle rotation. For this purpose, we discuss the changes upon the image that are caused by individual application of the "pure" motion components.

It is well-known that any rigid motion of an object in space between two points in time can be decomposed into a combination of translation and rotation. While many researchers have used a velocity-based formulation of the problem, the following treatment views motion in discrete time steps.

The viewing geometry involves a world coordinate system, (XYZ) as illustrated in FIG. 3. FIG. 3 shows the camera-centered coordinate system 130, lens center 126, image plane 128, and angles .phi., .theta. and .psi. of rotation. The origin 0 of coordinate system 130 is located at lens center 126. Focal length f is the distance between lens center 126 and image plane 128. Each 3-D point (XYZ) may be mapped onto image location (X,Y). Angles .phi., .theta. and .psi. specify angles of camera rotation about the X,Y and Z axes, respectively. Given the world coordinate system, a translation T=(U V W).sup.T applied to a point in 3-D X=(X Y Z).sup.T is accomplished through vector addition: ##EQU1## A 3-D rotation R about an arbitrary axis through the origin of coordinate system 130 can be described by successive rotations about its three axes: ##EQU2##

A general rigid motion in space consisting of translation and rotation is described by the transformation

Its six degrees of freedom are U, V, W, .phi., .theta., .psi..

This decomposition is not unique because the translation could be as well applied after the rotation. Also, since the multiplication of the rotation matrices is not commutative, a different order of rotations would result in different amounts of rotation for each axis. For a fixed order of application, however, this motion decomposition is unique.

To model the movements of the vehicle, the camera is considered as being stationary and the environment as being moving as one single rigid object relative to the camera. The origin 0 of coordinate system 130 is located in the lens center 126 of the camera.

The given task is to reconstruct the vehicle's or the camera's egomotion from visual information. It is therefore necessary to know the effects of different kings of vehicle or camera motion upon the camera image. Under perspective imaging, a point in space X=(X Y Z).sup.T is projected onto a location on the image plane X=(X Y).sup.T such that ##EQU3## where f is the focal length of the camera (see FIG. 3).

The effects of pure camera rotation are accounted for. Ignoring the boundary efforts when the camera is rotated around its lens center 126, the acquired image changes but no new views of the environment are obtained. Pure camera rotation merely maps the image into itself. The most intuitive effect results from pure rotation about the Z-axis of the camera-centered coordinate system 130, which is also the optical axis. Any point in the image moves along a circle centered at the image location x=(0 0 ). In practice, however, the amount of rotation .psi. of the vehicle about the Z-axis is small. Therefore, vehicle rotation is confined to the X- and Y-axis, where significant amounts of rotation occur.

The vehicle or camera undergoing rotation about the X-axis by an angle - .phi. and the Y-axis by an angle - .theta. moves each 3-D point X to point X' relative to the camera.

It is important to notice that this transformation contains no 3-D variables and is therefore a mapping of the image onto itself. This demonstrates that no additional information about the 3-D structure of the scene can be obtained under pure camera rotation.

An interesting property of this mapping should be mentioned at this point, which might not be obvious. Moving an image point on a diagonal passing through the center of the image at 45.degree. by only rotating the camera does not result in equal amounts of rotation about the X- and the Y-axis. This is again a consequence of the successive application of the two rotations R.sub..theta. and R.sub..phi., since the first rotation about the Y-axis also changes the orientation of the camera's X-axis in 3-D space. It also explains why the pair of equations in (7) is not symmetric with respect to .theta. and .phi..

In measuring the amount of camera rotation, the problem to be solved is the following: Given are two image locations x.sub.0 and x.sub.1, which are the observations of the same 3-D point at time t.sub.0 and time t.sub.1. The question here is the amount of rotation R.sub.100 and R.sub.74 which applied to the camera between time t.sub.0 and time t.sub.1, would move image point x.sub.0 onto x.sub.1 assuming that no camera translation occurred at the same time. If R.sub..phi. and R.sub.74 are applied to the camera separately, the points in the image move along hyperbolic paths. If pure horizontal rotation were applied to the camera, a given image point x.sub.0 would move on a path described by ##EQU8## Similarly pure vertical camera rotation would move an image point x.sub.1 along ##EQU9##

Since the 3-D rotation of the camera is modeled as being performed in two separate steps (R.sub..theta. followed by R.sub..phi.), the rotation mapping r.sub..phi. r.sub..theta. can also be separated into r.sub..theta. followed by r.sub..phi.. In the first step, applying pure (horizontal) rotation around the Y-axis r.sub..theta., point x.sub.0 is moved to an intermediate image location x.sub.c. The second step, applying pure (vertical) rotation around the X-axis r.sub..phi., takes point x.sub.c to the final image location x.sub.1. This can be expressed as ##EQU10##

FIG. 4 reveals a successive application of horizontal and vertical rotation. Image point x.sub.0 is to be moved to location x.sub.1 by pure horizontal and vertical camera rotation. Horizontal rotation (about the Y-axis) is applied first, moving x.sub.0 to x.sub.c, which is the intersection point of the two hyperbolic paths for horizontal and vertical rotation. In a second step, x.sub.c is taken to x.sub.1. Then the two rotation angles .theta. and .phi. are found directly. Further in FIG. 4, the image point x.sub.c =(x.sub.c y.sub.c) is the intersection point of the hyperbola passing through x.sub.0 resulting from horizontal camera rotation (10) with the hyperbola passing through x.sub.1 resulting from vertical camera rotation (11). Intersecting the two hyperbolae gives the image point x.sub.c, with ##EQU11##

The amount of camera rotation necessary to map x.sub.0 onto .sub.1 by applying R.sub..theta. followed by R.sub..phi. is finally obtained as ##EQU12##

When the vehicle or camera undergoes pure translation between time t and time t', every point on the vehicle is moved by the same 3-D vector T=(U V W).sup.T. Again, the same effect is achieved by keeping the camera fixed and moving every point X.sub.i in the environment to X.sub.i ' by applying -T.

Since every stationary point in the environment undergoes the same translation relative to the camera, the imaginary lines between corresponding points X.sub.i X.sub.i ' are parallel in 3-D space.

It is a fundamental result from perspective geometry that the images of parallel lines pass through a single point in the image plane called a "vanishing point." When the camera moves along a straight line, every (stationary) image point seems to expand from this vanishing point or contract towards it when the camera moves backwards. This particular image location is therefore commonly referred to as the focus of expansion (FOE) or the focus of contraction (FOC). Each displacement vector passes through the FOE creating the typical radial expansion pattern shown in FIG. 5. FIG. 5 reveals the location of the FOE. With pure vehicle translation, points in the environment (A,B) move along 3-D vectors parallel to the vector pointing from lens center 126 to the FOE in camera plane 128 (FIG. 3). These vectors form parallel lines in space which have a common vanishing point (the FOE) in the perspective image.

As can be seen in FIG 5, the straight line passing through the lens center of the camera and the FOE is also parallel to the 3-D displacement vectors. Therefore, the 3-D vector OF points in the direction of camera translation in space. Knowing the internal geometry of the camera (i.e., the focal length f), the direction of vehicle translation can be determined by locating the FOE in the image. The actual translation vector T applied to the camera is a multiple of the vector OF which supplies only the direction of camera translation but not its magnitude. Therefore,

Since most previous work incorporated a velocity-based model of a 3-D motion, the focus of expansion has commonly been interpreted as the direction of instantaneous heading, i.e., the direction of vehicle translation during an infinitely short period in time. When images are given as "snapshots" taken at discrete instances of time, the movements of the vehicle must be modeled accordingly as discrete movements from one position in space to the next. Therefore, the FOE cannot be interpreted as the momentary direction of translation at a certain point in time, but rather as the direction of accumulated vehicle translation over a period of time.

FIG. 6 illustrates the concept of FOE for discrete time steps. The motion of vehicle 132 between two points in time can be decomposed into a translation followed by a rotation. The image effects of pure translation (FOE.sub.a) are observed in image I.sub.o. FIG. 6 shows the top view of vehicle 132 traveling along a curved path at two instances in time t.sub.0 and t.sub.1. The position of the vehicle 132 in space is given by the position of a reference point on the vehicle P and the orientation of the vehicle is .OMEGA.. FIG. 6 also displays the adopted scheme of 3-D motion decomposition. First, the translation T is applied which shifts the vehicle's reference point (i.e., lens center 126 of the camera) from position P.sub.0 to position P.sub.1 without changing the vehicle's orientation .OMEGA.. The 3-D translation vector T intersects image plane 128 at FOE.sub.a. In the second step the vehicle is rotated b .omega. to the new orientation .OMEGA..sub.1. Translation T transforms image I.sub.0 into image I'.sub.1, which again is transformed into I.sub.1 by rotation .omega.. The important fact is that FOE.sub.a is observed at the transition from image I.sub.0 to image I'.sub.1, which is obtained by derotating image I.sub.1 by -.omega.. Throughout the present description , this scheme (FIG. 6) is used as a model for vehicle or camera motion.

The amount of camera translation can be measured. FIG. 7 shows the geometric relationships for the 2-D case. The amount of expansion from the FOE for discrete time steps is illustrated. FIG. 7 can be considered as a top view of the camera, i.e., a projection onto the X/Z-plane of the camera-centered coordinate system 130. The cross section of the image plane is shown as a straight line. The camera moves by a vector T in 3-D space, which passes through lens center 126 and the FOE in camera plane 128. The 3-D Z-axis is also the optical axis of the camera. The camera is translating from left to right in the direction given by T=(X.sub.f f).sup.T.

A stationary 3-D point is observed at two instances of time, which moves in space relative to the camera from X to X', resulting in two images x and x'. ##EQU13## Using the inverse perspective (8) transformation yields ##EQU14## From similar triangles (shaded in FIG. 7) ##EQU15## and therefore ##EQU16##

Thus, the rate of expansion of image points from the FOE contains direct information about the distance of the corresponding 3-D points from the camera. Consequently, if the vehicle is moving along a straight line and the FOE has been located, the 3-D structure of the scene can be determined from the expansion pattern in the image. However, the distance Z of a 3-D point from the camera can only be obtained up tot he scale factor .DELTA.Z, which is the distance that the vehicle advanced along the Z-axis during the elapsed time.

When the velocity of the vehicle (.DELTA.Z /t) in space is known, the absolute range of any stationary point can be computed. Alternatively, the velocity of the vehicle can be obtained if the actual range of a point in the scene is known (e.g., from laser range data). In practice, of course, any such technique requires that the FOE can be located in a small area, and the observed image points exhibit significant expansion away from the FOE. As shown below, imaging noise and camera distortion pose problems in the attempt to assure that both of the above requirements are met.

If a set of stationary 3-D points {(X.sub.i X.sub.i ')}is observed, then of course the translation in the Z-direction is the same for every point.

Therefore, the range of every point is proportional to the observed amount of expansion of its image away from the FOE ##EQU17## which renders the relative 3-D structure of the set of points.

The effects of camera translation T can be formulated as a mapping t of a set of image locations {x.sub.i } into another set of image locations {x'.sub.i }. Unlike in the case of pure camera rotation, this mapping not only depends upon the 3- D translation vector but also upon the actual 3-D location of each individual point observed. Therefore, in general, t is not simply a mapping of the image onto itself. However, one important property of t can be described exclusively in image plane 128, namely that each point must map onto a straight line passing through the original point and one unique location in the image (the FOE). This means that if vehicle 132 is undergoing pure translation , then there must exist an image location x.sub.f such that the mapping t satisfies the condition radial-mapping t (x.sub.f, I, I'): ##EQU18##

When the vehicle is not undergoing pure translation or rotation but combined 3-D motion (i.e., translation and rotation) of the form R.sub..phi. R.sub..theta. T, the effects in the image are described by a transformation d (for displacement) which is a combination f r.sub..phi. r.sub..theta. and t:

wherein I={x.sub.i }, I'={x'.sub.i } are the two sets of corresponding image points. FIG. 8a and 8b show a typical displacement field for a camera undergoing horizontal and vertical rotation as well as translation. The points x.sub.i .epsilon. I are marked with small circles 138. Rectangle 134 marks the area of search for the FOE. The derotated displacement is illustrated in FIG. 8b with the FOE marked by circle 136.

By decomposing a composite displacement field d into its three components r.sub..phi., r.sub..theta., and t, the vehicle's rotation and direction of translation in space can be computed from the information available in the image. This problem is addressed below. As discussed in the previous section, the 3-D motion M of the vehicle is modeled by a translation T followed by a rotation R.sub..theta. about the Y-axis and a rotation R.sub..phi. about the X-axis:

This results in a mapping d from the original image I.sub.O at time T.sub.0 into the new image I.sub.1 at time t.sub.1.

The intermediate image I'.sub.0 (26) is the result of the translation component of the vehicle's motion and has the property of being a radial mapping (23). Unlike the two images I.sub.0 and I.sub.1, which are actually given, the image I'.sub.0 is generally not observed, except when the camera rotation is zero. It serves as an intermediate result to be reached during the separation of translational and rotational motion components.

The question at this point is whether there exists more than one combination of rotation mappings r.sub..phi. and r.sub..theta. which would satisfy this requirement, i.e., if the solution is unique. It has been pointed out above that the decomposition of 3-D motion into R.sub..phi., R.sub..theta., R.sub..psi., and t is unique for a fixed order of application. This does not imply, however, that the effects of 3-D motion upon the perspective image are unique as well.

Related art has shown that seven points in two perspective views sufficient to obtain a unique interpretation in terms of rigid body motion and structure, except for a few cases where points are arranged in some very special configuration in space. Further art reports computer experiments which suggest that six points are sufficient in many cases and seven or eight points yield unique interpretations in most cases.

Due to its design and the application, however, the motion of a typical autonomous land vehicle (ALV) in space is quite restricted. The vehicle can only travel upright on a surface and its large wheelbase allows for only relatively small changes in orientation. It is also heavy and thus exhibits considerable inertia. Therefore, the final motion parameters must lie within a certain narrow range and it can be expected that a unique solution can be found even in cases when the number of points is near or above the minimum.

The fact that

suggests two different strategies for separating the motion components: (1) FOE from rotation--40 successively apply combinations of inverse rotation mappings

to the second image I.sub.1, until the resulting image I' is a radial mapping with respect to the original image I.sub.0. The locate the FOE x.sbsb.fk in I.sub.0 and (2) rotation from FOE --successively select FOE--locations (different directions of vehicle translation) X.sbsb.f1, X.sbsb.f2, . . . X.sbsb.f1 in the original image I.sub.0 and then determine the inverse rotation mapping r.sbsb..theta.1.sup.-1 r.sbsb..phi.1.sup.-1 that yields a radial mapping with respect to the given FOE x.sbsb.f1 in the original image I.sub.0.

Both alternatives were investigated under the assumption of restricted, but realistic vehicle motion, as stated earlier. It turned out that the major problem in the FOE-from-rotation approach is to determine if a mapping of image prints is (or is close to being) radial when the location of the FOE is unknown. Of course, in the presence of noise, this problem becomes even more difficult. The second approach was examined after it appeared that any method which extends the given set of displacement vectors backwards to find the FOE is inherently sensitive to image degradations.

Although there have been a number of suggestions for FOE-algorithms in the past, no results of implementations have been demonsstrated on real outdoor imagery. One reason for the absence of useful results might be that most researchers have tried to locate the FOE in terms of a single, distinct image location. In practice, however, the noise generated by merely digitizing a perfect translation displacement field may keep the resulting vectors from passing through a single pixel. Even for human observers it seems to be difficult to determine the exact direction of heading (i.e., the location of the FOE on the retina). Average deviation of human judgement from the real direction has been reported to be as large as 10.degree. and up to 20.degree. in the presence of large rotations. It was, therefore, an important premise in this work that the final algorithm should determine an area of potential FOE-locations (called the "fuzzy FOE") instead of a single (but probably incorrect) point.

The FOE may be obtained from rotation. In this method, the image motion is decomposed in two steps. First, the rotational components are estimated and their inverses are applied to the image, thus partially "derotating" the image. If the rotation estimate was accurate, the resulting displacement field after derotation would diverge from a single image location (the FOE). The second step verifies that the displacement field is actually radial and determines the location of the FOE. For this purpose, two problems have to be solved: (1) how to estimate the rotational motion components without knowing the exact location of the FOE, and (2) how to measure the "goodness of derotation" and locate the FOE.

The rotational components can be estimated. Each vector in the displacement field is the sum of vector components caused by camera rotation and camera translation. Since the displacement caused by translation depends on the depth of the corresponding points in 3-D space (equation 18), points located at a large distance from the camera are not significantly affected by camera translation. Therefore, one way of estimating vehicle rotation is to compute .theta. and .phi.from displacement vectors which are known to belong to points at far distance. Under the assumption that those displacement vectors are only caused by rotation, equations 14 and 15 can be applied to find the two angles. In some situations, distant points are selected easily. For example, points on the horizon are often located at a sufficient distance from the vehicle. Image points close to the axis of translation would be preferred because they expand from the FOE slower than other points at the same depth. However, points at far distances may not always be available or may not be known to exist in the image. In those cases, the following method for estimating the rotational components can be used. The design of the ALV (and most other mobile robots) does not allow rapid changes in the direction of vehicle heading. Therefore, it can be assumed that the motion of the camera between two frames is constrained, such that the FOE can change its location only within a certain range. If the FOE was located in one frame, the FOE in the subsequent frame must lie in a certain image region around the previous FOE location. FIG. 9a shows an image plane which illustrates this situation. The FOE of the previous frame was located at the center of square 140 which outlines the region of search for the current FOE, thus the FOE in the given frame must be inside square 140. Three displacement vectors are shown P1.fwdarw.P1', P2.fwdarw.P2', P3.fwdarw.P3'). The translational components (P1.fwdarw.Q1,P2.fwdarw.Q2, P3.fwdarw.Q3) of those displacement vectors and the FOE (inside square 140) are not known at this point in time but are marked in FIG. 9a.

The main idea of this technique is to determine the possible range of camera rotations which would be consistent with the FOE lying inside marked region 140. Since the camera rotates about two axes, the resulting range of rotations can be described as a region in a 2-D space. FIG. 9b shows this rotation space with the two axes .theta.and .phi.corresponding to the amount of camera rotation around the Y-axis and the X-axis, respectively. The initial rotation estimate is a range of .+-.10.degree. in both directions which is indicated by a square 142 in rotation space of FIG. 9b.

In general, the range of possible rotations is described by a closed, convex polygon in rotation space. A particular rotation (.theta.', .phi.') is possible if its application to every displacement vector (i.e., to its endpoint) yields a new vector which lies on a straight line passing through the maximal FOE-region. The region of possible rotations is successively constrained by applying the following steps for every displacement vector (FIG. 10a and 10b). First apply the rotation mapping defined by the vertices of the rotation polygon to the endpoint P' of the displacement vector P.fwdarw.P' 146. This yields a set of image points P.sub.i. Second, connect the points P.sub.i to a closed polygon 148 in the image. Polygon 148 is similar to the rotation polygon 144 but distorted by the nonlinear rotation mapping as shown in FIG. 10a. Third, intersect polygon 146 in the image with open triangle 150 formed by the starting point P of the displacement vector 146 and defined by two tangents 152 and 154 onto the maximal FOE-region. Rotations that would bring the endpoint of the displacement vector outside triangle 150 are not feasible. The result is a new (possibly empty) polygon 156 in the image plane. Fourth, new polygon 156 from the image plane back into the rotation space of FIG. 9b. Fifth, rotation polygon 158 is empty (number of vertices is zero), then stop. No camera rotation is possible that would make all displacement vectors intersect the given FOE-region. Repeat the process using a larger FOE-region.

FIGS. 11a, 11b and 11c show the changing shape of the rotation polygon during the application of this process to the three displacement vectors in FIG. 8.

Since the mapping from rotation space to the image plane is nonlinear (equation 9), the straight lines between vertices in the rotation polygon 160 do not correspond to straight lines in the image. They are, however, approximated as straight lines in order to simplify the intersection with the open triangle. The dotted lines in the image plane show the actual mapping of the rotation polygon onto the image. It can be seen that the deviations from straight lines are small and can be neglected. FIG. 11a shows rotation polygon after examining displacement vector P1.fwdarw.P1'. Any camera rotation inside the polygon would move the endpoint of the displacement vector (P1') into the open triangle formed by targets 164 and 166 through P1 to the maximal FOE-region given by square 68 in the image plane. The actual mapping of the rotation polygon into the image plane is shown with a dotted outline.

FIG. 11b reveals the rotation polygon 170 after examining displacement vectors P1.fwdarw.P1' and P2.fwdarw.P2'.

FIG. 11(c) shows the final rotation polygon after examining the three displacement vectors P1.fwdarw.P1', P2.fwdarw.P2'and P3.fwdarw.P3'. The amount of actual camera rotation (.theta.=-2.0.degree., .phi.=5.0.degree.) is marked with a small circle (arrow) 172.

Increasing the number of displacement vectors improves the rotation estimate. In practice, the amount of camera rotation can be constrained to a range of below 1.degree. in both directions. Rotation can be estimated more accurately when the displacement vectors are short, i.e., when the amount of camera translation is small. This is in contrast to estimating camera translation which is easier with long displacement vectors.

The situation when the rotation polygon becomes empty requires some additional considerations. As mentioned earlier, in such a case no camera rotation is possible that would make all displacement vectors pass through the given FOE-region. This could indicate one of the two alternatives. First, at least one of the displacement vectors belongs to a moving object. Second, the given FOE-region does not contain the actual location of the FOE, i.e., the region is not feasible. The latter case is of particular importance. If a region can be determined not to contain the FOE, then the FOE must necessarily lie outside this region. Therefore, the above method can not only be used to estimate the amount of camera rotation, but also to search for the location of the FOE. Unfortunately, if the rotation polygon does not become empty, this does not imply that the FOE is actually inside the given region. It only means that all displacement vectors would pass through this region, not that they have a common intersection inside this region. However, if not all vectors pass through a certain region, then this region cannot possibly contain the FOE. The following recursive algorithm searches a given region for the FOE by splitting it into smaller pieces (divide-and-conquer):

This algorithm searches for the smallest feasible FOE-region by systematically discarding sub regions from further consideration. For the case that the shape of the original region is a square, subregions can be obtained by splitting the region into four subsquares of equal size. The simple version shown here performs a depth-first search down to the smallest subregion (limited by the parameter "min-size"), which is neither the most elegant nor the most efficient approach. The algorithm can be significantly improved by applying a more sophisticated strategy, for example, by trying to discard subregions around the perimeter first before examining the interior of a region. Two major problems were encountered with the latter method. First, the algorithm is computationally expensive since the process of computing feasible rotations must be repeated for every subregion. Second, a small region is more likely to be discarded than a larger one. However, when the size of the region becomes too small, errors induced by noise, distortion or point-tracking may prohibit displacement vectors from passing though a region which actually contains the FOE.

Although this algorithm is not employed in the further treatment, it suggests an interesting alternative which departs significantly from traditional FOE-algorithms. Its main attractiveness is that it is inherently region-oriented in contrast to most other techniques which search for a single FOE-location. For the purpose of estimating the amount of rotation, the method using points at far distance mentioned earlier is probably more practical. Two other alternatives for locating the FOE once the rotation components have been estimated are discussed in the following.

Locating the FOE in a partially derotated image may be attempted. After applying a particular derotation mapping to the displacement field, the question is how close the new displacement field is to a radial mapping, where all vectors diverge from one image location. If the displacement field is really radial, then the image is completely derotated and only the components due to camera translation remain. Two different methods for measuring this property are noted. One method used the variance of intersection at imaginary horizontal and vertical lines. The second method computes the linear correlation coefficient to measure how "radial" the displacement field is. The variance of intersection in related art suggests to estimate the disturbance of the displacement field by computing the variance of intersections of one displacement vector with all other vectors. If the intersections lie in a small neighborhood, then the variance is small, which indicates that the displacement field is almost radial. The problem can be simplified by using an imaginary horizontal and vertical line instead, whose orientation is not affected by different camera rotations. FIG. 12 shows 5 displacement vectors P.sub.1 .fwdarw.P' . . . P.sub.5 ' intersecting a vertical line at x at y.sub.1 . . . Y.sub.5. Moving the vertical line from x towards x.sub.0 will bring the points of intersection closer together and will thus result in a smaller variance of intersection. The point of intersection of a displacement vector P.sub.i .fwdarw.P.sub.i ' with a vertical line at x is given by ##EQU19## The variance of intersection of all displacement vectors with the vertical line at position x is ##EQU20##

To find the vertical cross section with minimum intersection variance, the first derivative of (29) with respect to x is set to zero. The location x.sub.0 of minimum intersection variance is then obtained. Similarly, the position of a horizontal cross section with minimal intersection variance can be obtained.

The square root of the variance of intersection (standard deviation) at a vertical line was evaluated on the synthetic displacement field shown in FIG. 13. The actual FOE is located in the center of the image. Square 174 around the center (.+-.100 pixels in both directions) marks the region over which the error functions are evaluated.

FIGS. 14a, b, c and d show the distribution of the intersection standard deviation for increasing residual rotations in vertical direction in the absence of noise. The horizontal rotation is 1.degree. in all cases represented by these Figures. Locations of displacement vectors are presented by real numbers (not rounded to integer values). In FIG. 14a, no residual rotation exists, i.e., the displacement field is perfectly radial. The value of the horizontal position of the cross section varies .+-.100 pixels around the actual FOE. The standard deviation is zero for x=x.sub.f (the x-coordinate of the FOE) and increases linearly on both sides of the FOE. In FIGS. 14b-d, the residual vertical rotation is increased from 0.2.degree. to 1.0.degree.. The bold vertical bar marks the horizonal position of minimum standard deviation, the thin bar marks the location of the FOE (X.sub.f). It can be seen that the amount of minimum standard deviation rises with increasing disturbance by rotation, but that the location of minimum standard deviation does not necessarily move away from the FOE.

FIGS. 15-17 show the same function under the influence of noise. In FIG. 15a-d, noise was applied except that by merely rounding the locations of displacement vectors to their nearest integer values. These Figures show standard deviation of intersection (square root) at a vertical cross section at position x for different amounts of vertical rotation with no horizontal rotation. Uniform noise of .+-.1 pixels was added to the image locations in FIGS. 16a-d. In FIGS. 17a-d, uniform noise of .+-.2 pixels was applied to the image locations. It can be seen that the effects of noise are similar to the effects caused by residual rotation components. The purpose of this error function is to determine where the FOE is located, and how "radial" the current displacement field is.

If the displacement field is already perfectly derotated, then the location of minimum intersection standard deviation is the location of the FOE. Ideally, all vectors pass through the FOE, such that a cross section through the FOE yields zero standard deviation. The question is how well the FOE can be located in an image which is not perfectly derotated. FIG. 18a-d plot the location of minimum intersection standard deviation under varying horizontal rotation. The vertical rotation is kept fixed for each plot. Horizontal camera rotations from -1.degree. to +1.degree. are shown on the abscissa (rot). The ordinate (x0) gives the location of minimum standard deviation in the range of .+-.100 pixels around the FOE (marked xf). The location of minimum standard deviation depends strongly on the amount of horizontal rotation.

A problem is that the location of minimum standard deviation is not necessarily closer to the FOE when the amount of rotation is less. The function is only well behaved in a narrow range around zero rotation, which means that the estimate of the camera rotation must be very accurate to successfully locate the FOE. The second purpose of this error function is to measure how "radial" the displacement field is after partial derotation. This should be possible by computing the amount of minimum intersection standard deviation. Intuitively, a smaller amount of minimum intersection standard deviation should indicate that the displacement field is less disturbed by rotation. FIGS. 19a-d and 20a-d show that this is generally true by showing the amount of minimum intersection standard deviation under varying horizontal rotation. For the noise-free case in FIG. 19a, the amount of minimum intersection standard deviation becomes zero in the absence of horizontal and vertical rotations, indicating that the derotation is perfect. Unfortunately, the function is not well behaved even in this relatively small range of rotations (.+-.1.0.degree.). The curve exhibits some sharp local minima where an algorithm searching of an optimal derotation would get trapped easily. FIGS. 19a-d show the same function in the presence of .+-.2 pixels of uniform noise added to image locations.

The second method, utilizing linear correlation, of measuring how close a displacement field is to a radial pattern again uses the points y.sub.11 -y.sub.15 and y.sub.21 -y.sub.25 of intersection at vertical (or horizontal) lines x.sub.1 and x.sub.2 as illustrated in FIG. 21. The displacement vectors P.sub.1 .fwdarw.P'.sub.1 through P.sub.5 .fwdarw.P'.sub.5 are intersected by two vertical lines X.sub.1 and X.sub.2, both of which lie on the same side of the FOE that is within area 176. Since the location of the FOE is not known, the two lines X.sub.1 and X.sub.2 are simply located at a sufficient distance from any possible FOE-location. This results in two sets of intersection points {(x.sub.1, y.sub.1i)} and {(x.sub.2, y.sub.2i)}. If all displacement vectors emanate from one single image location, then the distances between corresponding intersection points in the two sets must be proportional, i.e, ##EQU21## Therefore, a linear relationship exists between the vertical coordinates of intersection points on these two lines. The "goodness" of this linear relationship is easily measured by computing the correlation coefficient for the y-coordinates of the two sets of points. The resulting coefficient is a real number in the range from -1.0 to +1.0. If both vertical lines are on the same side of the FOE, then the optimal value is +1.0. Otherwise, if the FOE lies between the two lines, the optimal coefficient is -1.0. The horizontal position of the two vertical lines is of no importance, as long as one of these conditions is satisfied. For example, the left and right border lines of the image can be used.

FIGS. 22a-d and 23a-d show plots for the correlation coefficient for intersection of displacement vectors at two vertical lines under varying horizontal rotations, under the same conditions as in FIGS. 19a-d and 20a-d. No noise was applied for FIG. 22a-d. In FIGS. 23a-d, a uniform noise of .+-.2 pixels was added to the image locations. The optimal coefficient is +1.0 (horizontal axis) in FIGS. 22a-d. The shapes of the curves of FIGS. 22a-d and 23a-d are similar, respectively, to FIGS. 19a-d and 20a-d for the minimum standard deviations shown above earlier, with peaks at the same locations. It is apparent, however, that each curve has several locations where the coefficient is close to the optimum value (+1.0), i.e., no distinct global optimum exists which is not only the case in the presence of noise (FIGS. 23a-d). This fact makes the method of maximizing the correlation coefficient useless for computing the FOE.

The main problem encountered in computing the FOE from rotation, just described above, is that none of the functions examined was well behaved, making the search for an optimal derotation and the location of the FOE difficult. Disturbances induced by noise and residual rotation components are amplified by extending short displacement to straight lines and computing their intersections. The method, rotation from FOE, described below avoids this problem by guessing an FOE-location first and estimating the optimal derotation for this particular FOE in the second step.

Given the two images I.sub.0 and I.sub.1 of corresponding points, the main algorithmic steps of this approach are: (1) Guess an FOE-location x.sub.f.sup.(i) in image I.sub.0 (for the current iteration i); (2) Determine the derotation mapping r.sub..theta..sup.-1, r.sub..phi..sup.-1 which would transform image I.sub.1 into an image I'.sub.1 such that the mapping (x.sub.f.sup.(i), I.sub.0, I'.sub.1) deviates from a radial mapping (equation 23) with minimum error E.sup.(i) ; and (3) Repeat steps (1) and (2) until an FOE-location x.sub.f.sup.(k) with the lowest minimum error E.sup.(k) is found.

An initial guess for the FOE-location is obtained from knowledge about the orientation of the camera with respect to the vehicle. For subsequent pairs of frames, the FOE-location computed from the previous pair can be used as a starting point.

Once a particular x.sub.f has been selected, the problem is to compute the rotation mappings r.sub..theta..sup.-1 and r.sub..phi..sup.-1 which, when applied to the image I.sub.1, will result in an optimal radial mapping with respect to I.sub.0 and x.sub.f.

To measure how close a given mapping is to a radial mapping, the perpendicular distances between points in the second image (x'.sub.i) and the "ideal" displacement vectors is measured. The "ideal" displacement vectors line on straight lines passing through the the FOE x.sub.f and the points in the first image x.sub.i (see FIG. 24) which illustrates measuring the perpendicular distance d.sub.i between lines from x.sub.f through points x.sub.i in the second image. The sum of the squared perpendicular distances d.sub.i is the final error measure. For each set of corresponding image points (x.sub.i .epsilon.I, x'.sub.i .epsilon.I'), the error measure is defined as ##EQU22##

In the following, it is assumed that the amount of residual image rotation in horizontal and vertical direction is moderately small (less than 4.degree.). In most practical cases, this condition is satisfied, provided that the time interval between frames is sufficiently small. However, should the amount of vehicle rotation be very large for some reason, a coarse estimate of the actual rotation can be found (as described above) and applied to the image before the FOE computation. With small amounts of rotation, the actual rotation mapping, where points move on horizontal and vertical hyperbolic paths, can be approximated by a horizontal and vertical shift with constant length over the entire image. Under this condition, the inverse rotation mapping r.sub..phi..sup.-1, r.sub..theta..sup.-1 can be approximated by a adding a constant vector s=(s.sub.x s.sub.y) which is independent of the image location:

Given two images I and I' the error measure (equation 31) becomes ##EQU23## where x.sub.i .epsilon. I and x'.sub.i .epsilon.I. For a given FOE-location x.sub.f, the problem is to minimize E with respect to the two unknowns s.sub.x and s.sub.y. To reduce this problem to a one-dimensional search, one point x.sub.g, called the "guiding point", is selected in image I which is forced to maintain zero error (see FIG. 25) wherein one vector x.sub.g is selected from the set of displacement vectors to determine the optimum 2-D shift to be applied to points x'.sub.i, given a FOE-location x.sub.f. First x'.sub.g is forced onto the line x.sub.f x.sub.g and then the entire image I'={x'.sub.1, x'.sub.2, . . . } is translated in the direction of this line until the error value reaches a minimum. Therefore, the corresponding point x'.sub.g must lie on a straight line passing through x.sub.f and x.sub.g. Any shift s applied to the image I' must keep x'.sub.g on this straight line, so

for all s, and thus,

For .lambda.=1, s=x.sub.g -x'.sub.g, which is the vector x'.sub.g .fwdarw.x.sub.g. This means that the image I' is shifted such that x.sub.g and x'.sub.g overlap. This leaves .lambda. as the only free variable and the error function (equation 33) is obtained as ##EQU24## Differentiating equation 35 with respect to .lambda. and forcing the resulting equation to zero yields the parameter for the optimal shift s.sub.opt as ##EQU25##

The optimal shift s.sub.opt and the resulting minimum error .SIGMA.(.lambda..sub.opt) for the given FOE-location x.sub.f is obtained by inserting .lambda..sub.opt into equations (34b) and (35) respectively, giving ##EQU26## The normalized error E.sub.n shown in the following results (shown in FIGS. 27-32) is defined as ##EQU27## where N is the number of displacement vectors used for computing the FOE.

Since in a displacement field caused by pure camera translation all vectors must point away from the FOE, this restriction must hold for any candidate FOE-location (as illustrated in FIG. 25). If after applying s.sub.opt (x.sub.f) to the second image I', the resulting displacement field contains vectors pointing towards the hypothesized x.sub.f, then this FOE-location is prohibited and can be discarded from further consideration such is the case at point x.sub.1 in FIG. 26. FIG. 26 shows a field of 5 displacement vectors. The optimal shift s.sub.opt for the given x.sub.f is shown as a vector in the lower right-hand corner. When s.sub.opt is applied to point x'.sub.1, the resulting displacement vector (shown fat) does not point away from the FOE. Since its projection onto the line x.sub.f x.sub.1 points towards the FOE, it is certainly not consistent with a radial expansion pattern.

The final algorithm for determining the direction of heading as well as horizontal and vertical camera rotations is the "find-FOE algorithm" which consists the steps: (1) Guess an initial FOE x.sub.f.sup.0, for example the FOE-location obtained from the previous pair of frames; (2) Starting from x.sub.f.sup.0, search for a location x.sub.f.sup.opt where E.sub.min (x.sub.f.sup.opt) is a minimum. A technique of steepest descent is used, where the search proceeds in the direction of least error; and (3) Determine a region around x.sub.f.sup.opt in which the error is below some threshold. The search for this FOE-area is conducted at FOE-locations lying on a grid of fixed width. In the examples shown, the grid spacing is 10 pixels on both x- and y-directions.

The error function E(x.sub.f) is computed in time proportional to the number of displacement vectors N. The final size of the FOE-area depends on the local shape of the error function and can be constrained not to exceed a certain maximum M. Therefore, the time complexity is O(MN).

The first set of experiments was conducted on synthetic imagery to investigate the behavior of the error measure under various conditions, namely, the average length of the displacement vectors (longer displacement vectors lead to a more accurate estimate of the FOE), the amount of residual rotation components in the image, and the amount of noise applied to the location of image points. FIGS. 27a-d shows the distribution of the normalized error E.sub.n (x.sub.f) for a sparse and relatively short displacement field containing 7 vectors. Residual rotation components of .+-.2.degree. in horizontal and vertical direction are present in FIGS. 27b-d to visualize their effects upon the image. This displacement field was used with different average vector lengths (indicated as length-factor) for the other experiments on synthetic data. The displacement vector through the guiding point is marked with a heavy line. The choice of this point is not critical, but it should be located at a considerable distance from the FOE to reduce the effects of noise upon the direction of the vector x.sub.f x.sub.g. FIGS. 27a-d, which show displacement field and minimum error at selected FOE-locations, the error function is sampled in a grid with a width of 10 pixels over an area of 200 by 200 pixels around the actual FOE, which is marked by small square 178. At each grid point, the amount of normalized error is (equation 41) indicated by the size of the circle 180. Heavy circles 180 indicate error values which are above a certain threshold. Those FOE-locations that would result in displacement vectors which point towards the FOE (as described above) are marked as prohibited (+). It can be seen that the shape of the 2-D error function changes smoothly with different residual rotations over a wide area and exhibits its minimum close to the actual location of the FOE. FIG. 27a represents no residual rotation, FIG. 27b represents 2.0.degree. of horizontal camera rotation (to the left), FIG. 27c represents 2.0.degree. of vertical rotation (upwards), and FIG. 27d represents -2.0.degree. vertical rotation (downwards).

FIGS. 28 to 33 show the effects of various conditions upon the behavior of this error function in the same 200 by 200 pixel square round the actual FOE as in FIGS. 27a-d.

FIGS. 28a-e show how the shape of the error function depends upon the average length (with length factors varying from 1 to 15) of the displacement vectors in the absence of any residual rotation or noise (except digitization noise). The minimum of the error function becomes more distinct with increasing amounts of displacement. FIGS. 29a-e show the effect of increasing residual rotation in horizontal direction upon the shape of the error function for relatively short vectors (length factor of 2.0) in absence of noise.

FIGS. 30a-e show the effect of residual rotation in vertical direction upon the shape of the error function for short vectors (length factor of 2.0) in absence of noise. Here, it is important to notice that the displacement field used is extremely nonsymmetric along the Y-axis of the image plane. This is motivated by the fact that in real ALV images, long displacement vectors are most likely to be found from points on the ground, which are located in the lower portion of the image. Therefore, positive and negative vertical rotations have been applied in FIGS. 30a-e.

In FIGS. 31a-j, residual rotations in both horizontal and vertical direction, respectively, are present, for short vectors with a length factor of 2.0. In FIGS. 31a-e, the error function is quite robust against rotational components in the image. FIGS. 31f-j show the amounts of optimal linear shift s.sub.opt under the same conditions. The result in FIG. 31e shows the effect of large combined rotation of 4.0.degree./4.0.degree. in both directions. Here, the minimum of the error function is considerably off the actual location of the FOE because of the error induced by using a linear shift to approximate the nonlinear derotation mapping. In such a case, it would be necessary to actually derotate the displacement field by the amount of rotation equivalent to s.sub.opt found at the minimum of this error function and repeat the process with the derotated displacement.

The effects of various amounts of uniform noise applied to image point coordinates for a constant average vector length of 5.0, are shown in FIGS. 32a-e. For this purpose, a random amount (with uniform distribution) of displacement was added to the original (continuous) image location and then rounded to integer pixel coordinates. Random displacement was applied in ranges from .+-.0.5 to .+-.4.0 pixels in both horizontal and vertical direction. The shape of the error function becomes flat around the local minimum of the FOE with increasing levels of the noise. The displacement field contains only 7 vectors. What is observed here is that the absolute minimum error increases with the amount of noise. FIGS. 32a-e thus serve as an indicator for the amount of noise present in the image and the reliability of the final result.

The length of the displacement vectors is an important factor. The shorter the displacement vectors are, the more difficult it is to locate the FOE correctly in the presence of noise. FIGS. 33a and 33b show the error functions for two displacement fields with different average vector lengths (length factor 2.0 and 5.0, respectively). For the shorter displacement field (length-factor 2.0) in FIG. 33, the shape of the error function changes dramatically under the same amount of noise (compare FIG. 31a). A search for the minimum error (i.e., local minimum ) inevitably converge towards an area 182 indicated by the small arrow, far off the actual FOE. For the image with length-factor 5.0 (FIG. 33b), the minimum of the error function coincides with the actual location of the FOE 178. The different result for the same constellation of points in the FIG. 32d is caused by the different random numbers (noise) obtained in each experiment. This experiment confirms that a sufficient amount of displacement between consecutive frames is essential for reliably determining the FOE and thus, the direction of vehicle translation.

The performance of this FOE algorithm is shown below on a sequence of real images taken from a moving ALV. Also, it is shown how the absolute velocity of the vehicle can be estimated after the location of the FOE has been determined. The essential measure used for this calculation is the absolute height of the camera above the ground which is constant ad known. Given the absolute velocity of the vehicle, the absolute distance from the camera of 3-D points in the scene can be estimated using equation 20.

Next, velocity over ground may be computed. After the FOE has been computed following the steps outlined above, the direction of vehicle translation and the amount of rotation are known. From the derotated displacement field and the location of the FOE, the 3-D layout of the scene can be obtained up to a common scale factor (equation 20). As pointed out above, this scale factor and, consequently, the velocity of the vehicle can be determined if the 3-D position of one point in space is known. Furthermore, it is easy to show that it is sufficient to know only one coordinate value of a point in space to reconstruct its position in space from its location in the image.

Since the ALV travels on a fairly flat surface, the road can be approximated as a plane which lies parallel to the vehicle's direction of translation (see FIG. 34). This approximation holds at least for a good part of the road in the field of view of the camera 184. FIG. 34 shows a side view of camera 184 traveling parallel to flat surface 186. Camera 184 advances in direction Z, such that a 3-D point on ground surface 186 moves relative to camera 184 from Z.sub.0 to Z.sub.1. Depression angle .phi. can be determined from the location of FOE 188 in image 128. Height of camera 184 above ground surface 186 is given.

Since the absolute height of camera 184 above the ground 186 is constant and known, it is possible to estimate the positions of points 192 and 194 on road surface 186 with respect to the vehicle 132 (of FIG. 6) in absolute terms. From the changing distances between points 192 and 194 and camera 184, the actual advancement and speed can be determined.

First, a coordinate system 130 is introduced which has its origin O in the lens center 126 of the camera 184. The Z-axis of coordinate system 130 passes through FOE 188 in the image plane 128 and aims in the direction of translation. The original camera-centered coordinate system (X Y Z) 130 is transformed into the new frame (X' Y' Z') merely by applying horizontal and vertical rotation until the Z-axis lines-up with FOE 188. The horizontal and vertical orientation in terms of "pan" and "tilt" are obtained by "rotating" FOE 188 (x.sub.f y.sub.f) into the center of image 128 (OO) using equations 14 and 15 in the following: ##EQU28##

The two angles .phi..sub.f and .phi..sub.f represent the orientation of camera 184 in 3-D with respect to the new coordinate system (X' Y' Z'). This allows determination of the 3-D orientation of the projecting rays passing through image points y.sub.0 and y.sub.1 by use of the inverse perspective transformation. A3-D point X in the environment whose image x=(xy) is given, lies on a ##EQU29## For points 192 and 194 on road surface 186 of FIG. 34, the Y-coordinate is -h which is the height of camera 184 above ground 186. Therefore, the value of .kappa..sub.s for a point on the road surface (x.sub.s y.sub.s) can be estimated as ##EQU30## and its 3-D distance is found by inserting .gamma..sub.s into equation 41 as ##EQU31##

If a point on the ground is observed at two instances of time, x.sub.s at time t and x'.sub.s at t', the resulting distances from the vehicle Z.sub.s at t and Z'.sub.s at t' yield the amount of advancement .DELTA.Z.sub.s (t, t') and estimated velocity V.sub.s (t, t') in this period as

Given a set of suitable displacement vectors S={x.sub.i -x'.sub.i }, the estimate of the distance traveled by the vehicle is taken as the weighed average of the measurements .DELTA.Z.sub.i on individual vectors ##EQU33## and the final estimate for the vehicle velocity is ##EQU34## This computation was applied to a sequence of real images which is described below.

In the following, results of the FOE-algorithm and computation of the vehicle's velocity over ground are shown on a real image sequence taken from the moving ALV. The original sequence was provided on standard video tape with a frame-rate of 30 per second. Out of this original sequence, images were taken in 0.5 second intervals, i.e., at a frame rate of 2 per second in order to reduce the amount of storage and computation. The images were digitized to a spatial resolution of 512.times.512, using only the Y-component (luminance) of the original color signal.

FIGS. 35a-i show the edge image of frames 182-190 of an actual test, with points 1-57 being tracked and labeled with ascending numbers. FIGS. 35a-i show the original image sequence taken from the moving ALV after edge detection and point detection. Selected points 1-57 are located at the lower-left corners of their marks. FIGS. 35j-p (frames 191-197), which include additional points 58-78, show the original image sequence taken after edge detection and point selection. An adaptive windowing technique was developed as an extension of relaxation labeling disparity analysis for the selection and matching of tracked points. The actual image location of each point is the lower left corner of the corresponding mark. The resulting data structure consists of a list of point observations for each image (time), e.g., time t.sub.o : ((P.sub.1 t.sub.0 x.sub.1 y.sub.1) (P.sub.2 t.sub.0 x.sub.2 y.sub.2) (P.sub.3 t.sub.0 x.sub.3 y.sub.3) . . . ) time t.sub.1 : ((P.sub.1 t.sub.1 x.sub.1 y.sub.1) (P.sub.2 t.sub.1 x.sub.2 y.sub.2) (P.sub.3 t.sub.1 x.sub.3 y.sub.3) . . . )

Points are given a unique label when they are encountered for the first time. After the tracking of a point has started, its label remains unchanged until this point is no longer tracked. When no correspondence is found in the subsequent frame for a point being tracked, either because of occlusion, or the feature left the field of view, or because it could not be identified, tracking of this point is discontinued. Should the same point reappear again, it is treated as a new item and given a new label. Approximately 25 points per image have been selected in the sequence shown in FIGS. 35a-i. In the search for the focus of expansion, the optical FOE-location from the previous pair of frames is taken as the initial guess. For the very first pair of frames (when no previous result as available), the location of the FOE is guessed from the known camera setup relative to the vehicle. The points which are tracked on the two cars (24 and 33) are assumed to be known as moving and are not used as reference points to compute the FOE, vehicle rotation, and velocity. This information is eventually supplied by the reasoning processes in conjunction with the qualitative scene model (FIG. 2).

FIGS. 36a-p illustrative the displacement vectors and estimates of vehicle motion for the image sequence shown in FIGS. 35a-p. shaded area 190 marks the possible FOE locations and circle 178 inside of area 190 is the FOE with the lowest error value. The FOE-ratio measures the flatness of the error function inside area 190.

Also, FIGS. 36a-p show the results of computing the vehicle's motion for the same sequence as in the previous Figure. Each frame t displays the motion estimates for the period between t and the previous frame t-1. Therefore, no estimate is available at the first frame (182). Starting from the given initial guess, the FOE-algorithm first searches for the image location, which is not prohibited and where the error function (equation 35) has a minimum.

The optimal horizontal and vertical shift resulting at this FOE-location is used to estimate the vehicle's rotations around the X- and Y-axis. This point, which is the initial guess for the subsequent frame, is marked as a small circle inside the shaded area. The equivalent rotation components are shown graphically on a .+-.1.degree. scale. They are relatively small throughout the sequence such that it was never necessary to apply intermediate derotation and iteration of the FOE-search. Along with the original displacement vectors (solid lines), the vectors obtained after derotation are shown with dashed lines.

After the location with minimum error has been found, it is used as the seed for growing a region of potential FOE-locations. The growth of the region is limited by two restrictions: (1) The ratio of maximum to minimum error inside the region is limited, i.e.,

(see equation 40 for the definition of the error function E.sub.n). No FOE-location for which the error ratio .rho..sup.i exceeds the limit .rho..sup.lim is joined to the region. Thus the final size of the region depends on the shape of the error function. In this example, the ratio .rho..sup.lim was set at 4.0. Similarly, no prohibited locations (FIG. 26) are considered; (2) The maximum size of the region M is given by the given FOE-region region regardless of their error values. The resulting error ratio .rho..sup.max =max(.rho..sup.i) for the points inside the region indicates the shape of the error function 72 for this area. A low value for the ratio .rho..sup.max indicates a flat error function. The value for .rho..sup.max is shown as FOE-RATIO in every image.

For the computation of absolute vehicle velocity, only a few prominent displacement vectors were selected in each frame pair. The criteria were that the vectors be located below the FOE and that their length be more than 20 pixels. The endpoints of the selected (derotated) vectors are marked with dark dots. The parameter used for the computation of absolute advancement is the height of the camera above the ground, which is 3.3 meters (11 feet).

FIG. 37 illustrates hardware embodiment 200 of the present invention. Camera 202 may be a Hitachi television camera having a 48.degree. vertical field of view and a 50.degree. horizontal field of view. Camera 202 may be set at a depression angle of 16.3.degree. below the horion. Image sequence 204 acquired by camera 202 are transmitted to tracking means 206 comprising an I.sup.2 S processor 208 and a VAX 11/750 computer 210. Tracking means 206 tracks tokens between frames of image sequences 204. The output of means 206 goes to means 212 for matching of tokens and corresponding images. For the two different computers 210 and 218, a VAX-Symbolics bidirectional network protocol means 214 is connected between means 212 and means 216 which includes Symbolics 3670 computer 218, though it is possible to use one computer thereby eliminating means 214. Computer 218 provides processing for obtaining fuzzy focus of expansion and rotations/translation estimates of motion. The language environment used with computer 218 is Common LISP.

The methods and processes in the embodiment of the present invention are implemented with the ensuing programming. ##SPC1##