Multi-Δt approach for peak-locking error correction and uncertainty quantification in PIV

A novel approach is devised for the quantification of systematic uncertainty due to peak locking in particle image velocimetry (PIV), which also leads to correction of the peak-locking errors. The approach, applicable to statistical flow properties such as time-averaged velocity and Reynolds stresses, relies on image recordings with multiple time separations Δt and a least-squares regression of the measured quantities. In presence of peak locking, the measured particle image displacement is a non-linear function of Δt due to the presence of measurement errors which vary non-linearly with the sub-pixel particle image displacement. Additionally, the measured displacement fluctuations are a combination of the actual flow fluctuations and the measurement error. When the image recordings are acquired with multiple Δt’s, a least-squares regression among the statistical results yields a correction where systematic errors due to peak locking are significantly diminished. The methodology is assessed for planar PIV measurements of the flow over a NACA0012 airfoil at 10 degrees angle of attack. Reference measurements with much larger Δt than the Δt’s of the actual measurements, such that relative peak-locking errors are negligible for the former, are used to assess the validity of the proposed approach.


Introduction
Peak-locking (also referred to as pixel-locking) is recognized as one of the major error sources in particle image velocimetry (PIV) measurements. Such an error source, mainly ascribed to particle image diameters small with respect to the sensor's pixel size, causes a bias of the measured particle image displacement towards the closest integer value (Westerweel Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI. 1997, Adrian and Westerweel 2011, Michaelis et al 2016, Raffel et al 2018. It is particularly relevant in high-speed PIV measurements with CMOS cameras, whose large pixel size (of the order of 10-20 µm) yields particle image diameters often smaller than one pixel. It is to be noted that the measurements with particle image diameters greater than one pixel are also affected by the peak-locking errors, although the peak-locking errors are less prominent than the random errors in this case (Westerweel 1997). Thus, the recent development of high-speed cameras has been influencing the study on peak-locking errors and their correction in PIV. Moreover, the peak-locking errors significantly affect the turbulence statistics extracted from the PIV measurements. When the flow velocity fluctuations encompass at least one pixel unit, the mean velocity is usually unaffected by the peak-locking errors in the instantaneous velocities. In such cases, the peak-locking errors appear instead in the fluctuations of velocity leading to inaccurate estimation of higher-order turbulence statistics, e.g. Reynolds stresses (Christensen 2004). Hence, the quantification and correction of the peak-locking errors in the PIV measurements is crucial for the evaluation of accurate turbulence statistics.
Several approaches have been proposed for quantification and correction of the peak-locking errors. A priori estimations and correction of the peak-locking errors based on theoretical models were presented by Angele and Muhammad-Klingmann (2005) and Cholemari (2007). Both the models have showed the effect of peak locking on the turbulence statistics following the work of Christensen (2004). The peak-locking errors have been corrected by assuming Gaussian distribution for the displacement and velocity probability density functions in the former model (Angele and Muhammad-Klingmann 2005), whereas the correction has been achieved by assuming sinusoidal variation of the peaklocking error with respect to the sub-pixel displacement in the latter model (Cholemari 2007). A number of works on the peak-locking correction at the processing (i.e. velocity estimation) and post-processing stages have been devised in the literature. Roth and Katz (2001) have applied a equalization to the sub-pixel particle image displacements for the correction. However, Hearst and Ganapathisubramani (2015) have demonstrated that pixel locking is non-uniform across an image. Therefore, identifying and adjusting for pixel locking with histograms computed based on the entire vector fields, as done by Roth and Katz (2001), may have been erroneous and the equalization process should be applied on a vectorby-vector basis (Hearst and Ganapathisubramani 2015). However, their approach is effective only in the absence of other error sources which might affect the histogram of the measured displacements and velocities. The peak-locking errors arising at the processing stage can be reduced by using state-ofthe-art processing algorithms (Scarano 2002, Roesgen 2003, Chen and Katz 2005, Liao and Cowen 2005. However, the peak-locking errors due to small particle image diameters are difficult to quantify or correct for. The concept of defocusing to increase the particle image diameter has become a common practice: a slight defocusing has been shown effective in reducing the peak-locking errors as reported by Overmars et al (2010). However, the optimal amount of defocusing is not possible to estimate and the excessive defocusing increases random errors in the detection of the particle images (Kislaya and Sciacchitano 2018). Furthermore, defocusing cannot be applied in tomographic measurements, where the same amount of defocusing cannot be imposed to the particles of the entire measurement volume. Michaelis et al (2016) proposed the use of optical diffusers to enlarge the particle image diameters by increasing the point spread function of the imaging system. A reduction of both systematic and random error components of the measured velocity by a factor of 3 was reported by Kislaya and Sciacchitano (2018). However, the effectiveness of the diffusers is limited for the CMOS cameras with large pixel size. Using the diffusers, a spread of about 10 µm for the incoming light can be achieved in the image plane (Michaelis et al 2016). Nevertheless, several CMOS cameras feature pixel sizes exceeding 10 µm. For those cameras, as discussed by Kislaya and Sciacchitano (2018), peak-locking errors remain present even when using two optical diffusers as proposed by Michaelis et al (2016); furthermore, the uncertainty associated with the peak locking errors remains unknown, as it cannot be easily evaluated via standard uncertainty quantification approaches (Sciacchitano et al 2015).
The discussion above shows that the problem of peaklocking errors due to small particle image diameters, especially in the case of CMOS cameras for high-speed PIV measurements, is still unsolved. In such a situation, a continuous development in the approaches based on multiple ∆t image acquisition (Nogueira et al 2009, Legrand et al 2012 has shown a high potential in peak-locking error quantification and correction, ∆t being a time separation between two frames in PIV. Using a set of different ∆t's for the same flow measurement allows for segregating the errors that scale with ∆t (e.g. peak-locking errors) from those that do not (Nogueira et al 2011). The recent work of Legrand et al (2018) has offered a 1-D analytical modeling of the peaklocking errors and has allowed for measurement correction. However, the method is iterative and computationally expensive to estimate the calibration coefficients mentioned in the algorithm. Also, selection of two ∆t's is not trivial in presence of turbulence in the flow and the ∆t values should be adjusted for different levels of turbulence. In the present work, a simple approach based on a least-squares regression of the measured time-averaged displacements and Reynolds stresses from multiple ∆t acquisitions is proposed to correct the peak-locking errors and quantify the uncertainty. The proposed methodology and its validation are explained in sections 2 and 3, respectively. The approach is then assessed for planar PIV measurements of the flow over a NACA0012 airfoil in a wind tunnel; the experimental setup and the results of such assessment are presented in sections 4 and 5, respectively.

Proposed methodology
The peak-locking uncertainty in PIV measurements is proposed to be quantified and reduced using regression analysis from acquisitions with multiple ∆t's. The approach, applicable for stationary processes (whose statistics do not change in time), is devised for the uncertainty in the time-averaged displacements (and velocities) as well as Reynolds stresses.

Mean displacement and velocity
Consider a stationary process, where the local actual timeaveraged flow velocity u true is constant in time. When PIV measurements are conducted with inter-frame time separation ∆t, the actual time-averaged particle image displacement in pixel units, not affected by any measurement error, is: where, M is the magnification factor relating the image plane to the object plane and s is pixel size of the camera sensor. From equation (1) it is evident that the actual particle image displacement increases linearly with the inter-frame time separation. In presence of the peak-locking errors, the measured particle image displacement is the sum of the true displacement and the measurement error: The peak-locking errors are often modeled as a sinusoidal (or close-to-sinusoidal, Cholemari 2007) function of the sub-pixel particle image displacements, hence they vary nonlinearly with ∆x true . Therefore, the measured particle image displacement becomes a non-linear function of ∆t, as illustrated in figure 1. In the proposed approach, the measurements are repeated with multiple ∆t's, provided that the ∆t's are selected such as to sample a sufficient portion of the peak-locking period (i.e. the variation of the time-averaged displacements is at least one pixel). Performing a linear regression of the measured displacements ∆x (∆t) allows to average out the systematic errors ε ∆x , thus yielding a regression displacement ∆x regr (∆t) which is a better estimate of ∆x true (see figure 1).

Quantification of peak-locking uncertainty.
Following Coleman and Steele (2009), the uncertainty of the time-averaged displacement is divided into systematic uncertainty and random uncertainty: where, U b,∆x and U r,∆x (∆t) are the systematic and random uncertainties in the time-averaged displacement ∆x, respectively. In this work, we assume that the systematic uncertainty U b,∆x is mainly ascribed to peak locking and not to the other bias error sources, such as modulation effects due to limited spatial resolution (Nogueira et al 2005) or bias errors due to velocity gradients (Keane andAdrian 1990, Cowen andMonismith 1997). The random uncertainty U r,∆x (∆t) can be attributed to several sources, including the actual flow fluctuations, noise in the recordings and out-of-plane particle motion (Raffel et al 2018, Sciacchitano 2019.
Considering the regression displacement as the best estimate of the true displacement, the systematic uncertainty can be estimated based on the difference between the measured and regression displacements: where, t C.I.,ν is the t-statistic describing the desired confidence interval and n is the number of ∆t acquisitions. In presence of severe peak locking, this systematic uncertainty constitutes a major part of the total uncertainty in the mean displacement, as can be seen in figure 1, where the measured timeaveraged displacements are underestimated or overestimated producing large systematic error and leading to large systematic uncertainty.
The random uncertainty in the time-averaged displacement is due to the finite number of measurement samples (Coleman and Steele 2009) and is given by: where N is the number of instantaneous measured displacements (∆x) at each ∆t acquisition. It is to be noted that the regression analysis does not affect the random uncertainty and random error in the measurements. The random uncertainty (and random error) in the time-averaged displacement at each ∆t acquisition can be reduced by acquiring large enough number of samples at that ∆t acquisition. Moreover, when time-resolved PIV measurements are conducted, the instantaneous measured displacements are correlated. As a consequence, an additional term appears in the calculation of the total uncertainty, namely the correlated systematic uncertainty U b,corr(∆x) (Coleman and Steele 2009): The value thus calculated in equation (6) should be added to the two terms under radical in equation (3) to calculate the total uncertainty in the time-averaged displacement.
The equation (6) can further be simplified considering that the covariance between the bias errors in the successive instantaneous displacements can be expressed as the product of the correlation coefficient (ρ) between them and the individual standard deviations (σ) assuming normal distribution for the bias errors (ε b ).
The estimation of the correlated systematic uncertainty U b,corr(∆x) is non-trivial as the bias errors in the instantaneous displacements are unknown and difficult to determine. Moreover, most of the measurements are followed by statistical analysis requiring uncorrelated samples in which case U b,corr(∆x) is negligible (or zero for completely uncorrelated samples). In this work, for sake of simplicity, only this case of uncorrelated samples is considered and the correlated systematic uncertainty is not included in the calculation of the total uncertainty. Finally, the total uncertainty in the measured time-averaged velocity u at each ∆t acquisition can be estimated using Taylor's approach from the uncertainty in the time-averaged displacement (calculated using equation (3)) at that particular ∆t acquisition: 2.1.2. Correction of peak-locking errors. The peak-locking errors in the time-averaged PIV measurements can be corrected for by removing the bias errors estimated via regression analysis. Thus, the measured displacements and velocities (∆x and u) are replaced with the displacements and velocity from regression (∆x regr and u regr ), respectively. The total standard uncertainty in the corrected displacement at each ∆t acquisition is then given by the uncertainty of the regression model (Coleman and Steele 2009): where, U b,∆x is the systematic peak-locking uncertainty in the measured mean displacement (calculated using equation (4)), n is the number of ∆t's and σ ∆t,∆t is the covariance in the ∆t's: The uncertainty in the corrected velocity is then derived from the uncertainty in the corrected displacement using Taylor's approach. It is to be noted that the uncertainty in the velocity is calculated at the mean of the ∆t's, that is where the regression error is the minimum (Montgomery et al 2011).

Reynolds stress
When the flow is turbulent and the velocity fluctuations encompass at least one pixel, the mean velocity is usually unaffected by the peak-locking errors in the instantaneous velocities. In such cases, the peak-locking errors appear instead in the velocity fluctuations leading to inaccurate estimation of higher-order turbulent statistic, e.g. Reynolds stress (Christensen 2004). The measured displacement fluctuations (∆x ′ ) are a combination of the actual flow fluctuations (∆x ′ true ), the random error of the displacement measurement (ε ∆x ) and the bias error fluctuations in the measured displacements (b ′ ∆x ), as also explained by Wilson and Smith (2013a), Sciacchitano and Wieneke (2016) and Scharnowski et al (2019): The equation (12) can be written in terms of total fluctuating measurement error (δ ∆x = ε ∆x + b ′ ∆x ) as: and in terms of velocity fluctuations as: The equation (14) can also be written in terms of the total fluctuating error in physical units δ ∆X = δ∆x ∆t · s M as: where, δ ∆X is the total fluctuating error in the measured displacement in physical units (e.g. m or cm). The measured variance of the velocity, or Reynolds stress, is then given by: where, R uu (or σ u 2 ) and R uu,true (or σ u,true 2 ) are the measured and 'true' Reynolds stresses, respectively. The relative importance of the different terms in the righthand-side of equation (16) is evaluated by means of Monte Carlo simulations. Two flow cases are considered, namely: (a) mean velocity of 10 m s −1 (same as the free stream velocity) and low flow fluctuations with turbulence intensity of 4%; (b) mean velocity of 2.5 m s −1 and high flow fluctuations with turbulence intensity of 31%.
The two flow cases are representative of measurements in a potential-like flow region and turbulent region, respectively. The turbulence intensity is calculated as the ratio of RMS velocity fluctuations to the free stream velocity. The maximum absolute values of the peak-locking bias error and the measurement random error are 0.25 pixel and 0.05 pixel, respectively and are kept the same for the two cases for ease of comparison. It is to be noted that the peak-locking bias error is modeled as a sinusoidal function of the sub-pixel displacement following Cholemari (2007). Also, the camera pixel size (s) of 20 µm and magnification factor (M) of 0.4 are considered while generating the synthetic data for ∆t's = 10, 11, …, 30 µs; those values are representative of a typical high-speed PIV experiment conducted with CMOS camera sensors. The measured Reynolds stresses for these ∆t's along with the true Reynolds stress, error variance over ∆t 2 and the covariance (between the velocity fluctuations and fluctuating error) over ∆t are plotted in figures 2(a) and (b) for the above mentioned synthetic data examples.
When the flow fluctuations are low compared to the peaklocking bias error, as in the case (i) of potential flow, the effect of peak-locking errors is visible in terms of the non-linear (sine-like function) variation in the measured Reynolds stress with respect to ∆t or the sub-pixel displacement. This is due to the variation in the variance of total fluctuating error (σ 2 δ,∆X ) and the covariance of the true flow fluctuations and fluctuating errors in the measurements [cov(u ′ true ,δ ∆X )] with different ∆t's. It is to be noted that the covariance can be positive or negative based on the amount of flow fluctuations and the magnitude of peak-locking errors which in turn depends on the measured sub-pixel displacement at a particular ∆t acquisition. In particular, when the average fractional displacement is close to zero, as it occurs for ∆t's = 10, 11, 14, 15 and 16 µs, positive flow fluctuations will yield negative peak-locking errors (and vice-versa), thus resulting in negative values of the covariance terms at those ∆t acquisitions (figure 2(a)). Instead, when the average fractional displacement is close to 0.5 pixels, positive (respectively negative) flow fluctuations will yield positive (negative) peak locking errors, thus resulting in positive values of the covariance terms (see for instance the results for ∆t's = 12 and 13 µs in figure 2(a)). It is to be noted that this behavior of the covariance [cov(u ′ true ,δ ∆X )] with respect to ∆t in the case (i) of potential flow is due to the flow fluctuations being smaller than the magnitudes of peaklocking errors. When the flow fluctuations are higher than the magnitudes of peak-locking errors, the behavior of covariance is not trivial as can be seen in the case (ii) of turbulent region in figure 2(b). However, the covariance terms in the case (ii) also depend on the flow fluctuations, amount of peak locking and the measured sub-pixel displacements as in the case (i). Therefore, similar overall trends can be seen for the measured Reynolds stresses, error variances (σ 2 δ,∆X ) and covariances [cov(u ′ true ,δ ∆X )] in both the synthetic data cases of potential flow region and turbulent region in figures 2(a) and (b), respectively. However, the non-linear variation with respect to ∆t is smother in the case (ii) compared to that in the case (i). It is due to the flow fluctuations being higher than the magnitude of peak-locking error in this case. Thus, in presence of peak locking, the Reynolds stresses may be overestimated or underestimated depending on the inherent flow fluctuations and the magnitude of peak-locking errors. In absence of peak locking, the Reynolds stresses are usually overestimated due to the random errors in the measurement as illustrated by Wilson and Smith (2013a) and also demonstrated by Scharnowski et al (2019) who estimated the turbulence levels of the flow facility using a multi-∆t approach.
It is clear from equation (16) that the Reynolds stresses measured at a particular spatial location are different when the PIV measurements are carried out with different inter-frame time separations (∆t's). Moreover, the error variance (σ 2 δ,∆X ) and the covariance [cov(u ′ true ,δ ∆X )] are functions of ∆t and in turn functions of sub-pixel displacement ∆x sub . Following Cholemari (2007), the fluctuating bias error due to peak locking is modeled as a sinusoidal function of measured sub-pixel displacement: Assuming that ε ∆x,i and b' ∆x,i are uncorrelated, the variance of the total errors in the displacement measurement at each ∆t acquisition is then given by: (18) Moreover, the synthetic data analysis shows that the covariance between the true velocity fluctuations and the fluctuating error varies non-linearly with ∆t or the sub-pixel displacement. Therefore, the covariance at each ∆t can be approximated as a cosine function of the mean measured sub-pixel displacement ∆x sub at that ∆t acquisition (figure 2): Substituting the equations (18) and (19) in the equation (16), the Reynolds stress is expressed as: It is to be noted that the variance of the random errors var(ε ∆x,i ) in the equation (20) is constant for different ∆t acquisitions. Moreover, s, M, A and B are constants. Therefore, the equation (20) can also be written as: Thus, it is proposed to perform a least-squares non-linear regression of the form of equation (21) to estimate the 'true' Reynolds stress from the measured Reynolds stresses at different ∆t's. The coefficients β 1 and β 2 are positive as they represent the variance and square of a constant, respectively.
For the ease of processing, the non-linear regression model in the equation (21) can be transformed to a linear regression model (Montgomery et al 2011) as: where, In principle, the accuracy of the measured Reynolds stresses is the maximum when a large ∆t is selected (∆t → ∞). However, this condition cannot be met in practice due to limitations caused by the out-of-plane particle motion and by the increasing truncation errors (Raffel et al 2018). Hence, the 'true' Reynolds stress is estimated by conducting measurements with different ∆t's, and then fitting a curve of the type mentioned in equation (22) by least-squares regression, where R uu,true represents the best estimate of the 'true' Reynolds stress.
The results of regression for the synthetic data can be seen in figures 3(a) and (b) for the potential flow region and turbulent region, respectively.
It is clear that the total Reynolds stress from regression is comparable to the measured Reynolds stress at each ∆t in both the cases. Also, depending on the ∆t, the measured Reynolds stress can differ significantly from the true Reynolds stress, as much as by 0.5 m 2 s −2 in the potential flow region and by 1 m 2 s −2 in the turbulent region. When the true Reynolds stress is estimated from the regression analysis, this difference is reduced to 0.07 m 2 s −2 and 0.05 m 2 s −2 , respectively. The quantification of the uncertainty of the estimated 'true' Reynolds stress is discussed in the following section 2.2.1.

Quantification of peak-locking uncertainty.
The estimate of the 'true' Reynolds stress obtained from regres-sionR uu,true represents a correction to the measured Reynolds stress, where the systematic error due to peak locking is removed, and can be used to quantify the systematic uncertainty in the measured Reynolds stress R uu at a given ∆t: where, t C.I.,ν is the t-statistic describing the desired confidence interval and n is the number of ∆t acquisitions.
Moreover, the random uncertainty in the measured Reynolds stress (R uu ) is given by (Sciacchitano and Wieneke 2016): where, N is the number of samples acquired at that ∆t.
The random errors of the measured instantaneous velocity vectors are known to yield overestimated Reynolds stress values (Wilson and Smith 2013b), consistently with the second term of the right-hand-side of equation (16). However, in presence of peak locking, the Reynolds stresses may be overestimated or underestimated depending on the covariance between the actual velocity fluctuations and the measurement error fluctuations. As a consequence, the uncertainty distribution of the Reynolds stresses is asymmetric: where, U Ruu,1 and U Ruu,2 are the total uncertainties in the measured Reynolds stresses, either of them representing the upper or lower uncertainty bounds depending on the shape of error distribution (Wilson and Smith 2013a). Moreover, the shape of error distribution is influenced by the flow fluctuations, affected by the peak-locking errors, and it is different in the flow regions with low fluctuations than those with high fluctuations compared to the magnitude of the peak-locking errors.
Therefore, in the case of overestimated Reynolds stresses U Ruu,1 and U Ruu,2 represent the lower and upper uncertainty bounds on the Reynolds stresses, respectively. Whereas, in the case of underestimated Reynolds stresses the lower and upper uncertainty bounds are given by U Ruu,2 and U Ruu,1 , respectively. The uncertainty of the estimated 'true' Reynolds stress from regressionR uu,true can be expressed in terms of the uncertainty of the total Reynolds stress from regression R uu,regr (Montgomery et al 2011): where, U Ruu,regr is the uncertainty of the regression model and is calculated as: and P is a n × 4 matrix of regression variables p 1 , p 2 , p 3 where, n is the number of ∆t acquisitions.

Validation of proposed methodology
It can be observed that the uncertainties in the measured mean velocity and Reynolds stress decrease with the increasing ∆t (or ∆x) as the relative peak-locking errors for larger displacements (at larger ∆t's) are less than those for the smaller displacements (at smaller ∆t's), as also reported by Wilson and Smith (2013b). It is to be noted that the considered range of ∆t's should be sufficiently small such that the effect of other error sources such as velocity gradients and out-ofplane particle motion can be considered negligible. This fact is used for validation of the proposed uncertainty quantification approach where reference measurements are carried out with a much larger ∆t than the actual measurement ∆t's. The errors in the measured quantities can then be calculated with respect to the reference quantities. Comparing the errors with the uncertainties estimated using the proposed methodology (explained in sections 2.1.1 and 2.2.1), it is possible to evaluate the effectiveness of the uncertainty quantification approach. Following Kislaya and Sciacchitano (2018), the reference measurements can be conducted with a time separation ∆t aux which is much larger than the ∆t's used in the actual measurements. Due to the large displacement ∆x aux in the reference measurements, the relative peak-locking error in ∆x aux is negligible with respect to that in the measured displacements at ∆t's. The reference displacement ∆x ref at each ∆t is thus calculated from the auxiliary displacement ∆x aux measured at ∆t aux : where, N is the number of instantaneous measured displacements (or velocities) in the reference measurements.
The errors in the measured quantities at each ∆t can then be calculated with respect to the reference quantities: These errors in the measured displacements, velocities and Reynolds stresses are used to estimate the reliability of the uncertainty quantification approach by means of uncertainty coverages. An uncertainty coverage is the percentage of measurements for which the error magnitudes are contained within the uncertainty limits (Timmins et al 2012): When U represents the standard uncertainty, the uncertainty coverage (C U ) should be equal to the error standard deviation coverage (C SD ); the latter is calculated as the percentage of measurements for which the errors are within the standard deviation of the errors distributions: In case of Gaussian error distribution, the above is analogous to comparing the one-sigma uncertainty coverage to 68%, as done in Sciacchitano et al (2015), among others. However, in the case of non-normal error distributions, it is necessary to first calculate the standard deviation of the errors and the corresponding standard deviation coverage to compare it with the uncertainty coverage. To assess the effectiveness of the uncertainty quantification approach, a term called relative coverage index (RCI) is proposed and is defined as the ratio between uncertainty coverage (C U ) and error standard deviation coverage (C SD ): Values of RCI close to one indicate accurate uncertainty estimations. Values of RCI less than one indicate that the uncertainty is underestimated, whereas values greater than one show that it is overestimated. Note that this parameter is applicable to any error distribution after estimation of the errors standard deviation.

Experimental setup
Planar PIV measurements of the flow over a NACA0012 airfoil at 10 degrees angle of attack were performed. Figure 4 illustrates a schematic of the experimental setup.
The experiment was conducted in the W-tunnel of Delft University of Technology, which is an open-jet open-return wind tunnel with an exit cross section of 0.4 × 0.4 m 2 and an area contraction ratio of 9. The maximum achievable freestream velocity is 30 m s −1 with 0.3% turbulence intensity (Tummers 1999). In this experiment, the free stream velocity was set to 10 m s −1 . The flow was seeded by a SAFEX seeding generator, which produces water-glycol droplets of 1 µm mean diameter. The particles were illuminated by a Quantronix Darwin-Duo laser (Nd:YLF, pulse energy of 25 mJ at 1 kHz, wavelength of 527 nm) and images were recorded with a LaVision High Speed Star 6 CMOS camera (12 bits, 20 µm pixel size, 1024 × 1024 pixels maximum resolution). The camera was equipped with a Nikon objective of 105 mm focal length and the optical aperture was set to f # = 5.6. A field of view of 45 mm × 45 mm was imaged with optical magnification of 0.46, which resulted in a theoretical particle image diameter of 0.5 pixel, following Raffel et al (2018). The small size of the particle image diameter caused peak locking in the measured displacements, influencing the estimation of mean velocity and Reynolds stresses. The image acquisition was conducted at a frequency of 200 Hz independently for multiple ∆t's from 10 µs to 24 µs in steps of 2 µs. The data sets at each ∆t consisted of 500 double-frame images. To validate the proposed multi-∆t approach, reference measurements were conducted with optical aperture f # = 11 and the time separation ∆t of 50 µs, which is larger than the ∆t's used in the actual measurements, following Kislaya and Sciacchitano (2018). It is to be noted that the effect of out-of-plane motions and gradient bias errors is negligible in the measurements with ∆t of 50 µs, as the free-stream particle image displacement is equal to 11.5 pixels. The processing was done via the LaVision DaVis 10 software, using Gaussian interrogation windows of 64 × 64 pixels with 75% overlap for the initial passes and 16 × 16 pixels with 75% overlap for the final passes.
The estimated time-averaged velocities (u) and RMS of the velocity fluctuations (u ′ rms ) in the reference measurements, normalized with respect to the free stream velocity (u ∞ ), are shown in figures 5(a) and (b), respectively. It is clear that the flow has varying degrees of fluctuations, e.g. low fluctuations in the potential flow region (0 < u ′ rms /u ∞ < 0.05) and relatively high fluctuations in the turbulent region (0.2 < u ′ rms /u ∞ < 0.6). As a consequence, the measured flow field is considered suitable to assess the effectiveness and limitations of the proposed approach in a range of flow conditions encountered in typical PIV measurements. A detailed analysis of the peak-locking errors correction and uncertainty quantification is conducted in regions I and II of figure 5(b), corresponding to the potential flow region and a turbulent region, respectively.

Results
The approach proposed in section 2.1 is applied to correct the peak-locking errors in the measured mean displacements and velocities, where the uncertainties in the corrected displacements and velocity obtained from regression are also provided. The analogous approach discussed in section 2.2 is then applied to correct the measured Reynolds stresses and quantify the uncertainty. The results are shown for the two selected regions in the flow as marked in figure 5(b).

Mean displacement and velocity
From the displacement plots of figure 6, it is evident that, while the reference displacement increases linearly with ∆t, the measured mean displacement does not, consistently with the discussion of section 2.1. It is clear from figures 6(a) and (b) that the regression displacement is much closer to the reference displacement, indicating a significant reduction in the peak-locking errors. The reduction in the peak-locking errors can further be seen in figures 7(a) and (b), where the errors in the measured mean displacement and the regression displacements at each ∆t, calculated with respect to the corresponding reference displacements, are shown. The maximum errors in the measured displacements are observed to be around ±0.15 pixels in both the regions I of potential flow (figure 7(a)) and II of high turbulence ( figure 7(b)). Whereas, the errors in the regression displacements are very close to zero in both the regions, being below 0.05 pixels. The peak-locking errors in the measured displacements are reduced up to 86% in the potential flow region I and up to 73% in the turbulent region II, when correcting the measured displacements by linear regression. The uncertainty bounds are also shown for both the uncorrected (measured) and corrected displacements (from regression) in figure 6. In the potential flow region, as shown in figure 6(a), the uncertainties in the measured displacements are around 0.07 pixels at different ∆t's, whereas those on the regression displacements are around 0.03-0.05 pixels.
Similarly, in the turbulent region, shown in figure 6(b), the uncertainties in the measured displacements are around 0.09 pixels at different ∆t's, and those on the regression displacements are around 0.03-0.05 pixels. Figures 8(a) and (b) show the distributions of the errors in the measured displacements and the regression displacements with respect to the reference displacements in the regions I and II, respectively. In the potential flow region ( figure 8(a)), the pdfs of error distributions for different ∆t acquisitions are narrow and separated from each other and from the pdf of regression error distribution. Whereas, in the turbulent region ( figure 8(b)), they are wide and close to each other and close to the regression one. The reason behind the pdfs being narrow in the potential flow region and wide in the turbulent region is the inherent fluctuations in the flow which are low in the potential flow region I and relatively high in the turbulent region II. Moreover, in the potential flow region, where the fluctuations are low, the peak-locking errors are significant and vary drastically (largely overestimated or underestimated) with varying ∆t. As a consequence, the error pdfs exhibit different mean values, from −0.13 pixels at ∆t = 12 µs to +0.16 pixels at ∆t = 10 µs. Although the errors at most ∆t's exhibit an approximately normal distribution, at some ∆t's the distribution is skewed positively (e.g. ∆t = 12 µs) or negatively (e.g. ∆t = 20 µs). Conversely, the error distribution stemming from the linear regression has a Gaussian shape, with 0.02 pixels mean and 0.001 pixels standard deviation. In contrast, in the turbulent region, the large velocity fluctuations induce particle image displacement variations larger than the peak-locking errors; as a consequence, the errors pdfs exhibit lower variations with varying ∆t. Nevertheless, the mean errors vary between −0.11 pixels and +0.08 pixels depending on the selected ∆t. The application of the linear regression reduces the mean bias error to −0.04 pixels. It should be remarked that, as it is clear from figures 8(a) and (b), the pdfs of displacement error distribution in different ∆t acquisitions lie on either side of the pdf of regression displacement error distribution. This is due to the fact that the measured mean displacements are overestimated or underestimated depending on the sub-pixel displacement which in turn depends on ∆t.
Analogous results in terms of time-averaged velocity are illustrated in figure 9. The errors in the measured velocities calculated with respect to the reference velocity are reduced up to 86% in the potential flow region I and up to 73% in the turbulent region II, after correcting the measured velocities by the regression analysis. In the potential flow region (figure 9(a)), the uncertainty in the corrected velocity estimated from the regression analysis is 0.07 m s −1 , which is much smaller than the uncertainties in the measured mean velocities (0.13-0.31 m s −1 ) at different ∆t's. Similarly, in the turbulent region ( figure 9(b)), the uncertainties of 0.17-0.36 m s −1 in the measured mean velocities are reduced to the uncertainty of 0.07 m s −1 in the corrected velocity obtained from regression. Thus, the overall reductions of the velocity uncertainties between 50% and 80% are achieved in both regions I and II.
A visual comparison of the reference, measured and regression velocities is conducted in figure 10. A particular zoomed area in the flow field is selected to make the comparison vivid. It is clear that, the mean measured velocities at different ∆t acquisitions are either overestimated or underestimated with respect to the velocities in the reference measurement depending on the sub-pixel displacement. For example, in the top region of the selected area, the measured velocities are overestimated at ∆t = 10 µs and 24 µs, whereas, they are underestimated at ∆t = 12 µs, consistent with figure 9(a) for the potential flow region. These measured velocities at different ∆t acquisitions are thus corrected by the regression analysis as can be seen in the third row of figure 10. It is clear that the corrected mean velocities from the regression are comparable to the corresponding mean velocities from the reference measurement.
The uncertainty coverages (C U ) for the uncertainties in the measured displacements and velocities are calculated following equation (39). The uncertainties exhibit a coverage of 52% and 62% considering the measurements at different ∆t acquisitions in the potential flow region I and turbulent region II, respectively. Whereas, the error standard deviation coverages (C SD ), calculated by the equation (40) are 61% and 70% in the regions I and II, respectively. Thus, the RCIs for the uncertainty quantification in the measured displacements and velocities are 0.85 and 0.89 in the potential flow region I and turbulent region II, respectively. This analysis proves the validity of the proposed approach for uncertainty quantification as the RCI's are close to one and the uncertainties are slightly underestimated in both the regions.

Reynolds stress
The results of the proposed methodology for correcting the Reynolds stresses are summarized in figures 11(a) and (b) in the potential flow region and turbulent region, respectively. It is clear from figure 11(a) that the effect of peak-locking errors is significant and can be seen by non-linear variation of the measured Reynolds stresses with respect to ∆t. Such result is associated with the magnitude of flow fluctuations being lower than the magnitude of peak-locking errors in this region I of potential flow, consistent with the results of the Monte Carlo simulations described in section 2.2. Conversely, in the turbulent region II, the measured Reynolds stresses show smoother variation with respect to the changing ∆t. In this region, the effect of the peak-locking errors is hidden due to the inherent flow fluctuations being higher. The reference Reynolds stress in the potential flow region is 0.08 m 2 s −2 which is comparatively low due to the low flow fluctuations in this region. However, the measured Reynolds stresses overestimate the reference value by up to 700%. The regression approach provides an estimate of the Reynolds stress that agrees within 50% with the true value. Conversely, in the turbulent region, higher flow fluctuations are encountered and the reference Reynolds stress is 2.0 m 2 s −2 .
The measured Reynolds stresses in this region are overestimated by 4% to 95%. When the regression approach is employed, the estimated Reynolds stress is slightly overestimated by 19% with respect to the reference value.
The errors in the measured Reynolds stresses are calculated with respect to the reference one using the equation (38)  fluctuations in the region I and high flow fluctuations in the region II, as also observed for the displacement error distributions. The errors in the regression Reynolds stress are comparatively lower than that in the measured Reynolds stresses at different ∆t's. Most of the errors in the potential flow region are positive, whereas those in the turbulent region are either positive or negative. Moreover, most of the pdfs for the measured Reynolds stresses lie on the positive side of the regression pdf in region I; in region II, instead, the pdfs of the measured Reynolds stresses lie on both positive and negative sides of the regression pdf. This fact is accounted for in the calculation of the total uncertainty in the Reynolds stresses at each ∆t using equations (28) and (29). For the Reynolds stresses overestimated with respect to the regression estimate, U Ruu,1 represents the lower uncertainty bound and U Ruu,2 is the upper uncertainty bound. Whereas, in case of the underestimated Reynolds stresses, U Ruu,1 and U Ruu,2 are the upper and lower uncertainty bounds, respectively, as can be seen in figures 11(a) and (b). In the potential flow region I, the largest uncertainty bound in the measured Reynolds stresses is 0.35 m 2 s −2 , whereas the uncertainty in the corrected Reynolds stress (i.e. the estimate of the 'true' Reynolds stress from regression) is 0.09 m 2 s −2 . At the point in the turbulent region II, the bigger uncertainty bound in the measured Reynolds stresses is around 0.95 m 2 s −2 , whereas the uncertainty in the corrected Reynolds stress is 0.61 m 2 s −2 . The uncertainties in the measured Reynolds stresses are reduced by around 75% in the potential flow region, and by around 35% in the turbulent region, after correction with the estimates of the 'true' Reynolds stresses from regression. The uncertainty coverages (C U ) are determined by comparing the calculated uncertainties with the errors (computed with respect to the reference measurements) in the measured Reynolds stresses as per the equation (39). The coverages for the errors standard deviation are also calculated using the equation (40). In the potential flow region I, the uncertainties exhibit a coverage (C U ) of 73% and the corresponding error standard deviation coverage (C SD ) is 61%. These values lead to the RCI of 1.2 indicating that the uncertainties in the measured Reynolds stresses are slightly overestimated in the potential flow region. In the turbulent region II, the uncertainty coverage is 64% and the error standard deviation coverage is 72%, resulting in RCI = 0.89. Thus the Reynolds stress uncertainties are slightly underestimated in the turbulent region. The results of the regression analysis are also summarized in figure 13, where the reference Reynolds stresses and the corrected Reynolds stresses are shown along with the measured Reynolds stresses at different ∆t acquisitions. It is to be noted that the measured values only at three ∆t's are plotted for sake of conciseness. Consistent with the findings of figure 11(b), the measured Reynolds stresses in the bottom region of the selected area are overestimated at ∆t's = 12 and 24 µs. The corrected Reynolds stresses by the regression analysis in this bottom region are slightly underestimated with respect to the reference values and are found to be comparable to the measured Reynolds stresses at ∆t = 10 µs acquisition. Similarly, in the potential flow region of this selected area in figure 13, the corrected Reynolds stresses are comparable to the reference ones, as also were seen in the figure 11(a). The results show that the peak-locking errors can be reduced by the regression analysis in the Reynolds stresses.

Selection of ∆t's
The proposed approach is based on the simple concept of leastsquares regression of the measured quantities from multiple ∆t acquisitions. The selection of ∆t's is an important step in the methodology. The first requirement for the selection of the ∆t's is that a minimum number of four ∆t's is considered. This requirement is a direct consequence of the fourdegree of freedom model (β 0 , β 1 , β 2 and β 3 ) employed for the Reynolds stress. Additionally, it is deemed likely that the accuracy of the regression result degrades when the range of the selected ∆t's covers a particle image displacement below one pixel, that is less than a period of the peak-locking error. Hence, the second requirement for selecting the ∆t's is that the particle image displacement range, defined as the difference between the largest time-averaged displacement (occurring in the acquisition with the largest ∆t) and the smallest time-averaged displacement (for the smallest ∆t), should be at least one pixel. In this section, the effect of the number of ∆t's selected for the regression analysis and the displacement range is investigated. Based on the two requirements, two sets of four values of ∆t's each are selected arbitrarily in the NACA0012 experiment. It is to be noted that, with four ∆t acquisitions, it is possible to perform the regression analysis for Reynolds stresses and to correct the peak-locking errors by the estimated 'true' Reynolds stress from the regression. However, it is not feasible to compute the uncertainty in the estimate of 'true' Reynolds stress, as the regression curve passes exactly through all the data points. The two sets of ∆t's are selected as follows: set I consists of measurements at ∆t's equal to 10, 12, 14, 16 µs, whereas in set II the measurements are performed at ∆t's of 10, 14, 18, 22 µs. The proposed approach of regression is applied for both the cases of ∆t's to correct the measured mean displacements (and velocities) and Reynolds stresses. It can be seen in figure 14(a) that the measured displacements cover more than one pixel range in both cases of ∆t's. Hence, the corrected displacements in both the cases are comparable to the reference displacements, and to each other. However, when the measured displacements do not cover one pixel range, the regression is not effective. This result is due to the majority of the measured displacements lying on one side of the reference displacements, resulting in the regression displacements to be on the same side and away from the reference one. Thus, for the effective displacements regression, it is necessary to have both the overestimated and underestimated measured displacements, which is possible when the range is greater than or equal to one pixel.
When analyzing the Reynolds stresses results ( figure 14(b)), it is clear that the corrected Reynolds stress in the case II is closer to the reference one than that in the case I. This is because of the larger range of ∆t's in the case II. Thus the regression in the case II performs slightly better than that in the case I. Nevertheless, the Reynolds stress estimated with the four ∆t's of set II is less accurate than the result obtained from the full set of eight ∆t's, which has been shown in figure 11(a). Hence, although four is the minimum number of time separations for the regression analysis, the use of more ∆t's can further increase the accuracy of the estimated Reynolds stresses.
Furthermore, to quantify the uncertainty in the estimate of 'true' Reynolds stress from regression, more than four data points or measured Reynolds stresses are required. Therefore, in the present work, eight ∆t's were selected with the gap of 2 µs between the consecutive values while acquiring the measurements, the results of which were presented in the sections 5.1 and 5.2. As a guideline for future applications of the proposed methodology, it is suggested that more than four ∆t's are selected, and that the particle image displacement range covers at least one pixel in the regions of interest, to ensure an effective evaluation and correction of the peaklocking errors.

Conclusions
A novel approach is proposed for the quantification of the peak-locking systematic uncertainty in PIV, which in turn also corrects the measured particle image displacement (and velocity) and Reynolds stress for peak-locking errors. The approach is based on the assumption that local flow statistics are constant in time, according to which the particle image displacement should vary linearly with the time separation (∆t) between two frames. However, in the presence of peak locking, the measured particle image displacement is a non-linear function of ∆t as the measurement error in the displacement varies non-linearly with sub-pixel particle image displacement. Similarly, in presence of peak locking the Reynolds stress may vary significantly depending on the selected ∆t. Hence, in the present approach, it is proposed to acquire the measurements at various ∆t's and perform a least-squares regression of the measured quantities (displacements and Reynolds stresses). The displacement (and velocity) and Reynolds stress from regression represent more accurate estimates of the 'true' values and thus used to quantify systematic errors and uncertainty due to peak locking in the measured displacement (and velocity) and Reynolds stresses, respectively. Guidelines for the selection of the ∆t's are also provided. The methodology is assessed for planar PIV measurements of the flow over a NACA0012 airfoil at 10 degrees angle of attack. The uncertainties in the measured velocity and Reynolds stresses are reduced by 50%-80% and 35%-75%, respectively, after correcting by the regression analysis. RCIs close to one are obtained for the mean velocity uncertainty and Reynolds stress uncertainty in both the potential and turbulent regions of the flow, indicating the ability of the methodology to quantify the uncertainty associated with peak-locking errors.