14th International Symposium on Particle Image Velocimetry – ISPIV 2021
August 1–5, 2021

Highly accurate optical flow method based on
volumetric segmentation for 3D PIV

Hua Yang1∗, Hang Shi1, Jin Lu1, Menggang Kang1, Zhouping Yin1

1 Huazhong University of Science and Technology, State Key Laboratory of Digital Manufacturing
Equipment and Technology, Wuhan, China

∗ huayang@hust.edu.cn

Abstract
In this study, we present a new three-dimensional optical flow method based on volumetric segmentation
for the velocity estimation of fluid flow. The proposed method uses a segmented smoothness term that
is designed on the assumption that the particle velocity varies continuously in each segmented volume
and discontinuously on the surfaces of the segmented volumes. Subsequently, the data term is proposed
on the basis of the segmented volumes and the fluid mass conservation equation, which is derived from
the Reynolds transport equation. In addition, the robust local level-set method is applied to segment the
particle volume according to the velocity distribution of fluid flow. The proposed method is evaluated
quantitatively on synthetic data and qualitatively on experimental data, and the velocity results are compared
to the advanced 3D velocity estimation methods. The results indicate that the proposed method can obtain
velocity fields with greater measurement accuracy for Tomo-PIV.

1 Introduction
As a non-contact velocity measurement technique, particle image velocimetry (PIV) has a very wide range of
applications such as aerodynamics, biological fluid mechanics, micro-scale complex flow, etc (Raffel et al.,
2018). With the development of laser lighting technology and image acquisition technology, PIV has been
extended from 2D velocity measurement to 3D velocity measurement (Lasinger et al., 2020). Tomographic
particle image velocimetry (Tomo-PIV), as a 3D three-component velocity measurement technology devel-
oped from 2D PIV, has been studied by many researchers (Scarano, 2013). In Tomo-PIV measurement, the
tracer particles are first seeded into the fluid flow. Then, a thick laser slice is used to illuminate a mea-
surement volume, which is observed by multiple cameras distributed in different views. Subsequently, the
particle distribution in the region of interest is reconstructed by the multiplicative algebraic reconstruction
technique (MART) from the multi-view image in each frame (Elsinga et al., 2006). Finally, the 3D velocity
field within a chosen interrogation volume is obtained by 3D cross-correlation method (Discetti and Astarita,
2012). Among the various technical problems of Tomo-PIV, improving the accuracy and resolution of the
velocity field estimated from the particle volume has received the most attention (Scarano, 2013).

After decades of development, the velocity estimation methods in Tomo-PIV can be divided into cross-
correlation-based methods and optical flow-based methods. In 3D cross-correlation methods, the particle
volumes are divided into several fixed-size interrogation volumes, the correlation coefficient of two interro-
gation volumes in two successive particle volumes is calculated through the fast Fourier transform, then the
velocity vector in each interrogation volume is determined by the position of the correlation peak (Lu et al.,
2021). Because of its strong robustness, the 3D cross-correlation methods have been widely used in Tomo-
PIV. Discetti and Astarita (2012) presented a 3D cross-correlation method based on voxel binning for quick
estimation of the predicted displacement field. A 3D cross-correlation method based on gradient iterative
volume deformation was proposed by Cheminet et al. (2014). Although these 3D cross-correlation methods
have achieved significant improvements in accuracy, they still cannot achieve accurate velocity field mea-
surement. The velocity vector obtained by the 3D cross-correlation methods can be interpreted as the spatial
average velocity in the interrogation volume. The resolution of the estimated velocity field is restricted by
the size of the interrogation volume, which makes it unsuitable for flow field measurement with a rapidly
changing velocity field.



Alternatively, the global optical flow method, which can provide velocity fields with pixel-level resolu-
tion, has been verified in 2D PIV and extended to Tomo-PIV (Ruhnau et al., 2005). The 3D global optical
flow methods use a global energy function composed of a data term and smoothness term to replace the in-
terrogation volume in the 3D cross-correlation methods, in which the data term assumes that the given points
retain the same intensity in the particle volume along their trajectories, and the smoothness term assumes
that all neighboring points have similar motions. Alvarez et al. (2009) first introduced the global optical flow
method to the research field of Tomo-PIV, which includes the incompressibility of the flow as a constraint to
the minimization problem. Lasinger et al. (2017) presented a 3D version of the global optical flow model,
augmented with a physically based smoothness term for incompressible fluids. Although these 3D global
optical flow methods are constantly improving the measurement accuracy of the velocity field in Tomo-PIV,
they still cannot obtain results with higher accuracy because of global energy function constraints (Lu et al.,
2019). Because the global energy function constraints will make the rapidly changing velocity field too
smooth, and the fine structure cannot be preserved. Therefore, improving the accuracy of 3D global optical
flow methods for rapidly changing velocity field measurement is an urgent problem to solve.

In this study, we propose a volumetric-segmentation-based optical flow (VS-OF) method for highly ac-
curate Tomo-PIV measurements. A segmented smoothness term is designed on the assumption that the
velocity varies continuously within each segmented volume and discontinuously on the surfaces of the seg-
mented volumes. Subsequently, the data term is proposed on the basis of the segmented volumes and the
fluid mass conservation equation, which is derived from the Reynolds transport equation. In addition, the
robust local level-set method is applied to segment the particle volume according to the velocity distribution
of the fluid flow. The rest of the present paper is organized as follows. In section 2, the proposed VS-OF
method is described in detail. Experimental analysis is given in section 3. Conclusions drawn from the
present work are summarized in section 4.

2 The proposed VS-OF method

2.1 Framework of the VS-OF method
In Tomo-PIV measurement, the accurate measurement of the velocity field is of great significance for the
subsequent in-depth analysis of the flow structure. In this study, a volumetric-segmentation-based optical
flow method is proposed. The main process of the VS-OF method is shown in Fig. 1. An initial velocity
field (u0,v0,w0) is first estimated by the 3D optical flow of the VS-OF method. Subsequently, the initial
velocity field is provided as an initialization for the next velocity field segmentation task. According to the
segmentation results of the velocity field, the next step of optical flow calculation is guided. After several
iterations, the final finer velocity field (u,v,w) is obtained.

In traditional 3D global optical flow methods, due to the global energy function constraints, the velocity
field is smoothed in the iterative calculation of optical flow. Therefore, in order to avoid smoothing, the VS-
OF method segments the 3D velocity field calculated in the previous step in each optical flow iteration. It
should be noted that the segmentation operation is implemented on the three components of the velocity field
rather than on the velocity magnitude, as shown in Fig. 1. Then, the particle volumes E (t) and E (t +dt)
are segmented according to the segmentation result of the velocity field, and the segmented velocity field
is taken as the initial value of the next optical flow iteration. The above process is repeated to obtain an
accurate 3D velocity field. The proposed VS-OF method has excellent comprehensive performance in terms
of flow velocity estimation due to the process of velocity field segmentation. The segmented smoothness
term and data term in the VS-OF method will be presented and discussed in the following subsections.

2.2 A segmented smoothness term
Inspired by the level set segmentation method (Li et al., 2007), a segmented smoothness term is designed
on the assumption that the velocity varies continuously within each segmented volume and discontinuously
on the surfaces of the segmented volumes. In addition, the segmented smoothness term can segment the
particle volume according to the 3D velocity distribution of the fluid flow. Without loss of generality, we
simplify the smoothness term to consider only the u component—that is, to minimize Eq. 1:∫

Ωu

Kσ|∇u|2dx+
∫

ΩC
u

Kσ|∇u|2dx+µ |Su| (1)



Figure 1: Concept of the proposed VS-OF method.

where u = (u,v,w)T is the velocity vector at location x, Ωu denotes the segmented volumes in the u compo-
nent, ΩC

u is the complement of Ωu, Su is the closed parameterized surface of the segmented volumes Ωu, |Su|
is the area of surface Su, µ is a parameter constraining the surface of a segmented volume to be as smooth as
possible, and Kσ denotes a Gaussian kernel with standard deviation σ.

Drawing lessons from the idea of (Li et al., 2007), the unknown evolution surface is replaced with the
level set function Φu(x) : Φu(x) > 0 if x is inside the closed surface Su, Φu(x) < 0 if x is outside Su, and
Φu(x) = 0 if x is on Su. We can approximate |Su| =

∫
Ω |∇H(Φu)|dx =

∫
Ω δ(Φu) |∇Φu|dx, and Eq. 1 can be

rewritten as follows:∫
Ω

(
H(Φu)Kσ|∇u|2 +(1−H(Φu))Kσ|∇u|2

)
dx+

∫
Ω

µδ(Φu) |∇Φu|dx

=
∫

Ω

(
Kσ|∇u|2 +µδ(Φu) |∇Φu|

)
dx (2)

where H(Φ) is the Heaviside function, which is equal to 1 if Φ ≥ 0 and to 0 if Φ < 0; δ(Φ) is the Dirac
function, and δ(Φ) = H ′(Φ). Generally, the Heaviside function H(Φ) in Eq. 2 is approximated by a smooth
function defined by:

H(Φ) =
1
2

[
1+

2
π

arctan(Φ)

]
(3)

The Dirac function δ(Φ), the derivative of H(Φ), is the following smooth function:

δ(Φ) = H ′(Φ) =
1
π

1
1+Φ2 (4)

We now extend the segmented smoothness term to three components, and Eq. 2 can be expressed as
follows:

EV S
smooth =

∫
Ω

(
Kσ|∇u|2 +µδ(Φu) |∇Φu|

)
dx+

∫
Ω

(
Kσ|∇v|2 +µδ(Φv) |∇Φv|

)
dx

+

∫
Ω

(
Kσ|∇w|2 +µδ(Φw) |∇Φw|

)
dx (5)



2.3 A segmented data term
Define a new segmentation surface S0 based on the segmented surface Su, Sv and Sw:S0 = Su ∪ Sv ∪ Sw =
{(x,y,z)|Φu(x,y,z) = 0,Φv(x,y,z) = 0,Φw(x,y,z) = 0}, and Ω0 denotes the segmented volume surrounded
by S0. Consistent with the smoothness term, we replace the segmented surface S0 with the level set function
Φ0. Therefore, we propose the segmented data term on the basis of the fluid mass conservation equation,
which is derived from the Reynolds transport equation:

∫
Ω0

(
∂E
∂t

+∇E ·u+E divu
)2

dx+
∫

ΩC
0

(
∂E
∂t

+∇E ·u+E divu
)2

dx (6)

In addition, for the segmented volumes, the material density in the fluid remains constant, and the
divergence of the velocity is zero: divu = ∇ · u = 0. In other words, the fluid is incompressible in each
segmented volume Ω0. We can derive a simplified segmented data term by constraining divu = 0:

EV S
data =

∫
Ω

H(Φ0)

(
∂E
∂t

+∇E ·u
)2

dx+
∫

Ω
(1−H(Φ0))

(
∂E
∂t

+∇E ·u
)2

dx

=
∫

Ω

(
∂E
∂t

+∇E ·u
)2

dx (7)

where Ω is the entire particle volume. Based on Eqs. 5 and 7, the model of the proposed VS-OF method is
as follows:

EV S=αEV S
data +EV S

smooth

=α
∫

Ω

(
∂E
∂t

+∇E ·u
)2

dx+
∫

Ω

(
Kσ|∇u|2+µδ(Φu) |∇Φu|

)
dx

+
∫

Ω

(
Kσ|∇v|2+µδ(Φv) |∇Φv|

)
dx+

∫
Ω

(
Kσ|∇w|2+µδ(Φw) |∇Φw|

)
dx (8)

where α is a parameter controlling the balance between the data term and smoothness term. Finally, the
multiscale technique and image warping operation are implemented in the segmented particle volumes to
avoid falling into local optimal solutions and obtain an accurate velocity field.

3 Experimental analysis

3.1 Experimental setup
In order to make a quantitative comparison with the proposed VS-OF method, we extend the multiple-pass
cross-correlation method provided by Thielicke and Stamhuis (2014) and the global optical flow method
proposed by Ruhnau et al. (2005) to three dimensions. In terms of numerical evaluation indicators, the
root mean square error (RMSE) and average angle error (AAE), which are often used to verify accuracy
in the field of optical flow (Alvarez et al., 2009; Lu et al., 2019), are used to quantitatively evaluate the
experimental results:

RMSE =

√
1
N

N

∑
i=1

|ut
i −ue

i |
2 (9)

AAE =
1
N

N

∑
i=1

arccos
(

ut
i ·ue

i

|ut
i| |ue

i |

)
(10)

where ut and ue denote the ground truth and the estimated velocity field, respectively. N is the total number
of voxels in the estimated particle volume, and the subscript i represents the voxel coordinates.



Method

Synthetic 

particle volume

Reconstructed 

particle volume

RMSE AAE RMSE AAE

3D  cross-correlation 0.46 9.77 0.62 13.38

3D global optical flow 0.24 6.68 0.47 12.57

The proposed VS-OF 0.21 5.76 0.39 10.49

Figure 2: Left: Velocity contours of the ground truth in the u component. Right: RMSE and AAE errors of
the velocity fields for the synthetic and reconstructed particle volumes estimated by different methods.

Figure 3: Velocity contours estimated from the synthetic particle volumes (top) and reconstructed particle
volumes (bottom) in the u component on the z = 0 plane. Left to right: Contours estimated by the 3D
cross-correlation method, the 3D global optical flow method, and the proposed VS-OF method, as well as
the ground truth.

3.2 Synthetic data
In the synthetic data, the ground truth velocity field from the Johns Hopkins Turbulence Database provid-
ed by Li et al. (2008) is used to quantitatively evaluate the proposed VS-OF method. The velocity field
comprises a numerically generated 3D isotropic turbulent flow with an average velocity of 1.3 voxels and a
maximum velocity of 2.9 voxels. A synthetic particle volume of 640×640×300 voxels is generated with a
random seed density of approximately 0.1 particles per pixel. The particles are projected onto four cameras
with 1024×1024 pixels at viewing angles of ±35◦ with respect to the x = 0 plane and ±18◦ with respect to
the y = 0 plane. Then, the particle volume is reconstructed using the reimplementation of MART proposed
by Scarano (2013) with a quality factor Q = 0.79.

To distinguish the influence of MART, the synthetic particle volume (without MART) and the recon-
structed particle volume (with MART) are used to evaluate the different methods. The results for the syn-
thetic and reconstructed particle volumes estimated by the different methods are presented in Fig. 2. Due to
the influence of MART, the RMSE and AAE errors of the reconstructed particle volume are larger than those
of the synthetic particle volume, but the errors of the proposed VS-OF method are smaller than those of the
other two methods on both the synthetic and the reconstructed particle volumes. The velocity contours in
the u component on the z = 0 plane estimated by the proposed VS-OF method exhibit behavior similar to
the ground truth, as shown in Fig. 3.



Figure 4: Left: Real experiment scene of cylinder wake (Michaelis et al., 2006). Right: Spectrum analysis
of the velocity fields for the experimental data.

Figure 5: Velocity contours in the xy-slice of the flow on the z = 203 plane. Top to bottom: Reference
flow field in u, v, and w components provided by Michaelis et al. (2006); estimated flow field in u, v, and w
components by the proposed VS-OF method.

3.3 Experimental data
In this subsection, we present qualitative results of the experimental data in the water flow (Michaelis et al.,
2006). The experimental data package contains of two particle volumes of 2107×1434×406 voxels. It is
a Karman vortex street behind a cylinder with a diameter D =12mm, which is positioned to the left of the
volume, with water flowing to the right, as shown in Fig. 4.

The velocity contours in the xy-slice of the flow on the z = 203 plane estimated by the proposed VS-OF
method are shown in Fig. 5. In addition to our own results, we also show the reference flow field provided in
the experimental data package, which is processed by the 3D cross-correlation method (Elsinga et al., 2006).
The results indicate that our method can better preserve the flow field structure with a rapidly changing
velocity field compared to the reference flow. To better investigate the performance of the proposed VS-OF
method, spectrum analysis is performed for the two processing methods, as shown in Fig. 4. The spectrum
obtained by the proposed VS-OF method is closer to the -5/3 spectrum slope. From this experiment, we can
conclude that the proposed VS-OF method is suitable for velocity field estimation in experimental data.

4 Conclusions
In this study, we propose a novel 3D optical flow method based on volumetric segmentation for the velocity
estimation of fluid flow. The proposed VS-OF method uses a segmented smoothness term designed on the
assumption that the particle velocity varies continuously in each segmented volume and discontinuously on
the surfaces of the segmented volumes. The data term is proposed on the basis of the segmented volumes and
the fluid mass conservation equation, which is derived from the Reynolds transport equation. In addition, the



robust local level-set method is applied to segment the particle volume according to the velocity distribution
of the fluid flow. To verify the performance of the proposed VS-OF method, we evaluate it quantitatively on
synthetic data and qualitatively on experimental data. It can be seen from the synthetic data that compared
with the 3D cross-correlation method and 3D global optical flow method, the proposed VS-OF method can
obtain velocity fields with greater accuracy. Finally, we demonstrate the good performance of the proposed
VS-OF method in a real cylinder wake flow experiment.

Acknowledgements
This work was supported by the National Natural Science Foundation of China (Grant No. 51875228), the
National Key R&D Program of China (Grant No. 2020YFA0405700), and the National Defense Science
and Technology Innovation Special Zone Project (Grant No. 193-A14-202-01-23).

References
Alvarez L, Castano C, and Garcı́a M (2009) A new energy-based method for 3D motion estimation of

incompressible PIV flows. Computer Vision and Image Understanding 113:802–810

Cheminet A, Leclaire B, and Champagnat F (2014) Accuracy assessment of a Lucas-Kanade based correla-
tion method for 3D PIV. in 17th International Symposium on Applications of Laser Techniques to Fluid
Mechanics, Lisbon, Portugal, July 7-10

Discetti S and Astarita T (2012) Fast 3D PIV with direct sparse cross-correlations. Experiments in fluids
53:1437–1451

Elsinga GE, Scarano F, and Wieneke B (2006) Tomographic particle image velocimetry. Experiments in
fluids 41:933–947

Lasinger K, Vogel C, and Pock T (2020) 3D fluid flow estimation with integrated particle reconstruction.
International Journal of Computer Vision 128:1012–1027

Lasinger K, Vogel C, and Schindler K (2017) Variational 3D-PIV for incompressible fluid flow estimation.
in 12th International Symposium on Particle Image Velocimetry, Busan, Korea, June 18-22

Li C, Kao CY, and Gore JC (2007) Implicit active contours driven by local binary fitting energy. in IEEE
Conference on Computer Vision and Pattern Recognition, Minneapolis, MN, USA, June 17-22

Li Y, Perlman E, and Wan M (2008) 3D fluid flow estimation with integrated particle reconstruction. Journal
of Turbulence 9:1–29

Lu J, Yang H, and Zhang Q (2019) A field-segmentation-based variational optical flow method for PIV
measurements of nonuniform flows. Experiments in Fluids 60:142

Lu J, Yang H, and Zhang Q (2021) An accurate optical flow estimation of PIV using fluid velocity decom-
positions. Experiments in Fluids 62:78

Michaelis D, Poelma C, and Scarano F (2006) A 3D time-resolved cylinder wake survey by tomographic
PIV. in 12th International Symposium on Flow Visualization, Göttingen, Germany, September 10-14

Raffel M, Willert CE, and Scarano F (2018) Particle image velocimetry: a practical guide. Springer

Ruhnau P, Kohlberger T, and Schnörr C (2005) Variational optical flow estimation for particle image ve-
locimetry. Experiments in Fluids 38:21–32

Scarano F (2013) Tomographic PIV: principles and practice. Measurement Science and Technology
24:012001

Thielicke W and Stamhuis E (2014) PIVlab–towards user-friendly, affordable and accurate digital particle
image velocimetry in MATLAB. Journal of open research software 2:e30