# Evidence of preferential sweeping during snow settling in atmospheric
turbulence

Jiaqi Li1,2  Aliza Abraham1,2  Michele Guala1,3  Jiarong Hong1,2
<EMAIL_ADDRESS>1Saint Anthony Falls Laboratory, University of Minnesota,
Minneapolis, MN, 55455, USA 2Department of Mechanical Engineering, University
of Minnesota, Minneapolis, MN, 55455, USA 3Department of Civil, Environmental
and Geo- Engineering, University of Minnesota, Minneapolis, MN, 55455, USA

###### Abstract

We present a field study of snow settling dynamics based on simultaneous
measurements of the atmospheric flow field and snow particle trajectories.
Specifically, a super-large-scale particle image velocimetry (SLPIV) system
using natural snow particles as tracers is deployed to quantify the velocity
field and identify vortex structures in a 22 m $\times$ 39 m field of view
centered 18 m above the ground. Simultaneously, we track individual snow
particles in a 3 m $\times$ 5 m sample area within the SLPIV using particle
tracking velocimetry (PTV). The results reveal the direct linkage among vortex
structures in atmospheric turbulence, the spatial distribution of snow
particle concentration, and their settling dynamics. In particular, with snow
turbulence interaction at near-critical Stokes number, the settling velocity
enhancement of snow particles is multifold, and larger than what has been
observed in previous field studies. SLPIV measurements show higher
concentration of snow particles preferentially located on the downward side of
the vortices identified in the atmospheric flow field. PTV, performed on high
resolution images around the reconstructed vortices, confirms the latter trend
and provides statistical evidence of the acceleration of snow particles, as
they move toward the downward side of vortices. Overall, the simultaneous
multi-scale particle imaging presented here enables us to directly quantify
the salient features of preferential sweeping, supporting it as an underlying
mechanism of snow settling enhancement in the atmospheric surface layer.

###### keywords:

atmospheric flows, particle/fluid flows, turbulent boundary layers

##  1 Introduction

Understanding the settling dynamics of inertial particles in turbulence is
important for predicting particle transport in the atmosphere, such as aeolian
transport of dust and sand (Durán et al., 2011), formation and growth of
droplets and particle aggregates in clouds (Shaw, 2003), and precipitation of
hydrometers, such as raindrops, graupels and snowflakes (Garrett et al., 2015;
Nemes et al., 2017; Zeugin et al., 2020; Li et al., 2021). Numerous laboratory
experiments and numerical simulations have been conducted to investigate the
effects of turbulence on the behavior of inertial particles. Two evident
manifestations of particle-turbulence interaction mechanisms are the formation
of particle clusters and the modulation of their settling velocity
(Balachandar and Eaton, 2010). These phenomena are observed in certain
conditions depending on the turbulence (through the Kolmogorov time scale,
$\tau_{\eta}$), and on the particle size (${\it D_{p}}$), density ($\rho_{p}$)
and aerodynamic properties, contributing to the definition of the particle
response time, $\tau_{p}$ (Maxey and Riley, 1983). The phenomena of clustering
and enhanced settling can be described as follows: as particles preferentially
concentrate in strain-dominated regions (e.g., in between vortices), they
settle along the downward side of swirling motions as clusters. As a result,
the fall speed of the particles on the downward side is increased. This
mechanism is known as preferential sweeping (Wang and Maxey, 1993). Studies
have shown that the average settling velocity of inertial particles in
turbulence can be enhanced significantly by the preferential sweeping
mechanism (Wang and Maxey, 1993; Yang and Lei, 1998; Aliseda et al., 2002;
Good et al., 2014; Falkinhoff et al., 2020), in particular under critical
conditions, i.e. when the Stokes number $\Stk=\tau_{p}/\tau_{\eta}\approx 1$
(Yang and Lei, 1998; Aliseda et al., 2002; Ferrante and Elghobashi, 2003).
There are also other mechanisms that have been described to hinder the
settling of inertial particles in turbulence such as loitering (Nielsen, 1993)
and vortex trapping (Tooby et al., 1977), but they usually tend to be
suppressed by preferential sweeping (Good et al., 2014; Rosa et al., 2016).

Despite the large number of laboratory experiments and simulations, field
measurements of inertial particles (e.g. snow particles, droplets, and dust)
settling in the atmospheric turbulence are scarce. The lack of field evidence
is mostly due to the fact that field measurements are experimentally
challenging (Shaw, 2003): local turbulent field conditions are difficult to
parameterize, and the effects of particle interaction and flow Reynolds
numbers on non-Stokesian particle kinematics is far from being clear (see
recent advancements by Petersen et al. (2019), Tom and Bragg (2019), and
Falkinhoff et al. (2020)). Moreover, the implementation of particle-turbulence
interaction mechanisms in predictive models of settling velocity at
geophysical scales is also limited, since the field conditions (e.g. wide
range of turbulence scales, complex particle shape) are often different from
those reproduced in laboratory experiments and simulations.

To enable spatially and temporally resolved flow measurements in the field, a
super-large-scale particle image velocimetry (SLPIV), using natural snow
particles as tracers, has been recently developed for studying the wake
structure downstream of a utility scale wind turbine in the atmospheric
boundary layer (Hong et al., 2014; Dasari et al., 2019; Abraham and Hong,
2020) and for the study of high Reynolds number wall turbulence (Toloui et
al., 2014; Heisel et al., 2018). Using the similar setup, Nemes et al. (2017)
quantified the settling trajectories of 87000 snow particles in a 4 m (width)
$\times$ 7 m (height) field of view using particle tracking velocimetry (PTV),
in parallel with a digital inline holography (DIH) system, to characterize the
size and morphology of snow particles. In the absence of direct estimates of
snow particle density, the acceleration probability density function (PDF)
obtained by PTV was used to estimate the Stokes number and the aerodynamic
particle response time of snow particles (Mordant et al., 2004; Bec et al.,
2006; Ayyalasomayajula et al., 2006). Nemes et al. (2017) found that the
settling velocity of snow particles measured using PTV showed multifold
enhancements in the atmospheric turbulence, in comparison to the still-air
terminal velocity ${\it W_{p}}=\tau_{p}\cdot{g}$ predicted using the
acceleration-based aerodynamic response time. Employing the same setup, Li et
al. (2021) investigated the settling and clustering of snow particles under
various turbulence and snow conditions. They observed intense clustering and
enhanced settling velocity during near-critical Stokes conditions, showing
statistical evidence of the correlation between enhanced settling velocity and
local preferential concentration, thus indirectly supporting the preferential
sweeping mechanism. Despite these major findings, both contributions could not
simultaneously provide flow measurements and trajectories of snow particles.

In the present study, we leverage the ability of SLPIV to measure large scale
flows and PTV to observe and track individual particles, allowing us to unveil
the direct, local, linkages between coherent vortex structures, snow
concentration distribution, and settling velocity at
$\Rey_{\lambda}\sim\textit{O}(10^{3})$. We quantify here both the preferential
concentration around vortices and enhanced settling velocity on the downward
side of vortices, thus highlighting the fundamental mechanism of preferential
sweeping. The experiment setup, atmospheric conditions, and turbulence
properties are introduced in section 2. In section 3, results are presented
for the quantification and visualization of preferential sweeping mechanism.
Conclusions and discussion follow in section 4.

Figure 1: (a) The experiment setup of simultaneous super-large-scale PIV
(SLPIV) and PTV measurements; (b) sample images showing the fields of view of
SLPIV (left) and PTV (right).

##  2 Methodology

###  2.1 Field experiment setup

The field deployment was conducted to acquire data on Jan. 13, 2020 between
18:00 and 21:00 local time, at the Eolos Wind Energy Research Field Station in
Rosemount, MN. The light sheet-based super-large-scale particle image
velocimetry (SLPIV) and particle tracking velocimetry (PTV), described in Hong
et al. (2014) and Nemes et al. (2017), respectively, have been applied to
capture the turbulent flow field, the trajectories and the concentration
distribution of snow particles. We used a 5-kW search light with a curved
mirror expanding the beam vertically into a light sheet to illuminate the snow
particles. For our current measurements, the light sheet thickness is
restricted to be 10 cm (different from our previous measurements with 30 cm
diameter light beam) at the ground and it increases to about 12 cm at 10 m
considering the divergence angle of our search light. The light sheet was
oriented to be parallel with the average wind direction and minimize the out-
of-plane motion. Throughout the deployment, the instantaneous wind direction
relative to the light sheet varied from -25 degrees to 15 degrees. An
11-minute duration dataset has been selected for the measurement presented in
this paper: within the selected period of time the wind direction was stable
and well-aligned with the light sheet direction with a deviation of less than
five degrees, and the snowfall intensity was steady, providing adequate
seeding density of around 150 snow particles per $\textnormal{m}^{2}$.

A Nikon D600 (Nikon Inc.) camera and a Sony A7RII (Sony Corp.) camera were
equipped to record the overall flow field at 30 fps and 1080 $\times$ 1920
$\textnormal{pixel}^{2}$, and the motion of snow particles at 120 fps and 720
$\times$ 1280 $\textnormal{pixel}^{2}$ respectively (referred as SLPIV and PTV
in the following sections). The start time of each recording for both cameras
is documented, and a large-scale turbulent structure visible as a void at the
beginning of both videos is used to further confirm the synchronization of the
two datasets. Both cameras were placed on tripods with measured tilt angles
from the horizontal direction (table 1). The relative locations of the two
cameras and the light sheet are illustrated in figure 1a with the defined
coordinate system (x, y, z) and corresponding velocity components (u, v, w).
The specifications (e.g. duration, field of view (FOV) elevation, size and
distance to the camera, etc.) for the two cameras are shown in table 1.

Figure 2: (a) Instantaneous flow field sample from SLPIV; (b) and snow
particle trajectory samples from PTV.

Sample images for SLPIV (Nikon camera) and PTV (Sony camera) are shown in
figure 1b, and sample results are shown in figure 2. PIV analysis is conducted
using LaVision Davis 8.2.0. A multi-pass setting was adopted with a final pass
of 32 $\times$ 32 $\textnormal{pixel}^{2}$ and 50$\%$ overlap. Around 700
vectors are obtained from each image pair. For PTV, we apply the learning-
based tracking method from Mallery et al. (2020) using a long short-term
memory (LSTM) network to acquire individual trajectories. Specifically, we
first implement tracking methods from Crocker and Grier (1996) and Ouellette
et al. (2006) to our PTV data. However, due to the lower quality of field data
and relatively high particle concentration, the conventional methods generate
substantially less tracks compared to those can be determined through manual
examination, potentially causing sampling bias (i.e., preferentially shorter
tracks, sampling only downward tracks). Therefore, good quality trajectories
generated by the conventional methods are manually selected as the training
set for the learning-based method. After iterations of training process, the
well-trained model generates significantly more tracks regardless of their
direction (upward or downward). In total, there are around 460,000
trajectories longer than $10\tau_{\eta}$ (where $\tau_{\eta}$ the Kolmogorov
time scale) being identified by the tracking algorithm. After the flow field
and trajectories are obtained, we further identify the distribution of
potential vortical structures in the flow field based on swirling strength
(Zhou et al., 1999) and calculate the Lagrangian velocity and acceleration
using the trajectory information. The swirling strength is defined as
$\lambda_{ci}$, the imaginary part of the complex eigenvalues of the velocity
gradient tensor ($\textbf{{D}}=\nabla\textbf{{u}}$). Under two-dimensional
measurement, D have either two real eigenvalues ($\lambda_{r}$) or a pair of
conjugate complex eigenvalues ($\lambda_{cr}\pm i\lambda_{ci}$), where
$\lambda_{cr}$ and $\lambda_{ci}$ are absolute values. Thus, the vortices can
be identified with finite $\lambda_{ci}$ (Adrian et al., 2000). Three
threshold values of the swirling strength (0.4 $\textnormal{s}^{-1}$, 0.5
$\textnormal{s}^{-1}$ and 0.65 $\textnormal{s}^{-1}$) are applied for
detecting the vortices, and the concentration and settling velocity are
analyzed using all three threshold values. In the results section, we will
show the preferential concentration and enhanced settling with the 0.4
$\textnormal{s}^{-1}$ threshold, and figures with the other thresholds will be
shown in the supplementary material.

Following the same procedures as in our previous studies (Toloui et al., 2014;
Hong et al., 2014), the traceability of snow particles for our SLPIV
measurement is analyzed. Specifically, the spatial resolution of the SLPIV is
usually limited by the smallest interrogation window size and light sheet
thickness (whichever is larger), i.e., $l=0.66$ m in our current SLPIV
measurements. Correspondingly, the flow time scale that our measurements can
resolve is estimated as $\tau_{f}=l/u_{\textnormal{rms}}=0.62$ s, where
$u_{\textnormal{rms}}=1.07$ m/s is the r.m.s. of the streamwise velocity
fluctuations. Thus, the particle Stokes number based on the particle response
time ($\tau_{p}=1.7-20$ ms, from acceleration PDF in section 3.1) and
$\tau_{f}$ is estimated to be $\Stk=\tau_{p}/\tau_{f}=0.0028-0.032$, much
smaller than the typical threshold for good traceability (i.e., 0.1 according
to Tropea et al. (2007)). As a result, turbulent flows above the spatial and
temporal resolution limits (i.e., 0.66 m and 0.62 s, respectively) are
reasonably captured in our measurements. Note also that the mean settling
velocity of the snow particles is subtracted from the SLPIV flow field.

A digital inline holography (DIH) (Nemes et al., 2017) system was deployed,
and around 1580 snow particles are captured during the 11 min of the SLPIV and
PTV data. The mean snow particle equivalent diameter is measured to be 0.39 mm
with a standard deviation of 0.29 mm, and the average aspect ratio of the
fitted minor and major ellipsoid axis is 0.62. The sample volume for the DIH
measurement is 42 $\textnormal{cm}^{3}$, thus the mean snow particle number
concentration is around 28,460 $\textnormal{m}^{-3}$, and the volume fraction
is $3.8\times 10^{-6}$.

SLPIV/PTV setup
---
| Deployment
---
datasets
| Duration
---
(min)
| Elevation
---
(m)
| FOV size
---
($\textnormal{m}^{2}$)
| Resolution
---
(mm/pixel)
| Camera angle
---
(deg.)
| Camera-to-light
---
distance (m)
SLPIV | 11 | 18.4 | 39.2 $\times$ 22.1 | 20.5 | 19.9 | 52.5
PTV | 11 | 7.9 | 5.3 $\times$ 3.0 | 4.15 | 18.4 | 19.1
Table 1: Summary of key parameters of the SLPIV and PTV measurement setups for the deployment dataset used in the present paper. U | $u_{\textnormal{rms}}$ | $w_{\textnormal{rms}}$ | $R_{b}$ | $L_{\textnormal{OB}}$ | $U_{\tau}$ | $z/L_{\textnormal{OB}}$ | L | $\tau_{L}$ | $\epsilon$ | $\eta$ | $\tau_{\eta}$ | $\lambda$ | $\Rey_{\lambda}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
m/s | m/s | m/s | - | m | m/s | - | m | s | $\textnormal{cm}^{2}/\textnormal{s}^{3}$ | mm | ms | mm | -
5.47 | 1.07 | 0.64 | 0.16 | 1651 | 0.48 | 0.0062 | 6.22 | 5.79 | 355 | 0.51 | 19.4 | 80.7 | 6478
Table 2: Estimated meteorological and turbulence conditions from the sonic
anemometer at $z=10\textnormal{m}$. See the text for the definition of the
symbols.

###  2.2 Atmospheric turbulence conditions

The atmospheric and turbulence conditions during the deployment are determined
using a meteorological tower instrumented with wind velocity, temperature and
humidity sensors at 50 m downstream of the light sheet. Four sonic anemometers
(20 Hz sampling rate, 5.8 cm horizontal measurement path length and 10 cm
vertical measurement path length, CSAT3, Campbell Scientific) are installed at
elevations of 10, 30, 80 and 129 m, and six cup-and-vane anemometers (1 Hz
sampling rate) are installed at elevations of 7, 27, 52, 77, 102 and 126 m.
Note that the measurement uncertainty of the sonic anemometer is $\pm$0.08 m/s
(Toloui et al., 2014), corresponding to 1.4$\%$ of the average wind speed.
Thus, we estimate that the uncertainties for the turbulence properties would
be less than 4$\%$. The key parameters are listed in table 2. The atmospheric
stability is estimated using the bulk Richardson number $R_{b}$ and the Monin-
Obukhov length $L_{\textnormal{OB}}$:

$R_{b}=-|g|\Delta\overline{\theta_{v}}\Delta
z/\left(\overline{\theta_{v}}\left[\left(\Delta V_{N}\right)^{2}+\left(\Delta
V_{W}\right)^{2}\right]\right)$ (1)

$L_{OB}=-U_{\tau}^{3}\overline{\theta_{v}}/\kappa
g\overline{w^{\prime}\theta_{v}^{\prime}}$ (2)

In the equations, g is the gravitational acceleration; $\theta_{v}$ is the
virtual potential temperature; $V_{N}$ and $V_{W}$ are the average wind
velocity components to the North and West respectively measured by the sonic
anemometers; $\kappa$ is the von Kármán constant; $U_{\tau}$ is the shear
velocity estimated from the Reynolds stresses (Stull, 1988), where
$U_{\tau}=\left(\left\langle
V_{N}^{\prime}V_{Z}^{\prime}\right\rangle^{2}+\left\langle
V_{W}^{\prime}V_{Z}^{\prime}\right\rangle^{2}\right)^{1/4}$. The average
velocity differences are calculated from the sonic anemometers at top (129 m)
and bottom (10 m) which yields a height difference $\Delta\textit{z}$ of 119
m. The Monin-Obukhov length and all other turbulence conditions are measured
with the data from the 20 Hz sonic sensor at 10 m. For the duration of the
analyzed dataset, the bulk Richardson number and the Monin-Obukhov length
indicate that the atmospheric boundary layer during the experiment is near-
neutrally stratified (typically for the near neutrally stratified atmospheric
boundary layer, $R_{b}$ ranges from 0 to 0.25 and stability parameter
($z/L_{\textnormal{OB}}$) ranges from 0 to 0.1 (Högström et al., 2002; Stull,
1988)).

The turbulence conditions are estimated using the methods described in Nemes
et al. (2017) and Li et al. (2021). Velocity data from the sonic anemometer at
10 m are used for the flow characterization, consistent with the sample area
elevation ranges of the SLPIV (from 3 m to around 40 m) and PTV (from 5.3 m to
10.6 m) measurements.The integral time scale $\tau_{L}$ and the length scale L
are estimated based on the temporal autocorrelation function $\rho_{uu}$:

$\rho_{uu}(\tau)=\left\langle
u^{\prime}(t)u^{\prime}(t+\tau)\right\rangle/u^{\prime 2}$ (3)

$\tau_{L}=\int_{0}^{T_{0}}\rho_{uu}(\tau)d\tau$ (4)

$L=u_{\textnormal{rms}}\tau_{L}$ (5)

In these equations, t is the variable time, $\tau$ is the time difference, and
$T_{0}$ is the first zero-crossing point the auto-correlation function. The
turbulence dissipation rate $\epsilon$ is estimated using the second-order
structure function of the streamwise velocity component, applying the Taylor
hypothesis to convert the measured time series into spatial velocity
variations:

$D_{11}(\tau)=\left\langle\left[u^{\prime}(t+\tau)-u^{\prime}(t)\right]^{2}\right\rangle$
(6)

$D_{11}(r)=C_{2}\epsilon^{2/3}r^{2/3}$ (7)

With the Kolmogorov prediction for the second-order structure function in the
inertial range (equation 7), where $C_{2}$ is a constant of around 2 in high-
Reynolds number turbulence (Saddoughi and Veeravalli, 1994), we can estimate
the turbulence dissipation rate
$\epsilon=\left(D_{11}/\left(C_{2}r^{2/3}\right)\right)^{3/2}$. Furthermore,
we calculated the Kolmogorov time and length scale,
$\tau_{\eta}=(v/\epsilon)^{1/2}$ and $\eta=\left(v^{3}/\epsilon\right)^{1/4}$,
the Taylor microscale, $\lambda=u_{\textnormal{rms}}(15v/\epsilon)^{1/2}$, and
the Reynolds number, $\Rey_{\lambda}=u_{\mathrm{rms}}\lambda/\nu$, where $\nu$
is the kinematic viscosity of air, and $u_{\mathrm{rms}}$ is the root mean
square (r.m.s.) of the velocity fluctuations $u^{\prime}$.

##  3 Results

###  3.1 Snow particle acceleration and Stokes number

The snow particle acceleration and vertical velocity obtained from PTV
analysis are evaluated in this section. Figure 3a shows the probability
density function (PDF) of the fluctuations of the snow particle acceleration
normalized by their r.m.s. value. The PDF is compared with data from previous
laboratory experiments and numerical simulations of tracers and inertial
particles in isotropic turbulence (Mordant et al., 2004; Bec et al., 2006;
Ayyalasomayajula et al., 2006). In figure 3a, the exponential tail of the in-
plane acceleration PDF curve of snow particles lies in between the curves with
Stokes numbers of 0.16 and 1.01 from Bec et al. (2006), while a comparison of
streamwise acceleration with Ayyalasomayajula et al. (2006), in a similar
boundary layer flow, seem to narrow the range to 0.09-0.15. As discussed in
Nemes et al. (2017), the acceleration kurtosis manifests the tendency of
inertial particles to experience only a portion of the high acceleration
events sustained by fluid parcels. The direct numerical simulation (DNS) by
Ireland et al. (2016) showed that the kurtosis of acceleration becomes
insensitive to the change of Reynolds number with $\Stk>0.1$ and
$\Rey_{\lambda}>398$ (e.g., as $\Rey_{\lambda}$ changes from 398 to 597 at
$\Stk=0.1$, the kurtosis of acceleration increases only 3$\%$, and the change
becomes smaller at higher $\Stk$). Therefore, following the reasoning in our
previous studies (Nemes et al., 2017; Li et al., 2021), we extend the
comparison of the acceleration PDFs to the atmospheric turbulence case with
high $\Rey_{\lambda}$ investigated here, and conservatively estimate the
Stokes number in the range of 0.09-1.01.

Figure 3: (a) PDFs of in-plane snow particle acceleration from PTV
(triangles), compared to $\Stk$ = 0 from Mordant et al. (2004) (dots), and Bec
et al. (2006) ($\Stk$ = 0.16, solid line; $\Stk$ = 0.37, dotted line; $\Stk$ =
1.01, dash-dotted line; $\Stk$ = 2.03, dashed line); (b) comparison of the
measured distribution (solid line) and average (dashed line) of settling
velocity with the estimated range of still-air terminal velocities (grey
region).

With the estimated $\Stk$, the aerodynamic particle response time of the
observed snow particles is predicted to be in the range from 1.7 ms to 19.6
ms, where $\tau_{p}=\Stk\cdot\tau_{\eta}$, leading to a still-air terminal
velocity defined by ${\it W_{p}}=\tau_{p}\cdot{g}$ and estimated between 0.02
m/s and 0.19 m/s. The estimated Stokes number indicates a near critical
condition ($\Stk\sim\textit{O}(1)$), anticipating the occurrence of
preferential concentration (clustering) and sweeping, as well as enhanced
settling velocity. In figure 3b, we compare the vertical velocity distribution
(solid line) from PTV (the average vertical velocity ($\langle w_{p}\rangle$)
of 0.73 m/s is indicated as a dashed line) with the estimated terminal
velocity range accounting for the uncertainty in Stokes number (grey region).
The increase is evident and multifold (around seven times larger on average,
$\langle w_{p}\rangle/\overline{W}_{p}$). This enhancement is consistent with
what has been observed in the previous study by Nemes et al. (2017) (around
three times enhancement on average). However, since our estimated range for
$\Stk$ is closer to the critical value, the observed enhancement here is
higher.

Note that the corresponding particle Reynolds number $\Rey={\langle
w_{p}\rangle}D/\nu$ based on the measured settling velocity and particle size
is $\sim$16.8, implying that a non-Stokesian drag correction is required. Due
to the disk-like shape of the observed snow particles, the Schiller-Neumann
approach based on the Reynolds-corrected sphere drag is not recommended, which
in part justifies the estimation of Stokes number from the acceleration PDF
with no explicit dependency on density, and size. The only option to account
for snow morphology is to use the semi-empirical $\chi$ number approach
proposed by Böhm (1989), corrected by Heymsfield and Westbrook (2010), and
also employed in Nemes et al. (2017). The resulting parameterization leads to
$\chi=997$ and a drag coefficient of $C_{\textnormal{D}}=3.53$, which is not
unusual given the relatively low particle Reynolds number (Westbrook and
Sephton, 2017). It is important to stress that the $\chi$ number accounts for
the snow morphology effects on drag (Garrett et al., 2015; Dunnavan et al.,
2019), not necessarily for the effect of ambient turbulence, which is the main
point of this work.

###  3.2 Preferential distribution of snow particle concentration

Snow particle concentration around vortices in the flow is first evaluated
using the SLPIV data. As shown in an instantaneous flow sample (figure 4a),
the vortices are identified using the swirling strength derived from the
velocity fields of SLPIV as described in detail in section 2.1. Subsequently,
the particle number concentration within and around the vortices is estimated
using the image intensity ($I(x,z,t)$) of the SLPIV data. The estimation of
concentration using image intensity is supported by Raffel et al. (2018) which
shows that the image intensity of PIV is proportional to the concentration of
particles with the same averaged diameters. However, factors such as
stochastic light attenuation by particles within the light sheet and in
between the light sheet and the cameras, as well as the power fluctuation of
the search light might lead to non-linear relationship between the light
intensity and local particle concentration (Kalt et al., 2007; Banko et al.,
2019). Nevertheless, due to relatively low volume fraction, the light
attenuation by particles between the light sheet and the cameras is not
inferred to be dominant as compared to the other two factors. Furthermore, to
minimize the spatial and temporal non-uniformity in background image intensity
due to the decay of light intensity with height and its fluctuation over time,
relative concentration $C^{*}=I(x,z,t)/I(x,z)_{\mathrm{10s,avg}}$ is defined
according to Li et al. (2021), where $I(x,z)_{\mathrm{10s,avg}}$ is an average
of the intensity of images recorded in a 10 s moving window. Figure 4b shows
the relative concentration of the flow sample corresponding to figure 4a. In
particular, we observe that snow concentration is low within the vortex cores.
The phenomenon is considered to be a result of the inertia bias (Maxey, 1987),
i.e., inertial particle trajectories are biased towards the region of low
vorticity. Remarkably, for all three strong vortices (i.e., vortices with high
swirling strength values) highlighted in the figure, including both prograde
and retrograde vortices, the particle concentration is preferentially higher
on the downward flow side of the vortices (referred to as preferential
concentration hereafter).

It is worth noting that the relative concentration map in figure 4b seems to
suggest the clustering of snow particles in the atmospheric surface layer
similar to those studied in our previous work (Li et al., 2021). However, the
clusters in figure 4b are on smaller scales and do not exhibit a clear
preferential orientation in comparison to those presented in Li et al. (2021)
with the January 2019 dataset. The difference is mainly caused by relatively
lower turbulence and snow concentration in the current cases which lead to a
weaker interaction between particles and turbulence. In addition, since the
current study focuses on elucidating the connection between the vortical
structures in the turbulent flow and settling of individual snow particle, the
PTV was designed to have a more focused field of view (3 m $\times$ 5 m,
smaller than the integral scale), limiting our ability to quantify the large-
scale clusters that extend beyond our region of interest.

Figure 4: (a) Swirling strength contour of an instantaneous sample (threshold
0.4 $\textnormal{s}^{-1}$ applied) with relative flow velocity vectors (white
arrows, subtracted by the convective velocity of 3.4 m/s of the prograde
vortex at the center left of the field of view). (b) Corresponding snow
particle concentration colormap ($C^{*}$ is the relative concentration) with
vortices detected from (a). Note that the white dashed lines represent vortex
boundaries, red arrows indicate prograde vortices and blue arrow indicates the
retrograde vortex, black circle and arrow define the coordinate system of the
local vortex.

To substantiate the observation of preferential concentration associated with
presence of vortices in the flow, the ensemble-averaged concentration contours
are calculated for prograde (28,700 prograde vortices identified) and
retrograde vortices (9,700 retrograde vortices), respectively. As shown in
figure 5, the averaged concentration is determined for both the central region
defined as the region within half effective radius
($R_{\mathrm{eff}}=\sqrt{A/\pi}$, where A is the area of the vortex region
detected through the swirling strength criterion) from the center of the
vortex, and in the rest of proximity defined as the region from
$0.7R_{\textnormal{eff}}$ to $1.4R_{\textnormal{eff}}$ from the center. The
latter is divided into twelve angular sectors with angle of $\pi/6$. For both
prograde (figure 5a) and retrograde vortices (figure 5b), it is observed that
the particle concentration is higher on the downward side, than that on the
upward side and in the center of vortices. For reference, the averaged
concentration distribution in the background (figure 5c) does not exhibit any
appreciable preferential concentration. For each prograde/retrograde vortex,
the reference background is defined as a circular region that has a radius
equal to the $R_{\textnormal{eff}}$ of the prograde/retrograde vortex with a
center at a randomly selected location in the SLPIV measurement field of view.
In addition, the averaged particle concentration is relatively higher at the
bottom of the downward side for retrograde vortices, possibly due to gravity
and downward fluid motion causing a stronger preferential concentration
effect. But for prograde vortices, relative concentration is more uniformly
distributed on the downward side. Such a discrepancy can be explained by the
different organizations of prograde and retrograde vortices in the atmospheric
surface layer. Specifically, the prograde vortices tend to form in packets
(Christensen and Adrian, 2001), predominantly located in the internal shear
layers in atmospheric surface layer, the interaction between snow particles
around a certain prograde vortex with adjacent vortices may smooth out the
particle concentration on the downward side.

It is well known that preferential concentration in turbulent flows is a
multi-scale phenomenon. As shown in (Baker et al., 2017), particles start
clustering at the Kolmogorov scale when they preferentially sample the high
strain regions in the flow. As the clusters yield larger response time than
that of individual particles, they can subsequently interact with larger scale
flow structures and grow in size up to the integral scale (see examples in Li
et al. (2021)). However, due to the limit of the SLPIV resolution, we can only
resolve vortices above a certain size ($\sim$ 66 cm) in our measurements.
Therefore, the preferential concentration statistics shown in figure 5 are
captured only by sampling large, energetic vortices leaving a signature
several times larger than the Taylor microscale in our coarsely resolved
turbulent flow. Nevertheless, the estimated Stokes number in section 3.1 (with
the upper range close to the critical condition) suggests strong interaction
between the particles and flow structure at the Kolmogorov scales, causing
preferential concentration and clustering that cannot be resolved with the
current SLPIV measurement. We thus acknowledge resolving a limited range of
scale contributing to particle clustering, but capturing the key mechanism
governing preferential sweeping at the resolved scale.

Figure 5: Comparison between the ensemble-averaged concentration from the
SLPIV dataset in the central region (dashed lines), defined as the region
within half $R_{\mathrm{eff}}$ from the vortex center where $R_{\mathrm{eff}}$
is the effective radius of the vortices, and in the proximity (solid lines),
defined as the region from $0.7R_{\mathrm{eff}}$ to $1.4R_{\mathrm{eff}}$ from
the center spanning a circumferential angle range of $\pi/6$, for (a) prograde
vortices, (b) retrograde vortices, and (c) reference background. For each
prograde/retrograde vortex, the reference background is defined as a circular
region that has a radius equal to the $R_{\mathrm{eff}}$ of the
prograde/retrograde vortex with a center at a randomly selected location in
the SLPIV measurement field of view.

Figure 6: (a) A schematic showing the box counting process to determine the
change of snow particle concentration due to the presence of a vortex (gray
solid and dashed contours) in the PTV sample area (black solid line). The
vortex is determined from the corresponding SLPIV measurement at the same time
instant as the PTV. A fixed counting box (black dashed line, 100 $\times$ 200
$\textnormal{pixel}^{2}$) is selected at the center of the PTV area. The
relative displacement (red arrow) between the vortex center
($x_{\mathrm{vor}}$, $z_{\mathrm{vor}}$) and box center ($x_{\mathrm{box}}$,
$z_{\mathrm{box}}$) is defined by equation 8. (b) The ensemble-averaged
particle concentration ($\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$) at different relative displacements
with respect to the vortex center. In total, 550 prograde (red line) and 380
retrograde (blue line) vortices are selected for the ensemble average,
respectively. The bin of relative displacement ($\Delta d$) used in ensemble
average has a width of $0.6R_{\mathrm{eff}}$ and spaced $0.4R_{\mathrm{eff}}$
from adjacent bins ($\sim$33$\%$ overlap). The bin size and spacing are
determined to ensure sufficient statistical convergence of the data.

Furthermore, we examine the images from the PTV measurements in which
individual snow particles can be counted and tracked within and around the
vortices (see movies 1 and 2 in the supplementary material), supporting the
observation of preferential concentration based on the change in the intensity
of images from SLPIV. We first select the vortices determined from the
corresponding SLPIV measurements at the same time instant as the PTV with
overlapping field of view. Most of these vortices have an equivalent diameter
$\gtrsim$1.3 m and are only partially inside the PTV images as they move
across the small (relative to SLPIV) sample area of PTV. As a result, the
previous method (that for SLPIV) for estimating particle concentration cannot
be applied. Instead, as illustrated in figure 6a, a box counting method is
used to determine the change of snow particle concentration due to the
presence of vortices in the PTV sample area (note that the uncertainty of
particle concentration associated with the changing light sheet thickness
within the PTV field is estimated to be less than 9$\%$). A fixed counting box
(100 $\times$ 200 $\textnormal{pixel}^{2}$) is selected at the center of the
PTV area. When the vortex appears in the sample area of PTV, the total number
of snow particles in the box is counted. The relative displacement (d) between
the vortex center ($x_{\mathrm{vor}}$, $z_{\mathrm{vor}}$) and box center
($x_{\mathrm{box}}$, $z_{\mathrm{box}}$) is calculated as:

$d=\frac{x_{\mathrm{box}}-x_{\mathrm{vor}}}{\left|x_{\mathrm{box}}-x_{\mathrm{vor}}\right|}\sqrt{\left(x_{\mathrm{box}}-x_{\mathrm{vor}}\right)^{2}+\left(z_{\mathrm{box}}-z_{\mathrm{vor}}\right)^{2}}$
(8)

Particularly, the sign of the relative displacement indicates which side the
particles are located with respect to the moving vortex. Similar to that in
the SLPIV data processing, to account for the variation of snow particle
concentration in the background, the relative concentration for PTV data
($C^{*}$) is estimated using total particle number counts in the box at each
time instant divided by the averaged total particle counts for the time
duration during which each vortex is present in the sample area of PTV.
Subsequently, the ensemble averaged $C^{*}$ in a bin of a width of
$0.6R_{\textnormal{eff}}$ ($\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$) at different relative displacements
is determined for prograde (550 in total) and retrograde (380 in total)
vortices, respectively (figure 6b). As the figure shows, for prograde
vortices, the $\left\langle C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$ is
evidently higher on the downward side of the vortices (positive $d$). For
retrograde vortices, the $\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$ also yields larger values on the
downward side (negative $d$), though the difference in $\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$ between two sides appears to be
smaller than that observed for prograde vortices. To determine the
significance of the difference in $\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$, we conduct a t-test to the
concentration distribution at the two sides of vortices: the particle
concentration on the downward side is in general 13$\%$-22$\%$ higher (at
95$\%$ confidence level) than that at the upward side for both prograde and
retrograde vortices. These trends (see also the supplementary material)
provide further support of the preferential concentration. Nevertheless, due
to the limited number of vortices that can be simultaneously detected by SLPIV
and PTV, the standard deviation of the $\left\langle
C^{*}\right\rangle_{0.6R_{\mathrm{eff}}}$ presented in the current analysis is
considerable. In addition, we would like to point out that the statistical
signature of preferential concentration observed in PTV is likely to be
underestimated potentially due to the relatively large counting box size used
in the analysis in comparison to the size of vortices.

###  3.3 Enhanced settling velocity due to preferential sweeping

In this subsection, we investigate how the settling velocity of snow particles
is influenced by the presence of the vortical structures in the flow.
Consistent with the method presented in the last section, the vortices are
first identified using SLPIV data. Subsequently, the particle trajectories
around these identified vortices are extracted from PTV for the following
analysis. For the prograde vortex (figure 7a), the settling velocity of
particles increases when moving toward the downward side (right side in the
sample) of the vortex. Similarly, for the retrograde vortex (figure 7c), the
settling of snow particles slows down; some particles are even lifted upward,
as they travel to the upward side (right side in the sample) of the vortex.
Both cases illustrate clearly higher settling velocities of the snow particles
situating on the downward side of vortices.

Figure 7: (a, c) Samples of snow particle trajectories around a prograde
vortex (a) and a retrograde vortex (c) (see also movie 1 and 2 uploaded in the
supplementary material). Black dashed lines represent vortex boundaries, and
trajectories are colored based on the dimensionless times $t^{*}$, defined as
the difference between the timestamps of snow particles and that of a selected
vortex at one time instant normalized by the Kolmogorov time scale. (b, d)
Histograms of the differences between the settling velocities on the downward
sides and upward sides for (b) prograde and (d) retrograde vortices (A total
of 4300 prograde and 1700 retrograde vortices are identified), normalized by
the average vertical velocity $\langle w_{p}\rangle$ of snow particles tracked
with zeros marked by the dash-dotted lines, and the dashed lines represent the
average velocity difference in each case. The bin size is chosen to be one
fifth of $\langle w_{p}\rangle$. The histograms are compared with the Gaussian
distributions marked by the solid lines.

To further substantiate these observations, the average vertical accelerations
conditioned on the downward ($\overline{a_{p,-}}$) and the upward
($\overline{a_{p,+}}$) sides of prograde and retrograde vortices are
calculated, and the histograms of settling velocity difference between the
downward ($\overline{w_{p,-}}$) and upward ($\overline{w_{p,+}}$) sides of
vortices, i.e., $\overline{w_{p,-}}-\overline{w_{p,+}}$, normalized by the
ensemble average snow particle vertical velocity ($\langle
w_{p}\rangle=0.73\mathrm{~{}m}/\mathrm{s}$) are presented for all prograde
(figure 7b) and retrograde (figure 7d) vortices, respectively. The
$\overline{a_{p,-}}$, $\overline{a_{p,+}}$, $\overline{w_{p,-}}$ and
$\overline{w_{p,+}}$ are calculated by averaging the vertical velocities of
particles within boxes ranging from $0.5R_{\mathrm{eff}}$ to
$1.5R_{\mathrm{eff}}$ to the center of vortices in the x direction and
covering the whole diameter ($2R_{\mathrm{eff}}$) in the vertical direction at
the two sides of vortices. Specifically, for prograde vortices, the two
conditionally averaged accelerations are $\overline{a_{p,-}}=-0.16\pm
2.20\mathrm{m/s^{2}}$ on the downward side and $\overline{a_{p,+}}=0.0065\pm
2.54\mathrm{m/s^{2}}$ on the upward side. While for retrograde vortices,
$\overline{a_{p,-}}=-0.33\pm 2.93\mathrm{m/s^{2}}$ and
$\overline{a_{p,+}}=0.12\pm 2.09\mathrm{m/s^{2}}$. These conditional averages
support the fact that particles on the downward side of vortices would
accelerate with the flow and particles falling on the upward side would
decelerate, or even be lifted up. Note that the variability in acceleration
(e.g. the standard deviation) is much higher than that for the settling
velocity due to the higher order derivative in the acceleration calculation.
Moreover, the settling velocity differences display a near Gaussian
distribution with the mean value above zero. As compared to Gaussian
distribution with the same mean value and standard deviation, the PDFs of
settling velocity difference exhibit higher probability near the mean values
and heavier tails on the right side for both prograde and retrograde vortices,
indicating statistically higher settling velocities on the downward side of
vortices. Specifically, about 78$\%$ of the prograde vortices yield a higher
settling velocity on the downward side with an average settling velocity
differences of 0.56$\langle w_{p}\rangle$, and the proportion is about 73$\%$
for the retrograde vortices with an average settling velocity difference of
0.48$\langle w_{p}\rangle$. Note that the total number of vortices identified
in figure 7 is larger than that for particle concentration in figure 6. It is
because we identify individual vortex from the PTV field for settling velocity
analysis, while vortices selected for preferential concentration in figure 6
are tracked across the region of interest.

However, not every single occasion is observed with higher settling velocity
at the downward side of vortices, due to the fact that the vortices detected
in the current study are planar projections of highly complex three-
dimensional vortices. In addition, the snow particles interacting with one
vortex can also be affected by adjacent vortices in atmosphere, which are
usually weaker and less appreciable in the SLPIV data as compared to the
vortices analyzed. Nevertheless, we observe an increasing percentage of events
showing such a trend when we sample vortices with higher swirling strength and
the percentage for prograde vortices is consistently higher than that for
retrograde vortices. Such discrepancies are likely due to the difference in
self-organization characteristics for prograde and retrograde vortices (i.e.,
prograde vortices are found to be predominantly located in the proximity of
internal shear layers as shown in the field PIV measurement in the atmospheric
surface layer by Heisel et al. (2018)). Therefore, prograde vortices could
have a cumulative and stronger effect on the enhanced settling of nearby snow
particles compare to retrograde vortices. Nevertheless, with the observed
preferential concentration from section 3.2 and statistically higher settling
velocity on the downward side of the vortices, we conclude that under near
critical conditions ($\Stk\sim\textit{O}(1)$) observed in our field
measurements, preferential sweeping plays a significant role in controlling
the settling velocity of hydrometeors in the atmospheric turbulence.

##  4 Conclusions and discussion

In this paper, we present the first field study of snow settling dynamics
based on simultaneous measurements of the atmospheric flow field using a
super-large-scale particle image velocimetry (SLPIV) and snow particle
trajectories using particle tracking velocimetry (PTV) within the SLPIV sample
area. Our results reveal the direct linkage among the coherent vortex
structures in the atmospheric turbulence, and the concentration distribution
and settling dynamics of snow particles in the field. Specifically, we observe
a settling velocity enhancement of around seven folds on average compared to
the estimated still-air terminal velocity. This value is larger than what has
been observed in our previous field study (Nemes et al., 2017), potentially
due to the fact that the Stokes number associated with the snow particles in
the present deployment is closer to the critical condition
($\Stk\sim\textit{O}(1)$). Using the SLPIV, we are able to detect the strong
vortices present in the atmospheric turbulent flow, and we show that the snow
concentration (represented by the variation in the particle image intensity in
SLPIV) is preferentially higher on the downward side for both prograde and
retrograde vortices in the flow. This observation is further substantiated by
counting individual snow particles around the vortices using PTV data. The
result indicates an average of 18$\%$ higher concentration on the downward
side of the detected strong vortices present in our study. In addition, the
samples of snow particle trajectories around vortices from PTV demonstrate
that snow particles accelerate as they move toward the downward side of
vortices and decelerate or even are lifted upward when traveling to the upward
side. Based on the histograms of snow settling velocity from PTV, we show that
the snow particles on the downward side of vortices yield a statistically
higher settling velocity than that at the upward side, with an average
difference of about 52$\%$ of the mean settling velocity.

Our results provide direct evidence and underlying mechanisms for the
preferential concentration and preferential sweeping during snow particles
settling in the atmospheric surface layer. While the presented results focus
on the quantification of particle-turbulence interaction mechanisms at the
scale of individual vortices, we recognize that atmospheric turbulence in the
near-neutrally stratified boundary layer in our study contains vortices over a
broad range of scales and intensity. Therefore, we hypothesize that the
preferential paths of highly concentrated particles along layers of vortices
are likely to produce a cumulative effect on enhancing the settling of the
snow particles. This conceptual framework has not been considered in current
snow settling models. Under this framework, our observed results of
preferential concentration and preferential sweeping may be critical to inform
a stochastic model to reproduce the observed fall speed under specific
micrometeorological conditions. By incorporating the framework, it may be
possible to model the settling of snow particles starting from the still-air
terminal velocity with specific drag coefficient and sampling a population of
vortices consistent with the inertial range in the air column. Under the
assumption of accelerated preferential paths, each snow particle has a large
probability to sample a number of vortices, and the effects of both the drag
force and vortex flow sweeping will progressively enhance the settling
velocity to a certain value. With this perspective, it is reasonable to
hypothesize that the multifold increase in observed settling velocity compared
with the still-air terminal velocity, as well as its large variability, is a
result of the cumulative effect of particles settling through the turbulent
air column occupied by many vortices. In addition, the disparity in settling
velocities observed on the downward and upward sides of vortices of 0.4 m/s
might also be directly affected by the azimuthal velocity of the flow around
vortices, manifested in vertical velocity fluctuations that are more readily
available from the measured flow field.

In the end, we acknowledge the uncertainties involved in our concentration
measurement using box counting from the PTV data. Such uncertainties are
largely caused by the limited data of synchronized SLPIV and PTV measurements
from our field deployment. Nevertheless, the main observations related to the
snow particle concentration and settling dynamics present in our study are
still statistically significant. A second relevant uncertainty is in the
estimate of the snow particle velocity in still air. While the Stokes number
range is conservatively defined, it would still be important to provide a
direct estimate of the snow particle density, combining measurements of single
particle volume with the weighing of particle ensemble in time. We expect that
more converged trends can be obtained with an increasing number of
deployments. In particular, we expect to extend our current field PTV to a
three-dimensional imaging system using multiple cameras. Such upgraded system
will allow us to quantify, more accurately, the particle concentration and
distribution, particle setting kinematics (e.g., curvature of the
trajectories) and dynamics (e.g., acceleration, inertial response) associated
with the presence of three-dimensional vortex structures in the atmospheric
turbulence.

## References

  * (1)
  * Abraham and Hong (2020) Abraham, A. and Hong, J.  2020 Dynamic wake modulation induced by utility-scale wind turbine operation, Appl. Energy, 257, 114003.
  * Adrian et al. (2000) Adrian, R. J., Christensen, K. T. and Liu, Z. C.  2000 Analysis and interpretation of instantaneous turbulent velocity fields, Exp. Fluids, 29 (3), 275-290.
  * Aliseda et al. (2002) Aliseda, A., Cartellier, A., Hainaux, F. and Lasheras, J. C.  2002 Effect of preferential concentration on the settling velocity of heavy particles in homogeneous isotropic turbulence, J. Fluid Mech., 468, 77-105.
  * Ayyalasomayajula et al. (2006) Ayyalasomayajula, S., Gylfason, A., Collins, L. R., Bodenschatz, E. and Warhaft, Z. 2006 Lagrangian measurements of inertial particle accelerations in grid generated wind tunnel turbulence, Phys. Rev. Lett., 97 (14), 144507.
  * Baker et al. (2017) Baker, L., Frankel, A., Mani, A. and Coletti, F.  2017 Coherent clusters of inertial particles in homogeneous turbulence, J. Fluid Mech., 833, 364-398.
  * Balachandar and Eaton (2010) Balachandar, S. and Eaton, J. K.  2010 Turbulent dispersed multiphase flow, Annu. Rev. Fluid Mech., 42, 111-133.
  * Banko et al. (2019) Banko, A. J., Villafañe, L., Kim, J. H., Esmaily, M. and Eaton, J. K. 2019 Stochastic modeling of direct radiation transmission in particle-laden turbulent flow, J. Quant. Spectrosc. Radiat. Transf., 226, 1-18.
  * Bec et al. (2006) Bec, J., Biferale, L., Boffetta, G., Celani, A., Cencini, M., Lanotte, A., Musacchio, S. and Toschi, F. 2006 Acceleration statistics of heavy particles in turbulence, J. Fluid Mech., 550, 349-358.
  * Böhm (1989) Böhm, H. P. 1989 A general equation for the terminal fall speed of solid hydrometeors, J. Atmos. Sci., 46 (15), 2419-2427.
  * Christensen and Adrian (2001) Christensen, K. T. and Adrian, R. J.  2001 Statistical evidence of hairpin vortex packets in wall turbulence, J. Fluid Mech., 431, 433-443.
  * Crocker and Grier (1996) Crocker, J. C. and Grier, D. G. 1996 Methods of digital video microscopy for colloidal studies, J. Colloid Interface Sci., 179 (1), 298-310.
  * Dasari et al. (2019) Dasari, T., Wu, Y., Liu, Y. and Hong, J. 2019 Near-wake behaviour of a utility-scale wind turbine, J. Fluid Mech., 859, 204-246.
  * Dunnavan et al. (2019) Dunnavan, E. L., Jiang, Z., Harrington, J. Y., Verlinde, J., Fitch, K. and Garrett, T. J. 2019 The shape and density evolution of snow aggregates, J. Atmos. Sci., 76 (12), 3919-3940.
  * Durán et al. (2011) Durán, O., Claudin, P. and Andreotti, B. 2011 On aeolian transport: Grain-scale interactions, dynamical mechanisms and scaling laws, Aeolian Res., 3 (3), 243-270.
  * Elghobashi and Truesdell (1992) Elghobashi, S. and Truesdell, G. C. 1992 Direct simulation of particle dispersion in a decaying isotropic turbulence, J. Fluid Mech., 242, 655-700.
  * Falkinhoff et al. (2020) Falkinhoff, F., Obligado, M., Bourgoin, M. and Mininni, P. D. 2020 Preferential concentration of free-falling heavy particles in turbulence, Phys. Rev. Lett., 125 (6), 064504.
  * Ferrante and Elghobashi (2003) Ferrante, A. and Elghobashi, S. 2003 On the physical mechanisms of two-way coupling in particle-laden isotropic turbulence, Phys. Fluids, 15 (2), 315-329.
  * Garrett et al. (2015) Garrett, T. J., Yuter, S. E., Fallgatter, C., Shkurko, K., Rhodes, S. R. and Endries, J. L. 2015 Orientations and aspect ratios of falling snow, Geophys. Res. Lett., 42 (11), 4617-4622.
  * Good et al. (2014) Good, G. H., Ireland, P. J., Bewley, G. P., Bodenschatz, E., Collins, L. R. and Warhaft, Z. 2014 Settling regimes of inertial particles in isotropic turbulence, J. Fluid Mech., 759 (R3).
  * Heisel et al. (2018) Heisel, M., Dasari, T., Liu, Y., Hong, J., Coletti, F. and Guala, M. 2018 The spatial structure of the logarithmic region in very-high-Reynolds-number rough wall turbulent boundary layers, J. Fluid Mech., 857, 704-747.
  * Heymsfield and Westbrook (2010) Heymsfield, A. J. and Westbrook, C. D. 2010 Advances in the estimation of ice particle fall speeds using laboratory and field measurements, J. Atmos. Sci., 67 (8), 2469-2482.
  * Högström et al. (2002) Högström, U., Hunt, J. C. R. and Smedman, A. S. 2002 Theory and measurements for turbulence spectra and variances in the atmospheric neutral surface layer, Bound.-Layer Meteorol., 103 (1), 101-124.
  * Hong et al. (2014) Hong, J., Toloui, M., Chamorro, L.P., Guala, M., Howard, K., Riley, S., Tucker, J. and Sotiropoulos, F. 2014 Natural snowfall reveals large-scale flow structures in the wake of a 2.5-MW wind turbine, Nat. Commun., 5, 4216.
  * Ireland et al. (2016) Ireland, P. J., Bragg, A. D. and Collins, L. R. 2016 The effect of Reynolds number on inertial particle dynamics in isotropic turbulence. Part 1. Simulations without gravitational effects, J. Fluid Mech., 796, 617-658.
  * Kalt et al. (2007) Kalt, P. A., Birzer, C. H. and Nathan, G. J. 2007 Corrections to facilitate planar imaging of particle concentration in particle-laden flows using Mie scattering, Part 1: Collimated laser sheets., Appl. Opt., 46 (23), 5823-5834.
  * Li et al. (2021) Li, C., Lim, K., Berk, T., Abraham, A., Heisel, M., Guala, M., Coletti, F. and Hong, J. 2021 Settling and clustering of snow particles in atmospheric turbulence, J. Fluid Mech., 912 (A49).
  * Mallery et al. (2020) Mallery, K., Shao, S. and Hong, J. 2020 Dense particle tracking using a learned predictive model, Exp. Fluids, 61 (10), 1-14.
  * Maxey (1987) Maxey, M. R. 1987 The gravitational settling of aerosol particles in homogeneous turbulence and random flow fields., J. Fluid Mech., 174, 441-465.
  * Maxey and Riley (1983) Maxey, M. R. and Riley, J. J. 2020 Equation of motion for a small rigid sphere in a nonuniform flow, Phys. Fluids, 26 (4), 883-889.
  * Mordant et al. (2004) Mordant, N., Crawford, A. M. and Bodenschatz, E. 2004 Experimental Lagrangian acceleration probability density function measurement, Physica D, 193 (1-4), 245-251.
  * Nemes et al. (2017) Nemes, A., Dasari, T., Hong, J., Guala, M. and Coletti, F. 2017 Snowflakes in the atmospheric surface layer: observation of particle–turbulence dynamics, J. Fluid Mech., 814, 592-613.
  * Nielsen (1993) Nielsen, P. 1993 Turbulence effects on the settling of suspended particles, J. Sediment. Res., 63 (5), 835-838.
  * Ouellette et al. (2006) Ouellette, N. T., Xu, H. and Bodenschatz, E.  2006 A quantitative study of three-dimensional Lagrangian particle tracking algorithms, Exp. Fluids, 40 (2), 301-313.
  * Petersen et al. (2019) Petersen, A. J., Baker, L. and Coletti, F. 2019 Experimental study of inertial particles clustering and settling in homogeneous turbulence, J. Fluid Mech., 864, 925-970.
  * Raffel et al. (2018) Raffel, M., Willert, C. E., Scarano, F., Kähler, C. J., Wereley, S. T. and Kompenhans, J. 2018 Particle Image Velocimetry: A Practical Guide. Springer.
  * Rosa et al. (2016) Rosa, B., Parishani, H., Ayala, O. and Wang, L. P. 2016 Settling velocity of small inertial particles in homogeneous isotropic turbulence from high-resolution DNS, Int. J. Multiph. Flow, 83, 217-231.
  * Saddoughi and Veeravalli (1994) Saddoughi, S. G. and Veeravalli, S. V. 1994 Local isotropy in turbulent boundary layers at high Reynolds number, J. Fluid Mech., 268, 333-372.
  * Shaw (2003) Shaw, R. A. 2003 Particle-turbulence interactions in atmospheric clouds, Annu. Rev. Fluid Mech., 35 (1), 183-227.
  * Stull (1988) Stull, R. B. 1988 An introduction to boundary layer meteorology. Kluwer Academic Publishers.
  * Toloui et al. (2014) Toloui, M., Riley, S., Hong, J., Howard, K., Chamorro, L. P., Guala, M. and Tucker, J. 2014 Measurement of atmospheric boundary layer based on super-large-scale particle image velocimetry using natural snowfall, Exp. Fluids, 55 (5), 1-14.
  * Tom and Bragg (2019) Tom, J. and Bragg, A. 2019 Multiscale preferential sweeping of particles settling in turbulence, J. Fluid Mech., 872, 995-995.
  * Tooby et al. (1977) Tooby, P. F., Wick, G. L. and Isaacs, J. D. 1977 The motion of a small sphere in a rotating velocity field: a possible mechanism for suspending particles in turbulence, J. Geophys. Res., 82 (15), 2096-2100.
  * Tropea et al. (2007) Tropea, C., Yarin, A. L. and Foss, J. F. 2007 Springer handbook of experimental fluid mechanics. Berlin: Springer.
  * Vaillancourt and Yau (2000) Vaillancourt, P. A. and Yau, M. K. 2000 Review of particle-turbulence interactions and consequences for cloud physics, Bull. Am. Meteorol. Soc., 81 (2), 285-298.
  * Wang and Maxey (1993) Wang, L. P. and Maxey, M. R. 1993 Settling velocity and concentration distribution of heavy particles in homogeneous isotropic turbulence, J. Fluid Mech., 256, 27-68.
  * Westbrook and Sephton (2017) Westbrook, C. D. and Sephton, E. K. 2017 Using 3‐D‐printed analogues to investigate the fall speeds and orientations of complex ice particles, Geophys. Res. Lett., 44 (15), 7994-8001.
  * Yang and Lei (1998) Yang, C. Y. and Lei, U. 1998 The role of the turbulent scales in the settling velocity of heavy particles in homogeneous isotropic turbulence, J. Fluid Mech., 371, 179-205.
  * Zeugin et al. (2020) Zeugin, T., Krol, Q., Fouxon, I. and Holzner, M. 2020 Sedimentation of snow particles in still air in stokes regime, Geophys. Res. Lett., 47 (15), e2020GL087832.
  * Zhou et al. (1999) Zhou, J., Adrian, R. J., Balachandar, S. and Kendall, T. M. 1999 Mechanisms for generating coherent packets of hairpin vortices in channel flow, J. Fluid Mech., 387, 353-396.


# A Data-Driven Analysis of
Robust Automatic Piano Transcription

Drew Edwards, Simon Dixon, Emmanouil Benetos, Akira Maezawa, Yuta Kusaka DE,
SD, and EB are with Queen Mary University of London (QMUL), London, E1 4NS, UK
(a.c.edwards@qmul.ac.uk); AM and YK are with Yamaha Corporation, Hamamatsu,
Shizuoka 430-8650, Japan. DE is a PhD student at the UKRI CDT in AI and Music,
supported by UKRI (EP/S022694/1), QMUL, and Yamaha. EB is supported by
RAEng/Leverhulme Trust Research Fellowship (LTRF2223-19-106). This project
used the Tier 2 HPC facility JADE2, funded by EPSRC (EP/T022205/1). For the
purpose of open access, the authors have applied a Creative Commons
Attribution (CC BY) license to any Accepted Manuscript version arising.

###### Abstract

Algorithms for automatic piano transcription have improved dramatically in
recent years due to new datasets and modeling techniques. Recent developments
have focused primarily on adapting new neural network architectures, such as
the Transformer and Perceiver, in order to yield more accurate systems. In
this work, we study transcription systems from the perspective of their
training data. By measuring their performance on out-of-distribution annotated
piano data, we show how these models can severely overfit to acoustic
properties of the training data. We create a new set of audio for the MAESTRO
dataset, captured automatically in a professional studio recording environment
via Yamaha Disklavier playback. Using various data augmentation techniques
when training with the original and re-performed versions of the MAESTRO
dataset, we achieve state-of-the-art note-onset accuracy of 88.4 F1-score on
the MAPS dataset, without seeing any of its training data. We subsequently
analyze these data augmentation techniques in a series of ablation studies to
better understand their influence on the resulting models.

###### Index Terms:

piano transcription, data augmentation

##  I Introduction

Automatic music transcription (AMT) is the signal processing task of
converting an audio recording of a musical performance into symbolic form. The
target output could be a low-level physical description of each note played by
any instrument in the recording, such as MIDI, or sheet music with
representations of musical features such as meter, time signature, rhythm,
instrumentation, and dynamics.

In this work, we focus on AMT applied to piano, or automatic piano
transcription (APT). APT has improved dramatically over the past decade. The
two largest relative improvements in state-of-the-art accuracy were from a new
architecture to jointly predict onsets and frames [1] and from the release of
the MAESTRO dataset [2]. The Onsets and Frames architecture yielded a more
than fifty percent increase in the note-onset F1-score, when trained and
evaluated on the MAPS dataset [3]. The MAESTRO dataset yielded an additional
five percent relative improvement on MAPS test evaluation in a zero-shot
setting, along with becoming the de-facto standard in training APT systems.
Since then, the research has largely focused on new model architectures and
adapting these systems for multi-instrument transcription [4, 5, 6, 7].

However, most of this research ignores out-of-distribution111In the context of
this research, “out-of-distribution” refers to solo piano recordings drawn
from a different distribution than the training data, e.g. with different
pianos, recording conditions, pieces and musical styles. performance of their
trained models. The emphasis has been on state-of-the-art results on held out
test sets, despite the well-documented phenomenon of overfitting to acoustic
conditions shared between test and training sets [8]. This lack of focus on
generalization can lead to models which exhibit shortcut learning [9], and
thus the reported evaluation metrics can overstate the true performance of
these systems. Furthermore, publications may forgo techniques to make their
models more robust and generalizable, as this has been demonstrated to
slightly hurt their test set results. This is in contrast to computer vision
literature where data augmentation has become standard practice in training
deep learning models [10].

In this research, we analyze the robustness of APT with a data-driven
methodology. We demonstrate how existing systems can overfit to the acoustic
conditions of their training data, which can hamper generalization. We then
perform a series of experiments to analyze the impact of data augmentation on
out-of-distribution transcription performance. Our core contributions are the
following:

  * •

a detailed analysis of data augmentation techniques applied to the task of
APT;

  * •

an automatically re-performed version of all $\approx 200$ hours of the
MAESTRO dataset captured in a studio recording environment via Yamaha
Disklavier technology;

  * •

showing how a combination of data augmentation and increased timbral diversity
results in state-of-the-art (SOTA) performance on the MAPS dataset, and making
the model weights publicly available.

Figure 1: Illustration of the data augmentation pipeline employed. The process
flows from left to right, starting with a Random 7-Band EQ, where each band’s
gain is randomly adjusted within a range of -10dB to 5dB. Next, random
background noise is added, selected from a set of 65 different 10-second
clips, with the signal-to-noise ratio meticulously maintained between 17.5 and
25dB. Then we apply a random pitch shift, varying between -10 and 10 cents, to
account for slight tuning discrepancies in production data. Another Random
7-Band EQ is then applied. The final stage of the pipeline introduces reverb,
utilizing one of 14 unique impulse responses, to simulate different spatial
acoustic characteristics. Each component of this pipeline is applied
independently with a probability of 0.5, ensuring a rich and diverse set of
augmentations.

##  II Related Work

Deep learning models have become the current SOTA in AMT. The first end-to-end
neural network approach for piano transcription was by Sigtia et al. [11].
This was soon improved by Hawthorne et al. [1] by formulating a dual-objective
of predicting both note onsets and frame activity. Recognizing that their
approach was limited by data, the same authors introduced the MAESTRO dataset
[2]: 200 hours of virtuoso solo piano performances with finely-aligned MIDI
labels. Before MAESTRO, the largest available dataset for piano transcription
was MAPS [3], which has 18 hours of audio. Soon after the MAESTRO publication,
Kong et al. [4] introduced a novel technique to regress on the onset times,
allowing higher resolution than the spectrogram frame width. Subsequent
attempts have introduced the Transformer [12, 5] and Perceiver [6]
architectures for piano and multi-instrument transcription. The current SOTA
is from the recent work of Toyama et al., which relies upon the self-attention
mechanism of Transformers to capture long-term dependencies in the frequency
and time axes, yielding a note-level transcription F1 score of 97.4 evaluated
on the MAESTRO test set [7].

The main technique we explore for regularizing deep transcription models is
data augmentation. There has been much work on using data augmentation for
audio based learning tasks [13, 14, 2, 15]. Thickstun et al. [16] explore
label-preserving pitch-shift transformations. Simon et al. [17] use mixtures
of monophonic transcription examples to pre-train their multi-instrument
polyphonic transcription transformer. Lu et al. [6] use cross-dataset mixtures
for training, along with label-preserving pitch shifts. The audiomentations
library222https://github.com/iver56/audiomentations provides a suite of audio-
based data augmentations to apply on the fly during machine learning training,
which we make use of in this work.

###  II-A Robustness of Piano Transcription Systems

Despite the great progress on this task, nearly all of the results mentioned
above fail to report out-of-distribution evaluation metrics. The one exception
is the MAESTRO paper [2], which emphasises the importance of using data
augmentation when evaluating on the MAPS test set. A very promising direction
in robust transcription research comes from Maman and Bermano [18]: the
authors describe a technique to automatically align scores with live
performances. Using only synthetic data and automatically aligned scores, they
achieve the current SOTA on the MAPS test set and report very encouraging
cross-dataset performance.

Many factors seem to contribute to a more robust piano transcription system.
First, certain data augmentation techniques can improve generalization. For
example, models can overfit to the tuning of pianos in their training data,
and small amounts of random pitch shifting can mitigate this effect. Data
augmentation is the core technique explored in this publication and is
explored in detail in Section IV. Second, timbral and acoustic diversity seems
to improve generalization. This is likely the explanation for the very robust
performance of Maman and Bermano, since they train on MusicNet [19], which
includes audio from a wide variety of sources and includes multi-instrumental
performances. Robustness to noisy training data is another factor: Kong et al.
[4] demonstrate that their triangulation technique permits learning when note-
onset labels are off by 50ms, while the architecture from Onsets and Frames is
susceptible to much lower accuracy in this regime.

##  III Studio MAESTRO Dataset

The audio in the MAESTRO dataset contains recordings from ten years of piano
competitions. All pieces are performed by a human on a Yamaha Disklavier piano
in a concert hall with an audience present. While the instrument used in the
various years may be the same, the tuning will certainly vary year to year.
Furthermore the exact placement of microphones will not be identical. On one
hand, this is desirable: the data is very realistic and contains useful noise
and variations for better generalization. On the other hand, we may seek to
remove these differences to better understand the influence of data
augmentation and the extent of overfitting that state-of-the-art transcription
models exhibit.

To that end, we create a data collection apparatus to automate the recording
of acoustic piano performances. The equipment used consists of a Yamaha C6X
grand piano manufactured in 2008 with a Disklavier ENSPIRE PRO player piano
controller unit, two Neumann U87Ai large dual-diaphragm microphones, two
Schoeps MK 6 capsule microphones, an iPhone 12 Pro, a Fireface UFX audio
interface, and a MacBook Pro. A Python tool was developed to playback MIDI
through the Yamaha Disklavier while simultaneously maintaining an audio input
stream. The timing of the MIDI playback is synchronized to the clock of the
audio interface to prevent drift. The software used for data collection is
made available on GitHub 333https://github.com/almostimplemented/piano-
capture. The apparatus was used to re-record all 200 hours of the MAESTRO
dataset, which is released as part of this publication and is available on
Zenodo 444https://zenodo.org/records/10082144.

##  IV Training Experiments

All experiments in this paper use Kong et al.’s [4] model architecture and
implementation. They use a deep neural network architecture that combines
convolutional layers for feature extraction from input spectrograms with gated
recurrent units to capture temporal dependencies. Precise onset and offset
time predictions are achieved via regression heads. However, this choice is
somewhat arbitrary: other powerful transcription models could have been used.
Instead of focusing on modifications to the model architecture, we analyze
various training conditions characterized by the training data.

###  IV-A Data Augmentation Experiment

Training with the data collected from the studio recording setup exhibited
severe overfitting. Our initial training experiment used only the studio audio
with no data augmentation. We evaluate on the test split of MAESTRO, first
using the Studio MAESTRO audio and then using the original audio. On the
Studio MAESTRO data, the note-onset metrics [20] are
($\mathit{Precision}=99.75$, $\mathit{Recall}=95.09$, $\mathit{F1}=97.32$). On
the original MAESTRO data, the performance falls to
($\mathit{Precision}=77.08$, $\mathit{Recall}=85.69$, $\mathit{F1}=80.77$). We
believe the lack of acoustic and timbral variety limits the ability of the
network to generalize to new pianos. To address these shortcomings, we deploy
a novel data augmentation pipeline and supplement the training data with the
original MAESTRO audio plus six additional audio renderings (varying in
timbre) of the MIDI by Modartt Pianoteq version 7
555https://www.modartt.com/pianoteq_overview.

Figure 1 shows the general scheme of our data augmentation pipeline, which is
implemented with the audiomentations package and is inspired by the pipeline
of Hawthorne et al. [2] with a few modifications. We apply two random seven
band parametric equalizers, additive background noise with a variable signal-
to-noise (SNR) ratio, random pitch shifting (between $\pm 0.1$ semitones), and
reverb. The background noise is inspired by the pub noise of the Audio
Degradation Toolbox (ADT) [14]: the audio files consist of recordings of noisy
public environments, typically at cafes and bars. We use restaurant08.wav from
ADT plus four recordings available on freesound.org, cut into 10-second
segments to apply at training time (each training example is also 10 seconds).
The reverb is similarly generated from 14 impulse responses from
echothief.com. For each training example, there is a 50% chance that each
stage of the pipeline will be applied.

The neural network model is taken from Kong et al. [4]. During training we use
a batch size of 32, a learning rate of $5\times 10^{-4}$, and run for 200000
steps. Since we have 8 different sources of audio (original MAESTRO, Studio
MAESTRO, and the six synthesized versions), we apply a sampling scheme of
$1/4,1/4,1/12,...,1/12$ respectively. We forgo training the pedal predictor
and focus on note-onset prediction, which is independent of the model’s pedal
predictions.

###  IV-B Results

We present results of evaluation on the MAPS test set configuration [11].
Note-onset precision, recall, and F1 are reported in Table I. In all cases,
the models were not trained on the MAPS dataset. Maman and Bermano [18] was
the existing SOTA on the MAPS test set, which can be attributed to the
diversity of their training data. Toyama et al. [7] hold the SOTA on the
MAESTRO test set ($Precision=99.64$, $Recall=95.44$, $F1=97.44$), but their
model fails to reach comparable performance in the out-of-distribution setting
on MAPS. Our model achieves a new SOTA for MAPS, without ever training on the
MAPS train split. The model architecture is identical to Kong et al. [4], so
the 6% improvement owes exclusively to our data-driven methodology. It is
worth noting that all models (except Maman and Bermano) suffer a considerable
drop in performance when evaluating on the MAPS test set compared to that of
MAESTRO. This motivates further investigation into regularization techniques
of these models. In the next section, we look at the effect of our chosen data
augmentations in finer detail.

TABLE I: Out-of-distribution performance of state-of-the-art piano transcription models evaluated on the MAPS test set.  | Note-level metrics
---|---
Model | Precision | Recall | F1
Hawthorne et al. [2] | 87.5 | 85.6 | 86.4
Kong et al. [4] | 78.3 | 87.2 | 82.4
Maman and Bermano [18] | 88.2 | 86.5 | 87.3
Toyama et al. [7] | 84.6 | 85.7 | 85.1
Ours | 89.5 | 87.4 | 88.4

##  V Data Augmentation Analysis

We present several experiments to explore the effect of data augmentation on
training piano transcription systems. The first experiment (V-A) applies our
data augmentations independently at test time to demonstrate the sensitivity
of a model trained without augmentation. The second experiment (V-B) applies
each augmentation at training time and measures the impact on out-of-
distribution evaluation. The final experiment (V-C) is an ablation study,
again at training time, that removes one component of the augmentation
pipeline at a time.

[h]Data Degradation Results, showing the effect of data augmentation on
performance. Kong et al.’s model is trained only on MAESTRO’s training data.
Our model benefits from a diversified training set, including original,
studio, and synthesized MAESTRO recordings, all augmented. Evaluation spans
test splits from MAESTRO, Studio MAESTRO, and MAPS, focusing on the note-onset
F1 score. MAESTRO Studio MAPS Kong [4] Ours Kong [4] Ours Kong [4] Ours No
Augmentation 96.8 96.6 95.7 97.7 82.4 88.4 Background Noise 93.4 95.4 91.4
97.1 82.6 88.1 EQ 96.6 96.5 95.1 97.7 81.8 88.4 Pitch Shift 85.2 94.1 76.5
96.0 72.4 87.3 Reverb 89.1 95.2 87.0 96.9 72.3 87.1

###  V-A Data Degradation Study

As an initial investigation of the robustness of published models, we show in
Table V how simple augmentations to test data can cause a dramatic drop in
accuracy. We decompose our data augmentation pipeline and apply each
augmentation in isolation during test set evaluation. We perform this on three
datasets: MAPS, MAESTRO, and Studio MAESTRO. We evaluate Kong et al.’s and our
best model trained with data augmentation. Kong’s model suffers a significant
drop in F1 score across all conditions, as much as 19.2 percentage points. Our
model is mostly invariant to the perturbations, with the largest drop of only
2.5.

###  V-B Single Augmentations

To better understand the influence of data augmentation during training, we
conduct a series of training experiments during which we apply only one
augmentation. Furthermore, we look at the influence of selection of training
data, using either the original MAESTRO or the Studio MAESTRO re-recording.
Each augmentation is applied during training time, again with a 50%
probability of being applied per example. Due to resource constraints, we
limit the training to 28,000 steps (roughly 10 epochs over the ground truth
labels). All other training parameters match those described in Section IV-B.
The results are shown in Table II. The importance of pitch shifting and reverb
as data augmentation techniques is evident, with out-of-distribution
improvements in F1 score of 3.1% and 2.8% respectively. In these experiments,
the additive background noise augmentation seems insignificant or even
detrimental to performance. This is likely due to the clean recording
conditions of the MAPS test set. The EQ augmentation also does not appear to
help. This could indicate poor parameterization of the equalizer, or perhaps
too weak boosts and cuts.

TABLE II: Single Augmentation Results:
Each value corresponds to the note-onset F1 score when evaluating on the MAPS test set. Each experimental condition used a single data augmentation and either the original MAESTRO data or the Studio MAESTRO data. |  | Training data
---|---|---
|  | MAESTRO | Studio MAESTRO
Augmentation | No Augmentation |  82.4 666For this condition, we report the evaluation of the published model. | 75.2
Background | 82.7 (+0.3) | 73.3 (-1.9)
Pitch Shift |  85.5 (+3.1) |  78.3 (+3.1)
Reverb |  85.2 (+2.8) | 76.6 (+1.4)
EQ | 82.1 (-0.3) | 74.4 (-0.8)

###  V-C Ablation Study

In the next series of experiments, we again run training for 28,000 steps, but
this time we remove only one component from the full data augmentation
pipeline. The results are reported in Table III, with the first row being the
full data augmentation pipeline. Once again, the strong influence of pitch
shifting and reverb are demonstrated by the relatively large drop in
performance when skipping these augmentations. Some interaction effects seem
to be present between the various augmentations. In particular, the EQ
augmentation appears inconsequential or even detrimental when applied as a
single augmentation, but removing it from the pipeline in the ablation study
hurts results.

TABLE III: Ablation Experiments:
In these experiments, we train the model with versions of the data augmentation pipeline with one augmentation removed. Note-onset F1 score on the MAPS test set is reported. |  | Training data
---|---|---
|  | MAESTRO | Studio MAESTRO
Ablation | Full Augmentation | 86.4 | 79.0
Skip Background | 85.4 (-1.0) | 78.8 (-0.2)
Skip Pitch Shift |  82.9 (-3.5) |  75.6 (-3.4)
Skip Reverb |  82.8 (-3.6) | 77.1 (-1.9)
Skip EQ | 86.4 (-0.0) | 77.8 (-1.2)

##  VI Conclusion

This work explored the role of data and data augmentation techniques in
enhancing the robustness of APT systems, particularly regarding their
performance on out-of-distribution data. Utilizing a new studio recording of
the MAESTRO dataset and simple data augmentation, our research yields a new
SOTA on the MAPS test set and new insights into the influence of data
augmentation on piano transcription models.

Despite the strong performance on the MAPS test set, the findings underscore
the nuanced and occasionally divergent effects of different data
augmentations. Particularly, pitch shifting and reverb emerged as substantial
in improving model generalization, while others like EQ augmentation did not
yield anticipated benefits. A limitation of our approach is that the
combinatorial nature of exploring mixtures of data augmentations made an
exhaustive study prohibitively expensive. Furthermore, we were unable to
conduct precise tests of significance, which would require running the
training task for every condition multiple times and measuring variance across
runs. Additional evaluation datasets are necessary to further measure the
robustness of our model; in particular, live performances in a noisy
environment would provide a better indication of whether the background noise
augmentation is useful.

Moving forward, one avenue for future work lies in leveraging the advanced
physical modeling features of Pianoteq to synthesize piano audio with diverse
tunings and string inharmonicities, which could further enrich the timbral
diversity of training data. Additionally, the studio recording of MAESTRO may
present a valuable resource for other research areas, such as audio style
transfer, offering a unique opportunity to delve into style translations
between different recordings. The gap between performance on MAESTRO test and
MAPS test motivates further research into the differences between these
datasets, to ensure accuracy of their labels and identify any inconsistencies.

In light of our findings, we advocate for a heightened emphasis on out-of-
distribution evaluation in future music transcription research, not only to
benchmark model performance but also to delve deeper into understanding and
enhancing model robustness and generalization across varied and potentially
challenging acoustic environments.

## References

  * [1] Curtis Hawthorne, Erich Elsen, Jialin Song, Adam Roberts, Ian Simon, Colin Raffel, Jesse Engel, Sageev Oore, and Douglas Eck,  “Onsets and frames: Dual-objective piano transcription,”  in Proceedings of the 19th International Society for Music Information Retrieval Conference, New Orleans, Louisiana, USA, 2018, pp. 50–57.
  * [2] Curtis Hawthorne, Andriy Stasyuk, Adam Roberts, Ian Simon, Cheng-Zhi Anna Huang, Sander Dieleman, Erich Elsen, Jesse Engel, and Douglas Eck,  “Enabling factorized piano music modeling and generation with the MAESTRO dataset,”  in International Conference on Learning Representations, 2019.
  * [3] Valentin Emiya, Nancy Bertin, Bertrand David, and Roland Badeau,  “MAPS - A piano database for multipitch estimation and automatic transcription of music,”  Research report inria-00544155, Telecom ParisTech, July 2010.
  * [4] Qiuqiang Kong, Bochen Li, Xuchen Song, Yuan Wan, and Yuxan Wang,  “High-resolution piano transcription with pedals by regressing onset and offset times,”  IEEE/ACM Transactions on Audio, Speech and Language Processing, vol. 29, pp. 3707–3717, 2021.
  * [5] Joshua P Gardner, Ian Simon, Ethan Manilow, Curtis Hawthorne, and Jesse Engel,  “MT3: Multi-task multitrack music transcription,”  in International Conference on Learning Representations, Online, 2022.
  * [6] Wei Tsung Lu, Ju-Chiang Wang, and Yun Ning Hung,  “Multitrack music transcription with a time-frequency perceiver,”  in ICASSP 2023 - 2023 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP), Rhodes Island, Greece, June 2023, pp. 11771–11775, IEEE.
  * [7] Keisuke Toyama, Taketo Akama, Yukara Ikemiya, Yuhta Takida, Wei-Hsiang Liao, and Yuki Mitsufuji,  “Automatic piano transcription with hierarchical frequency-time transformer,”  in Proceedings of the 24th International Society for Music Information Retrieval Conference, Milan, Italy, 2023, pp. 215–222.
  * [8] Arthur Flexer and Dominik Schnitzer,  “Effects of album and artist filters in audio similarity computed for very large music databases,”  Computer Music Journal, vol. 34, no. 3, pp. 20–28, 2010.
  * [9] Robert Geirhos, Jörn-Henrik Jacobsen, Claudio Michaelis, Richard Zemel, Wieland Brendel, Matthias Bethge, and Felix A. Wichmann,  “Shortcut learning in deep neural networks,”  Nature Machine Intelligence, vol. 2, no. 11, pp. 665–673, Nov 2020.
  * [10] Connor Shorten and Taghi M. Khoshgoftaar,  “A survey on image data augmentation for deep learning,”  Journal of Big Data, vol. 6, no. 1, July 2019.
  * [11] Siddharth Sigtia, Emmanouil Benetos, and Simon Dixon,  “An end-to-end neural network for polyphonic piano music transcription,”  IEEE/ACM Transactions on Audio, Speech and Language Processing, vol. 24, no. 5, pp. 927–939, May 2016.
  * [12] Curtis Hawthorne, Ian Simon, Rigel Swavely, Ethan Manilow, and Jesse Engel,  “Sequence-to-sequence piano transcription with transformers,”  in Proceedings of the 22nd International Society for Music Information Retrieval Conference, Online, 2021, pp. 246–253.
  * [13] Justin Salamon and Juan Pablo Bello,  “Deep convolutional neural networks and data augmentation for environmental sound classification,”  IEEE Signal Processing Letters, vol. 24, no. 3, pp. 279–283, 2017.
  * [14] Matthias Mauch and Sebastian Ewert,  “The audio degradation toolbox and its application to robustness evaluation,”  in Proceedings of the 14th International Society for Music Information Retrieval Conference, Curitiba, Brazil, 2013, pp. 83–88.
  * [15] Brian McFee, Eric J. Humphrey, and Juan Pablo Bello,  “A software framework for musical data augmentation,”  in Proceedings of the 16th International Society for Music Information Retrieval Conference, Málaga, Spain, 2015, pp. 248–254.
  * [16] John Thickstun, Zaid Harchaoui, Dean P. Foster, and Sham M. Kakade,  “Invariances and data augmentation for supervised music transcription,”  in International Conference on Acoustics, Speech, and Signal Processing (ICASSP), Calgary, Alberta, Canada, 2018, pp. 2241–2245.
  * [17] Ian Simon, Josh Gardner, Curtis Glenn-Macway Hawthorne, Ethan Manilow, and Jesse Engel,  “Scaling polyphonic transcription with mixtures of monophonic transcriptions,”  in Proceedings of the 23rd International Society for Music Information Retrieval Conference, Bengaluru, India, 2022, pp. 44–51.
  * [18] Ben Maman and Amit H. Bermano,  “Unaligned supervision for automatic music transcription in the wild,”  in International Conference on Machine Learning, ICML 2022, 17-23 July 2022, Baltimore, Maryland, USA. 2022, vol. 162 of Proceedings of Machine Learning Research, pp. 14918–14934, PMLR.
  * [19] John Thickstun, Zaid Harchaoui, and Sham M. Kakade,  “Learning features of music from scratch,”  in International Conference on Learning Representations (ICLR), Toulon, France, 2017, pp. 1378–1391.
  * [20] Mert Bay, Andreas F. Ehmann, and J. S. Downie,  “Evaluation of multiple-f0 estimation and tracking systems,”  in Proceedings of the 10th International Society for Music Information Retrieval Conference, Kobe, Japan, 2009, pp. 315–320.