Dense velocity reconstruction with VIC-based time-segment assimilation

The vortex-in-cell time-segment assimilation (VIC-TSA) method is introduced. A particle track is obtained from a finite number of successive time samples of the tracer’s position and velocity can be used for reconstruction on a Cartesian grid. Similar to the VIC + technique, the method makes use of the vortex-in-cell paradigm to produce estimates of the flow state at locations and times other than the measured ones. The working principle requires time-resolved measurements of the particles’ velocity during a finite time interval. The work investigates the effects of the assimilated length on the spatial resolution of the velocity field reconstruction. The working hypotheses of the VIC-TSA method are presented here along with the numerical algorithm for its application to particle tracks datasets. The novel parameter governing the reconstruction is the length of the time-segment chosen for the data assimilation. Three regimes of operation are identified, based on the track length and the geometrical distance between neighbouring tracks. The regime of adjacent tracks arguably provides the optimal trade-off between spatial resolution and computational effort. The VIC-TSA spatial resolution is evaluated first by a numerical exercise; a 3D sine wave lattice is reconstructed at different values of the particles concentration. The modulation appears to reduce (cut-off delay) when the time-segment length is increased. Large-scale PIV experiments in the wake of a circular cylinder at Red = 27,000 are used to evaluate the method’s suitability to real data, including noise and data outliers. Both primary vortex structures in the Kármán wake as well as interconnecting ribs are present in this complex flow field, with a typical diameter close to the average inter-particle distance. When the time-segment is increased to adjacent tracks and beyond, a more regular time dependence of local and Lagrangian properties is observed, confirming the suitability of the time-segment assimilation for accurate reconstruction of sparse velocity data.


Introduction
Three-dimensional velocity measurements based on PIV have been transitioning from domain cross-correlation (Tomographic PIV, Elsinga et al. 2006) towards more efficient particle-based analysis techniques, notably the Shake-the-Box algorithm (Schanz et al. 2016a, b). The main advantages offered by particle tracking velocimetry (PTV) techniques (Nishino et al. 1989;Malik et al. 1993;Maas et al. 1993) are that the velocity is defined pointwise, opposed to the spatial averaging over the interrogation window (cross-correlation analysis). A second important advantage of particle-based analysis is the significant level of data compression, where only the particle positions and their velocity need to be stored, opposed to reconstructed volumes often counting billions of discrete volumetric elements (voxels). One of the widely reported shortcomings of the PTV analysis has been the low reliability in particle image pairing. The latter problem grows exponentially with the tracers' concentration (Schanz et al. 2016a, b). The introduction of iterative algorithms for 3D particle field reconstruction (IPR, Wieneke 2012) has paved the way to a more robust particle tracking analysis (Schanz et al. 2016a, b). As a result, the tractable seeding density for particle-based methods has attained those practised for tomographic PIV (Raffel et al. 2018).
One of the main problems arising from the direct treatment of spatially sparse data is the difficulty in performing a number of data reduction operations, ranging from simple visualisations, statistical analysis, spatial differentiation and evaluating more complex operators such as for pressure evaluation (van Oudheusden 2013). Reducing first the data onto a Cartesian grid (Cartesian grid reduction, or CGR in the remainder) is currently the norm. Several works have aimed at optimising such process in terms of accuracy, reliability, resilience to noise and computational cost. Early works on CGR have proposed interpolation techniques, mostly using a polynomial basis for either interpolation or multi-dimensional regression. Specific basis functions for PTV data have also been proposed, notably the adaptive Gaussian windowing technique (AGW, Agüí and Jimenez 1987). The main shortcoming of interpolation techniques is the mismatch between the solutions admitted by the basis functions (especially for low-order polynomials) and those representing the physical flow behaviour. For instance, a linear interpolation of scattered data produces an unphysical discontinuous distribution of the velocity gradient and the related flow properties (vorticity, rate of shear). Methods that enforce the solution to obey physical laws have been initially proposed for PIV data (already on a Cartesian grid) as a mean to reduce the measurement errors. Solenoidal filters that impose mass conservation (Vedula andAdrian 2005, Schiavazzi et al. 2014) fall in the above category.
More complex methods based on data assimilation have been introduced for CGR seeking to enhance the spatial resolution. By enforcing compliance of the instantaneous flow measurement (velocity and its material derivative) with the governing equations of fluid motion (conservation of mass and momentum), data assimilation techniques have shown to enable resolving the flow at scales smaller than those attained by interpolation or cross-correlation analysis. The VIC + technique (Schneiders and Scarano 2016) is based on the assimilation of the data to the dynamical equation of vorticity (Vortex-in-Cell, Christiansen 1973) and it has shown to provide spatial resolution significantly higher than that of cross-correlation analysis (typically used for tomographic PIV data). Physical processes, such as vorticity fluctuations dynamics and turbulent kinetic energy dissipation, taking place at the right-side of the spectrum of flow scales were captured to a fuller extent by the VIC + data assimilation, compared to the cross-correlation analysis . Further studies have confirmed the benefits of imposing mass and momentum conservation for CGR applied to particle trajectories measured with the STB technique (FlowFit, Schanz et al. 2016a, b). Furthermore, the same group (Ehlers et al. 2020) has investigated ways to improve the temporal consistency of FlowFit CGR by introduction of artificial Lagrangian tracers (virtual particles).
In the work of Jeon et al. (2018), significant advancements to the VIC + technique have been reported, both in terms of computational efficiency (multi-grid computations) and for the treatment of domain boundaries.
In summary, CGR techniques can be qualified by their accuracy and spatial resolution. Figure 1 compares such methods, stating the relatively poor resolution and accuracy of spatial averaging methods, compared to interpolators. Further advantages to the latter are obtained with the use of physical constraints as done by data assimilation techniques.
Most of the proposed approaches perform the assimilation on the basis of the instantaneous flow properties (velocity, vorticity, acceleration). The potential of 4D data assimilation has been already explored in the framework of variational methods in computational fluid dynamics problems (Chandramouli et al. 2020, among others).
The application to time-resolved sparse 3D data is recently being considered. The assimilation of measurements comprising a finite time interval (time-segment assimilation or TSA) has been introduced recently by Schneiders and Scarano (2016) in the framework of vortex-in-cell technique. In the latter work, the principle was introduced along with the analysis of synthetic data indicating the potential for delaying the cut-off wavenumber for spatial reconstruction.
A simpler, but computationally more efficient 4D (i.e. space-time), approach to the problem has been developed by Jeon et al. (2019) where the time-segment is spanned using Taylor approximation of the flow field velocity. The latter can be considered as an advection-based time-segment assimilation (ADV-TSA) that models the flow behaviour along a finite time interval. The study has brought evidence that an improved dense reconstruction can be achieved using different values of the time-segment. Yet, no conclusive Fig. 1 Summary of CGR methods and qualitative indication of their performance in terms of measurement accuracy and spatial resolution. The height of the bars is purely indicative criterion has been stated that defines an optimum for the time-segment length.
The above 4D techniques can be categorised as pouring time into space viz. enhancing the spatial resolution of the reconstruction from spatially sparse and temporally well resolved data. The latter conditions are encountered in a number of 3D flow measurements, in particular after the recent development of large-scale 3D PIV experiments by means of helium-filled soap bubbles (Bosbach et al. 2009) where the flow velocity is sampled at the sparse tracers' locations.
The present work discusses the VIC-TSA data assimilation principle. A first-principle analysis based on geometry and kinematics lays the foundations of the operation regimes for this technique. Furthermore, the mathematical background and the numerical algorithm are presented and discussed. The features and performances of the VIC-TSA technique are based on the analysis of a simplified synthetic PTV experiment, using a steady 3D lattice of sinusoids, following the approach of the 4th PIV challenge (Kahler et al. 2016). The main parameter governing the behaviour of this method is the time interval (time-segment) considered for the data assimilation. In the present work, the effect of the time-segment length on the spatial resolution is also examined from the analysis of experimental data obtained with large-scale PIV measurements in the wake of a circular cylinder at Re D = 27,000.

Working principle of VIC-TSA
The working principle of VIC-TSA may be seen as based upon two concepts: 1. Increasing the amount of data used for CGR, by including measured velocity data over a finite time interval (particles tracks) 2. Requiring a solution that satisfies the governing equations for vorticity over a finite time interval The first concept is elaborated in the section below. The second one is the subject of Sect. 3.

Effective inter-particle distance for multi-frame data
Particle tracers immersed in the fluid flow are assumed to be homogeneously distributed in space. Their concentration C is defined as the number of particles present in a unit volume. The concentration of particle tracers directly determines the length-scales of the flow that can be instantaneously resolved.
Consider the flow schematically depicted in Fig. 2: three tracer particles [ p 1 , p 2 , p 3 ] travel along their respective trajectories [ Γ 1 , Γ 2 , Γ 3 ] (solid blue lines). Their position at a given time instant t 0 is indicated by a dark red circle. Provided that the measurement of the particles position is obtained in time-resolved mode (e.g. by means of high-repetition rate PIV), the sampled positions at Nt + 1 consecutive times [ t −Nt∕2 , … , t −1 , t 0 , t 1 , … , t Nt∕2 ] separated by the sampling time interval ∆t s are shown by hollow circles. We denote the time-segment length encompassed by the Nt + 1 samples as T = Nt • Δt s . The array of Nt + 1 positions for each particle is denoted as a track ℑ (represented by solid black lines in Fig. 2-right). Furthermore, let us denote the average distance among neighbouring particles by r and the minimum distance among neighbouring tracks by r . The condition r ≤ r is always satisfied. Moreover, the ratio r∕r decreases when the track length L (viz time-segment length T ) increases.
For the CGR of the flow velocity from sparse tracers, the spacing of the Cartesian grid depends fundamentally upon the mutual distance among neighbouring tracer particles, in turn depending upon the seeding concentration. A second parameter is the scheme used to interpolate or reconstruct the velocity field between the scattered data points. The VIC + data assimilation method (Schneiders and Scarano 2016) makes use of the instantaneous velocity and acceleration for the assimilation of the instantaneous volumetric Other samples along the time-segment are marked by hollow circles. Right: particle tracks ℑ (solid black lines overlay-ing the trajectories) and distance between neighbouring particles r at time instance t 0 (dashed red) and minimum distance between tracks r (dashed black) velocity field, to retrieve smaller flow scales than the interparticle distance. In the above work, the criterion h = ¼ r is proposed, where h is the grid spacing. In the present work, the attention focuses on the distance between particle tracers during a finite time interval of length T , that translates into the distance between tracks ∼ r . The latter depends not only upon the inter-particle distance, but also on the length of the path L travelled during T .
The track's distance has the following property with time: if a track length is considered such that L ≪r (defined here as impulsive regime), then r ≈ r and little-to-no improvement can be expected in the CGR process compared with instantaneous reconstruction methods. If L ≈r (adjacent tracks regime), then r < r , and an improvement is expected in terms of spatial resolution that should be proportional tor∕r . It is conjectured here that the equivalent increase in volumetric concentration may scale as r∕r 3 . Finally, any further increase in L beyond r (stringy tracks regime) is expected to yield diminishing benefits since no substantial decrease of r will occur. The working hypothesis of the time-segment assimilation (TSA) stems from the above considerations. The dependency of spatial resolution is shifted from the particles' distance to the tracks' distance. As a consequence, the duration of the assimilated time-segment T becomes the parameter that governs such shift. The three resulting regimes are illustrated in Fig. 3.
With a short time-segment (for instance given by only three time samples, as in Fig. 3-left), the average tracks' distance differs very little from the particles' distance and the spatial resolution of CGR is expected to remain similar to that of methods based on the instantaneous flow properties. On the other side for extremely long tracks as in the stringy tracks regime , the distance between neighbouring tracks has already reached the minimum limit of distance between trajectories. The spatial resolution of CGR is expected to approach that obtained with an experiment where the tracers concentration corresponds to particles' distance of r . The intermediate condition of adjacent tracks regime is identified as a potential optimum considering the computational efficiency of the VIC calculations, where the track length is increased until most benefits of reducing the tracks distance are included. The criterion for selecting the time-segment length that corresponds to the adjacent tracks regime reads is as follows: where U is a reference velocity of the flow, T is the length of the time-segment, and r is the average inter-particle distance between the particles in the domain. The following sections aim to prove the stated rationale: the accuracy of the flow reconstruction with TSA depends upon the above introduced regimes.
The above discussion established a criterion for selecting the optimal time-segment length for reducing the effective inter-particle distance r . Increasing the time-segment length further does not further reduce r , but does require the solution resulting from VIC-TSA to satisfy the flow governing equations over a longer time duration.
Here, the optimum time duration for VIC-TSA is associated with the flow timescales. When taking an integration time similar or larger than the expected flow timescale, the solution is consequently required to behave coherently over this time duration.

VIC-TSA algorithm
The aim of the VIC-TSA technique is to achieve a velocity field at the middle of the time-segment that minimises the difference between the measured particles velocity and the VIC computation along the entire time-segment. The current implementation was introduced by Schneiders and Scarano (2018) and is recalled here for completeness. The input to the algorithm is based on time-resolved 3D tracks data, as issued for instance from tomographic PTV (Novara and Scarano 2013) or Shake-the-Box analysis (aka Lagrangian particle tracking, Schanz et al. 2016a, b). The set of particle tracks is referred to as ℑ and it entails spatially scattered measurements of particles position x p and velocity V p sampled Nt times along their where ω = ∇ × U is the vorticity. The initial vorticity field, ω(x,t 0 ), is numerically integrated in time using a vortexparticle discretization (Vortex-in-Cell, Christiansen 1973).
Within the vortex-in-cell (VIC) framework, velocity is calculated from vorticity at each integration time instant by solution of a Poisson equation, In the solution of the above equation, it is assumed that boundary conditions can be selected that make the problem well posed. In case no boundary conditions can be defined a-priori, they can be considered as additional unknowns and increase the degrees of freedom in the optimization procedure. The above assumptions have been discussed in Schneiders and Scarano (2016) and Jeon et al. (2019).
Solving (2) and (3) for the time-segment T yields the velocity temporal evolution in this measurement time-segment. At every measurement time instant t n , the computed velocity field U h can be evaluated at the instantaneous particle location x p to yield the instantaneous disparity J p relative to a single particle: When the operation is performed for the entire particle ensemble and across the time-segment T, one obtains the total disparity: that is referred in the remainder as the cost function J, which indicates the quality of the numerical simulation, considering the measurements as a reference. Therefore, J is at the basis of the optimization algorithm.
In the present work, the boundary conditions are kept fixed, based on an a-priori estimation and only the initial vorticity field is taken as degree of freedom. The optimization problem is defined accordingly as follows: The sampling frequency adopted in experiments is usually high enough to temporally resolve the relevant timescales of the flow problem. An illustration of the numerical process is given in Fig. 4: the process begins by initialising the velocity field at the time t 0 . This is achieved by interpolating the tracks data V p onto a Cartesian grid, yielding U h . The vorticity field ω h is thus estimated by finite-difference at time t 0 . The vorticity field along the time-segment is computed using the Vortex-in-Cell framework, similarly to the time super-sampling technique (Schneiders et al. 2014). In Fig. 4, the sparse and gridded variables are typed in black and red, respectively.
An iterative problem is built upon these steps to minimize the cost function J (grey line arrow in Fig. 4) by tuning the optimization variable, the gridded vorticity, in every iteration. The problem is solved using the limited-memory Broyden-Fletcher-Goldfarb-Shanno method (L-BFGS, Liu and Nocedal 1989), an approximation of the BFGS, based on a quasi-Newton method, with a limited amount of memory.
The resulting method remains computationally intensive: VIC-TSA requires approximately 500 and 20 times more CPU time compared to linear interpolation and the The cost function J is evaluated from the differences between the computed velocity and the tracks data. The cost function is minimised iteratively (grey line arrow) yielding the optimized velocity field at t 0 VIC + algorithm, respectively. The computational effort increases linearly with the length of the time-segment and with the third power of grid resolution.

Numerical assessment
A synthetic flow field is generated for the numerical assessment of the VIC-TSA technique. The field is described by the steady 3D sine wave lattice shown below where x, y and z are the Cartesian coordinates, A = 2 m/s is the velocity amplitude of the waves along x and y and = 0.2 m is their wavelength. The velocity is considered uniform along z.
Sine-based velocity distributions have been extensively used in previous works for steady numerical assessments in 2D (Willert and Gharib 1991;Scarano and Riethmuller. 2000) and in 3D (4th PIVchallenge, Kahler et al. 2016;Schneiders and Scarano 2018). The analysis is based on the sinusoid reconstruction, with attention to the modulation of the peak value, as the resolution of the measurement technique deteriorates. The current approach considers a simulated volume of size 0.4 × 0.4 × 0.8 m 3 , where randomly distributed particle tracers at a uniform concentration are generated and propagated according to the underlying velocity field. No particle reconstruction errors are simulated. The particle motion is calculated using a fourth-order Runge-Kutta method. The concentration is maintained constant by introducing at the inlet side as many particles as those exiting.
A spatial response analysis is carried out for the mentioned velocity field varying the particles concentration C, and the number of samples Nt composing the assimilated time-segment. The VIC-TSA method is applied with the shortest duration of observation, three snapshots (central scheme in time), a longer segment, Nt = 11, and a significantly longer segment, Nt = 31. The non-dimensional time length of these segments gives an estimation of the regime of the cases explained in Sect. 2. For a time separation between samples Δt S = 0.001 s, a convection velocity of the particles A = 2 m/s and an average inter- respectively. For the current analysis, the grid spacing used for the TSA reconstruction is kept constant for the different cases,h = 0.01m. A first qualitative analysis of the spatial response is given in Fig. 5, where the u velocity component is displayed. The analysis compares the results when the average inter-particle distance is varied ( r = 0.08, 0.05, 0.03 ). The upper row illustrates the tracers captured in a slice of the domain of 0.2 m height. The reference velocity field is shown at the right of the figure and it presents the sine waves pattern of Eq. 7. A general trend is observed, whereby more peaks are retrieved and with higher values as the assimilation is performed over a longer time. It is remarked here that the technique reconstructs a coherent velocity pattern a low densely seeded flows by using long assimilation times ( r = 0.08 for 31 samples). Such feature is particularly interesting for experiments conducted at low seeding density, as it is the case when helium-filled soap bubbles (HFSB) are used as flow tracers in wind tunnels. Figure 6-left shows the amplitude modulation evaluated with three values of the time-segment duration. The spatial response of a top-hat moving average filter (solid black line) over a cubic element containing 5 tracer particles is also shown for reference. The latter can be representative of the cross-correlation analysis as used in tomographic PIV. The horizontal axis represents the normalized wavenumber r*, computed as the ratio between the average inter-particle distance and the flow wavelength r * = r∕ (lower axis), or between the linear size of the cubic element used as linear filter and the flow wavelength (l* on the upper axis).
The effect of spatial modulation by decreasing the seeding density is evident for all curves. The most important observation is that the − 3 dB roll-off point (or critical inter-particle distance r * c ) is progressively delayed when a longer sequence is considered. For a sequence with only three frames, the peak value is halved at r * c = 0.2. In comparison, for a sequence with 31 frames, one obtains r * c = 0.35. Such shift corresponds to a fivefold increase of particle concentration.
The present results should, however, be taken with a word of caution, considering that for a steady velocity field, the temporal integration of information yields additional spatial information. Further work should be devoted to defining appropriate analytical test cases that also include the unsteady development of the flow. The experimental section in the present work deals with such unsteady conditions, in absence, however of a reference velocity field as in the present case.
Given the affinity of the TSA-VIC technique with the VIC + method, results from the study of Schneiders et al. (2017, appendix) are reported here for comparison. To some surprise, VIC + data feature a flatter behaviour compared to TSA-VIC for sequences of 3 to 11 frames. This effect is ascribed to a more accurate optimisation by the use of radial basis functions in the VIC + algorithm. Radial basis functions have not been implemented in the present work for computational affordability and the TSA-VIC algorithm is currently based on the discrete finite-differences method. The TSA-VIC algorithm requires the numerical evaluation of Eq. 2, whereby time marching is required. To this purpose, the effect of the internal time integration on the reconstructed velocity field is evaluated for the chosen values of the time-segment at a seeding particles concentration of C = 3000 particles/m 3 and for a normalized inter-particle distance r * = 0.2. Figure 6-right shows the amplitude modulation of the reconstructions for the normalized time marching frequency, dt * = Δt march ∕Δt s . The CFL number computed from the reference velocity A = 2 m/s and grid spacing h = 0.01 m is given as top x-axis (Ferziger et al, 2020). In the inquired range (CFL < 0.5), the result appears to be independent of the choice on the internal time separation.
Finally, the spatial response to the choice of the computational grid is examined. Figure 7 compares the results obtained with a fine grid ( h = 0.01 m) with those returned using a too coarse grid (h = 0.05 m). The latter results in the amplification of numerical errors, especially when a long-time interval is considered. Such effects can also be observed in Fig. 6-left for small values of r*. The grid spacing considered in Fig. 7, however, well beyond the ¼ criterion suggested in Schneiders and Scarano (2016). Nevertheless, the present result needs to be retained in mind as caveat, when exploring ways to assimilate data with VIC-TSA over a coarse grid.

Cylinder wake experiment
The experimental case considered for VIC-TSA analysis is the Karman wake behind a circular cylinder. This easily reproducible experiment is a common benchmark for measurement techniques dealing with unsteady flows. Yet, the large range of length-and timescales, from the Karman vortices to the shear layer turbulence, poses a challenge to the spatial resolution of field measurement techniques.

Experimental apparatus and procedures
The experiments are conducted in a large low-speed wind tunnel, the open jet facility (OJF) at the laboratories of the aerospace engineering faculty of TU Delft. The OJF features an exit cross section of 2.85 × 2.85 m 2 . A 2-m-long cylinder diameter D = 10 cm is installed vertically and immersed in the free-stream flow at 4 m/s. The resulting value of the Reynolds number is 27,000. The domain of interest is the near wake of the cylinder. The measurement covers a domain of 5D(X) × 6D(Y) × 3. 8D(Z). For the measurements over such an extended region, sub-millimetre helium-filled soap bubbles are used as flow tracers . A HFSB seeder hosting 204 generators with cross-sectional area of 0.6 × 0.9 m 2 is placed in the wind tunnel settling chamber. After the wind tunnel contraction, a seeded streamtube of approximately 50 cm diameter is produced with particles concentration of about 1 bubble/cm 3 . Figure 8(top-left) illustrates the HFSB generator used for the experiment.
Illumination is provided by two LaVision LED-Flashlight 300 (72 elements, white light) directed along the span of the cylinder. The pulsed width is 50 ms (one tenth of the pulse separation time Dt s = 0.5 ms, at a repetition rate of 2 kHz). A high-speed tomographic PIV system is employed to record the light scattered by the HFSB tracers with three LaVision HighSpeedStar 6 CMOS cameras (1024 × 1024 px, 5400 frames/s, 12-bit, 20 μm pixel pitch) equipped with Nikon 50 mm focal length objectives, set at numerical aperture f # = 22. The cameras subtend a tomographic angle β = 65 degrees. The main components of the PIV system and measurement parameters are summarised in Table 1.
Three sequences of 2000 frames are acquired at a rate of 2 kHz. The particle's motion between subsequent exposures is approximately ∆X = 2 mm in the free stream. Considering a shedding frequency based on a Strouhal number of 0.2, each sequence captures approximately eight shedding cycles. Synchronisation of illumination and image recording is performed with a LaVision programmable timing unit (PTU 10). Data recording, storage, and Shake-the-Box analysis operations are performed with the DaVis 10 software.

Data processing and TSA analysis
The images are pre-processed with a time-frequency domain high-pass filter (Sciacchitano and Scarano 2014) to reduce background illumination. Raw and pre-processed recordings are illustrated in Fig. 9 and further detailed in the videos. An image-to-object geometrical mapping function is obtained recording the image of a calibration plate with controlled translation across the measurement domain. The residual errors are reduced to below 0.1 pixels using the volume self-calibration algorithm (Wieneke 2008). Subsequently, Fig. 7 Amplitude response convergence history as a function of the reconstruction grid an optical transfer function calibration (Schanz et al. 2013) is performed prior to the interrogation of the recordings sequence.
The tracers motion is evaluated with the Shake-the-Box algorithm (Schanz et al. 2016a, b). Approximately 40,000 tracers are detected in every recording. The use of three cameras and the tracks acceptance criterion yields on average 12,000 measurements at every time instant. The resulting tracks concentration is C tracks = 0.12/cm 3 and the corresponding average distance among neighbouring tracers on the detected tracks is 1.25 cm (discussed as r in Sect. 2). Table 2 summarises the parameters used for the particles motion analysis. An illustration of particle tracks with colorcoded instantaneous velocity is given in Fig. 10 (left). The instantaneous vorticity field reconstructed with VIC# (Jeon et al. 2019) is shown in Fig. 10 (right). Recalling the above discussion, the condition of adjacent tracks ( T AT ) is therefore attained when the time-segment length for assimilation is T ~ 6 ms (13 exposures). The VIC-TSA analysis is performed The internal time marching for the VIC algorithm is chosen as Δt march = 0.5 ms, which complies with the CFL conditions (CFL = 0.25) and corresponds to the time separation between recorded samples. The relevant parameter of the VIC-TSA analysis is reported in Table 3.

Vorticity analysis
The ability to describe the range of flow scales can be examined from the measurement of the instantaneous vorticity distribution. According to the discussion presented in Fig. 6right, the limit condition where spatial scales are resolved is when u* is just above 0.5. First, a planar visualisation of instantaneous velocity and primary vorticity (along the cylinder axis) illustrates the behaviour of the reconstruction when τ * varies from the impulsive to the stringy regime. In Fig. 11 (top-left), the time-segment of 1 ms corresponds to τ * = 0.15 (impulsive regime). This condition approximates that of the velocity-acceleration assimilation method (VIC + Schneiders et al. 2016). Although the overall oscillatory motion is captured with the distribution of the velocity vectors, the blobs of fluctuating vorticity are comparatively less intense, when compared to the analysis performed with a longer time-segment (e.g. bottom-right corresponds to τ * = 2.8).
A quantitative analysis of the effect of τ * on the spatial response of the reconstruction method is presented in Fig. 12. The standard deviation of the axial vorticity (ω y ) dominates over the other two components. This is consistent with the alignment of the mean shear with the cylinder axis leading to the formation of spanwise vortices. The values of streamwise vorticity τ x exceed those of ω z given the formation of secondary vortices along the streamwise direction (ribs or fingers, Williamson 1996). The analysis performed with VIC-TSA by varying τ * indicates that the amplitude of vorticity fluctuations increases monotonically with the length of the time-segment (Fig. 12-left). However, diminishing benefits are found as τ * > 2. It may be concluded that the condition τ * = 2 corresponds to a practical optimum (also considering the computational effort required for the VIC time marching and the numerical optimisation process) for the choice of the time-segment in the VIC-TSA technique. Moreover, the comparison between the impulsive (τ * < 0.1) and stringy (τ * > 2) regimes yields an overall increase in approximately 30% for the vorticity fluctuations in the present experiment. Compared to the velocity reconstruction based on linear interpolation (shown in blue in Fig. 12), the VIC-TSA reconstructs vorticity from 25 to 50% higher than the linear method. An additional comparison is made with tomographic PIV. In this case, the recordings are reconstructed with the sequential MART method (SMTE, Lynch and Scarano 2015) and the crosscorrelation analysis is performed with interrogation boxes of 52 × 52 × 52 voxels with 75% overlap, yielding the same grid spacing of 8 mm. The amplitude of the vorticity fluctuations is higher than that obtained with the linear interpolation, which apparently contradicts the analysis of Schneiders and Scarano  (2018), where for a given set of particles in 3D space, linear interpolation yields a higher spatial resolution than the crosscorrelation analysis. However, in the present case, the linear interpolation is performed using only the velocity vectors of the accepted tracks, whereas the cross-correlation interrogation may also include particles not reconstructed with IPR. Overall, the level of vorticity fluctuations is similar to that of VIC-TSA with very short time-segment.
A question that may arise is whether the increase in vorticity fluctuations amplitude is not to be ascribed to an overall increase in the reconstruction noise by larger values of τ *. An indication can be found in Fig. 12-right, where the time derivative of the vorticity is shown. Here, the time derivative of vorticity appears little affected by the variation in τ *, indicating that the increase in vorticity fluctuations amplitude with the longer time-segment shall not be ascribed to noise. Results from tomographic PIV exhibit comparatively larger values. The measurement noise associated with twopulse systems can be lowered with multi-frame methods as discussed in Raffel et al. (2018).
The three-dimensional instantaneous flow organization is presented by iso-surfaces of the vorticity components in Fig. 13. The illustration shows the unsteady flow organisation dominated by spanwise rollers forming the Kármán street, producing vorticity fluctuations of highest amplitude. These vortices are approximately aligned with the cylinder axis and are shed with alternating sign of vorticity (red for counterclockwise and blue for clockwise). Additional to the Kármán vortices, secondary elongated structures are formed that interconnect the primary rollers from the outer side of the wake where low-speed fluid is ejected under the effect induced by the primary vortices. These structures develop along the streamwise (X) and the vertical (Z) direction and are produced in the region of flow subject to positive strain between two Kármán rollers. Such organisation is most often observed at flow regimes with lower Reynolds number (Re D < 1000, Prasad and Williamson 1997;Kanaris et al. 2011, among others) where these structures are reported as fingers or ribs. They become less coherent in the high-Reynolds regimes, yet the mechanism of ribs interconnecting the Kármán rollers is reported in works conducted at Re D = 3900 (Parnaudeau et al. 2008). In Fig. 13, these structures are identified by considering the vorticity components along X and Z and they are colour-coded using the sign of ω x (yellow for counterclockwise and green for clockwise). In the present context, Fig. 13 is used to elucidate the overall effect of the time-segment length used for the VIC-TSA data reduction. The first column shows the reconstruction obtained with T = 1 ms (τ * = 0.3), where the main rollers are identified along with the secondary structures. The overall vorticity distribution is visibly scattered and the secondary structures do appear, but mostly interrupted. When entering the adjacent tracks regime (τ * = 2) in the middle column of Fig. 13, the vorticity distribution appears more coherent with a clear distinction between the primary rollers and the ribs. Moving further to T = 19 ms (τ * = 6), one can observe an overall increase in vorticity coherence. Although, the analysis shown in Fig. 12 does not yield an important increase in vorticity fluctuations beyond τ * = 2, the flow reconstruction shows well-connected ribs and vorticity with peak fluctuations at a scale of 1-2 cm.

Conclusions
A data assimilation technique has been introduced that considers data distributed over a finite time-segment to reconstruct the optimal 3D vorticity and velocity fields correspond to values obtained with linear interpolation on the same data mesh. Red-black labels correspond to results obtained with tomographic PIV. Impulsive, adjacent and stringy regimes indicated by red, green and blue areas, respectively. Right: temporal derivative of vorticity from VIC-TSA and by tomographic PIV (red-black labels) described within such segment. The technique is coupled to the Vortex-in-Cell paradigm to produce a dense reconstruction of the velocity field based on spatially sparse measurements of tracers' trajectories. A first-principle analysis indicates the potential of VIC-TSA in terms of spatial resolution, by identifying three main data reduction regimes, namely impulsive, adjacent tracks and stringy.
A numerical analysis that follows the approach presented in the literature is carried in a 3D lattice of sinusoids, where the velocity amplitude modulation is inspected for varying conditions of tracers concentration, sinusoids wavelength and, most importantly, time-segment length. The cut-off of amplitude reconstruction is delayed when the time-segment increases, consistently with the hypotheses formulated in the first-principle analysis.
Experiments are conducted with large-scale PIV, where the wake of a circular cylinder at Re D = 27,000 is measured using helium-filled soap bubbles and a tomographic imaging system comprising three high-speed cameras. The tracers motion analysis is delivered using the Shake-the-Box algorithm, and the database of particle trajectories is examined using VIC-TSA. The statistical analysis of vorticity fluctuations reveals that the modulation of vortex structures reconstructed by VIC-TSA decreases as the time-segment is made longer. An optimum trade-off between spatial resolution and computational burden is identified with the criterion τ * = 2, i.e. twice the value of adjacent tracks regime. The topology analysis of the vortex structure confirms the salient features of the Kármán wake with sub-structures being resolved at a scale comparable to the average inter-particle distance.