IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

1

Orientation in Manhattan: Equiprojective Classes
and Sequential Estimation
André T. Martins,

Pedro M. Q. Aguiar, Member, IEEE, and Mário A. T. Figueiredo, Senior Member, IEEE

Abstract— The problem of inferring 3D orientation of a camera from video sequences has been mostly addressed by firstcomputing correspondences of image features. This intermediatestep is now seen as the main bottleneck of those approaches. In
this paper, we propose a new 3D orientation estimation method
for urban (indoor and outdoor) environments, which avoids
correspondences between frames. The scene property exploitedby our method is that many edges are oriented along threeorthogonal directions; this is the recently introduced Manhattan
world (MW) assumption.
The main contributions of this paper are: the definition ofequivalence classes of equiprojective orientations; the introduction of a new small rotation model, formalizing the fact that thecamera moves smoothly; and the decoupling of elevation andtwist angle estimation from that of the compass angle. We build
a probabilistic sequential orientation estimation method, basedon an MW likelihood model, with the above listed contributions
allowing a drastic reduction of the search space for eachorientation estimate. We demonstrate the performance of our
method using real video sequences.
Index Terms— Camera orientation, sequential estimation,
Manhattan world assumption, camera calibration.

I. I NTRODUCTION
Applications in areas such as digital video, virtual reality,
mobile robotics, and visual aids for blind people, require
efficient methods to estimate the 3D pose of a video camera
from the images it captures.
The most popular approaches to 3D pose estimation arefeature-based. In the multi-view case, this requires finding
correspondences between features [2], [3], [4]. In the singleimage case, typical methods involve feature grouping [5], [6],[7]. Naturally, in both cases, feature detection (e.g., corners,
edges) is an indispensable first step. However, it is widely
accepted that automatic feature matching or grouping areserious bottlenecks. Moreover, by basing all inference ona usually small feature set (relative to the whole image),potentially useful information may be prematurely discarded.
In the multi-view case, methods that estimate the 3D structure directly from the image intensity values, i.e., without
involving feature detection and matching, have been proposed[8], [9]. These approaches lead to complex time-consuming
algorithms and strongly rely on the assumption that the
The authors are with the Department of Electrical and Computer Engineering, Instituto Superior Técnico, Technical University of Lisbon, 1049-001
Lisboa, Portugal. P. Aguiar is also with the Institute for Systems and Robotics.
M. Figueiredo is also with the Institute for Telecommunications. Email addresses: aftm@mega.ist.utl.pt, aguiar@isr.ist.utl.pt, andmtf@lx.it.pt.
An early version of this work has appeared in [1].

brightness pattern remains (approximately) constant from view
to view.
Recently, a very different approach has been proposed
which avoids dealing with features in the single-image case,by using prior knowledge about the structure of the scene.
Specifically, in typical indoor and outdoor urban scenes, manyedges are aligned with one of the three directions defining an
orthogonal coordinate system. Under this so-called Manhattan
world (MW) assumption, Coughlan and Yuille [10], [11] used
Bayesian inference to estimate the rotational component of
the 3D pose (i.e., 3D orientation) of the camera, with respect
to this coordinate system, from a single image. The MW
assumption was also used by [12] for camera calibration andextended by [13] to more general urban environments.
In this paper, we propose a new method for 3D orientation
estimation from image sequences in MW environments. The
novelties in our method are the following:• while in [10], [11], the MW prior is used to perform
3D orientation estimation from a single image, we extend
its use for sequences of images;• we introduce a new small rotation (SR) model that
expresses the fact that the video camera undergoes asmooth 3D motion;• by defining the 3D orientation in terms of the equivalenceclasses of equiprojective orientations, we reduce the spacein which the solution has to be searched;• we show how the estimate of the elevation and twist
angles can be computed independently of the compass
angle, thus reducing the computational load.
The paper is organized as follows. In Section II, we review the geometry of camera orientation. The concept ofequiprojective orientations and the small rotation (SR) model
are introduced in Section III and IV, respectively. Section
V describes the sequential estimation method. Experimental
results are shown in Section VI and Section VII concludes
the paper.
II. C AMERA O RIENTATION AND VANISHING P OINTSLet (x, y, z) and (n, h, v) be the Cartesian coordinate
systems of the MW and the camera, respectively. These are
related through the equation (n, h, v)T = O · (x, y, z)T ,where O ∈ SO(3) is the orientation matrix, i.e., the camera
orientation. In the following text, we often denote orientation
as O ≡ O(α, β, γ), expressing the fact that it is parameterizedwith three angles: α, the compass (azimuth) angle, corresponding to rotation about the z axis; β, the elevation angle above
the xy plane; and γ, the twist about the principal axis (seeFig. 1).

IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

z

γ

v0 v
y

0

x

α

β

v0
P

h0

h0
h
image planenFig. 1. Parameterization of the camera orientation. Left: compass angle α
and elevation angle β. Right: twist angle γ represented on the image plane.
(Note: the image plane is placed in front of the optical center.)

The principal point P lies on the sphere with center at the
optical center 0 (chosen as the origin of the MW reference
frame) and radius equal to the focal length f . Its 3D coordinates are related with the compass and elevation angles via
P = (Px , Py , Pz )T = f (cos α cos β, sin α cos β, sin β)T .
(1)
The orientation O(α, β, γ) can be determined by finding wherethe vanishing points (VPs) of the MW axes project on theimage plane [2], [3]. In fact, let (h, v) be the reference frameof this plane and let the 2D principal point p be its origin,
i.e., p = (0, 0)T . Assuming a pinhole and radial-distortion-freecamera, the 2D coordinates, vx , vy , vz , of the VP projections
are related with O(α, β, γ) via³´Tαvx = f Rγ − tan
,
−
tan
β
,
cos β
³
´T
α
vy = f Rγ cot,(2)cos β , − tan β
vz

=

T

f Rγ (0, cot β) ,where Rγ is the twist matrix,
·
cos γ
Rγ =
− sin γ

sin γ
cos γ

¸
.

(3)

In the above, Cartesian coordinates are used only for simplicity; vanishing points at infinity can be handled by using
homogeneous coordinates.
III. E QUIPROJECTIVE O RIENTATIONSConsider the problem of determining the camera orientationfrom the set of three VPs on a single image. Since it is notknown which VP corresponds to which MW axis, the problem
has multiple solutions. This ambiguity motivates the conceptof equiprojectivity.
byDefinition 1 (equiprojective orientations): Denote
V(O) = {vx , vy , vz } the set of VPs determined by anorientation O. Two orientations O and O∗ are termed
equiprojective iff they have identical sets of VPs, i.e., iff
V(O) = V(O∗ ).

2

Equiprojectivity, as just defined, is reflexive, symmetric, andtransitive; therefore, it is an equivalence relation. The following result provides a way to find the complete equivalenceclass of a given orientation, i.e., the set of all orientations
which are equiprojective with it.
Proposition 2: Let O be an orientation and P =(Px , Py , Pz )T the corresponding principal point. The equivalence class of O always has 24 elements. Each O(n) =O(αn , βn , γn ), for n = 1, ..., 24, corresponds to a principal
point P(n) related to P through P(n) = Mn P, where Mn is a
3×3 signed permutation matrix (i.e., entries in {−1, 0, 1}, with
one nonzero entry per row and per column) with det Mn = 1.
The angles αn and βn are obtainable from P(n) according
to (1); the twist angles γn depend on O(α, β, γ) and P(n) as
follows:(n)⇐ MTn z = (0, 0, 1)T (Pz = Pz )
γ

(n)

T
T

γ±π
⇐ Mn z = (0, 0, −1) (Pz = −Pz ) γ + atan tan α ± π ⇐ MTn z = (1, 0, 0)T (Pz(n) = Px )
sin β
γn =
(n)
α
γ + atan tan⇐ MTn z = (−1, 0, 0)T (Pz = −Px )


sin β


(n)
α

γ − atan cot± π ⇐ MTn z = (0, 1, 0)T (Pz = Py )

sin β



(n)
T
T
cot α
γ − atan sin β⇐ Mn z = (0, −1, 0) (Pz = −Py ).

Proof: Given an orientation O, the corresponding image
plane can be seen as the plane that is tangent to the sphere{w : ||w|| = f } in P. The intersection of each MW axis
x, y, and z with the image plane defines its respective VP.
Hence, a necessary condition for an orientation O(n) to beequiprojective with O is that their corresponding principal
points (respectively P(n) and P) have the same coordinates
up to permutations and/or sign changes, which is equivalent
to the existence of a signed permutation matrix Mn satisfying
P(n) = Mn P. Any permutation matrix satisfies det Mn =±1; however, not all matrices of this kind yield a solution.
Particularly, if P and P(n) differ by a single permutation or bya single sign change, the triangles formed by the VPs at eachcase have opposite orientations, i.e., they are “reflected”. Since
the composition of two reflections is the identity, the numberof permutations plus the number of sign changes defined byany matrix Mn must be even; this is equivalent to imposing
det Mn = 1. Because the number of possible permutations
in a 3-vector is 3! = 6, and the number of sign changes is
23 = 8, we can combine permutations and sign changes in
48 different ways; since half of these correspond to “mirror
images”, the cardinality of the set {Mn } is 24 (see illustrationin Fig. 2).
For each Mn , we are able to know which VP in V(O)
corresponds to which VP in V(O(n) ). Namely, for every i, j ∈(n){x, y, z}, the VP vi and vj correspond iff jT Mn i = ±1, i.e.,(n)
iff Pi = ±Pj . Taking j = z, we have:((n)][vz pvi ]⇐ zT Mn i = 1(Pz = Pi )γn −γ =
(n)][vz pvi ] ± π ⇐ zT Mn i = −1 (Pz = −Pi ).
Finally, from (2)-(3), we obtain (2).
The concept of equiprojectivity is useful in any problem oforientation estimation, or VP location, since it allows reducing
the search spaces. This was also pointed out in [12], where analgorithm was proposed to round a quaternion to a canonicalIEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

z

=f

z = 0

z =  f
x = f

x = 0

x =  f

x = 0

x = f

y = 0

y = f

y = 0

y =  fy = 0

Fig. 2. 3D locations of the principal points of equiprojective orientations, on
the octants of a sphere with radius f . Here, we have two equivalence classes:
the white and the black points. Black points correspond to “mirror images”
of white points.

3

consistent with the SR(ξ) model iff there exists a small fixed
angle ξ such that |ρk | ≤ ξ for any k.
In our experiments, we have used a SR(5◦ ) model, which
implies that for a sampling rate of 12.5 Hz the rotation angle
is always less than 62.5◦ in each second; this is an intuitively
reasonable assumption.
The following proposition expresses how the variations of
the compass, elevation and twist angles between consecutive
frames are bounded due to the SR model.
Proposition 5: If the camera motion is consistent with theSR(ξ) model, then, at any frame k, the following bounds hold:• The elevation variation, ∆β = βk − βk−1 , satisfies|∆β| ≤ ξ.

value in SO(3)/C, where C is the octohedral group of cube
symmetries. We formalize this search space reduction in thefollowing proposition (proved in the Appendix).
Proposition 3: Every orientation O has an equiprojectiveO∗ = O(α∗ , β ∗ , γ ∗ ) such that:i π πii π πiα∗ ∈ − ,, β∗ ∈ − ,, and γ ∗ ∈ ]−ϕ, ϕ] , (4)
4 4
4 4
√
where ϕ = atan 2 ≈ 54.7◦ . An equivalent statement is: forany camera orientation O, there exists at least one VP inside
the region of the image plane shown in Fig. 3.

•(5)The compass variation, ∆α = αk − αk−1 , satisfies|∆α| ≤ aξ (βk , βk−1 ) ≡(³´cos |∆β|−cos ξ
acos 1− cos⇐ |βk−1 +βk | ≤ π−ξβk−1 cos βk(6)π⇐ otherwise.
2
If Ok−1 is in the region defined by (4), then, independently of βk and βk−1 :|∆α| ≤ acos(2 cos ξ − 1).

•(7)

The twist variation, ∆γ = γk − γk−1 , satisfies|∆γ| ≤ gξ (βk−1 ),(8)where gξ is an even function that increases in the subdomain [0, π2 ] from gξ (0) = ξ to gξ ( π2 ) = π. If Ok−1 is in
the region defined by (4), then |βk−1 | ≤ π4 and³π´|∆γ| ≤ gξ,(9)4
Fig. 4 plots gξ in the subdomain [0, π4 ], for ξ = 5◦ ; thisvalue of ξ leads to |∆γ| ≤ 7.08◦ .

ϕ ∼ 54.7
0ffilm

Fig. 3. Representation of the image plane. It is guaranteed that there exists
at least one vanishing point in the shaded region.

7.5

Maximum gamma variation for ξ = 5°
7
6.5

g (β

) = max |∆γ | [°]ξ k–1

IV. S MALL ROTATIONS M ODEL6

Let us now assume that the camera is moving and acquiring a sequence of frames {I1 , . . . , IN }. We denote byOk (αk , βk , γk ) the orientation at the k-th frame. The sequenceof orientations {O1 , ..., ON } depends only on the rotational
component of the motion. In typical video sequences, thecamera orientation evolves in a smooth continuous way. We
formalize this property by introducing the small rotations (SR)model, next described.

5.5

Definition 4: Let Rk (ρk , ek ) be the rotational componentof the camera motion between the (k−1)-th and k-th frames,
where ρk and ek denote the angle and the axis of rotation,
respectively. Independently of ek , we say that the camera is

5

β

k–1

4.5
0

5

10

15

20

25

30

35

40

[°]
45

Fig. 4. Maximum variation for the twist angle as function of the initial
elevation angle, using a SR(5◦ ) model.

Proof: Rk (ρk , ek ) is the composition of two rotations:
Rk1 (ρk1 , ek1 ) transforming the principal point Pk−1 in Pk ,followed by Rk2 (ρk2 , ek2 ) that twists the camera through theprincipal axis. Composing these two rotations, and taking intoρρaccount that ek1 ⊥ek2 , we obtain cos ρ2k = cos k21 cos k22 .

IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

Therefore, the SR(ξ) condition |ρk | ≤ ξ implies bothρρ
cos k21 ≥ cos 2ξ and cos k22 ≥ cos 2ξ , i.e., |ρk1 | ≤ ξ and|ρk2 | ≤ ξ. Since cos ρk1 = f −2 PTk Pk−1 , from (1) we obtaincos ρk1 = cos βk cos βk−1 cos ∆α+sin βk sin βk−1 ≤ cos ∆β.
This suffices to prove (5). Now, rewriting the latest inequality
for cos ∆α and simplifying leads to (6). If Ok−1 is in theregion defined by (4), then |βk + βk+1 | ≤ π/4 + π/4 + ξ ≤π −ξ. The maximum value of ∆α occurs for βk = βk−1 = π4 ,which leads to (7).
For ∆γ, we couldn’t find an simple closed-form expression
for gξ (βk−1 ). Instead, since ρk is a function of βk−1 , βk , ∆αand ∆γ, we can study gξ assuming that αk−1 = γk−1 = 0.
Spherical symmetry implies that gξ is an even function; also,a simple geometric argument shows that gξ (βk−1 ) increases
with |βk−1 |. Writing Rk as a composition of the three individual compass, elevation and twist rotations, and using theformula for the product of quaternions, yields√
AB − C B 2 + C 2 − A2|∆γ| = 2 acos,(10)
B2 + C 2∆βwhere ³A = cos ρ2k , B = cos ∆αand C =
2´ cos 2∆β∆αsin 2 cos 2 sin βk − cos βk sin ∆β . Numerical maximization of (10) w.r.t. ∆α and βk (for ρk = ξ) approximatesgξ .

If the orientation Ok−1 lies in the minimal region definedby (4), the search space for Ok is significantly reduced by thebounds imposed by Proposition 5. In particular, with ξ = 5◦ ,we have |∆α| ≤ 7.08◦ , |∆β| ≤ 5◦ and |∆γ| ≤ 7.08◦ . If Ok−1
does not lie in this minimal region, there is an equiprojectiveorientation that does. This shows how the SR model and theequiprojective orientations can be used together to reduce thesearch space.
V. S EQUENTIAL O RIENTATION E STIMATION
A. Estimation Criterion
To estimate the sequence of camera orientations{O1 , . . . , ON } from the observed image sequence{I1 , . . . , IN }, we adopt a probabilistic sequential estimation
framework, making use of the MW and SR assumptions.
The MW assumption states that the images contain manyedges consistent with the x, y and z axes; hence, the statisticsof the image intensity gradient ∇Ik of each image carryinformation about the corresponding camera orientation Ok
via a likelihood function P (∇Ik |Ok ) [10], [11]. In this paper,
we embed this idea in a sequential estimation framework,
using a maximum a posteriori (MAP) criterion:nob k−1 ) , (11)b k = arg max log P (∇Ik |Ok ) + log P (Ok |OOOkb k−1 ) penalizes large changes between
where the prior P (Ok |O
consecutive orientation estimates.
A fully Bayesian sequential estimation approach would
require computationally expensive Monte Carlo methods [14],[15]. Our results show that the simplified criterion in (11) leads
to good results and, by exploiting the equiprojectivity results
and the SR assumption introduced in the previous section, can4

be implemented in near real time. An alternative scheme wasproposed in [13], in which Ok is estimated via an iterativeb k−1 .
(EM) algorithm initialized with O

B. Likelihood Function
In this subsection, to simplify the notation, we will omit thetime index k, and derive the likelihood function P (∇I|O) fora generic image. Let Eu = (Eu , φu ) denote the element of
the image gradient ∇I at pixel u, where Eu is the gradientmagnitude and φu the gradient direction. As in [10], [11], thelikelihood function is derived as follows:
Each pixel u has a class label mu ∈ {1, 2, 3, 4, 5}. Pixels
in classes 1, 2, 3 belong to edges consistent with the x,y, z axes, respectively. Pixels in class 4 are on edges notconsistent with those axes. Non-edge pixels are in class
5. These classes have prior probabilities {P (mu )} (weadopt the values used in [10], [11]).
• The gradient magnitude and direction are conditionallyindependent, given the class label. Naturally, the gradient magnitude is also conditionally independent of thecamera orientation and of the pixel location. Thus,•

P (Eu |mu , O, u) = P (Eu |mu ) P (φu |mu , O, u), (12)where½
P (Eu |mu ) =Pon (Eu ) , if mu 6= 5
Poff (Eu ) , if mu = 5,

(13)

and Pon (Eu ) and Poff (Eu ) are the probability mass functions of the quantized gradient magnitude, conditioned onwhether pixel u is on or off an edge, respectively. These
probabilities are learned off-line.
• Let θx (O, u), θy (O, u), θz (O, u) be the gradient directions that would be ideally observed at location u
if mu = 1, 2, 3, respectively. The gradient direction
probability function is
Pang(φu −θx (O, u)) ⇐ mu = 1Pang(φu −θy (O, u)) ⇐ mu = 2
P (φu |mu , O, u) =Pang(φu −θz (O, u)) ⇐ mu = 3U(φu )⇐ mu = 4, 5,
(14)
where
½ 1−²
⇐ t ∈ [−τ, τ ]Pang (t) = 2τ²⇐t ∈ ] − π/2, −τ [ ∪ ]τ, π/2],
π−2τ
and U (·) is the uniform pdf on ]− π2 , π2 ]. In our experiments, we use ² = 0.1, and τ = 4◦ .
• Finally, the joint likelihood is obtained by marginalizing(summing) over all possible models at each pixel, andassuming independence among different pixels:
P (∇I|O) = P ({Eu } |O) =5
Y X
P(Eu |mu ) P(φu |mu , O, u) P (mu ) . (15)u mu =1

IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

5

C. Locating the Estimates
The maximization in (11), with the likelihood function(15), is a 3-dimensional optimization problem with respect
to α, β, and γ. We propose an approximate solution which
decouples the problem into two simpler steps: a 2-dimensional
optimization w.r.t. β and γ, followed by a one-dimensional
search w.r.t. α. This approximation is supported on the factthat the vanishing point vz does not depend on the compass
angle α, as is clear from (2).
In the first step, we estimate β and γ, for frame k, accordingto³´
nβbk , γbk
= arg max log P ({Eu }k |β, γ) +β,γolog P (β, γ|βbk−1 , γbk−1 ) , (16)P (φu |mu , β, γ, u) =½
Pang (φu − θz (β, γ, u)) ⇐ mu = 3(17)
U (φu )⇐ mu = 1, 2, 4, 5.
Notice that the use of a uniform distribution is simply away of ignoring angle information from all pixels but thosecorresponding to the z axis (mu = 3), when estimating βk andγk ; it doesn’t mean that those angles are actually uniformly
distributed.
P (β, γ|βbk−1 , γbk−1 ) is a truncated bivariate Gaussian withmean [βbk−1 , γbk−1 ]T , defined over the region β ∈ ]βbk−1 −ξ, βbk−1 +ξ] and γ ∈ ]b
γk−1−gξ (βbk−1 ), γ
bk−1 +gξ (βbk−1 )]. Thisprior formalizes the SR assumption (see (5) and (8)) as wellas angle variation smoothness. The variance of this Gaussian
controls the tradeoff between the smoothness of the estimatedsequence of angles and the accuracy of this estimates. In thefirst frame, the prior is flat over the entire domain (β, γ) ∈]−45◦ , 45◦ ] × ]−54.7◦ , 54.7◦ ], according to (4).
Given βbk and γ
bk , we then estimate the compass angle αk
using
n
α
bk = arg max log P ({Eu } |α, βbk , γ
bk )+
α
o
log P (α|bαk−1 , βbk−1 , βbk ) , (18)where the prior P (α|bαk−1 , βbk−1 , βbk ) is a truncated Gaussian with mean αbk−1 , defined over the interval ]b
αk−1 −aξ (βbk , βbk−1 ), αbk−1 + aξ (βbk , βbk−1 )] (see (6)). For the first
frame, the prior is flat over ]−45◦ , 45◦ ]. The maximizationsin (16) and (18) are carried out by exhaustive search.
b k (b
If a given estimate O
αk , βbk , γbk ) is located outside of theminimal region defined in (4), we replace it by an equiprojective orientation inside that region. As explained in the lastparagraph of subsection IV, this allows aξ (βbk , βbk−1 ) to be lessthan 7.1◦ , hence keeping a small search space. As a final step,at each frame k, we select an orientation from the equivalenceb k , such that the resulting sequence satisfies the SR
class of O
model.

20

Degrees

where the likelihood P ({Eu }k |β, γ) is a version of (15)which only models direction information of those edges consistent with the z axis. More specifically, instead of (14), we
use here

α
β
γ

40

0

–20
0

20

40
60
Frame number

80

Fig. 6. Top: frames 20, 30, 40, and 50 of another video sequence. Bottom:
camera angle estimates.

VI. E XPERIMENTS
The algorithm was tested with outdoor MPEG-4 video
sequences, acquired with a hand-held camera. Although thesequences are of low quality due to radial distortion and several over- and under-exposed frames, our algorithm was able
to successfully estimate the camera orientation, as illustrated
in Fig. 5.
The images in Figs. 6 and 7 show frames from two othersequences. Notice that the algorithm is able to estimate thecorrect orientation, despite the many edges not aligned with
the MW axes (e.g., people in Fig. 7). The plots in the same
figures represent the estimates of the orientation angles, forthese two sequences. Note that the estimates on the plot ofFig. 7 are slightly noisier than those in Fig. 6, due to the lower
image quality. The smoothness of these estimates is controlled
by the prior variances referred in Subsection V-C; here, thesevariances are the same for both sequences and the three angles.
Of course, there is a tradeoff between smoothness and ability
to accurately follow fast camera rotations.
Typical processing time for each (288×360)-pixels frame
is below one second, on a 3.0 GHz Pentium IV, using a MATLAB implementation. The only effort made to speed up thecomputation was the exclusion of non-relevant pixels by nonmaxima suppression followed by thresholding of the gradientmagnitude. We are currently working on a C implementationto achieve frame-rate.
VII. C ONCLUSION
We have proposed a probabilistic approach to estimatingcamera orientation from video sequences of urban scenes.

IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)

Fig. 5.

6

Orientation estimates (superimposed cubes represent the estimated MW axes) for the first and several other frames of a video sequence.

which gives us both the Euclidean distance di = (viT vi )1/2
between points vi and p = (0, 0)T , and the angle θij =vT v
acos dii djj formed by the two lines [pvi ] and [p vj ], with i 6= j.
Consider now the disk D with radius f centered at p, i.e.,
D = {(u, v) ∈ R2 : u2 + v 2 ≤ f 2 }. We have vi ∈ D iff
di ≤ f , which, by (19), is equivalent to Pi2 ≥ f 2 /2. Since
Px2 + Py2 + Pz2 = f 2 , the condition Pi2 ≥ f 2 /2 implies that
Pj2 ≤ f 2 /2 for any j 6= i, which means that there cannot exist
more than one VP in the interior of disk D. Furthermore, thethree VPs are all in the exterior or at the boundary of D iff
Pi2 ≤ f 2 /2, for i ∈ {x, y, z}.
To complete our proof, we need the following intermediate
result:80αβγ

Proposition 6: Any two VPs vi and vj , with i 6= j, verifycos θij ≤ 0. Furthermore, if vk ∈ D, with k 6= i and k 6= j,then cos θij ≥ − 31 .

Degrees

60
40

The method avoids standard intermediate steps such as featuredetection and correspondence, or edge detection and linking.
Experimental results show that the method is able to handle
low-quality video sequences, even with many spurious edges.

Proof: The first statement comes directly from (19).
To prove the second statement, we obtain min cos θij =2− minfdi dj , w.r.t. the principal point coordinates Pi and Pj ,over the domain defined by Pi2 + Pj2 ≤ f 2 /2. The minimum
occurs for |Pi | = |Pj | = f2 with value −1/3.
√
Since 12 acos(− 13 ) = atan 2 ≈ 54.7◦ , the shaded area in
Figure 3 is a simple consequence of Proposition 6. To show(4), consider an orientation O and let vi be a VP in thisshaded area. Proposition 2 then guarantees the existence ofan equiprojective orientation O∗ , satisfying: (i) vz∗ = vi ,and (ii) dx ≤ dy . From (2) –¤ (3) we√ have,√due¤ to (i), thatβ ∗ ∈ ]−π/2, π/2] and γ ∗ ∈ −atan 2, atan 2 , and due to(ii), that α∗ ∈ ]−π/2, π/2].

A PPENDIXR EFERENCES

Here we prove Proposition 3. From (1) – (2), we have (fori, j ∈ {x, y, z})
³ 2
´
(
⇐i=j
f 2 Pf 2 − 1T(19)vi vj =i−f 2
⇐ i 6= j,[1] A. Martins, P. Aguiar, and M. Figueiredo. Navigating in Manhattan:
3D orientation from video without correspondences. In Proc. IEEE Int.
Conf. on Image Processing, Barcelona, Spain, 2003.
[2] O. Faugeras. Three-Dimensional Computer Vision. MIT Press, Cambridge, MA, 1993.
[3] R. Hartley and A. Zisserman. Multiple View Geometry. Cambridge
University Press, 2000.

20
0
0

50

100
150
Frame number

200

250

Fig. 7. Top: frames 110, 130, 150, and 170 of a third video sequence.
Bottom: camera angle estimates.

IEEE TRANSACTIONS ON PATTERN ANALYSIS AND MACHINE INTELLIGENCE, 2005 (ACCEPTED)[4] J. Mundy and A. Zisserman, editors. Geometric Invariants in Computer
Vision. MIT Press, Cambridge, MA, 1992.
[5] E. Lutton, H. Maı̂tre, and J. Lopez-Krahe. Contribution to the determination of vanishing points using Hough transform. IEEE Trans. on
PAMI, 16(4):430–438, 1994.
[6] S. Utcke. Grouping based on projective geometry constraints anduncertainty. In Proc. IEEE Int. Conf. on Computer Vision, Bombay,
India, 1998.
[7] J. Kosecka and W. Zhang. Video compass. In European Conf. on Comp.
Vision, Springer Verlag, LNCS 2350, 2002.
[8] B. Horn and E. Weldon Jr. Direct methods for recovering motion. Int.
Jour. Computer Vision, 2(1):51–76, 1988.
[9] G. P. Stein and A. Shashua. Model-based brightness constraints: On
direct estimation of structure and motion. IEEE Trans. on PAMI,
22(9):992–1015, 2000.
[10] J. Coughlan and A. Yuille. Manhattan world: Compass direction froma single image by Bayesian inference. In Proc. IEEE Int. Conf. on
Computer Vision, Corfu, Greece, 1999.
[11] J. Coughlan and A. Yuille. The Manhattan world assumption: Regularities in scene statistics which enable Bayesian inference. In Proc. Neural
Information Processing Systems, Denver, CO, USA, 2000.
[12] J. Deutscher, M. Isard, and J. MacCormick. Automatic camera calibration from a single Manhattan image. In Proc. European Conf. on
Computer Vision, Springer Verlag, LNCS 2350, 2002.
[13] G. Schindler and F. Dellaert.
Atlanta world: An expectationmaximization framework for simultaneous low-level edge grouping andcamera calibration in complex man-made environments. In Proc.
IEEE Comp. Soc. Conf. on Computer Vision and Pattern Recognition,
Washington, DC, 2004.
[14] N. Gordon A. Doucet, N. Freitas. Sequential Monte Carlo Methods inPractice. Springer-Verlag, New York, 2001.
[15] M. Isard and A. Blake. CONDENSATION – conditional density
propagation for visual tracking. Int. Jour. of Computer Vision, 29(3):5–
28, 1998.

7