14th International Symposium on Particle Image Velocimetry – ISPIV 2021
Chicago, USA, August 1-4, 2021

Particle position prediction based on Lagrangian
coherency for flow over a cylinder in 4D-PTV

Ali Rahimi Khojasteh1∗, Dominique Heitz1, Yin Yang1, Lionel Fiabane1

1 INRAE, OPAALE, F-35044, Rennes, France
∗ correspondent author: ali.rahimi-khojasteh@inrae.fr

Abstract
Recent developments in time-resolved Particle Tracking Velocimetry (4D-PTV) consistently improved track-
ing accuracy and robustness. We propose a novel technique named ”Lagrangian coherent predictor” to esti-
mate particle positions within the 4D-PTV algorithm. We add spatial and temporal coherency information
of neighbour particles to predict a single trajectory using Lagrangian Coherent Structures (LCS). We found
that even a weak signal from coherent neighbour motions improves particle prediction accuracy in complex
flow regions. We applied Finite Time Lyapunov Exponent (FTLE) to quantify local boundaries (i.e. ridges)
of coherent motions. Synthetic analysis of the wake behind a smooth cylinder at Reynolds number equal
to 3900 showed enhanced estimation compared with the recent predictor functions employed in 4D-PTV.
Results of the experimental study of the same flow configuration are reported. We compared predicted posi-
tions with the optimised final positions of Shake The Box (STB). It was found that the Lagrangian coherent
predictor succeeded in estimating particle positions with minimum deviation to the optimised positions.

1 Introduction
This paper discusses a novel approach in time-resolved Particle Tracking Velocimetry (4D-PTV) to predict
particle positions over space and time based on Lagrangian Coherent Structures (LCS). In the classic 4D-
PTV algorithm, predicted particle positions are given to the optimisation process for further corrections. The
optimisation can deal with slight deviations between the predicted position and the true position. However,
the optimisation fails to find the true position if the deviation is large enough to have multi-candidates for a
single particle at tn+1. This shows the importance of having an appropriate prediction in dense and complex
motions. The proposed idea is initiated by arguing that predictions in PTV techniques focus on a single
particle individually, while this single particle is not acting alone [1]. We propose this approach to locally
differentiate coherent and non-coherent motions of neighbour particles around a single particle position to
improve prediction accuracy, as shown in Figure 1. On the other hand, the information of coherent particles
should be shared with each neighbour particle to predict their behaviour accurately and avoid misprediction.
The present study was designed as a complementary function for 4D-PTV algorithms such as STB [2] and
KLPT [3].

Briefly, it is assumed that particle positions are known for n time steps (four to five). Afterwards, a
mathematical prediction function is implemented to estimate particle positions for time-step n+1 followed
by ”shaking” and position refinement. It should be noted that the ”shaking” process tries to look for a can-
didate true position very close to the predicted position. If misprediction happens, no matter how many
times we perform ”shaking”, the true position is not achievable. This implies the importance of producing
accurate predictions. This paper seeks to investigate the possibilities of improvements in motion estimation
by adding meaningful physics into the prediction function. A simple prediction approach is a Polynomial
function, suggested by Schanz et al. [4], resulting in reasonable predictions and 3D particle position recon-
struction in simple flows [4, 5, 6]. However, significant off prediction occurs in case of flow associated with
complexities such as high turbulence level, high Reynolds number, and mixing flows [7]. In such conditions,
even by increasing the order of the Polynomial predictor functions from 3 to 10 [7], off prediction stays re-
mained. The solution for this challenge is implementing optimal temporal filtering such as the Wiener filter,
which has been first examined in 4D-PTV experiments by Schröder et al. [5]. Since then, this concept
became consistent in the STB studies due to its high robustness and accurate motion estimations [8, 7]. As



Figure 1: Particle prediction scenario from tn to tn+1: a), Known particle positions from history starting from
tn−4 up to tn (particle size is increasing gradually to show time-step differences); b), The trajectory (golden
line) obtained from filtered curve fitting of known particle positions; c), Prediction based on extrapolating
of the fitted trajectory (red dash line) from tn to tn+1; d), Modified prediction (grey dash line) using velocity
and acceleration information of coherent particles in neighbourhood of the target particle at tn.

mentioned, the Wiener filter showed robust behaviour in prediction with complex flows but still suffers in
high motion gradients. This implies the fact that the prediction function suffers from a lack of information
to find true positions. Worth mentioning that these prediction-based techniques rely on one particle individ-
ually, excluding it from surroundings. All the information we know from an individual particle is its history.
Even if we implement filtering and smoothing schemes such as STB using Wiener filter [8], our information
is limited by the history of the target particle, ignoring that every particle is spatially and temporally coherent
with a specific group of other particles following the same behaviour. This motivated us to take into account
a group of coherent motions for predicting a single particle. We propose to locally determine information
of coherent and non-coherent particles during the trajectory procedure by using the Finite-Time Lyapunov
Exponent (FTLE). More details of coherent motion detection are discussed in Section 2. After that, we
address the prediction function with the minimisation approach in Section 3. In the following Sections 4
and 5, we study and evaluate our proposed technique using synthetic and experimental case studies of the
wake over and behind a smooth cylinder at Reynolds number equal to 3900.

2 Coherent and non-coherent neighbours
Every single particle is spatially and temporally coherent with a specific cluster of other particles following
the same behaviour [1]. Figure 2 shows a schematic of coherent and non-coherent groups of neighbour
particles in different colours evolving by time (t1 − t5). A particle can spatially meet a group of other
particles in which there is no coherency link between them. There are many available concepts to identify
Lagrangian Coherent Structures (LCS) from looking for separatrix lines or surfaces which divide structures
into different coherent regions [9]. We introduced such a concept for the track initialisation in 4D-PTV
studies, named Lagrangian Coherent Track Initialisation (LCTI), in synthetic flow configurations of the
recent LPT Challenge [10] and for a real jet impingement experiment [1]. LCTI showed that FTLE is
locally applicable for surrounding sparse trajectories to determine coherent neighbours.

To this end, we defined a local Eulerian frame around each particle, while all neighbourhood particles
inside this area must be classified as coherent or non-coherent with the target particle. This frame is fixed
during a series of time-steps that can provide an Eulerian view of the neighbourhood behaviour. Velocity
values of the target particle are used to quantify the Eulerian frame size in each direction. If 2D / 3D
velocity values are equal in each direction, then the shape would be a circle / sphere around the target
particle; however, more directional gradients would adjust the shape (see Figure 1.d). All particles inside
the Eulerian frame in the same phase or with phase delay are considered neighbourhoods. By tracking each



Figure 2: Schematic of particle trajectories in 2D pair vortices starting from t1 to t5: t1), The target particle
with coherent neighbour particles located in a clockwise vortex (golden cluster), non-coherent particles
belong to different clusters; t2), The target particle trajectory (golden line) approaching to particles in the
red cluster; t3), The target particle separation with non-coherent particles in the red cluster and approaching
to particles in the grey cluster; t4), Separation of non-coherent particles in the grey cluster with the target
particle; t1−t5), Full trajectory view of the target particle alongside with coherent particles in golden cluster.

Lagrangian particle in the flow over a finite time (see Figure 2), we can compute the FTLE value based on
neighbour trajectories as,

σ
t
t0 =

1
|T |

√
λmax (∆) =

1
|T |

log(
δx(t)
δx(t0)

) (1)

where σt
t0 is named the scalar FTLE value showing the amount of stretching over the interval time

T = t − t0, and δx is the displacement of neighbour particles. λmax is the maximum eigenvalue of the
right Cauchy-Green deformation tensor [11] obtained from the displacement map. A lower FTLE value
means a neighbouring particle is coherent and acts in the same behaviour with the target particle spatially
over a specific temporal scale. As a result of this classification, the mean values, including velocity and
acceleration of the coherent neighbour particles, can be superimposed (weighted averaged) on the target
particle (see Figure 1). This physics-based information is added to the resulting prediction function. The
prediction function of this study is a 3rd order polynomial predictor that must satisfy particle history and
additional coherent dynamics. Details of the prediction function is addressed in Section 3.

3 Prediction Function
Polynomial function is the most simple predictor that can be used in the time resolved particle tracking
techniques. The polynomial coefficients must be determined optimally by minimising mean square error
such that the corresponding polynomial curve with order of n best fits the given positions. This can be
formalised as,

n

∑
j=n−m

n+1

∑
i=1

ai.t i−1
j = yn−m,n, (2)

where ai is unknown coefficients of the predictor function. yn−m,n is known positions of the last finite frames
(in this study m = 4). Therefore, the least square cost function is,

J =
1
n

n

∑
i=n−m

(Xi− yi)
2 , (3)



Figure 3: Lagrangian particle trajectories obtained from transport of synthetic particles through DNS Eule-
rian fields.

The coherent prediction function adds information obtained from temporal / local spatial Lagrangian
coherent particles to come up with additional constraints in the polynomial cost function (Equation 3).
In the worst-case scenario where there is no coherent neighbour information, the prediction function is
just a simple polynomial function without additional constraints. Each particle carries sets of information,
including position, first and higher-order derivative values. Assuming positions of at least three time steps
m are known. In this function, we impose first and second-order derivatives of each coherent particle into
the prediction function as below,

n+1

∑
i=2

ai,n+1.(i−1)∗ t i−2
n = ÿ LCS

n

n+1

∑
i=3

ai,n+1.(i−1)∗ (i−2)∗ t i−3
n = ÿ LCS

n .

(4)

Therefore, each particle will end up with the weighted averaged of local coherent velocity and coherent
acceleration values (ẏLCS , ÿLCS). In the present study, we take four time step histories of particles to minimise
the cost function and then predict the next step. The first and second derivatives of all coherent particles are
weighted averaged based on their FTLE level and distance to the target particle. Two weighted averaged
values, velocity and acceleration, of each target particle can be obtained for the estimation. Therefore, the
modified cost function (i.e. coherent predictor) can be written as,

J =
1
n

n

∑
i=1

(Xi− yi)
2 +(Ẋn− ẏn)

2 +(Ẍn− ÿn)
2. (5)

Velocity constraint controls the direction of the prediction function, while in the case of having high
acceleration gradients, a second-order constraint is required to control the acceleration of the prediction.
The solution for the cost function in Equation 5 is not only smooth on the history of the target particle but
also satisfies local coherent first and second-order derivatives. We compared the performance of the coherent
predictor with three other prediction functions, as listed in Table 1. DNS predictor was defined as a reference
using the Euler equation to transport particle positions by the ground truth DNS velocity. Predicted particle
positions are followed by shaking or other optimisation techniques. Therefore, all particles are either tracked
or untracked except for the inlet and outlet trajectories. In every time-step, untracked particles are like
additive noises and might gradually cause to collapse of the whole trajectory process. Due to this, untracked
particles must be fed by other complementary treatments. To this end, new information can be extracted
if any groups of tracked particles are found to be located in the neighbourhood of untracked particles with
a time step phase delay (i.e. tn+1). This phase delay means that new tracked particles at time-step tn+1 are



Method Fit parameters Cost function Prediction

a) DNS predictor - - Xn+1 = ẊDNS · tn+1

b) Polynomial predictor ∑
n
j=n−` ∑

n+1
i=1 ai, j · t i−1

j = yn−`,n J = 1
n ∑

n
i=n−m (Xi− yi)

2 yn+1 = ∑
n+1
i=1 ai,n+1 · t i−1

n+1

c) Wiener filter ∑
`
i=1 wi.un = yn J = (Xn−uT

n w)2 yn+1 = ∑
`+1
i=2 wi.un

d) Coherent predictor ∑
n
j=n−` ∑

n+1
i=1 ai, j.t i−1

j = yn−`,n J = 1
n ∑

n
i=1 (Xi− yi)

2 +(Ẋn− ẏn)
2 +(Ẍn− ÿn)

2 yn+1 = ∑
n+1
i=1 ai,n+1 · t i−1

n+1

∑
n+1
i=2 ai,n+1.(i−1)∗ t i−2

n = ẏ LCS
n

∑
n+1
i=3 ai,n+1.(i−1)∗ (i−2)∗ t i−3

n = ÿ LCS
n

Table 1: Prediction function formulation

locally coherent with one specific untracked particle at time-step tn. Another technique to reduce the number
of untracked particles is to use backward prediction, which is well established in classic schemes such as
nearest neighbour trajectory. Similarly, we implemented the backward predictor to search for additional
information from the coherent particles to estimate in reverse pace followed by backward shaking. On the
one hand, surrounding information of an untracked particle can provide the least information to predict in
backward pace. In addition, this treatment can also connect spilt tracklets for reconstructing longer particle
trajectories. This process is iterative, meaning that every forward step is embedded with a series of backward
estimations from the current time-step up to the first step. In the present study, we evaluated and compared
the performance of coherent predictor only in forward prediction.

Figure 4: Normal pdf of particle position error in x direction of four predictor functions.

4 Synthetic Evaluation
To evaluate our novel particle position prediction scheme, we used a DNS simulation of the wake behind a
smooth cylinder at Reynolds number equal to 3900 computed by an open-access code named Incompact3d
[12]. Particle trajectories around the cylinder are shown in Figure 3. Particles are transported by every 10
DNS time step using the fourth order Runge Kutta temporal and trilinear spatial schemes. The synthetic
dataset is available to the public for interested readers [13, 14]. In the synthetic case, particle trajectories are
smooth and predictable when the synthetic temporal scale is with the same order of the DNS time step due
to the small travelling distance between two time steps (less than the Kolmogorov timescale). However, the
travelling distance is comparably large for a real PTV experiment. To mimic the real experiment, we created



Figure 5: Position estimation error averaged in z direction: a), DNS predictor; b), Polynomial Predictor; c),
Wiener filter; d), Coherent predictor.

around 150,000 ground truth trajectories associated with noise for every 20 DNS time step. By increasing
the temporal scale, less particle position temporal information is available, and then the prediction would be
more challenging. Therefore, even a weak signal of coherent motion would lead to a better prediction.

We compared position prediction of four schemes with the ground truth particle trajectories (see Table
1). First, we predicted particle positions based on known ground truth DNS velocity with a linear Euler
transport function. In such a scenario, we can estimate the minimum uncertainty level that can be achieved
in this sparse temporal scale. Both Wiener filter and Polynomial predictors are also selected to be compared
with the LCS based predictor (i.e. coherent predictor). Figure 4 shows normal pdf of the predicted position
error in x direction of four schemes. Position error in x direction shows that deviations for coherent predictor
stay virtually below 0.05 ε/D, where D is the cylinder diameter. On the contrary, a significant number of
particles are mispredicted in both Polynomial and Wiener filter techniques. Similar significant improvements
of using coherent predictor are observed in y and z directions. Figure 5 shows the projected distribution of the
position error on xy plane for each predictor function. Interestingly, the prediction error is highly correlated
with the flow behaviour. Although the DNS predictor (see Figure 5.a) uses known ground truth velocity
information, the travelling distance is large enough to introduce the prediction errors, particularly inside
the wake region. As shown in Figure 5.b, third order polynomial has the worst prediction error, which can
be up to 0.2 ε/D around the cylinder leading edge and inside the wake region. The polynomial prediction
error distribution is fully shaped by the flow motion meaning that any gradients inside the flow create huge
estimation error. Overall and local performance of the Wiener filter is better than the Polynomial predictor.
Wiener filter succeeded to reduce the prediction error in most of the peak regions (see Figure 5.b.c). Error
distribution reveals that coherent predictor has the best performance locally and globally compared to Wiener
and Polynomial predictors. Worth mentioning that a small prediction error reduces the probability of picking
a wrong particle from the surroundings in the optimisation process of 4D-PTV.

5 Experimental Evaluation
An experimental study of the cylinder wake flow at Reynolds number equal to 3900 (same value as the
synthetic data) was performed. We designed an experimental setup with four cameras as shown in Figure
6.c. Four CMOS SpeedSense DANTEC cameras with a resolution of 1280×800 pixels and the maximum
frequency of 3 kHz are empowered. Cameras are equipped with Nikon 105 mm lenses. The first two cameras
are positioned in backward light scattering, while the second two cameras received maximum intensity
signal in forward scattering. The calibration error was lower than 0.06 pixel and reduced to 0.04 after the



a) b)

c)

100 200 300 400 500 600 700 800 900 1000

100

200

300

400

500

600

700

800

900

1000

2

3

4

5

5.5

V
el
o
ci
ty

(m
/
s)

Figure 6: The cylinder wake flow at Reynolds 3900 with four cameras: a), Snapshot of the 4D-PTV ex-
periment; b), Side view of particle trajectories superimposed by vorticity iso-surfaces; c), Schematic of the
experimental setup.

volume self calibration. The volume of interest was 200 mm× 150 mm× 46 mm starting from roughly
4D downstream of the cylinder, knowing that the vortex formation zone ends at 4D. The aperture was set at
11 to achieve 46 mm depth of focus. We used an LED system to illuminate this large volume. The seeding
particles were Helium Filled Soap Bubbles (HFSB) [15] resulting desired intensity signal with appropriate
particle size. However, bubbles are limited by three main factors in the wind tunnel experiments, including
generation rate, lifetime, and image glare points. We placed 50 bubble generator nozzles with airfoil-shaped
structures inside the wind tunnel chamber. The nozzles were far upstream of the measurement section to
ensure a sufficient number of bubbles are created, and the main flow field is not disturbed by the existence of
nozzles. The bubble lifetime is very short (less than 2−3 minutes) inside the wind tunnel, mainly because
they explode by passing through honeycomb layers. To overcome this issue, we injected bubbles inside the
chamber for up to 5 minutes when the wind tunnel is off before starting the acquisition. We found that
particles larger than three pixels create two glare points on two sides of the bubble. This requires more



-5 -2.5 0 2.5 5

Position error (pixels)

0

0.2

0.4

0.6

0.8

1

n
or

m
al

p
d
f

Coherent predictor

Polynomial predictor

Wiener -lter

Figure 7: Experiment normal pdf results of particle position error in x direction of three predictors. Each
predictor is compared with final optimised positions of STB Davis.

image treatments before running the 4D-PTV algorithm to avoid false particle reconstruction. However, the
intensity of two glare points can diffuse and merge if the particle size is around two pixels. Therefore, we
adjusted the camera magnification to reach two particle pixel sizes on average to surpass the glare point
issue. One snapshot of the experiment is shown in Figure 6.a. Trajectory results of the current experiment
with superimposed vorticity iso-surfaces are shown in Figure 6.b.

To quantify the results of different schemes, we compared predictions with optimised positions obtained
from STB Davis. As a result of the experiment, STB managed to successfully build nearly 12000 particles as
shown in Figure 5.b. Noisy particle reconstruction of four time steps was used as an input of the prediction
functions. We compared three techniques, Polynomial, Wiener filter, and coherent predictors, with final
optimised positions. The deviation of position estimated of each technique is shown in Figure 7. The
distribution shows that the coherent predictor has more accurate estimations within 1 pixel deviation from
the optimised positions. Position estimations of Wiener filter and coherent predictors stay below 2.5 pixels
deviation for nearly all particles. On the contrary, the Polynomial predictor has maximum deviation with
STB Davis.

6 conclusion
We proposed a robust technique to predict particle positions based on their local temporal and spatial coher-
ent motions. LCS can classify and divide the coherent neighbour motions. We imposed first and second-
order derivatives of the neighbour coherent motions into the predictor function in addition to the particle
history. To assess the proposed method named coherent predictor, we performed the synthetic analysis
of the wake behind a smooth cylinder at Reynolds number equal to 3900. We compared three predictor
functions. Polynomial predictor showed maximum deviation with the ground truth data. Whereas coherent
predictor provided the most accurate position estimation. We found that the flow regions highly impact the
estimation error. Inside the wake region, particularly the vortex formation zone and the two sideward shear
layers, cause more challenges in prediction. These mentioned regions are featured by high acceleration and
3D directional motions. We also performed the 4D-PTV experiment of the wake flow behind a cylinder at
the same Reynolds number. It was found that the coherent predictor is reliable to estimate particle positions
very close to the optimised positions.



References
[1] Ali Rahimi Khojasteh, Yin Yang, Dominique Heitz, and Sylvain Laizet. Lagrangian Coherent Track

Initialisation (LCTI). In arXiv, volume 9, 2021.

[2] Daniel Schanz, Sebastian Gesemann, and Andreas Schröder. Shake-The-Box: Lagrangian particle
tracking at high particle image densities. Experiments in Fluids, 57:1–27, 5 2016.

[3] Yin Yang and Dominique Heitz. Kernelized Lagrangian Particle Tracking. 4 2021.

[4] Daniel Schanz, Andreas Schröder, Sebastian Gesemann, Dirk Michaelis, and Bernhard Wieneke.
’Shake The Box’: A highly efficient and accurate Tomographic Particle Tracking Velocimetry (TOMO-
PTV) method using prediction of particle positions. In 10th International Symposium on Particle Image
Velocimetry, 2013.

[5] Andreas Schröder, Daniel Schanz, Dirk Michaelis, Christian Cierpka, Sven Scharnowski, and Chris-
tian J. Kähler. Advances of PIV and 4D-PTV ”shake-The-Box” for Turbulent Flow Analysis -the Flow
over Periodic Hills. Flow, Turbulence and Combustion, 95(2-3):193–209, 2015.

[6] Daniel Schanz, Andreas Schröder, and Sebastian Gesemann. ’Shake The Box’ - a 4D PTV algorithm:
Accurate and ghostless reconstruction of Lagrangian tracks in densely seeded flows. 17th International
Symposium on Applications of Laser Techniques to Fluid Mechanics, pages 7–10, 2014.

[7] Shiyong Tan, Ashwanth Salibindla, Ashik Ullah, and Mohammad Masuk. An open-source Shake-
the-Box method and its performance evaluation. In 13th International Symposium on Particle Image
Velocimetry, volume i, 2019.

[8] Andreas Schröder, Daniel Schanz, Reinhard Geisler, Sebastian Gesemann, and Christian Willert. Near-
wall turbulence characterization using 4D-PTV Shake-The-Box. In 11th International Symposium on
Particle Image Velocimetry, 2015.

[9] George Haller. Lagrangian coherent structures. Annual Review of Fluid Mechanics, 47:137–162, 2015.

[10] Ali Rahimi Khojasteh, Dominique Heitz, Yin Yang, and Sylvain Laizet. Lagrangian Coherent Track
Initialisation. In 3rd Workshop and 1st Challenge on Data Assimilation & CFD Processing for PIV
and Lagrangian Particle Tracking, 2020.

[11] Shawn C. Shadden, Francois Lekien, and Jerrold E. Marsden. Definition and properties of Lagrangian
coherent structures from finite-time Lyapunov exponents in two-dimensional aperiodic flows. Physica
D: Nonlinear Phenomena, 212(3-4):271–304, 12 2005.

[12] Sylvain Laizet and Ning Li. Incompact3d: A powerful tool to tackle turbulence problems with up to
O(105) computational cores. International Journal for Numerical Methods in Fluids, 67:1735–1757,
12 2011.

[13] Ali Rahimi Khojasteh, Sylvain Laizet, Dominique Heitz, and Yin Yang. Lagrangian and Eulerian
dataset of flow over a circular cylinder at Reynolds number 3900, 2021.

[14] Ali Rahimi Khojasteh, Sylvain Laizet, Dominique Heitz, and Yin Yang. Lagrangian and Eulerian
dataset of the wake over a smooth cylinder at a Reynolds number equal to 3900. Data in Brief, 2021.

[15] Fulvio Scarano, Sina Ghaemi, Giuseppe Carlo Alp Caridi, Johannes Bosbach, Uwe Dierksheide, and
Andrea Sciacchitano. On the use of helium-filled soap bubbles for large-scale tomographic PIV in
wind tunnel experiments. Experiments in Fluids, 56:1–12, 2 2015.


	Introduction
	Coherent and non-coherent neighbours
	Prediction Function
	Synthetic Evaluation
	Experimental Evaluation
	conclusion