User:Dan.mankowitz

Range Flow is a topic in computer vision  which is discussed in CVonline

Introduction
Range Flow can be seen as an extension to the concept of Optical Flow. Optical Flow involves determining whether objects in an image are moving as well as the velocity with which an object is moving by analysing sequences of images. This is computed and displayed as a 2D velocity vector field.

Range Flow is an extension of this concept to 3 dimensions. Here, the motion of deformable surfaces are determined from analysing sequences of depth maps. The movement of objects as a result of the depth map sequence will produce a 3D velocity vector field describing the relevant object's movement. The 3D velocity is a vector of 3 components and is represented as the vector $$f = [U, V, W]^{T}$$. This vector will be referred to for the remainder of the article as the 'Range Flow'.

Range Flow Constraint Equation
When analysing the motion of objects in a sequence of image frames, we can represent the time-varying motion of the object as a depth function $$\textstyle Z(X,Y,t)$$. An assumption is made that the depth of an object in consecutive images is approximately constant. This means that a small change in depth represented by $$\Delta X, \Delta Y, \Delta t$$ should have a negligable effect on the image. That is,


 * $$ Z(X+ \Delta X, Y+\Delta Y, t+\Delta t) - Z(X,Y,t) = 0$$

In addition to this assumption, it is also assumed that the object is made up of local planar patches. This means that the depth function can be represented by a first order Taylor series yielding:


 * $$Z(X,Y,t) = Z_{0}(t) + Z_{X}X(t) + Z_{Y} Y(t)$$

Here, $$Z_{X}, Z_{Y}$$ represent partial derivatives and $$X(t), Y(t)$$ represent the local world coordinates around the point of interest. This equation is then derived with respect to time and results in the following equation:


 * $$\textstyle \frac{dZ}{dt} = \frac{d Z_0}{d t} + Z_X \frac{dX}{dt} + Z_Y \frac{dY}{dt} + X \frac{dZ_X}{dt} + Y \frac{dZ_Y}{dt}$$

The slope does not change as it is assumed that the object undergoes a pure translation. This is because the object moves an extremely small distance between consecutive depth map images. The movement can therefore be modelled as a translation. The last two terms equal zero due to the translation assumption resulting in the equation below.


 * $$Z_{x}U + Z_{y}V - W + Z_{t} = 0 $$

This is referred to as the range flow constraint equation and is analogous to the brightness change constraint equation that is used in the calculation of Optical flow. The equation is derived based on the assumption that an object observed between two consecutive frames will have a similar depth value as stated previously. The parameters $$U, V, W$$ denote the range flow and are equal to $$\frac{d}{dt}[X Y Z]^{T} $$ respectively. This equation represents a constraint plane with parameters $$U, V, W$$ defined in the velocity solution space. The normal to the plane is $$[Z_{X}, Z_{Y}, -1]^{T}$$.

Estimating the Partial Derivatives
Once the range flow constraint equation has been derived, the next logical step would be to calculate the range flow variables $$ U, V, W $$. However, there are a number of unknowns. These include the partial derivatives $$Z_{x}, Z_{y}, Z_{t}$$ as well as the range flow variables $$U, V, W$$. It is possible to estimate the partial derivatives using some simple matrix manipulation techniques.

The following calculation is based on the assumption that the derivatives $$Z_{x}, Z_{y}, Z_{t}$$ are approximately constant in a local image neighborhood. By considering an infinitesimally small movement represented by $$\Delta X, \Delta Y, \Delta t $$ a first order approximation of the function can be computed yielding the equation below:

$$Z(X+ \Delta X, Y+\Delta Y, t+\Delta t) = Z(X,Y,t) + \frac{dZ}{dX}\Delta X + \frac{dZ}{dX}\Delta Y + \frac{dZ}{dt}\Delta t$$

This equation can be extended to a linear system of equations of the form $$Az = b$$.



(A) \begin{align} \begin{bmatrix} \Delta X_1 & \Delta Y_1 & \Delta t_1 \\ \Delta X_2 & \Delta Y_2 & \Delta t_2 \\ &. &  \\ &.  &  \\ &.  &  \\ \Delta X_N & \Delta Y_N & \Delta t_N \\ \end{bmatrix}

(z)\begin{bmatrix} Z_X\\ Z_Y\\ Z_t\\ \end{bmatrix} &= (b)\begin{bmatrix} b_1 \\ .\\ .\\ .\\ b_N\\ \end{bmatrix}. \end{align} $$

The system of equations are representative of a group of pixels in a local neighborhood where the partial derivatives are constant and the column vector $$ b = Z(X+ \Delta X, Y+\Delta Y, t+\Delta t) - Z(X,Y,t)$$. Each pixel is offset by a value $$\Delta X_{i}, \Delta Y_{i}, \Delta t_{i}$$ respectively which forms the above-mentioned set of linear equations. This will enable the partial derivatives to be estimated and calculated using a pseudo inverse of the form:

$$ \mathbf{z} = (\mathbf{A}^{T}\mathbf{A})^{-1}\mathbf{A}^{T}\mathbf{b}$$

The aperture problem
The aperture problem is an effect that causes ambiguities in determining the direction of motion of objects in a scene. For example, consider looking through a small aperture at a plane with a uniform texture. It is impossible to visually determine the motion of the plane simply by looking through the aperture. This is extended to range flow and is a significant problem when determining the motion of the range flow components.



Once the partial derivatives have been calculated using the pseudo-inverse, the range flow variables $$ U, V, W $$ need to be determined. However, in calculating these variables, a number of problems exist. Firstly, this is an underconstrained problem. Only one constraint equation exists but three unknowns need to be calculated. This results in ambiguities in the velocity solution space. One constraint equation does not provide enough information to unambiguously determine the motion of the range flow components.

The range flow constraint equation defines a constraint plane in $$(U, V, W)$$ space with a vector $$ \overrightarrow{n} = [Z_X, Z_Y, -1]^{T} $$ normal to the plane. Due to the aperture problem, it is only possible to observe motion that is perpendicular to the plane. The motion of the velocity components that lie on the plane will not be visible and thus motion ambiguities arise.

The only vector that can capture motion of the plane is the normal vector $$ [Z_X, Z_Y, -1]^{T} $$. The best normal vector describing this motion is the minimal vector between the velocity origin and the constraint plane and is referred to as raw normal flow.



$$ f_{raw} = \frac{-Z_{t}}{Z_{X}^{2} + Z_{Y}^2 + 1} \begin{bmatrix} Z_{X} \\ Z_{Y} \\ 1 \\ \end{bmatrix} $$

This is one instance of the aperture problem.

In order to prevent these range flow ambiguities, more constraint planes are required to intersect the current constraint plane. These planes will slice up the velocity solution space, ultimately constraining the number of ambiguous range flow solutions. There are three different kinds of image neighborhoods that can be encountered in three dimensional range data. These include corner or point-like structures, lines and planes. The range flow associated with each of these image neighborhoods are termed Full flow, Line flow and Plane flow respectively.



This flow representation is possible on corners or points of objects. It is possible to determine all three components, $$ U, V, W $$, of range flow motion at these locations. This is due to at least three independent constraint planes intersecting at a single point producing a unique solution for range flow.
 * Full Flow

When two constraint planes intersect with one another, the range flow component along the constraint line resulting from the intersection cannot be determined. However, the range flow can be resolved in all other possible directions. A vector normal to this constraint line defines the flow of the line. This vector is called the line flow and is perpendicular to the linear structure. The point closest to the velocity origin which lies on the constraint line, defines the appropriate line flow vector.
 * Line Flow

On a planar surface, it is impossible to detect any changes in velocity within the plane. The only change in velocity that can be detected is any movement that occurs perpendicular to the plane. This vector component is termed plane flow. This vector is constructed from a point on the constraint plane that lies closest to the velocity origin.
 * Plane Flow



Local Total Least Squares (TLS)
This method is an extension of the structure tensor algorithm used in estimating optical flow. This algorithm assumes that a local image neighborhood has constant range flow. The image neighborhood consists of n pixels and each pixel has its own constraint equation. The normal to each pixel's constraint plane in velocity space is $$\overrightarrow{d} = [Z_{X}, Z_{Y}, 1, Z_{T}]^{T}$$ and the range flow velocity components are defined as $$\overrightarrow{u} = [U, V, W, 1]^{T}$$ respectively. Placing each of the pixel normals in a matrix $$D = [d_1 d_2 ... d_n]^{T}$$, it is possible to estimate the range flow as:

$$ argmin ||Du||_{2}$$ which is subject to the normalising constraint $$\overrightarrow{u}^{T}\overrightarrow{u} = 1$$

This is structured as a total least squares formula. That is, the pixel with the minimum norm to the velocity origin defines the range flow for that image neighborhood.

Calculating these normal flows can be achieved by analysing the structure tensor matrix. This matrix will determine the dominant direction in which the depth gradient is changing. The matrix is calculated using the following formula:

$$ S = D^{T}D =

\begin{bmatrix} Z_{XX}^{avg} & Z_{XY}^{avg} & Z_{X}^{avg} & Z_{XT}^{avg}\\ Z_{YX}^{avg} & Z_{YY}^{avg} & Z_{Y}^{avg} & Z_{YT}^{avg}\\ Z_{X}^{avg} & Z_{Y}^{avg} & 1^{avg} & Z_{T}^{avg}\\ Z_{TX}^{avg} & Z_{TY}^{avg} & Z_{T}^{avg} & Z_{TT}^{avg}\\ \end{bmatrix} $$

This creates a 4 by 4 matrix where $$Z_{AB}^{avg}$$ represents the average of the partial derivative product of $$Z_{A}$$ and $$Z_{B}$$ respectively. The average is computed locally using a Box or Binomial Filter. By analysing the eigenvalues and eigenvectors of this matrix, the desired range flow can be calculated. It turns out that the smallest eigenvalue of this matrix is $$\lambda_{4}$$ corresponding to the eigenvector $$\overrightarrow{e}_{4}$$. This enables the calculation of the full flow $$f_{full}$$:

$$f_{full} = \frac{1}{a_{44}} \begin{bmatrix} a_{14}\\ a_{24}\\ a_{34}\\ \end{bmatrix} $$

where $$a_{ij}$$ corresponds to element $$(i,j)$$ in the eigenvector $$\overrightarrow{e}_{j}$$.

In order to calculate the line flow, the two most significant eigenvalues and therefore eigenvectors $$e_{1}, e_{2}$$ define the minimum norm solution. This solution, $$f_{line}$$ is presented below:

$$f_{line} = \frac{1}{1 - a_{41}^{2} - a_{42}^{2}} \begin{bmatrix} a_{41} \begin{bmatrix} a_{11}\\ a_{21}\\ a_{31}\\ \end{bmatrix} + a_{42} \begin{bmatrix} a_{12}\\ a_{22}\\ a_{32}\\ \end{bmatrix} \\ \end{bmatrix} $$

Finally, calculating the plane flow is achieved from the eigenvector $$e_{1}$$ corresponding to the largest eigenvalue $$\lambda_{1}$$. This results in the plane flow norm $$f_{plane}$$.

$$f_{plane} = \frac{a_{41}}{a_{11}^2 + a_{21}^{2} + a_{31}^{2}} \begin{bmatrix} a_{11}\\ a_{21}\\ a_{31}\\ \end{bmatrix} $$

It is important to note that if the trace of the structure tensor matrix is below a certain threshold, say $$\rho_{1}$$, then the range flow for that image neighborhood will not be calculated. This is because the image neighborhood does not display significant movement and therefore the range flow does not need to be calculated saving processing time.

Global Smoothing
Another method that can be used to estimate the range flow is global smoothing. This algorithm tries to minimise distortions in the flow over the entire image, hence the use of the word smoothness. It is largely based on the optical flow smoothness algorithm called the Horn-Schunck Method. In order to perform this procedure, a smoothness term is defined and this term needs to be minimised imposing a smoothness constraint. The smoothness term is shown below:

$$argmin(\Sigma_{i=1}^{3} (\nabla f_{i})^2)$$

Here, $$(\nabla f_{i})$$ represents the partial derivative with respect to the $$i^{th}$$ velocity component. This is then added to a global energy functional which needs to be minimised. This functional is of the form shown below.

$$E = \int_A (Z_{X}U + Z_{Y}V + W + Z_{t})^2 + \alpha^2 \Sigma_{i=1}^{3} (\nabla f_{i})^2 dr$$

Minimisation of this functional can be solved iteratively using the Euler-Lagrange equations yielding a far smoother flow field. The variable $$\alpha$$ determines the importance of the smoothness term. Thus the larger the value of $$\alpha$$, the smoother the flow field becomes. However, $$\alpha$$ should not be increased without bound as important flow field detail may be lost if the image is over-smoothed. Therefore a trade-off results between producing a smooth flow field without much detail (large values of $$\alpha$$) or maintaining a noisy flow field (small values of $$\alpha$$) with a significant amount of flow field detail.

The global energy functional is determined by integrating over an area $$A$$. This area is representative of a segmented part of an image. Therefore, it is important to note that this method only makes sense if the image is segmented into different moving objects. Minimising the energy functional and determining a suitable value for $$\alpha$$ ultimately results in an estimation of the range flow.

The Addition of Intensity
There are certain instances where intensity data may be available in addition to the depth data. If this is the case, it is possible to construct another constraint equation of the form:

$$I_{X}U + I_{Y}V + I_{t} = 0$$

This has the same structural form as that of the optical flow constraint equation. This results in two constraint equations which can be utilised in the methods mentioned previously to estimate range flow. For Local TLS, the new constraint equation results in a new structure tensor matrix. Summing the two tensor matrices together yields a new structure tensor which can then be used to estimate the range flow.

This equation can also be added to the smoothness algorithm with a regularisation parameter $$\beta$$. The parameter acts much like a tuning-knob as it is used to control the influence of the intensity data and the depth data on the specified image. Low values of $$\beta$$ give more influence to the intensity data whereas higher values of $$\beta$$ creates a bias towards the range data.

The utilisation of colour information as an additional constraint to estimate the range flow in an image has also been implemented. The addition of this data creates an additional constraint which ultimately enables increased accuracy in estimating the range flow. This is because it has removed some of the motion ambiguities resulting from the aperture problem mentioned previously.

Examples
An example of range flow can be seen in a synthetic depth map in figure a below. The rendered version of the figure is shown in figure b. It should be noted that these images are artificial but they do illustrate the three different types of range flow. Figures c,e and g illustrate the X-Y range flow components for full flow, line flow and plane flow respectively. Figures d,f and h illustrate the X-Z range flow components for full flow, line flow and plane flow respectively.



A toy tiger was used as a test subject in order to determine if there is a significant increase in the accuracy of the range flow calculations, when intensity information is added to the system. The data was gathered using a laser range finder. As can be seen in the figure, intensity data provides a significant amount of additional information yielding a more detailed range flow diagram. The diagram is more dense and the directional error as well as the magnitude error of the individual range flow components decrease. Both the local TLS method and the global smoothing method experience a significant increase in detail on their respective range flow diagrams.





Skinning
One notable example of range flow is its use in analysing soft tissue deformation. By attaching speckled garments to a subjects arm, it is possible to determine the optical flows that show the deformation of the human skin. Using depth data gathered from these observations, corresponding range flows can be calculated and can then be used in analysing the deformations.

Analysing this deformation leads to one application of optical flow called skinning. This is used in 3D animation to animate the deformation of skin around articulations as characters perform various kinds of movements.



Plant Leaf Motion
Another application is determining the motion of living castor oil leaves. This can be used to determine the diurnal motion patterns of these leaves.





It should also be noted that there are many other applications for range flow such as object tracking, motion estimation and gesture recognition.