Probabilistic Analysis of Projected Features in Binocular Stereo

1. Introduction

Some geometrical relationships between projected primitives in binocular stereo systems will be analysed in the next sections with the aim of providing a characterization from a probabilistic point of view. To this end, we will consider the parallel stereo system model and the well known pinhole camera model [1].

The characterizations that will be derived will be readily usable as valuable sources of information to solve the correspondence problem in stereo systems [2] and their nature will be that of a priori information sources in Bayesian models.

To begin with, we will introduce the stereo system model that will be used for the analysis together with the notation that will be employed and the parameters that will be necessary for the calculations. Afterwards, we will use this model to derive the joint probability density function (pdf) of the orientation of the projections on the image planes of arbitrary small edges. In this case, we will find a cumbersome expression so, then, we will focus on the derivation of a tractable pdf of a convenient function of the orientation of the projections.

Later, we will turn our attention to the so called disparity gradient, which defines important relationships between projections in stereo systems. We will find three different usable pdfs of the disparity gradient that can be used to solve the correspondence problem in parallel stereo systems. Finally, a brief summary will be drawn.

2. Geometric relationships in the parallel stereo system model

In order to perform our analysis, we consider a common model for stereo image acquisition systems. The two cameras of the stereo system are considered to be identical. These cameras are modelled using the well known pinhole camera model with focal length $f$, parallel optical axes and image planes defined on the same geometric plane [3], [1]. This description defines the so called parallel stereo system model. An illustration of the geometry and the projection process with this model is represented in Fig. ▭.

For simplicity, the centre of the real world coordinate system is considered to be equidistant to the optical centres of the two cameras of the system ($C_l$ and $C_r$). The optical centers of the cameras are separated a distance $b$: the baseline. As shown in Fig. ▭, the $X$ axis is parallel to the linebase $b$ and the $Z$ axis is perpendicular to the image planes.

Figure 1.

Parallel stereo system model.

In Fig. ▭, $A$ and $B$ represent the edges of a straight segment $\overline{AB}$ of length $\delta$. $A$ is located at $(X,Y,Z)$ in the world coordinate system. The segment has an arbitrary orientation described by the angles $\alpha$ and $\beta$ defined with respect to the $XY$ and $XZ$ planes, respectively.

The edge points and the segment are projected onto the left and right image planes of our parallel stereo system. Thus, we find the projected points $A_l$ and $B_l$ on the left image and the projected segment ${\delta }_l$ on the same image. Also, the angle between ${\delta }_l$ and the horizontal on the left image is denoted ${\theta }_l$. Similarly, on the right image plane we find $A_r$, $B_r$, ${\delta }_r$ and the angle ${\theta }_r$.

Recall that the optical axes of the two cameras are parallel in our stereo model. Also, we consider that equally numbered horizontal lines on the two image planes comply with the epipolar constraint [4].

The segments on the image planes that correspond to the projection of the same segment in the real world are partially characterized and related by their respective orientations on the left and right images. This orientation can be analysed to be used to solve the correspondence problem in stereo systems.

Using the model selected, we will focus in the next sections on the orientation of the projection of small straight edges (${\delta }_l$ and ${\delta }_r$). Then, we will also consider a well known feature: the disparity gradient [5], and we will show how to develop probability characterizations of this feature under different conditions [6].

3. Joint probability density function of the orientation of projected edges

Making use of the geometrical relationships established in the previous section and in Fig. ▭, we will derive a relationship between the location and orientation of the edgel $\delta$ [7] in the real world, and the orientations of its projections described by the angles ${\theta }_l$ and ${\theta }_r$ (Fig. ▭) in the corresponding image planes. Then, under appropriate hypotheses, we will find the description of the joint probabilistic behaviour of the projected angles.

Consider the definitions and the geometry shown in Fig. ▭ where the length of the segment $\delta$ is arbitrarily small. We can write the location of the projected points in the left and right images using their coordinates on the corresponding image planes [1]. Let $B_l=(B_{lx},B_{ly})$, then, using the geometry involved and using $B$ and its projections as starting reference, we can write $A_l=(A_{lx},A_{ly})=(B_{lx}+{\delta }_{l}cos{\theta }_l,B_{ly}+{\delta }_{l}sin{\theta }_l)$.

Now, let's look at the right ($r$) image. Under the hypotheses described previously, and using the length of the projected edgel on the right image, ${\delta }_r$, making use of the fact that the $y$ coordinates must be the same in the two images, it is simple to observe that ${\delta }_{l}sin{\theta }_l={\delta }_{r}sin{\theta }_r$ and, so, ${\delta }_r=\frac{{\delta }_{l}sin{\theta }_l}{sin{\theta }_r}$.

After these observations, the coordinates of the projections of $A$ and $B$ can be written as follows:

$$\begin{array}{ccc}& & A_r=(A_{rx},A_{ry})=(B_{rx}+\frac{{\delta }_{l}sin{\theta }_l}{sin{\theta }_r}cos{\theta }_r,B_{ry}+{\delta }_{l}sin{\theta }_l)\hfill \\ & & B_r=(B_{rx},B_{ry})\hfill \end{array}$$uid2

But our objective must be to find the relation between the projections and the orientation of the edgel in the real world, such orientation is described by the angles $\alpha$ and $\beta$ in Fig. ▭. Working in this direction, the following relations can be observed:

$$\left\{\begin{array}{c}\alpha =arctan\frac{B_z-A_z}{B_x-A_x}\hfill \\ \beta =arctan\frac{B_y-A_y}{\sqrt{{(B_x-A_x)}^2+{(B_z-A_z)}^2}}=arcsin\frac{B_y-A_y}{\sqrt{{(B_x-A_x)}^2+{(B_y-A_y)}^2+{(B_z-A_z)}^2}}\hfill \end{array}\right.$$uid3

On the other hand, using the projection equations of the pinhole camera model [1], the following relations can be found:

$$\begin{array}{ccc}\hfill X& =& -\frac{x_{l}b}{x_r-x_l}-\frac{b}{2}\hfill \\ \hfill Y& =& -\frac{y_{l}b}{x_r-x_l}\hfill \\ \hfill Z& =& \frac{fb}{x_r-x_l}\hfill \end{array}$$uid4

where $(X,Y,Z)$ correspond to the coordinates of a generic point in the real world and $(x_l,y_l)$, $(x_r,y_r)$ correspond to its projections on the left and right images, respectively.

Now, using eqs. (▭) to (▭) together with eqs. (▭) and (▭), it is possible to find the expressions of the following terms involved in the calculation of the projected angles:

$$\begin{array}{ccc}\hfill & & B_x-A_x\hfill \\ & & B_y-A_y\hfill \\ & & B_z-A_z\hfill \end{array}$$uid5

Then, using these expressions in eq. (▭) and writing all the terms as functions of the real world coordinates of $A$, the coordinates of $B_r$, the camera parameters $f$ and $b$ and the orientation of the projections of the edgel (${\theta }_l$ and ${\theta }_r$), we find the equations that lead us from $(\alpha ,\beta )$ to $({\theta }_l,{\theta }_r)$:

$$\left\{\begin{array}{c}\alpha =arctan\left[\frac{Zsin({\theta }_r-{\theta }_l)}{Xsin({\theta }_r-{\theta }_l)-\frac{b}{2}sin({\theta }_r+{\theta }_l)}\right]\hfill \\ \beta =arctan\left[\frac{bsin{\theta }_{l}sin{\theta }_r-Ysin({\theta }_r-{\theta }_l)}{\sqrt{{[Xsin({\theta }_r-{\theta }_l)-bsin({\theta }_r+{\theta }_l)]}^2+{[Zsin({\theta }_r-{\theta }_l)]}^2}}\right]\hfill \end{array}\right.$$uid6

After these operations, we are ready to derive the joint pdf of the orientation of the projections of the segment: $f_{{\theta }_l,{\theta }_r}({\theta }_l,{\theta }_r)$. To this end, only the pdf of $(\alpha ,\beta )$ is required at this stage.

Since there is no reason to think differently, we will assume that these two parameters are independent uniform random variables (rv's) ranging from 0 to $\pi$ [8]. Under these hypotheses, it is evident that the joint pdf of $(\alpha ,\beta )$ is $f_{\alpha \beta }(\alpha ,\beta )=\frac{1}{{\pi }^2}$. So, in order to derive the desired expression, we only need to calculate the modulus of the Jacobian of the transformation [9]:

$$|J_d|=\left|\begin{array}{cc}\frac{\partial \alpha }{\partial {\theta }_l}& \frac{\partial \alpha }{\partial {\theta }_r}\\ \frac{\partial \beta }{\partial {\theta }_l}& \frac{\partial \beta }{\partial {\theta }_r}\end{array}\right|$$uid7

Thus, we must find the partial derivatives of $\alpha$ and $\beta$ with respect to ${\theta }_l$ y ${\theta }_r$. These are not simple expressions because of the functions involved. As an example, observe the result obtained for the last element of $J_d$:

$$\begin{array}{c}\frac{\partial \beta }{\partial {\theta }_r}=\frac{[bsin{\theta }_{l}cos{\theta }_r-Ycos({\theta }_r-{\theta }_l)]\{{[Xsin({\theta }_r-{\theta }_l)-bsin({\theta }_r+{\theta }_l)]}^2+Z^2{sin}^2({\theta }_r-{\theta }_l)\}-...}{\{{[Xsin({\theta }_r-{\theta }_l)-bsin({\theta }_r+{\theta }_l)]}^2+...}\cdots \hfill \\ \cdots \frac{...[bsin{\theta }_{l}sin{\theta }_r-Ysin({\theta }_r-{\theta }_l)]\{[Xsin({\theta }_r-{\theta }_l)-bsin({\theta }_r+{\theta }_l)]...}{...Z^2{sin}^2({\theta }_r-{\theta }_l)+{[bsin{\theta }_{l}sin{\theta }_r-Ysin({\theta }_r-{\theta }_l)]}^2\}...}\cdots \\ \hfill \cdots \frac{...[Xcos({\theta }_r-{\theta }_l)-bcos({\theta }_r+{\theta }_l)]+Z^{2}sin({\theta }_r-{\theta }_l)cos({\theta }_r-{\theta }_l)\}}{\cdots \sqrt{{[Xsin({\theta }_r-{\theta }_l)-bsin({\theta }_r+{\theta }_l)]}^2+Z^2{sin}^2({\theta }_r-{\theta }_l)}}\end{array}$$uid8

Since analytical expressions for all the required terms can be found by direct calculations, it is possible to obtain the desired pdf operating in the usual way [9]:

$$f_{{\theta }_l,{\theta }_r}({\theta }_l,{\theta }_r)=\frac{1}{{\pi }^2}\left|J_d\right|$$uid9

Unfortunately, this expression far from being simple because of the complexity of the terms involved. This fact should encourage us to search for a more usable expression capable of statistically describing a certain relation between the orientation of the projected segments. In the next section, we find such expression by using a function of cot ${\theta }_l$ and cot ${\theta }_r$.

4. Probability density function of the difference of the cot of the orientation of projected segments

A tractable expression to relate the orientation of projected segments can be found by defining a suitable function of the projected angles shown in Fig. ▭. Let $f_K\left(k\right)$, with $k$ a function of $\{{\theta }_l,{\theta }_r\}$ denote such function.

More specifically, the pdf of the modulus of the difference of the cot of the projected angles in the selected binocular stereo system will be derived.

Taking into account the scene depicted in Fig. ▭, let $\overline{AB}$ define, again, a straight segment with arbitrary length $\delta$. The orientation of this segment is described by the angles $\alpha$ y $\beta$ as shown in the figure.

Now, the location of the edges of the segment in the real world coordinate system will be written as follows:

$$A:(A_x,A_y,A_z)=(X,Y,Z)$$uid10

$$B:(B_x,B_y,B_z)=(X+\delta cos\beta cos\alpha ,Y-\delta sin\beta ,Z-\delta cos\beta sin\alpha )$$uid11

And taking into account the geometry selected, the coordinates of the projections of the edges of the segment can be written as:

$$\begin{array}{ccc}\hfill A_{rx}=-\frac{f}{A_z}(A_x-\frac{b}{2})& & A_{lx}=-\frac{f}{A_z}(A_x+\frac{b}{2})\hfill \\ \hfill A_{ry}=-\frac{f}{A_z}{A}_y& & A_{ly}=-\frac{f}{A_z}{A}_y\hfill \end{array}$$uid12

$$\begin{array}{ccc}\hfill B_{rx}=-\frac{f}{B_z}(B_x-\frac{b}{2})& & B_{lx}=-\frac{f}{B_z}(B_x+\frac{b}{2})\hfill \\ \hfill B_{ry}=-\frac{f}{B_z}{B}_y& & B_{ly}=-\frac{f}{B_z}{B}_y\hfill \end{array}$$uid13

Now, let

$$k=|cot\left({\theta }_l\right)-cot\left({\theta }_r\right)|$$uid14

Substituting the cot functions by the corresponding expressions in terms of the projections of the edges of the segment, using the projection equations (▭) to (▭), multiplying by $A_z{B}_z$, substituting $B_i$ as a function of the coordinates of $A$ and dividing by $cos\beta$, the following expression is found:

$$k=\left|\frac{-bsin\alpha }{Ztan\beta -Ysin\alpha }\right|$$uid15

This expression will be used to derive the pdf of $k$.

To begin with, the joint pdf of $k$ and $\alpha$ will be derived. To this end, the following transformation equations will be used:

$$\left\{\begin{array}{c}k=\left|\frac{-bsin\alpha }{Ztan\beta -Ysin\alpha }\right|\hfill \\ \alpha =\alpha \hfill \end{array}\right.$$uid16

The modulus of the Jacobian of the transformation can be easily determined:

$$\left|J\right|=\left|\begin{array}{cc}\frac{\partial k}{\partial \alpha }& \frac{\partial k}{\partial \beta }\\ \frac{\partial \alpha }{\partial \alpha }& \frac{\partial \alpha }{\partial \beta }\end{array}\right|=\frac{bsin\alpha Z{sec}^2\beta }{{(Ztan\beta -Ysin\alpha )}^2}$$uid17

With all this, the joint pdf of $k$ and $\alpha$ can be readily obtained [10], [9]:

$$f_{k,\alpha }(k,\alpha )=\sum _{r}f(\alpha (k_r,{\alpha }_r),\beta (k_r,{\alpha }_r))\frac{1}{|J_r|}$$uid18

where $r$ represents the set of roots of the transformation of $(\alpha ,\beta )$ as a function of $(k,\alpha )$. Two different solutions can be found for this transformation because of the modulus operation in equation (▭):

$$\left\{\begin{array}{c}\left\{\begin{array}{cc}\beta =arctan\left[sin\alpha \left(\frac{kY+b}{kZ}\right)\right],\hfill & withk=\frac{bsin\alpha }{Ztan\beta -Ysinalpha}\hfill \\ \beta =arctan\left[sin\alpha \left(\frac{kY-b}{kZ}\right)\right],\hfill & withk=\frac{-bsin\alpha }{Ztan\beta -Ysinalpha}\hfill \end{array}\right.\hfill \\ \alpha =\alpha \hfill \end{array}\right.$$uid19

Assuming, that the orientation angles $\alpha$ and $\beta$ behave as uniform random variables [8] with range $(0,\pi )$ and assuming independence, it is clear that $f(\alpha ,\beta )=\frac{1}{{\pi }^2}$ [9]. Then, equation (▭) can be written, after substitution of the terms involved as:

$$f_{k,\alpha }(k,\alpha )=\frac{1}{{\pi }^2}{\left.\frac{{(Ztan\beta -Ysin\alpha )}^2}{bsin\alpha Z{sec}^2\beta }\right|}_{\beta =arctan\left[sin\alpha \frac{kY+b}{kZ}\right]}+\frac{1}{{\pi }^2}{\left.\frac{{(Ztan\beta -Ysin\alpha )}^2}{bsin\alpha Z{sec}^2\beta }\right|}_{\beta =arctan\left[sin\alpha \frac{kY-b}{kZ}\right]}$$uid20

Now, $\alpha$ and $\beta$ can be expressed in terms of $\alpha$ and $k$, making use of the following identity: $sec[arctana]=\sqrt{1+a^2}$. Thus, the following expression is found after some simplifications:

$$f_{k,\alpha }(k,\alpha )=\frac{1}{{\pi }^2}\frac{bsin\alpha }{k^{2}Z\left[1+{sin}^2\alpha {\left(\frac{kY+b}{kZ}\right)}^2\right]}+\frac{1}{{\pi }^2}\frac{bsin\alpha }{k^{2}Z\left[1+{sin}^2\alpha {\left(\frac{kY-b}{kZ}\right)}^2\right]}$$uid21

Now, the last step to reach our objective is to integrate with respect to $\alpha$. The two terms of the previous fdp can be integrated similarly. It will be shown how the first one is handled:

$$\begin{array}{c}{I}_1={\int }_{\alpha =0}^{\pi }\frac{b}{{\pi }^2{k}^{2}Z}\frac{sin\alpha }{\left[1+{sin}^2\alpha {\left(\frac{kY+b}{kZ}\right)}^2\right]}d\alpha =\left\{\begin{array}{c}cos\alpha =x\\ -sin\alpha d\alpha =dx\end{array}\right\}\Rightarrow \hfill \\ \frac{b}{{\pi }^2{k}^{2}Z}{\int }_{x(\alpha =0)}^{x(\alpha =\pi )}\frac{-dx}{1+(1-x^2){\left(\frac{kY+b}{kZ}\right)}^2}=\\ \frac{b}{{\pi }^2{k}^{2}Z\left[1+{\left(\frac{kY+b}{kZ}\right)}^2\right]}{\int }_{x(\alpha =0)}^{x(\alpha =\pi )}\frac{-dx}{1-x^2\frac{{\left(\frac{kY+b}{kZ}\right)}^2}{1+{\left(\frac{kY+b}{kZ}\right)}^2}}=\left\{\begin{array}{c}x\frac{\frac{kY+b}{kZ}}{\sqrt{1+{\left(\frac{kY+b}{kZ}\right)}^2}}=y\\ dx=\frac{\sqrt{1+{\left(\frac{kY+b}{kZ}\right)}^2}}{\frac{kY+b}{kZ}}dy\end{array}\right\}\Rightarrow \\ \frac{b}{{\pi }^2{k}^{2}Z\sqrt{1+{\left(\frac{kY+b}{kZ}\right)}^2}\frac{kY+b}{kZ}}{\int }_{y\left(x\right(\alpha =0\left)\right)}^{y\left(x\right(\alpha =\pi \left)\right)}\frac{-dy}{1-y^2}=\\ \hfill \frac{2b}{{\pi }^{2}k\sqrt{1+{\left(\frac{kY+b}{kZ}\right)}^2}(kY+b)}arctanh\left(\frac{kY+b}{\sqrt{k^2{Z}^2+{(kY+b)}^2}}\right)\end{array}$$uid22

The second term can be integrated likewise.

Finally, the target pdf, $f_k\left(k\right)$, can be written:

$$\begin{array}{c}{f}_k\left(k\right)=\frac{2b}{{\pi }^{2}k\sqrt{1+{\left(\frac{kY+b}{kZ}\right)}^2}(kY+b)}arctanh\left(\frac{kY+b}{\sqrt{k^2{Z}^2+{(kY+b)}^2}}\right)+\hfill \\ \hfill \frac{2b}{{\pi }^{2}k\sqrt{1+{\left(\frac{kY-b}{kZ}\right)}^2}(kY-b)}arctanh\left(\frac{kY-b}{\sqrt{k^2{Z}^2+{(kY-b)}^2}}\right),k>0\end{array}$$uid23

This is the expression we were looking for. The behaviour of this function is represented in Fig. ▭.

Figure 2.

Probability density function of the modulus of the difference of the cot of the orientation of projected segments (Y=0).

5. The disparity gradient

The disparity gradient has been successfully used in the process of establishment of the correspondence relationships in stereo vision systems. Although the probabilistic behaviour of this feature has been used previously [11], [12], the process to derive some of the pdfs related to the disparity gradient has not been detailed. In this section, we will focus on the specific procedure to find different approximations of the probabilistic characterization of the disparity gradient. Thus, we will derive several expressions of the pdf of the disparity gradient$DG$ represents the random variable whereas $dg$ represents a realization of $DG$.:

$$f_{DG}\left(dg\right)$$uid26

We will pay attention to the assumptions required to derive the pdfs and to the approximations used in the different cases considered.

5.1. Comments on the disparity gradient

The disparity gradient has been successfully used as a source of information to solve the correspondence problem in stereo systems [13], [14], [5], [15], [16], [17], [12].

Generally speaking, the disparity gradient provides a priori information regarding how the real world scene is projected onto the image planes of a stereo system and, consequently, how different matching points in the projected images must be related in terms of geometrical (disparity related) relationships

The disparity refers to the difference between the coordinates of the projections of a certain point of the 3D world onto the image planes of a stereo system. Obviously, the disparity gradient refers to the rate of change of the disparity between nearby or related points [5].

Furthermore, it has been confirmed that the human visual system shows certain limitations related to the disparity gradient when matching stereo images [18]. More specifically, it was proved that 1 represents the limit of the disparity gradient for most of the subjects evaluated. On the other hand, other experiments were performed by other authors that showed that, under certain conditions, the disparity gradient can be over that threshold but with low probability. In fact, Pollard [19] derived a probability function for the disparity gradient in a stereo system with fixation point.

Additionally, the disparity gradient is able to consider other important constraints often employed for the analysis of three dimensional scenes such as figural continuity, ordering of projected features or continuity of the disparity gradient itself [17], [5].

5.2. Stereo system for the probabilistic analysis of the disparity gradient

In the following sections devoted to the probabilistic analysis of the disparity gradient in a parallel binocular stereo system, the specific geometry that will be considered is shown in Fig. ▭. According to this figure, the locations in the real world of the points $A$ and $B$, that define a straight segment with its mid-point at $(X_0,Y_0,Z_0)$ and length $2\delta$, are given by the following expressions:

$$\begin{array}{ccc}\hfill A=\left(X_0+\delta cos\beta cos\alpha ,Y_0+\delta cos\beta sin\alpha ,Z_0-\delta sin\beta \right)& & \\ \hfill B=\left(X_0-\delta cos\beta cos\alpha ,Y_0-\delta cos\beta sin\alpha ,Z_0+\delta sin\beta \right)& \end{array}$$uid29

Figure 3.

Parallel binocular stereo system for the analysis of the disparity gradient.

Then, the projections of the edge points of the segment onto the right and left image planes are given by:

$$\begin{array}{ccc}\hfill A_r=\left(-\frac{f}{A_z}(A_x-\frac{b}{2}),-\frac{f}{A_z}{A}_y\right)& & \\ \hfill B_r=\left(-\frac{f}{B_z}(B_x-\frac{b}{2}),-\frac{f}{B_z}{B}_y\right)& & \\ \hfill A_l=\left(-\frac{f}{A_z}(A_x+\frac{b}{2}),-\frac{f}{A_z}{A}_y\right)& & \\ \hfill B_l=\left(-\frac{f}{B_z}(B_x+\frac{b}{2}),-\frac{f}{B_z}{B}_y\right)& \end{array}$$uid31

In this scenario, the disparity gradient is defined as the quotient between the difference of disparity between the two points observed and their Cyclopean separation [19]:

$$dg=\frac{Differenceofdisparity}{Cyclopeanseparation}$$uid32

Taking into account that the Cyclopean projections of $A$ and $B$ are given by the following equation:

$$\frac{A_r+A_l}{2}and\frac{B_r+B_l}{2}$$uid33

and using the disparity vectors associated to the points $A$ and $B$ given by

$$(A_l-A_r)and(B_l-B_r)$$uid34

respectively. Then the disparity gradient can be written as follows:

$$dg=2\frac{\left|\right|(A_r-B_r)-(A_l-B_l)\left|\right|}{\left|\right|(A_r-B_r)+(A_l-B_l)\left|\right|}$$uid35

Now, by substitution of the expressions of $A_l$, $B_l$, $A_r$ and $B_r$, multiplying by $A_z{B}_z$, substituting by their expressions in terms of $\delta$, $\beta$ and $Z_0$, after some simplifications and reordering all the terms, the following expression is found:

$$dg=\frac{|bsin\beta |}{\left|\right|\left(-X_{0}sin\beta -Z_{0}cos\beta cos\alpha ,-Y_{0}sin\beta -Z_{0}cos\beta sin\alpha \right)\left|\right|}$$uid36

This is the main equation that will be used to derive different expressions of the disparity gradient in different scenarios.

The following sections describe the scenarios and the procedures issued to derive the different probability density functions.

5.3. Primitives centred in the world reference system

In our first scenario, we will be able to derive an exact analytical expression of the pdf of the disparity gradient This expression can be considered to be illustrative of the behaviour of $dg$. Moreover, in the next subsection, we will show how the same expression is found under different conditions and assumptions.

In this first scenario, we will assume that $X_0=0$, $Y_0=0$ and $\alpha =0$ (see Fig. ▭). Then, the expression of the disparity gradient (eq. (▭)) is readily simplified to give:

$$dg=\frac{b}{Z_0}|tan\beta |$$uid38

We will assume that the angle of orientation $\beta$ behaves as a uniform random variable in the range $(0,\pi )$.

Paying attention to the symmetry of $dg$, it is possible to pose the problem in a more convenient way. Without loss of generality, the modulus of $tan\beta$ in eq. (▭) can be removed by simply allowing the random variable $\beta$ to be defined as a uniform random variable in $(0,\frac{\pi }{2})$. The application of this and other symmetry conditions that will be considered later will allow us to avoid some expressions that involve the calculation of the modulus of certain functions and thus the analysis and some of the expressions involved will remain conveniently more simple.

According to equation (▭), it is quite simple to obtain the derivative of the disparity gradient with respect to $\beta$. Let $g\left(\beta \right)=dg$, then $g^{\text{'}}\left(\beta \right)=\frac{b}{Z_0{cos}^2\beta }$. On the other hand, it is possible to obtain $\beta$ as $g^{-1}\left(dg\right)=arctan\left(\frac{Z_0}{b}dg\right)$. Thus, finally, the pdf of $DG$ is directly obtained:

$$f_{DG}\left(dg\right)={\left.\frac{\frac{\pi }{2}}{\frac{b}{Z_0}\left|\frac{1}{{cos}^2\beta }\right|}\right|}_{\beta =g^{-1}\left(dg\right)}={\left.\frac{\frac{2Z_0}{\pi b}}{1+{tan}^2\beta }\right|}_{\beta =g^{-1}\left(dg\right)}=\frac{\frac{2Z_0}{\pi b}}{{\left[tan\left(arctan\left(\frac{Z_0}{b}dg\right)\right)\right]}^2+1}=$$

$$f_{DG}\left(dg\right)=\frac{\frac{2}{\pi }\frac{b}{Z_0}}{d{g}^2+{\left(\frac{b}{Z_0}\right)}^2},dg\in (0,\infty )$$uid39

Figure 4.

Probability density function of the disparity gradient when the primitives projected are centred in the world reference system.

Figure 5.

Distribution function of the disparity gradient when the primitives projected are centred in the world reference system.

In this expression (eq. (▭) and Fig. ▭), a unilateral Cauchy probability density function should be identified. In our scenario, this Cauchy function is tuned by the parameters 0 and $\frac{b}{Z_0}$ [20]. The distribution function can be easily found (See Fig. ▭):

$$F_{DG}\left(dg\right)=\frac{2}{\pi }arctan\left(\frac{Z_0}{b}dg\right),dg\in (0,\infty )$$uid42

5.4. Narrow field of view cameras

In this section, another step in the analysis of the behaviour of the disparity gradient will be done. We will consider a binocular stereo system with cameras of narrow field of view satisfying the epipolar constraint. This is a scenario that can be applied in numerous cases. Moreover, we can consider this scenario as a basic model for the analysis of stereo systems and suitable for practical applications.

In this scenario, the disparity gradient is given by:

$$\begin{array}{c}dg=2\frac{\left|\right|\left(\frac{fb}{A_z},0\right)-\left(\frac{fb}{B_z},0\right)\left|\right|}{\left|\right|\left(-\frac{2f}{A_z}(X_0+\delta cos\beta cos\alpha ),-\frac{2f}{A_z}(Y_0+\delta cos\beta sin\alpha )\right)...}\cdots \hfill \\ \\ \hfill \cdots \frac{}{...-\left(-\frac{2f}{B_z}(X_0-\delta cos\beta cos\alpha ),-\frac{2f}{B_z}(Y_0-\delta cos\beta sin\alpha )\right)\left|\right|}\end{array}$$uid44

After the substitution of $A_z$ and $B_z$ by their respective expressions in terms of $X_0$, $Y_0$, $Z_0$, $\alpha$, $\beta$ and $\delta$ and reordering all the terms the following expression can be found:

$$dg=\frac{\sqrt{b^2{sin}^2\beta }}{\sqrt{(X_0^2+Y_0^2){sin}^2\beta +Z_0^2{cos}^2\beta +2Z_{0}sin\beta cos\beta (X_{0}cos\alpha +Y_{0}sin\alpha )}}$$uid45

We will derive the desired pdf making use of this equation.

The fact that the cameras of the stereo system have a narrow field of view implies that the coordinates in the real world of the projected objects should satisfy the following condition: $Z_0\gg X_0,Y_0$. On the other hand, the angle $\beta$ should not be equal to $\frac{\pi }{2}$ (as a matter of fact, being $\beta$ a continuous random variable, this conditions represents and event with zero probability).

Under the hypotheses described, removing $X_0$ and $Y_0$ from the expression of the disparity gradient, because of the narrow field approximation, and assuming that $Z_0\ll Z_0^2$, the following simplified expression is found:

$$dg\approx \frac{\sqrt{b^2{sin}^2\beta }}{\sqrt{Z_0^2{cos}^2\beta }}=\frac{bsin\beta }{Z_0|cos\beta |}$$uid46

In this scenario, the symmetry of the geometry and the behaviour of the random variables $\alpha$ and $\beta$ allows us to consider the following range for the uniform random variables $\alpha$ and $\beta$: $(-\frac{\pi }{2},\frac{\pi }{2})$ and $(0,\frac{\pi }{2})$, respectively. And then, the expression of the disparity gradient can be written as:

$$dg=\frac{bsin\beta }{Z_{0}cos\beta }$$uid47

Now, in order to derive the behaviour of the disparity gradient, we will observe the region in which the random variable $DG$ is smaller than a certain value $dg$. Then, $Prob\{DG<dg\}$ is given by the probability that the random variables $\alpha$ and $\beta$ are such that $DG<dg$. Let $C_{dg}$ denote the region in the $\alpha$-$\beta$ plane that complies with this condition:

$$Prob\{DG<dg\}=Prob\{(\alpha ,\beta )\in C_{dg}\}$$uid48

This probability can be easily found by integrating the joint pdf of $\alpha$ and $\beta$ in the region $C_{dg}$:

$$F_{DG}\left(dg\right)=\int {\int }_{C_{dg}}{f}_{\alpha ,\beta }(\alpha ,\beta )d\alpha d\beta$$uid49

where, according to the selected hypotheses, the joint pdf required is given by $f_{\alpha ,\beta }(\alpha ,\beta )=\frac{2}{{\pi }^2}$.

In order to define the region $C_{dg}$, eq. (▭) must be used in order to obtain the solutions of $\beta$:

$$\beta =arctan\left(\frac{dg{Z}_0}{b}\right)$$uid50

So, the region in the $\alpha$-$\beta$ plane that defines $C_{dg}$ is given by the following relations:

$$\left\{\begin{array}{c}\alpha \in \left(-\frac{\pi }{2},\frac{\pi }{2}\right)\hfill \\ \beta \in \left(0,arctan\left(\frac{dg{Z}_0}{b}\right)\right)\hfill \end{array}\right.$$uid51

Thus, it is possible to derive the probability distribution function of the disparity gradient solving the following integral:

$$F_{DG}\left(dg\right)={\int }_{\alpha =-\frac{\pi }{2}}^{\alpha =\frac{\pi }{2}}{\int }_{\beta =0}^{\beta =arctan\left(\frac{dg{Z}_0}{b}\right)}\frac{2}{{\pi }^2}d\beta d\alpha$$uid52

which is given by:

$$F_{DG}\left(dg\right)=\frac{2}{\pi }arctan\left(\frac{Z_0}{b}dg\right)$$uid53

Then, the probability density function can be readily obtained:

$$f_{DG}\left(dg\right)=\frac{\frac{2}{\pi }\frac{b}{Z_0}}{d{g}^2+{\left(\frac{b}{Z_0}\right)}^2}$$uid54

Observe that, under different conditions and hypotheses, the same expressions for the behaviour of the disparity gradient as in the case of primitives centred in the world coordinate system (Sec. ▭ ) have been obtained. Of course, this fact comes from the assumption that $Z_0\gg X_0,Y_0$ which asymptotically leads to the more specific case in which $X_0=0$ and $Y_0=0$.

5.5. General case. Approximate expression

Under general conditions, a close analytic solution for the probability density function or the probability distribution function of the disparity gradient has not been found. So, we will face the derivation of an approximate solution.

To this end, consider the following approximate expression of the disparity gradient in our stereo system (Fig. ▭):

$$dg=\frac{b}{\sqrt{X_0^2+Y_0^2+Z_0^2{cot}^2\beta +2Z_{0}cot\beta K(X_0+Y_0)}}$$uid56

In this expression, obtained after eq. (▭), the terms $(X_{0}cos\alpha +Y_{0}sin\alpha )$ have been substituted by $K(X_0+Y_0)$. Note that $K$ should not modify the region in which the disparity gradient is properly defined: $DG\in [0,\infty )$. Using this idea, it is possible to arrive at the desired goal. Now the procedure is described.

We know that if $\beta \to 0$, then $dg\to 0$. So, we can find a condition to impose on $K$ so that $max\left\{DG\right\}\to \infty$. To this end, the minimum of the denominator in eq. (▭) can be found in the usual way, deriving the expression in the square root with respect to $\beta$ and finding the roots:

$$\frac{\partial }{\partial \beta }\left[X_0^2+Y_0^2+Z_0^2{cot}^2\beta +2Z_{0}cot\beta K(X_0+Y_0)\right]=0$$uid57

$$-2Z_0^{2}cot\beta {csc}^2\beta -2Z_0{csc}^2\beta K(X_0+Y_0)=0$$uid58

Now, since $csc\beta \ne 0$ $\forall \beta$, the following must be fulfilled:

$$Z_{0}cot\beta +K(X_0+Y_0)=0$$uid59

Thus, the following relation is found:

$$cot\beta =-\frac{K(X_0+Y_0)}{Z_0}$$uid60

Recall that in the minimum the denominator in eq. (▭) must be zero. Substituting $cot\beta$ according to the previous expression in the denominator of eq. (▭), the following must be fulfilled:

$$X_0^2+Y_0^2+Z_0^2{\left[-\frac{K(X_0+Y_0)}{Z_0}\right]}^2+2Z_0\left[-\frac{K(X_0+Y_0)}{Z_0}\right]K(X_0+Y_0)=0$$uid61

which leads to the following expression:

$$K=\sqrt{\frac{X_0^2+Y_0^2}{{(X_0+Y_0)}^2}}$$uid62

Thus, the approximation of the disparity gradient that will be used is given by:

$$dg\approx \frac{b}{\sqrt{X_0^2+Y_0^2+Z_0^2{cot}^2\beta +2Z_0\sqrt{X_0^2+Y_0^2}cot\beta }}$$uid63

Now, the probability distribution function will be found. Consider $C_{dg}$ as the region in which $DG<dg$ and let $C_{dg}(\alpha ,\beta )$ denote the region in the $\alpha$-$\beta$ plane such that $DG<dg$. Then, again:

$$F_{DG}\left(dg\right)=\int {\int }_{C_{dg}(\alpha ,\beta )}{f}_{\alpha ,\beta }(\alpha ,\beta )d\alpha d\beta$$uid64

(a)

(b)

(c)

(d)

Figure 6.

Probability distribution function of the disparity gradient {1}. General case: simulation results (solid line) and analytic approximation (dashed line).

Since $DG$ does not depend on $\alpha$ (eq. (▭)), the region $C_{dg}(\alpha ,\beta )$ can be defined as a function of $\beta$, exclusively:

$$F_{DG}\left(dg\right)={\int }_{C_{dg}\left(\beta \right)}{\int }_{\alpha }{f}_{\alpha ,\beta }(\alpha ,\beta )d\alpha d\beta ={\int }_{C_{dg}\left(\beta \right)}\frac{1}{\pi }d\beta$$uid66

In order to define $C_{dg}\left(\beta \right)$, $dg$ must also be written as a function of $\beta$; the following result if easily obtained:

$$cot\beta =-\frac{\sqrt{X_0^2+Y_0^2}}{Z_0}\pm \frac{b}{dg{Z}_0}$$uid67

Let ${\beta }_1$ and ${\beta }_2$ represent the two solutions of this equation, then the region $C_{dg}\left(\beta \right)$ is defined by the following intervals:

$$C_{dg}\left(\beta \right)=\left\{\begin{array}{c}(-\frac{\pi }{2},min({\beta }_1,{\beta }_2))\hfill \\ \bigcup \hfill \\ (max({\beta }_1,{\beta }_2),\frac{\pi }{2})\hfill \end{array}\right.$$uid68

With all this, the desired solution, the probability distribution function of the disparity gradient, is given by (Figs. ▭ and ▭):

$$F_{DG}\left(dg\right)=1-\frac{1}{\pi }\left[arccot\left(-\frac{\sqrt{X_0^2+Y_0^2}}{Z_0}-\frac{b}{dg{Z}_0}\right)-arccot\left(-\frac{\sqrt{X_0^2+Y_0^2}}{Z_0}+\frac{b}{dg{Z}_0}\right)\right]$$uid69

(e)

(f)

(g)

(h)

Figure 7.

Probability distribution function of the disparity gradient {2}. General case: simulation results (solid line) and analytic approximation (dashed line).

Note that this solution is mathematically correct, however some considerations must be taken into account so that $F_{DG}\left(dg\right)$ behaves as a proper probability distribution function sec. 2.2[10]. Specifically, the function $arccot$ returns an angular value which, ultimately, can be seen as a periodic function with period $\pi$. This means that there is an infinite number of solutions of $arccot$, although the main solution is often considered to be in the interval $(-\frac{\pi }{2},\frac{\pi }{2})$. In our specific development, the function derived behaves properly if the solutions of the function $arccot$ are selected in the range $(-\pi ,0)$.

After the probability distribution function (eq. (▭)), the probability density function (pdf) of the disparity gradient is readily found [9]:

$$\begin{array}{c}{f}_{DG}\left(dg\right)=\frac{1}{\pi }\frac{2bZ\left[d{g}^2{Z}_0^2+b^2+d{g}^2(X_0^2+Y_0^2)\right]}{d{g}^4{Z}_0^4+b^4+d{g}^4{(X_0^2+Y_0^2)}^2+2d{g}^2{Z}_0^2{b}^2...}\cdots \hfill \\ \hfill \cdots \frac{}{...+2d{g}^4{Z}_0^2(X_0^2+Y_0^2)-2b^{2}d{g}^2(X_0^2+Y_0^2)}\end{array}$$uid71

which is a usable expression of the pdf of the disparity gradient that completes the analysis of the probabilistic behaviour of this parameter under the conditions and hypotheses selected.

6. Concluding summary

In this chapter, we have dealt with the probabilistic behaviour of certain relations established between the projection of features onto the image planes of a parallel stereo system. Specifically, we have considered relations between the orientation of projected edgels and the disparity gradient.

The projected edgels are simple features that can be considered in a matching stage [7]. The relation between their orientations constitutes an a priori source of information that, using the models proposed, can be used in the matching processes [21] of stereo systems. The formulae of the relation between the orientation of the projections derived are perfectly suited for application in Bayesian models for stereo matching [22].

The disparity gradient is an important parameter for stereo matching systems [14]. In this chapter, it has been analysed under different conditions to find proper probability density functions usable in a probabilistic context.

The functions derived can be used alone to match random dot stereo pairs [23], [24], [17], [25]. Also, these functions can contribute and collaborate with other matching models in the solution of the correspondence problem in stereo systems. Specifically, Bayesian approaches can be employed to solve the correspondence problem [26] using the proposed models of the disparity gradient [12].

Acknowledgements

This work was supported by the Ministerio de Economía y Competitividad of the Spanish Government under Project No. TIN2010-21089-C03-02 and Project No. IPT-2011-0885-430000.