Simultaneous Estimation of Optical Flow and Its Boundaries based on
the Dynamical System Model
Yuya Michishita, Noboru Sebe, Shuichi Enokida and Eitaku Nobuyama
Kyushu Institute of Technology, Iizuka, Fukuoka, Japan
Keywords:
Optical Flow, Optical Flow Boundaries, Dynamical System Model.
Abstract:
Optical flow is a velocity vector which represents the motion of objects in video images. Optical flow esti-
mation is difficult in the neighborhood of flow boundary. To resolve this problem, Sasagawa (2014) proposes
a modified dynamical system model in which one assumes that, in the neighborhood of flow boundaries,
the brightness flows in the perpendicular direction, and considers the resulting corrections to the brightness
constancy constraint. However, in that model, the correction is occurred even in place where the flow is contin-
uous. We propose a new model, which switches the conventional model and the proposed model in Sasagawa
(2014). As a result, we expect improvement of the estimate accuracy in place where the flow is continuous.
We conduct numerical experiments to investigate the improvements that the proposed model yields in the
estimation accuracy of optical flows.
1 INTRODUCTION
Optical flow is a velocity vector which represents the
motion of objects in video images. The estimation of
optical flow is a fundamental tool for the object mo-
tion measurements using a visual sensor, and is uti-
lized in various fields. Recently, high-accurate esti-
mation method of optical flow has been expected to
improve the performance of various video image anal-
ysis systems. Actually, the accurate boundaries infor-
mation of optical flow will contribute to establish ac-
curate estimation of the optical flow, and vice versa.
However, it is a chicken and egg dilemma. In this pa-
per, a new estimation method which has the capability
of simultaneous estimation of the optical flow and its
boundary is proposed.
The optical flow equation describes the brightness
constancy constraint and has been used to determine
optical flows in Horn and Schunck (1981) and Lu-
cas and Kanade (1981). Sebe et al. (2009) is one
of the optical flow estimation methods which utilize
the optical flow equation. Sebe et al. (2009) regard
the optical flow equation as a dynamical system and
apply Kalman filter to that dynamical system to esti-
mate the optical flow. Based on the method in Sebe et
al. (2009), an estimation method was proposed to es-
timate dense optical flow (Fukami et al., 2011). There
are mainly two advantages of the use of Kalman fil-
ter for dense optical flow estimation. One is that it
allows us to obtain the covariance matrix of the es-
timation error. The covariance matrix of the estima-
tion error provides a confidence of estimated optical
flow. This enables us to assess the results of optical
flow estimation, even if the actual values are not avail-
able in practice. The other advantage is that the mea-
surement residual, which can be obtained as the dif-
ference between the estimated intensity and the mea-
sured intensity, can detect the optical flow boundary.
At the boundary of optical flow, the brightness con-
stancy constraint does not hold. In other words, the
dynamical system model used for estimation is not
correct. This causes a large measurement residual.
Accordingly, measurement residual enables us to de-
tect the boundary of optical flow. However, in Fukami
et al. (2011), the boundary information of optical flow
is not used to improve the accuracy of estimation.
One difficulty in optical flow estimation is that
the estimation accuracy deteriorates near flow bound-
aries, such as occlusion. To resolve this problem,
Sasagawa (2014) proposes a modified dynamical sys-
tem model in which the brightness conserved quan-
tity flows in the perpendicular direction in the neigh-
borhood of flow boundaries. He also considers the
resulting corrections to the brightness constancy con-
straint. The model proposed in Sasagawa (2014) im-
proved estimation accuracy in the neighborhood of
flow boundaries. However, as Sasagawa (2014) did
not use any information about the measurement resid-
Michishita Y., Sebe N., Enokida S. and Nobuyama E.
Simultaneous Estimation of Optical Flow and Its Boundaries based on the Dynamical System Model.
DOI: 10.5220/0006101303770385
In Proceedings of the 12th International Joint Conference on Computer Vision, Imaging and Computer Graphics Theory and Applications (VISIGRAPP 2017), pages 377-385
ISBN: 978-989-758-227-1
Copyright
c
 2017 by SCITEPRESS – Science and Technology Publications, Lda. All rights reserved
377

ual, the brightness corrections arose even in places
where the flow was continuous. This deteriorate the
estimation accuracy in regions not near image bound-
aries.
To address this difficulty, we propose a new model
in which measurement residuals are used to switch
between the model proposed in Sasagawa (2014) and
conventional models. In our proposed model, we as-
sume that boundaries are present in regions with large
measurement residuals, and thus apply the corrections
in such regions; in regions where measurement resid-
uals are small we do not apply the corrections. Conse-
quently, no corrections arise in regions where bound-
aries do not exist, whereupon we expect improved
estimation accuracy. We conduct numerical experi-
ments to investigate the improvements that the pro-
posed model yields in the estimation accuracy of op-
tical flows.
2 DYNAMICAL-SYSTEM-BASED
ESTIMATION MODEL
In this section, following the conventional works in
Fukami et al. (2011) and Sasagawa (2014), we review
the dynamical system based optical flow estimation.
2.1 The Optical Flow Equation
The gradient method for optical flow estimation is a
technique in which estimations are made based on the
assumption that brightness values do not vary while
moving over infinitesimal time intervals (Horn and
Schunck, 1981; Lucas and Kanade, 1981). Denoting
the brightness value at time t and coordinates (x, y)
by I(x, y,t), they derive the following equation by as-
suming that brightness values remain unchanged after
moving:
dI
dt
=
∂I
∂x
u +
∂I
∂y
v +
∂I
∂t
= 0, (1)
where u and v are the optical flows in the x and y direc-
tions, this relation is known as the optical flow equa-
tion. Sebe et al. (2009) regard the optical flow equa-
tion (1) as a dynamical system and apply Kalman fil-
ter to estimate the state of dynamical system, i.e., the
optical flow.
2.2 Quasi-dense Estimation
Here we review the quasi-dense estimation method
of Fukami et al. (2011). The quasi-dense estima-
tion refers to an estimation method in which, rather
than applying optical flow estimation to all pixels, it
decimates the set of pixels and performs optical flow
estimation on just the remaining subset. After opti-
cal flow estimation of the subset of pixels, an inter-
polation operator is applied to interpolate optical flow
values, which are used to estimate brightness of all
pixels. The quasi-dense estimation considers, optical
flow estimation for a video image of height M and
width N with L frames. The vectors q
u
, q
v
, whose
components are the optical flows of each pixel [u,v]
T
,
are defined as follows.
q
u(k)
=

u
(1,1,k)
u
(1,2,k)
... u
(M,N,k)

T
, (2)
q
u(k)
=

v
(1,1,k)
v
(1,2,k)
... v
(M,N,k)

T
. (3)
The decimated flows ˜q
u
, ˜q
v
are given by
q
u
= M ˜q
u
, (4)
q
v
= M ˜q
v
, (5)
where M is an interpolating operator that interpolates
˜q
u
and ˜q
v
to obtain the flows q
u
and q
v
at all the pixels.
Also, using the brightness I of each pixel, q
I
is defined
as follows:
q
I(k)
=

I
(1,1,k)
I
(1,2,k)
... I
(M,N,k)

T
. (6)
Using these quantities, the state of dynamical system
q
(k)
is defined as follows:
q
(k)
=
h
q
T
I(k)
˜q
T
u(k)
˜q
T
v(k)
i
T
. (7)
Also, f (q
(k)
) is defined by
f (q
(k)
) =




q
I(k)
− [diag{M ˜q
u(k)
}D
x
+ diag{M ˜q
v(k)
}D
y
]q
I(k)
˜q
u(k)
˜q
v(k)




,
(8)
where D
x
, D
y
are partial differential operators with re-
spect to x, y. Using these equations, Fukami et al.
(2011) regard optical flow equations of all pixels as a
nonlinear dynamical system of the form
q
(k+1)
= f (q
(k)
), (9)
I
(k)
= Hq
(k)
, (10)
where I
(k)
, H are given by
I
(k)
=

I
(1,1,k)
I
(1,2,k)
... I
(M,N,k)

T
, (11)
H =

E
MN×MN
O
MN×2L

. (12)
The matrices E
m×n
, O
m×n
are the m×n identity matrix
and zero matrix, respectively.
2.3 Kalman Filter
Kalman filter is one of the most effective method for
estimating the state of dynamical systems. Fukami et
al. (2011) apply Kalman filter to the dynamical sys-
tem (9), (10) to estimate the optical flow. For this pur-
pose, they consider the following noise-added version
of the model:
VISAPP 2017 - International Conference on Computer Vision Theory and Applications
378



q
I(k+1)
˜q
u(k+1)
˜q
v(k+1)


= f (q
(k)
) +


η
I(k)
η
u(k)
η
v(k)


, (13)
I
(k)
= Hq
(k)
+ ζ
(k)
, (14)
where η
(k)
denotes system noise, while ζ
(k)
represents
measurement noise; they take these to be Gaussian
noise processes with covariance matrices Q
(k)
, R
(k)
.
Also, the Jacobian F
(k)
of f (q
(k)
) in equation (8) is
F
(k)
=


E + N(q
(k)
) J
u
(q
(k)
) J
v
(q
(k)
)
O E O
O O E


, (15)
where N(q
(k)
), J
u
(q
(k)
), J
v
(q
(k)
) are given by
N(q
(k)
) = −[diag{M ˜q
u(k)
}D
x
+diag{M ˜q
v(k)
}D
y
], (16)
J
u
(q
(k)
) = −diag{D
x
q
I(k)
}M , (17)
J
v
(q
(k)
) = −diag{D
y
q
I(k)
}M . (18)
The procedures of extended Kalman filter time up-
dates and measurement updates for the model of
equations (13), (14), and (15) are as follows.
Time updates:
ˆq
(k+1|k)
= f ( ˆq
(k|k)
), (19)
ˆ
F
(k)
=


E + N( ˆq
(k|k)
) J
u
( ˆq
(k|k)
) J
v
( ˆq
(k|k)
)
O E O
O O E


,
(20)
P
(k+1|k)
=
ˆ
F
(k)
P
(k|k)
ˆ
F
T
(k)
+ Q
(k)
. (21)
Measurement updates:
ˆq
(k+1|k+1)
= ˆq
(k+1|k)
+ K
(k+1)
{I
(k+1)
− H ˆq
(k+1|k)
},
(22)
P
(k+1|k+1)
= P
(k+1|k)
− K
(k+1)
HP
(k+1|k)
, (23)
K
(k+1)
= P
(k+1|k)
H
T
{HP
(k+1|k)
H
T
+ R
(k)
}
−1
, (24)
where ˆq
(m|k)
and P
(m|k)
are the estimate and the es-
timation error covariance matrix at time m from the
state at time k, respectively. The estimation process
proceeds by repeating the time update and measure-
ment update for each time point.
2.4 Correcting Estimated Values using
Measurement Residuals
In this subsection we review the method proposed in
Sasagawa et al. (2013) for correcting estimated val-
ues by using measurement residuals. The estimation
via extended Kalman filter may fail due to factors
such as nonlinearities or flow boundaries. If the esti-
mation fails, it does not subsequently recover, and for
this reason the estimated values should be corrected.
For cases in which they conclude that the estimation
has failed, based on measurement residuals, they ap-
ply a reset procedure of the following form.
• The estimated brightness are reset to that of mea-
sured values.
• Because the flow value is unknown, the estimated
flow values are reset to 0.
• For the covariance matrix, the diagonal elements
are reset to large initial values to facilitate effec-
tive estimation; the off-diagonal elements are re-
set to 0 because they have no information about
multivariable correlations.
2.5 Introduction of Artificial Diffusion
Term
In this subsection, we review the diffusion term intro-
duced by Sasagawa et al. (2013) . In methods (Sebe
et al., 2009; Fukami et al., 2011), the dynamics of the
optical flows are modeled as
˜q
u(k+1)
= ˜q
u(k)
, (25)
˜q
v(k+1)
= ˜q
v(k)
. (26)
These are augmented by a diffusion term capturing
the smoothing effect of the estimation, as follows:
˜q
u(k+1)
= D(D
2
x
+ D
2
y
) ˜q
u(k)
, (27)
˜q
v(k+1)
= D(D
2
x
+ D
2
y
) ˜q
v(k)
, (28)
where D is a diffusion coefficient that determines the
weight of the diffusion term. With the introduction of
the diffusion term, equation (8) takes the form
f (q
(k)
) =




q
I(k)
− [diag{M ˜q
u(k)
}D
x
+ diag{M ˜q
v(k)
}D
y
]q
I(k)
D(D
2
x
+ D
2
y
) ˜q
u(k)
D(D
2
x
+ D
2
y
) ˜q
v(k)




.
(29)
Also, the Jacobian of equation (29) becomes
F
(k)
=


E + N(q
(k)
) J
u
(q
(k)
) J
v
(q
(k)
)
O D(D
2
x
+ D
2
y
) O
O O D(D
2
x
+ D
2
y
)


.
(30)
2.6 Estimation Model with Brightness
Constancy Corrections
In this subsection, we review an optical flow estima-
tion model that takes into account the brightness con-
stancy constraint in the neighborhood of flow bound-
aries (Sasagawa, 2014). When the flow becomes dis-
continuous due to occlusions or other factors, the
Simultaneous Estimation of Optical Flow and Its Boundaries based on the Dynamical System Model
379

brightness of objects that had previously been visible
disappears, because it is obscured by some covering
object. However, the optical flow equation used in the
estimation model is based on the brightness constancy
constraint and does not account for the possibility that
brightness may be extinguished. For this reason, opti-
cal flow estimation is difficult in the neighborhood of
flow boundaries. To address this difficulty, Sasagawa
et al. (2013) reason that, in the event of an occlu-
sion, brightness - being a constant quantity - flows in
the perpendicular direction, and they add corrections
to the optical flow equation describing the volume of
the perpendicular brightness flow. These corrections
should be applied only at flow boundaries, thus they
construct the corrections to be proportional to the in-
ner product of the flows boundary normal vector and
the flow difference vector. However, because the nor-
mal vector to the flow boundary is unknown, they
replace it with the brightness gradient vector, which
they assume to point in the direction perpendicular to
the boundary. Then the equation to be satisfied for
each pixel takes the form shown in (33). The bright-
ness gradient
−→
∆I, which they use in place of the nor-
mal vector, and the difference vector are defined as
follows.
−→
∆I =
"
∂I
∂x
∂I
∂y
#
, (31)

∆u
∆v

=
"
∂u
∂x
∂u
∂y
∂v
∂x
∂v
∂y
#
−→
∆I
q
I
2
x
+ I
2
y
. (32)
They reason that the correction should be proportional
to the inner product of the vectors (31) and (32), so
the brightness constancy constraint is modified as fol-
lows.
∂I
∂t
= −

∂I
∂x
u +
∂I
∂y
v

−
−→
∆I
T

∆u
∆v

z. (33)
These corrections are large when large movements
are made in the direction of the increasing brightness
gradient. In the second term on the right-hand side
of equation (33), they estimate only the quantity z by
taking all other factors as proportionality constants.
Using this model, they expect improved accuracy in
the flow estimation near boundaries.
The above equation is a model for the equation
that should hold for each pixel. They next present
their model-description method for the image as a
whole. They begin by adding z to the state of the esti-
mation model.
q
z(k)
=

z
(1,1,k)
z
(1,2,k)
... z
(M,N,k)

T
, (34)
where they take ˜q
z
to be
q
z
= M ˜q
z
. (35)
As for the optical flows ˜q
u
, ˜q
v
, they apply decimation
to q
z
to guarantee observability. They apply the es-
timation procedure to the decimated vector ˜q
z
. Then
the state q
(k)
may be defined as follows:
q
(k)
=
h
q
T
I(k)
˜q
T
u(k)
˜q
T
v(k)
˜q
T
z(k)
i
T
. (36)
Similar to the optical flows ˜q
u
, ˜q
v
, the dynamics of the
perpendicular flow is modeled as
˜q
z(k+1)
= ˜q
z(k)
. (37)
Using equations (36) and (37), f (q
(k)
) in equation
(29) is modified as follows.
f (q
(k)
) =














q
I(k)
− {diag{M ˜q
u(k)
}D
x
+diag{M ˜q
v(k)
}D
y
}q
I(k)
−diag{diag{I
x(k)
}∆u
(k)
+diag{I
y(k)
}∆v
(k)
}M ˜q
z(k)




D(D
2
x
+ D
2
y
) ˜q
u(k)
D(D
2
x
+ D
2
y
) ˜q
v(k)
˜q
z(k)










.
(38)
The quantities I
x(k)
, I
y(k)
, ∆u, ∆v here are defined as
follows.
I
x(k)
=
h
∂
∂x
I
(1,1,k)
∂
∂x
I
(1,2,k)
...
∂
∂x
I
(M,N,k)
i
T
=D
x
q
I(k)
, (39)
I
y(k)
=
h
∂
∂y
I
(1,1,k)
∂
∂y
I
(1,2,k)
...
∂
∂y
I
(M,N,k)
i
T
=D
y
q
I(k)
, (40)
∆u
(k)
=

∆u
(1,1,k)
∆u
(1,2,k)
... ∆u
(M,N,k)

T
,(41)
∆v
(k)
=

∆v
(1,1,k)
∆v
(1,2,k)
... ∆v
(M,N,k)

T
. (42)
The elements ∆u
(x,y,k)
, ∆v
(x,y,k)
in equations (41) and
(42) are the quantities ∆u, ∆v of equation (32) for
pixel (x, y) at time k. The calculation of the Jacobian
is given in APPENDIX.
3 NEW PROPOSED METHOD
In this section we describe a new model that improves
the model proposed in Sasagawa (2014).
3.1 Switching between Estimation
Models based on Measurement
Residuals
Using the model proposed in Sasagawa (2014), we
successfully confirmed improved estimation accuracy
VISAPP 2017 - International Conference on Computer Vision Theory and Applications
380

in cases where the flow becomes discontinuous due
to occlusions or other factors. However, their model
also applied corrections to regions in which the flow
was continuous; this caused an outflow of brightness
and yielded the problematic result that the estimation
accuracy of the model was inferior to that without the
correction. To address this difficulty, we here propose
a new model in which measurement residuals are used
to switch between the model accounting for bright-
ness constancy corrections and existing models. Our
model modifies equation (33) as follows.
∂I
∂t
= −

∂I
∂x
u +
∂I
∂y
v

−
−→
∆I
T

∆u
∆v

Az. (43)
The quantity A is defined by
A = diag(a
(1,1,k)
, a
(1,2,k)
, ... , a
(M,N,k)
), (44)
where,
a
(i, j,k)
=

1 |(I
(i, j,k)
− q
I(k|k−1)
| ≥ τ),
0 |(I
(i, j,k)
− q
I(k|k−1)
| < τ).
(45)
In equation (45), I
(i, j,k)
− q
I(k|k−1)
is the measure-
ment residual and τ is a threshold value. The func-
tion a(i, j, k) is a step function taking the value 1 if
the measurement residual exceeds the threshold, and
0 otherwise. This ensures that, in a region of flow
discontinuities, the correction represented by the sec-
ond term of equation (43) is present, while this term
is removed in regions of continuous flow, and ensures
that no correction is applied. Equation (33) has the
flaw that the correction term is present even in re-
gions where there is no occlusion and the flow is con-
tinuous, and thus results in an outflow of brightness
that degrades the estimation accuracy. In contrast,
equation (43) uses measurement residuals and a step
function to ensure that the model that attempts to ap-
ply brightness constancy corrections is not used in re-
gions where the flow is continuous, and thus promises
improved estimation accuracy. Because we modified
equation (33) according to (43), f (q
(k)
) in equation
(38) is modified to
f
A
(q
(k)
) =














q
I(k)
− {diag{M ˜q
u(k)
}D
x
+diag{M ˜q
v(k)
}D
y
}q
I(k)
−diag{diag{I
x(k)
}∆u
(k)
+diag{I
y(k)
}∆v
(k)
}AM ˜q
z(k)




D(D
2
x
+ D
2
y
) ˜q
u(k)
D(D
2
x
+ D
2
y
) ˜q
v(k)
˜q
z(k)










.
(46)
The Jacobian of (46) is given in APPENDIX.
4 NUMERICAL EXPERIMENTS
In this section we conduct numerical experiments to
test the efficacy of the model proposed in Section 3.1.
4.1 Description of Experiments
The video images used in these numerical experi-
ments are taken from the sleeping2 sample in the MPI
Sintel benchmark (Butler et al., 2012) . For our tests,
we extract the regions delineated by white frames in
Figures 1 and 2. We performed estimation on a 60
60 pixel region. To track flow boundaries through
all frames, we displace the test region by 1 pixel per
frame in the imaging direction (Brox et al., 2004).
The covariance matrices are set to
Q
I(k)
= 0.001E, (47)
Q
u(k)
= Q
v(k)
(48)
= 0.4E, (49)
R
(k)
= 0.1E. (50)
We set the measurement residual τ in equation (45) to
10. For numerical integration we used the 4th-order
Runge-Kutta method with 4 time steps per frame. We
used the endpoint error to assess the optical flow es-
timation accuracy for each pixel. The endpoint error
is the distance between the tips of the estimated flow
vector
˜
U and the actual flow vector U and may be de-
termined by using equation (51).
q
(
˜
U −U)
T
(
˜
U −U). (51)
Figure 1: Frame 1 of video image containing a boundary.
Figure 2: Frame 1 of video image not containing a bound-
ary.
Simultaneous Estimation of Optical Flow and Its Boundaries based on the Dynamical System Model
381

4.2 Investigating the Difference of
Estimate Accuracy Among Three
Models
We investigate the difference of the estimate accu-
racy among the model proposed in Section 3, Fukami
model and the model in Sasagawa (2014) by numeri-
cal experimants. Here, Fukami model is that Fukami
et al. (2011) includes correcting estimated values us-
ing measurement residuals and diffusion term. Fig-
ures 3 and 4 show the mean endpoint error for im-
age frames containing boundaries and not containing
boundaries, respectively. The mean endpoint errors
for all frames are listed in Table 1.
0 5 10 15 20 25 30 35 40 45 50
0
0.1
0.2
0.3
0.4
0.5
Number of Frame
Average of EndPointError
Fukami model
Sasagawa (2014)
Proposed model
Figure 3: Mean endpoint error for each frame of the bound-
ary containing images.
0 5 10 15 20 25 30 35 40 45 50
0
0.1
0.2
0.3
0.4
0.5
Number of Frame
Average of EndPointError
Fukami model
Sasagawa (2014)
Proposed model
Figure 4: Mean endpoint error for each frame of the non
boundary containing images.
Table 1: Mean endpoint errors for all frames.
Images
containing
boundary
Images not
containing
boundary
Fukami model 0.3875 0.1642
Sasagawa (2014) 0.3549 0.1674
Proposed model 0.3454 0.1636
Figure 5 shows the estimation frame containing
a boundary. For this image, Figures 6 and 7 show
the endpoint error of frame 30 for Sasagawa (2014)
and the proposed model, respectively. Next, Figure
8 shows the estimation frame that does not contain a
boundary. For this image, the endpoint error of frame
40 obtained with Sasagawa (2014) and the proposed
model are shown in Figures 9 and 10, respectively.
Figure 5: Estimation region for frame 30 of the image con-
taining a boundary.
Figure 6: Endpoint er-
ror for image of Figure
5 (frame 30, Sasagawa
(2014)).
Figure 7: Endpoint er-
ror for image of Figure
5 (frame 30, proposed
model).
Figure 8: Estimation region for frame 40 of the image not
containing a boundary.
Figure 9: Endpoint er-
ror for image of Figure
8 (frame 40, Sasagawa
(2014)).
Figure 10: Endpoint er-
ror for image of Figure
8 (frame 40, proposed
model).
For the images containing a boundary, Figures 11
and 12 show values of the correction term z deter-
mined by Sasagawa (2014) and the proposed model,
respectively. For the images not containing a bound-
ary, Figures 13 and 14 show values of the correction
term z determined by Sasagawa (2014) and the pro-
posed model, respectively.
VISAPP 2017 - International Conference on Computer Vision Theory and Applications
382

10
20
30
40
50
60
10
20
30
40
50
60
0
10
20
30
40
50
60
Figure 11: Correction term z determined by
Sasagawa (2014) for frame 30 of the image con-
taining a boundary.
10
20
30
40
50
60
10
20
30
40
50
60
0
10
20
30
40
50
60
Figure 12: Correction term z determined by the
proposed model for frame 30 of the image con-
taining a boundary.
10
20
30
40
50
60
10
20
30
40
50
60
0
10
20
30
40
50
60
Figure 13: Correction term z determined by
Sasagawa (2014) for frame 40 of the image not
containing a boundary.
10
20
30
40
50
60
10
20
30
40
50
60
0
10
20
30
40
50
60
Figure 14: Correction term z determined by the
proposed model for frame 40 of the image not con-
taining a boundary.
5 DISCUSSION
Our numerical experiments to compare the estimation
accuracy of Sasagawa (2014) and the proposed model
confirm the efficacy of our proposed model. For es-
timation of the image containing a boundary, Figure
3 shows that the estimation accuracy of the proposed
model is superior to that of Sasagawa (2014). As we
see in Figures 6 and 7, this is because the measure-
ment residuals are small, whereas the model proposed
in this paper does not apply corrections in regions of
continuous flow, and so yield improved estimation ac-
curacy. From Figure 4 we see that the proposed model
also achieves superior estimation accuracy when esti-
mating regions that do not contain boundaries. From
Figures 12 and 14 we see that the proposed model
yields a value of 0 for the correction term z in re-
gions where the flow is continuous. This confirms
that, in regions not containing boundaries, our pro-
posed model correctly switches to the model that ig-
nores corrections.
At the end of this section, we briefly give a com-
ment on the computational cost. The examples are
carried out by MATLAB (Release 2011a) on a PC
(core i7-4790K, 4.00GHz with 16GB RAM). We use
4 cores to carry out the example with the proposed
method, and it takes about 40.3 [sec]. As the extended
Kalman filter calculates a huge error covariance ma-
trix, enormous computational cost is required. To re-
duce the computational cost, mathematical tools used
in data assimilation, e.g. the ensemble Kalman filter,
might be applicable. Reducing the computational cost
is one of our important future works.
6 CONCLUSIONS
Our numerical experiments confirm the effectiveness
of our proposed model. In the proposal model, we
switched between two models: one that accounts for
Simultaneous Estimation of Optical Flow and Its Boundaries based on the Dynamical System Model
383

corrections to the brightness constancy constraint and
one that does not use measurement residuals. This en-
sures that corrections are not applied in regions where
the flow is continuous, and so the estimation accu-
racy is superior to that of Sasagawa (2014). In addi-
tion, a step function is used to switch between models
based on measurement residuals in this model. In fu-
ture work, we plan to replace this with a sigmoid or
other nonlinear function to achieve further improve-
ments in estimation accuracy. Furthermore, the infor-
mation of measurement residual can be used to mod-
ify the diffusion coefficient. This modification real-
izes the nonlinear diffusion used in Brox and Weicket
(2002). This is also our future work.
REFERENCES
Brox, T., Bruhn, A., Pepenberg, N., and Weickert, J. (2004).
High accuracy optical flow estimation based on a the-
ory for warping. European Conference on Computer
Vision (ECCV), pages 25–36.
Butler, D. J., Wulff, J., Stanley, G. B., and Black, M. J.
(2012). A naturalistic open source movie for optical
flow evaluation. Conf. on Computer Vision (ECCV),
pages 611–625.
Fukami, T., Enokida, S., Sebe, N., and Nobuyama, E.
(2011). Quasi-dense optical flow estimation based on
dynamical system models. SICE Annual Conference,
pages 802–811.
Horn, B. K. P. and Schunck, B. G. (1981). Determining
optical flow. Artificial Intelligence, pages 185–203.
Lucas, K. and Kanade, T. (1981). An interactive image reg-
istration techniques with an application to stereo vi-
sion. In Proceedings of 7th IJCAL, pages 674–679.
Sasagawa, T. (2014). Optical flow estimation accuracy im-
provement in the flow near the boundary. Kyushu In-
stitute of Technology Graduate School of master’s the-
sis.
Sasagawa, T., Sebe, N., Enokida, S., Nobuyama, E., and
Fukami, T. (2013). Improvement of the accuracy of
the optical flow estimation by simultaneous estimation
of flow boundary. The 32nd SICE Kyushu Branch An-
nual Conference, pages 151–154.
Sebe, N., Nobuyama, E., and Enokida, S. (2009). Optical
flow estimation method in the moving image by the
dynamic systems theory approach. The 38th Sympo-
sium on Control Theory, pages 357–362.
Thomas Brox, J. W. (2002). Nonlinear matrix diffusion
for optic flow estimation. Pattern Recognition (Proc.
DAGM), pages 446–453.
APPENDIX
A.1. Jacobian of Sasagawa’s Model
The Jacobian F
(k)
of (38) is given.
F
(k)
=






E + N(q
(k)
) J
u
(q
(k)
)
−diag{g
I(k)
} −diag{g
u(k)
}M
O D(D
2
x
+ D
2
y
)
O O
O O
J
v
(q
(k)
) − diag{g
v(k)
}M J
z
(q
(k)
)
O O
D(D
2
x
+ D
2
y
) O
O E




, (52)
where
J
z
(q
(k)
) = −

diag{diag{I
x(k)
}∆u
(k)
+ diag{I
y(k)
}∆v
(k)
}

M . (53)
The correction term g is defined by
g =
−→
∆I
T

∆u
∆v

z. (54)
The notation g(x, y, k) is used to denote the value of
equation (54) at pixel (x, y) and time k. The quantities
g
I(k)
, g
u(k)
, g
v(k)
in equation (52) may be expressed as
follows in terms of partial differentials of q
I
, q
u
, q
v
.
g
I(k)
=
h
∂
∂I
g
(1,1,k)
∂
∂I
g
(1,2,k)
...
∂
∂I
g
(M,N,k)
i
T
, (55)
g
u(k)
=
h
∂
∂u
g
(1,1,k)
∂
∂u
g
(1,2,k)
...
∂
∂u
g
(M,N,k)
i
T
, (56)
g
v(k)
=
h
∂
∂v
g
(1,1,k)
∂
∂v
g
(1,2,k)
...
∂
∂v
g
(M,N,k)
i
T
.
(57)
Also, g
I
x
(k)
, g
I
y
(k)
, g
u
x
(k)
, g
u
y
(k)
, g
v
x
(k)
, g
v
y
(k)
are given
by
g
I
x
(k)
=
h
∂
∂I
x
g
(1,1,k)
∂
∂I
x
g
(1,2,k)
...
∂
∂I
x
g
(M,N,k)
i
T
, (58)
g
I
y
(k)
=
h
∂
∂I
y
g
(1,1,k)
∂
∂I
y
g
(1,2,k)
...
∂
∂I
y
g
(M,N,k)
i
T
, (59)
g
u
x
(k)
=
h
∂
∂u
x
g
(1,1,k)
∂
∂u
x
g
(1,2,k)
...
∂
∂u
x
g
(M,N,k)
i
T
, (60)
g
u
y
(k)
=
h
∂
∂u
y
g
(1,1,k)
∂
∂u
y
g
(1,2,k)
...
∂
∂u
y
g
(M,N,k)
i
T
, (61)
g
v
x
(k)
=
h
∂
∂v
x
g
(1,1,k)
∂
∂v
x
g
(1,2,k)
...
∂
∂v
x
g
(M,N,k)
i
T
, (62)
g
v
y
(k)
=
h
∂
∂v
y
g
(1,1,k)
∂
∂v
y
g
(1,2,k)
...
∂
∂v
y
g
(M,N,k)
i
T
. (63)
VISAPP 2017 - International Conference on Computer Vision Theory and Applications
384

The partial differentials of g in equation (54) are given
by
∂g
∂I
=
∂g
∂I
x
∂I
x
∂I
+
∂g
∂I
y
∂I
y
∂I
, (64)
∂g
∂u
=
∂g
∂u
x
∂u
x
∂u
+
∂g
∂u
y
∂u
y
∂u
, (65)
∂g
∂v
=
∂g
∂v
x
∂v
x
∂v
+
∂g
∂v
y
∂v
y
∂v
. (66)
Also,
∂g
∂I
x
= {I
3
x
u
x
+ I
3
y
(u
y
+ v
x
)
+I
x
I
2
y
(2u
x
− v
y
)}(I
2
x
+ I
2
y
)
−
3
2
z, (67)
∂g
∂I
y
= {I
3
y
v
y
+ I
3
x
(u
y
+ v
x
)
+I
2
x
I
y
(2v
y
− u
x
)}(I
2
x
+ I
2
y
)
−
3
2
z, (68)
∂g
∂u
x
= I
2
x
(I
2
x
+ I
2
y
)
−
1
2
z, (69)
∂g
∂u
y
= I
x
I
y
(I
2
x
+ I
2
y
)
−
1
2
z, (70)
∂g
∂v
x
= I
x
I
y
(I
2
x
+ I
2
y
)
−
1
2
z, (71)
∂g
∂v
y
= I
2
y
(I
2
x
+ I
2
y
)
−
1
2
z. (72)
From equations (64), (65), (66), the quantities
g
I
(k)
, g
u
(k)
, g
v
(k)
in equations (55), (56), (57) may be
computed as follows:
g
I(k)
=
h
diag{g
I
x
(k)
}D
x
+ diag{g
I
y
(k)
}D
y
i
q
I(k)
,
(73)
g
u(k)
=
h
diag{g
u
x
(k)
}D
x
+ diag{g
u
y
(k)
}D
y
i
q
u(k)
,
(74)
g
v(k)
=
h
diag{g
v
x
(k)
}D
x
+ diag{g
v
y
(k)
}D
y
i
q
v(k)
.
(75)
A.2. Jacobian of the Proposed Model
The Jacobian F
(k)
of (46) is given as follows:
F
A(k)
=






E + N(q
(k)
) J
u
(q
(k)
)
−diag{g
AI(k)
} −diag{g
Au(k)
}M
O D(D
2
x
+ D
2
y
)
O O
O O
J
v
(q
(k)
) − diag{g
Av(k)
}M J
Az
(q
(k)
)
O O
D(D
2
x
+ D
2
y
) O
O E




, (76)
where
J
Az
(q
(k)
) = −

diag{diag{I
x(k)
}∆u
(k)
+ diag{I
y(k)
}∆v
(k)
}

AM . (77)
Also, the correction term g
A
is defined as follows:
g
A
=
−→
∆I
T

∆u
∆v

Az. (78)
We denote by g
A(x,y,k)
the value of equation
(78) for pixel (x, y) at time k; The quantities
g
AI(k)
, g
Au(k)
, g
Av(k)
in equation (76) are defined as fol-
lows:
g
AI(k)
= Ag
I(k)
, (79)
g
Au(k)
= Ag
u(k)
, (80)
g
Av(k)
= Ag
v(k)
. (81)
Simultaneous Estimation of Optical Flow and Its Boundaries based on the Dynamical System Model
385