B Dipierro

M Abid

Rayleigh-Taylor instability in variable density swirling flows

Using linear instability theory and nonlinear dynamics, the Rayleigh-Taylor instability of vari-able density swirling flows is studied. It is found that the flow topology could be predicted, when the instability sets in, using a function χ dependent on density and axial and azimuthal velocities. It is shown that even when the inner axial-flow is heavier than the outer one (a favorable case for the development of the Rayleigh-Taylor instability thanks to the centrifugal force) the instability is not necessarily Rayleigh-Taylor-dominated. It is also shown that when the Rayleigh-Taylor instability develops, it is helical.

Introduction

Swirling flows are important for many application devices and as a fundamental problem considering their relevance to aircraft trailing vortices, vortical transport of momentum and energy in meteorology, and vortex bursting (or breakdown) [1][2][3][4]. Hence, they are widely studied in the literature ( [1][2][3][4][5][6][7], among many others) but the physical mechanisms of their instabilities, when density variations are present, are not generally discussed.

In a variable density swirling flow (for example a columnar vortex with an axial flow, a configuration that will receive a particular attention in this study) there are ingredients for the development of two fundamental instabilities: Rayleigh-Taylor [8][START_REF] Taylor | Proc. R. Soc. Lond. A[END_REF][START_REF] Sharp | [END_REF] (due to density variations) and Taylor-Couette [11][12][13] (promoted by differential rotation). Intuitively, the instability of the flow is a mixing of Rayleigh-Taylor and Taylor-Couette instabilities. Theref o r e ,an a t u r a lq u e s t i o n( f r o map h y s i c a lpo i n to fv i e w )i s how to formulate the instability problem to account for the contribution of each of these latter instabilities to the general, three-dimensional, instability of a variable density swirling flow? Note that the key instability dynamics of Rayleigh-Taylor and Taylor-Couette can both be represented with a minimum of two dimensions. The former is characterized by "bubbles" for positively buoyant fluid and "spikes" for negatively buoyant fluid, the latter being axisymmetric with the organization of the flow as toroidal patterns.

The objective of the present study is to answer the above question based on the recent asymptotic analysis by [14]. To this end, since Taylor experiments [11]o fa fluid between two rotating cylinders have indicated that the disturbance could be axisymmetric when the instabila e-mail: abid@irphe.univ-mrs.fr ity occurs, let us define two particular instabilities: an instability with an axisymmetric centrifugal-force that could be identified with an axisymmetric Taylor-Couette instability, and a two-dimensional "orthogonal state" to the previous one, with no axisymmetry of the centrifugal force, that could be identified with a two-dimensional Rayleigh-Taylor instability (see below for more precision about the orthogonality). We show that the nonlinear evolution of the three-dimensional instability is a mixing of these two (pure) states with a coupling term of axial and azimuthal velocities. This coupling term is responsible for the threedimensionality of the flow, giving it a helical form.

Linear asymptotics

Consider an inviscid and incompressible swirling flow with a basic velocity field (in a cylindrical coordinate system (r, θ, x) with the basis (e r , e θ , e x ); (x, y, z)beingtheusual cartesian coordinate system),

U 0 (x)=(0,V(r),W(r)) (er,e θ ,ex) ,r = y 2 + z 2 (1)
and a density field, ρ b (r). Let f (r) be the amplitude, k the axial wavenumber, m the azimuthal wavenumber and ω the complex angular frequency of a perturbation to the flow. For a helical perturbation of the form f (r)exp(i(mθ + kx -ωt)), in the limit of large wavenumbers, it is shown in [14] that the leading order (in ǫ = O(1/k)=O(1/m)) of the perturbation growth rate σ ≡ℑ (ω)i se q u a lt oχ(m, k, r 0 ) where the function χ is defined as: and r 0 is the location where χ is extremal. In the χ expression, Ω = V/r is the local angular velocity, H is defined as:

χ 2 (m, k, r)= m 2 G 2 + k 2 r 2 H 2 +2mkrΩW ′ m 2 + k 2 r 2 , (2)
H 2 = G 2 -Φ, G (the analogue of the Brunt- Väisäla frequency) is defined as: G 2 = -(V 2 /r)ρ ′ b /ρ b and Φ = r -3 d(r 2 V 2 )/dr =2 ΩΘ (the Rayleigh discriminant).
The prime, ′ , is the derivative with respect to the r variable.

It is also shown in [14] that for the particular case (m ≫ 1,k = 0), the leading order (in ǫ = O(1/m)) of the growth rate is equal to G(r 0 )( r 0 is the location where G is extremal), and for the particular case (m =0 ,k ≫ 1) the leading order (in ǫ = O(1/k)) of the growth rate is equal to H(r 0 )( r 0 is the location where H is extremal). In the present paper we note that these two particular cases could be related to two classically well known instabilities. In the case with (m =0 ,k = 0), the perturbation grows only in the azimuthal direction leading eventually to the formation of finger-like structures of heavy fluid penetrating the light one, as explained and sketched in Figure 1a. This is the Rayleigh-Taylor instability where the centrifugal acceleration replaces the gravity. In the case with (m =0 ,k = 0), the perturbation grows axisymmetrically, and is wavy in the axial direction, leading to the organization of the flow in toroidal structures characteristic of the Taylor-Couette instability (Fig. 1b). Therefore, we will call the two orthogonal eigenstates (m =0 ,k =0 )a n d( m =0 ,k = 0) Rayleigh-Taylor and Taylor-Couette, respectively. By inspection of the χ expression (Eq. ( 2)), clearly the three-dimensional instability of a variable density swirling flow is a mixing of both Rayleigh-Taylor instability (with a weight equal to m 2 /(m 2 + k 2 r 2 0 )), Taylor-Couette instability (with a weight equal to k 2 r 2 0 /(m 2 + k 2 r 2 0 )) and an instability due to the coupling of rotation and the axial shear (with a weight equal to 2mkr 0 /(m 2 + k 2 r 2 0 )). Note that the coupling term is responsible for the three-dimensionalization of the flow instability. Note also that for a fixed k and m →∞the function χ → G and the instability should be Rayleigh-Taylor-dominated. For, m fixed and k →∞ the function χ → H and the instability should be Taylor-Couette-dominated.

The linear asymptotical-results of Di Pierro and Abid [14] are valid for large azimuthal and axial wavenumbers. It will be shown here, numerically, that they work well even for low values of m and k and that they are also valid in the nonlinear regime.

Nonlinear study

We now study the instability behavior beyond the point that linearizing assumptions are valid, and use numerical simulation to study this non-linear regime. The governing equations are:

∂ρ ∂t + u •∇ρ =0, (3)
∂u ∂t +(u •∇)u + ∇P ρ - 1 Reρ ∆u =0, (4)
∇•u =0, ( 5
)
where ρ is the variable (in space and time) density field, Re =ρU a/µ is the Reynolds number with µ the dynamic viscosity, ρ a density scale (chosen here as the density at infinity), U a characteristic velocity (chosen as the velocity on the axis) and a the characteristic length scale (chosen as the vortex-core length). These equations are the full variable density Navier-Stokes equations and there is no Boussinesq approximation [15,16]. Therefore large density variations are allowed. To appreciate the analogy with a Rayleigh-Taylor system note that for Re = ∞ any basic velocity (Eq. ( 1)), U 0 (x), and density, ρ b (r), are equilibrium solutions of the previous equations (see the appendix) provided that ∇P = ρ b (V 2 /r)e r (to be compared with an equilibrium of a planar density-stratification in the gravity field, g = -ge x ,f o rw h i c h∇P = -ρ b ge x ).

To study this possible competition, we solve the variable density, constant dynamical viscosity, Navier-Stokes equations using numerical simulations. The time is made nondimensional using the time scale a/U .W eu s ep e r iodic boundary conditions and the initial condition is a superposition of the basic flow (Eq. ( 1)) and a perturbation flow obtained as an eigenmode of the linearized Navier-Stokes equations (about the basic flow) [17]. We seek the nonlinear behavior of the solution (after the time of the saturation of the linear instability, noted t s and defined in Fig. 2). The solution is obtained numerically using a pseudo-spectral (Fourier) method for the spatial discretization and a second order Runge-Kutta scheme for the time discretization. Aliased modes are removed using the "2/3" rule [18,19]. The computation is solved in a cartesian domain with a volume equal to (4πa) 3 . The use of a pseudo-spectral algorithm in a cubic (or parallelipipedic) box, for a flow that is initially axisymmetric, is possible for flows with an exponential decay of their vorticity Fig. 2. Plot against time of the energy, E mk , of the perturbation modes a) (m =1,k = 2) for the experimental profile (Exp) and (b) (m =3 ,k = 1) for the Batchelor-q-Vortex (BqV). Solid-line: the linear instability-theory prediction. Dot-dashed line: the energy growth obtained using the numerical solution of the variable density Navier-Stokes equations. The saturation time is defined as the time when the exponential growth, predicted by linear instability theory, is no longer valid due to nonlinearities of Navier-Stokes equations: for t<t i s (i denoting one of the two profiles) the dynamics is linear, and for t>t i s it is nonlinear. Note that for t<t i s the growth rate of the perturbation predicted using linear instability theory (solid line) is well recovered in the numerical simulation (dot-dashed line): this is a validation of the numerical procedure.

with the radial distance. The advantages of this method (compared to a method using a cylindrical coordinate system) are: (a) this pseudo-spectral method is easy to implement and parallelize; (b) there is no singularity at the origin of the coordinate system, hence there is no need to implement the pole conditions (or equations) [20]. In the literature, many studies have been done with a pseudo-spectral algorithm in a parallelipipedic box, for a flow that is initially axisymmetric, principally for jets, swirling jets and wakes [21][22][23][24][25]. Note that these numerical studies are done in the constant density case. In the present work, we extended the numerical method to the variable density case, as described below.

Using standard notations, at time t n = nδt (n ∈ N), where δt is the time-step, equations (3)-( 5)a r es o l v e da s follows. Let

P ρ =(I d -(1/ρ)∇L -1 ρ ∇•)whereI d is the iden- tity operator, L -1
ρ is the inverse (in the space of zero mean functions) of L ρ ≡∇•(1/ρ∇)and∇• is the divergence operator. Note that P ρ is an exact projector, i.e., P 2 ρ = P ρ , and that ∇•P ρ v = 0, for any differentiable vector field v and density profile ρ (i.e., P ρ is a projector onto the space of divergence-free vector fields). The equation

∂u ∂t n = -(u n •∇)u n + 1 Reρ n ∆u n - ∇P n ρ n ≡ F(u n ,ρ n ) - ∇P n ρ n (6)
is advanced in time using a modified, fully explicit, second order Runge-Kutta scheme with two intermediate projection steps:

K * 1 = F(u n ,ρ n ), (7)
K n 1 = P ρ K * 1 , (8)
K * 2 = F(u n + δtK n 1 ,ρ n ), (9)
K n 2 = P ρ K * 2 ( 10
)
u n+1 = u n + 1 2 δt(K n 1 + K n 2 ). ( 11
)
Note that this time integration is chosen to be fully explicit because all the terms in F(.) are nonlinear (especially the diffusive term: (∆u n /(ρ n Re))). Note also that the projection step needs the inversion of L ρ that will be detailed when the pressure is calculated. The density is calculated with a classical second-order Runge-Kutta scheme,

∂ρ ∂t n = -u n •∇ρ n (12)
≡ G(u n ,ρ n ),

ρ n+1 = ρ n + 1 2 δt(G(u n ,ρ n ) (13) + G(u n ,ρ n + δtG(u n ,ρ n )). ( 14
)
Consistency of the numerical scheme with variable density Navier-Stokes equations leads to the following equation for the pressure:

L ρ P n = ∇• ∇P n ρ n = 1 2 ∇•(K * 1 + K * 2 ), (15)
with K * 1 and K * 2 expressed above. Noting that the modified Runge-Kutta scheme, for the velocity, needs the computation of P ρ K * 1 and P ρ K * 2 separately, equation ( 15)i s written as Poisson equations,

∆P n i = ∇ρ n •∇P n i ρ n + ρ n ∇•K * i ≡ h n i ,i =1, 2, (16)
and, if needed, the pressure is computed as P n =( P n 1 + P n

2 )/2. The solution of equation ( 16), for the "partial " pressure P i ,giv esL -1 ρ ∇•K * i that is necessary in the computation of P ρ K * i . To this end, equation ( 16)i ss o l v e d iteratively. Let k (∈ N) be the current iteration index. Hence, P n i,k is computed as:

P n i,k = ∆ -1 ∇ρ n •∇P n i,k-1 ρ n + ρ n ∇•K * i ≡ ∆ -1 h n i,k-1 , (17)
and the iterative process is stopped when the residual,

||∆P n i,k -h n i,k-1 ||
, is smaller than a desired accuracy (for a given norm). Note that ∆ -1 is the laplacian inverse for zero mean functions. The number of iterations needed to converge is very dependent on the initial value P n i,0 . Therefore, to accelerate slightly this process, the previous converged value of P i is used to initialize the next one: P n+1 i,0 = P n i,c ,w h e r eP n i,c is the converged "partial" pressure at the previous time. In practice, we find that six iterations are sufficient to obtain a velocity divergence O(10 -12 ).

The overall numerical procedure is validated as follow: (a) we check the spectral accuracy by verifying the exponential decay of the power spectrum density (at high wavenumbers) of each computed field (velocity, pressure and density); (b) we check that the numerical time-error, ǫ t , is of order two when decreasing the time-step δt, i.e., ǫ t = O((δt) 2 ) in the limit δt → 0; (c) we compute the growth rates of flow perturbations at early times (when nonlinearities are small: t ≪ t i s , i denoting one of the two profiles used, see below and Fig. 2) and compare them to growth rates obtained using diagonalization of linearized Navier-Stokes equations [17]: we check that differences between growth rates are less than 5%. For example, for the radial velocity (which is initially zero), r , computed using the Navier-Stokes equations, let us define the energy in a helical structure as:

E mk (t)= D dxu 2 r (x,t)e -i(kx+mθ) , ( 18
)
where D is the volume of the computational box. In Figure 2, the energy E mk o ft h em o d e( m =1 ,k =2 )i s plotted as a function of time (for s =2a n du s i n gt h e experimental profile as a basic flow, see below for its definition). The energy growth-rate predicted using the linear instability theory (the continuous line) is 0.67. The en-ergy growth-rate obtained using the numerical solution of Navier-Stokes equations (the dot-dashed line) is: 0.66 for t<t Exp s . We study the competition between Rayleigh-Taylor and Taylor-Couette instabilities for two particular flows. The first one is the well known Batchelor q vortex (Bqv) [26], that we extend to the variable density case:

V (r)/W (0) = q 1 -exp(-(r/a) 2 ) (r/a) , W (r)/W (0) = exp(-(r/a) 2 ), ρ b (r)/ρ ∞ =1+(s -1) exp(-(r/a) 2 ), ( 19
)
where q is the swirl number and s is the density ratio between the vortex center and the ambient fluid. The second one is a fit of experimentally obtained vortices with axial flows [27], that we extend also to the variable density case:

V (r)/W (0) = q(r/a)exp(-(r/a) 4 ), W (r)/W (0) = 1 1 + (exp((r/a) 2 log(2)) -1) 3 , ρ b (r)/ρ ∞ =1+ s -1 1 + (exp ((r/a) 2 log(2)) -1) 3 . (20)
In what follows the length scale is chosen as the vortexcore length a, the velocity scale as the axial velocity at r =0(U = W (0)) and the density scale as the density at infinity ρ = ρ ∞ ≡ ρ b (r →∞). Two, new, nondimensional numbers could be introduced: the Atwood number At = (ρ b (0)-ρ ∞ )/(ρ b (0)+ρ ∞ )=(s-1)/(s+1), and the Taylor number Ta = Ω 2 a 4 /(µ/ρ ∞ ) 2 = q 2 Re 2 (with Ω = qU/a the angular-velocity scale).

The azimuthal velocity of the Bqv has an algebraic decay. To accommodate this velocity profile in the periodic domain of computation we remark that it is diverg e n c ef r e ea n di t sv o r t i c i t y( ω 0 = ∇×U 0 ) has an exponential decay. Therefore, we compute the velocity as: U 0 = -∆ -1 ∇×ω 0 , which is automatically periodic if the operator ∇ is the gradient on the space of periodic functions.

In practice, if L y and L z are the lateral dimensions of the computational box, we choose L y ≫ a and L z ≫ a (see Fig. 3) to ensure that the vortex neighbors, introduced by periodicity, do not affect the growth rates of the linear instability. We used L x = L y = L z ≃ 14a. Two resolutions (192 3 and 256 3 ) are used to ensure convergence. For azimuthal modes, the characteristic scale is 2πa (since they grow near the vortex core, r = a,w h e r e the gradients are high). Therefore, with lateral dimensions L y = L z ≃ 14a,w eh a v e( 2 πa/14a)256 ≃ 114 points. If we consider, for example, that 5 points per wavelength are sufficient to accurately represent a wave [START_REF] Threfethen | Spectral Methods in Matlab[END_REF], we are able to accurately represent azimuthal modes up to m ≃ 22. This is enough for the physics we are considering (instability, transition to turbulence and not turbulence).

Computational box

L y L z Vortex 2a
Fig. 3. Geometry of the computational box in the lateral plane (y, z). The lateral lengths are chosen such that Ly ≫ a and Lz ≫ a to ensure that the vortex neighbors, introduced by the periodization process, do not affect the growth rates of the linear instability.

The reynolds number is chosen as Re = 1000, that is (initially) the self-advection is 1000 times more important than the viscous diffusion. This value of Re is selected as the highest possible value that retains numerical stability, at the used resolutions. It is interesting to mention that this value of Reynolds number is close to ones used in a variable-density jet experiment [START_REF] Boujemaa | [END_REF], and a constantdensity swirling-jet experiement [27].

Note that these profiles are exact solutions of Euler equations. When these profiles are used as initial conditions for the Navier-Stokes equations, the corresponding velocity profiles diffuse (see the appendix). For example, when s = 1, the axial velocity of the Bqv diffuses such that a = a 0 1+4t/Re,w h e r ea 0 is the initial vortex core length (indeed the axisymmetric diffusion equation has a solution of the form exp(-(r/a) 2 )/a provided that a(da/dt)=2 /Re). Note also that the azimuthal velocity of the Bqv is stable to axisymmetric perturbations, however that of experimental profiles is unstable to this kind of perturbations [START_REF] Landau | Fluid Mechanics[END_REF][START_REF] Bayly | [END_REF].

For the Bqv, in the linear framework (using the function χ), we find that χ is dominated by G if s>1 for all m, k and q. Furthermore, the coupling term is not identically zero. Therefore, we expect the development of a three-dimensional Rayleigh-Taylor instability and the formation of helical spikes in the nonlinear regime. This is confirmed by the results presented, for example, in Fig- ure 4. The evolution of the Bqv( f o rs =2andq = -0.5), initially disturbed with the mode k =1a n dm =3 ,i s shown. At the beginning, three helical fingers grew radially. The heavy fluid, initially localized in the vortex core, is gradually transported far from the axis, owing to these helical structures. When nonlinearities become important, fingers roll-up into spikes and at the end of the simulation, the heavy fluid is localized around r ∼ 2. The simulation is repeated with other unstable (m, k)m o d e s , and the same flow topology is observed. Hence, the helical Rayleigh-Taylor instability is dominant for the Bqvw i t h s>1. It is interesting to mention that for higher values of Re, vorticity would be greater and one would expect tighter roll-up of spikes. Therefore, for the late evolution of this roll-up a time-varying Reynolds number, based on the evolving length-scale of the spikes, would be more appropriate.

For s<1, for the Bqv, we find that χ is dominated by the coupling term ΩW ′ . The topology of the flow in the nonlinear regime is different from the previously described Rayleigh-Taylor instability. It is shown in Figure 5 for s = 0.5a n dq = -0.5 and for the same mode (k =1a n d m = 3), for example. Helical structures appear and there is no spikes formation.

For the experimental profile we find that the function χ is principally dominated by H for all q, s, m and k values. We therefore expect that in the nonlinear regime the flow develops a topology similar to that resulting from the Taylor-Couette instability. This is confirmed in Figure 6. We see that disturbances, even in the nonlinear stage, are concentrated near r = 1, as if the flow was bounded (note that a confined dynamics may be also observed when important density stratifications are present in Rayleigh-Taylor systems [32,33]). Furthermore, although the fluid topology is helical it presents parts of toroidal structures similar to Taylor vortices. Note that even if s>1( af avorable case for the development of the Rayleigh-Taylor instability), the instability is not Rayleigh-Taylor dominated.

4Conclusion

We now summarize our main findings. The three-dimensional instability of swirling flows, submitted to density variations, is formulated as a competition among three distinct instabilities: Rayleigh-Taylor, Taylor-Couette and an instability due to the coupling between the rotation and the axial shear. We have presented the results of a linear asymptotic analysis, valid for incompressible inviscid fluids in the limit of large wavenumbers. We have demonstrated here that the analysis is also valid for nonlinear development at low wavenumbers, as confirmed by numerical simulations. It is found that the flow topology could be predicted, when the instability sets in, using the function χ dependent on density and axial and azimuthal velocities. As an application it is shown that even when the inner axial-flow is heavier than the outer one (a favorable case for the development of the Rayleigh-Taylor instability thanks to the centrifugal force) the instability is not necessarily Rayleigh-Taylor-dominated.

Appendix A

In a cylindrical coordinate system (r, θ, x)t h ev a r iable density Navier-Stokes equations for velocity, u = (u r ,u θ ,u x ), and density, ρ,a r e : For any stationary (axisymmetric) velocity, U =( u r = 0,u θ = V (r),u x = W (r)), and density, ρ = ρ b (r), profiles the following terms vanish: In the constant density case (choose for example ρ =1),if the swrling-flow velocity U is used as an initial condition for the Navier-Stokes equations, the velocity components u θ and u x will diffuse according to: The pressure is given by:

∂ρ
∂u θ ∂t = 1
u 2 θ r = ∂P ∂r ,
and there is no generation of radial velocity.

In the numerical simulations, we checked that when a swirling jet is used as an initial condition there is no generation of azimuthal modes m =4andm = 8 by periodicity. That is to say the swirling jet only diffuses and there is no generation of radial velocity by periodicity. This is also a validation that the values of lateral lengths of the computational box are sufficient and the ghost swirling-flows have no effects on the flow of the computational box.

We would like to acknowledge valuable comments of anonymous referees that ameliorated the presentation of the present work. The "Centre Informatique National de l'enseignement supérieur" (CINES) is also acknowledged for computational facility.

Fig. 1 .

 1 Fig. 1. (a) A schematic of a cylindrical Rayleigh-Taylor instability. The flow acceleration, Γ, is oriented outward from the vortex core, and if the vortex is heavier than the ambient fluid (i.e. ∇ρ is oriented from the outside to the vortex core), this leads to a Rayleigh-Taylor unstable configuration. (b) Axisymmetric Taylor-Couette instability.

Fig. 4 .=1 . 5 .

 45 Fig. 4. Nonlinear evolution of a variable density Bqv, with q = -0.5( Ta =2 .5 × 10 5 )a n ds =2( At =0 .33). The initial "interface" disturbance has k =1,m = 3 and an initial amplitude O(10 -4 ) times that of the initial velocity. (a)-(d) Contour plot of the density at x =0fort/t Bqv s =0, 0.6, 1.1a n d1 .8 respectively (t Bqv s is the linear-instability time-saturation). (d) Isosurface of the density at a medium level for t/t Bqv s =1 .5. Three helical spikes are ejected radially. Only the central part (near the vortex) of the computational box is shown. Note that spatial scales are in unit of a.

Fig. 5 .

 5 Fig. 5. Same conditions as in Figure 4 but with s =0 .5( At = -0.33). Helical unstable structures due to the coupling of rotation and axial shear appear. Only the central part (near the vortex) of the computational box is shown.

Fig. 6 .=1 . 5 .

 65 Fig. 6. Nonlinear evolution of a variable density experimental profile with q = -0.5( Ta =2 .5 × 10 5 )a n ds =2( At =0 .33). The initial "interface" disturbance has k =2 ,m = 1 and an initial amplitude O(10 -4 ) times that of the initial velocity. (a)-(d) Contour plot of the density at y =0 ,f o rt/t Exp s

  •∇)u r , (u •∇)u θ , (u •∇)u x , u r u θ r ,and the continuity equation (A.5) is automatically verified. Therefore (U,ρ b ) are exact solutions of Euler equations (Re = ∞) provided that:

  (r, t =0)=W (r).


Ramana Sundararaman
email: ramanasubramanyam.sundararaman@polytechnique.edu

Cédric De

Almeida Braga

Eric Marchand
email: eric.marchand@irisa.fr

Julien Pettré
email: julien.pettre@inria.fr

Tracking Pedestrian Heads in Dense Crowd

Tracking humans in crowded video sequences is an important constituent of visual scene understanding. Increasing crowd density challenges visibility of humans, limiting the scalability of existing pedestrian trackers to higher crowd densities. For that reason, we propose to revitalize head tracking with Crowd of Heads Dataset (CroHD), consisting of 9 sequences of 11,463 frames with over 2,276,838 heads and 5,230 tracks annotated in diverse scenes. For evaluation, we proposed a new metric, IDEucl, to measure an algorithm's efficacy in preserving a unique identity for the longest stretch in image coordinate space, thus building a correspondence between pedestrian crowd motion and the performance of a tracking algorithm. Moreover, we also propose a new head detector, HeadHunter, which is designed for small head detection in crowded scenes. We extend HeadHunter with a Particle Filter and a color histogram based re-identification module for head tracking. To establish this as a strong baseline, we compare our tracker with existing state-of-the-art pedestrian trackers on CroHD and demonstrate superiority, especially in identity preserving tracking metrics. With a light-weight head detector and a tracker which is efficient at identity preservation, we believe our contributions will serve useful in advancement of pedestrian tracking in dense crowds. We make our dataset, code and models publicly available at https://project.inria.fr/crowdscience/ project/dense-crowd-head-tracking/.

Introduction

Tracking multiple objects, especially humans, is a central problem in visual scene understanding. The intricacy of this task grows with increasing targets to be tracked and remains an open area of research. Alike other subfields in Computer Vision, with the advent of Deep Learning, the task of Multiple Object Tracking (MOT) has remarkably advanced its benchmarks [START_REF] Dave | Tao: A large-scale benchmark for tracking any object[END_REF][START_REF] Geiger | Vision meets robotics: The kitti dataset[END_REF][START_REF] Geiger | Are we ready for autonomous driving? the kitti vision benchmark suite[END_REF][START_REF] Leal-Taixé | MOTChallenge 2015: Towards a benchmark for multitarget tracking[END_REF][START_REF] Milan | MOT16: A benchmark for multi-object tracking[END_REF][START_REF] Voigtlaender | MOTS: multi-object tracking and segmentation[END_REF] since its Figure 1. Comparison between head detection and full body detection in a crowded scene from CroHD. HeadHunter detects 36 heads whereas Faster-RCNN [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF] can detect only 23 pedestrians out of 37 present in this scene.

inception [START_REF] Ferryman | Pets2009: Dataset and challenge[END_REF]. In the recent past, the focus of MOTChallenge benchmark [START_REF] Dendorfer | CVPR19 tracking and detection challenge: How crowded can it get?[END_REF] has shifted towards tracking pedestrians in crowds of higher density. This has several applications in fields such as activity recognition, anomaly detection, robot navigation, visual surveillance, safety planning etc. Yet, the performances of trackers on these benchmark suggests a trend of saturation 1 . Majority of online tracking algorithms today follow the tracking-by-detection paradigm and several research works have well-established object detector's performance to be crucial in tracker's performance [START_REF] Bergmann | Tracking without bells and whistles[END_REF][START_REF] Bewley | Simple online and realtime tracking[END_REF][START_REF] Ciaparrone | Deep learning in video multi-object tracking: A survey[END_REF]. As the pedestrian density in a scene increases, pedestrian visibility reduces with increasing mutual occlusions, leading to reduced pedestrian detection as visualized in Figure 1. To tackle these challenges yet track humans efficiently in densely crowded environments, we rekindle the task of MOT with tracking humans by their distinctly visible part -heads. To that end, we propose a new dataset, CroHD, Crowd of Heads Dataset, comprising 9 sequences of 11,463 frames with head bounding boxes annotated for tracking. We hope that this new dataset opens up opportunities for promising future research to better understand global pedestrian motion in dense crowds.

Supplementing this, we develop two new baseline methods on CroHD, a head detector, HeadHunter and a head tracker, HeadHunter-T. We design HeadHunter peculiar for head detection in crowded environments, distinct from standard pedestrian detectors and demonstrate state-of-theart performance on an existing head detection dataset. HeadHunter-T extends HeadHunter with a Particle Filter framework and a light-weight re-identification module for head-tracking. To validate HeadHunter-T to be a strong baseline tracker, we compare it with three published top performing pedestrian trackers on the crowded MOTChallenge benchmark, evaluated on CroHD. We further perform comparisons between tracking by head detection and tracking by body detection to illustrate the value of our contribution.

To establish correspondence between a tracking algorithm and pedestrian motion, it is necessary to understand the adequacy of various trackers in successfully representing ground truth pedestrian trajectories. We thus propose a new metric, IDEucl to evaluate tracking algorithms based on their consistency in maintaining the same identity for the longest length of a ground truth trajectory in the image coordinate space. IDEucl is compatible with our dataset and can be extended to any tracking benchmark, recorded with a static camera. In summary, this paper makes the following contributions (i) We present a new dataset, CroHD, with annotated pedestrian heads for tracking in dense crowd, (ii) We propose a baseline head detector for CroHD, HeadHunter, (iii) We develop HeadHunter-T, by extending HeadHunter as the baseline head tracker for CroHD, (iv) We propose a new metric, IDEucl, to evaluate the efficiency of trackers in representing a ground truth trajectory and finally, (v) We demonstrate HeadHunter-T to be a strong baseline by comparing with three existing state-of-the-art trackers on CroHD.

Related Work

Head Detection Benchmarks: The earliest benchmarks in head detection are [START_REF] Hoai | Talking heads: Detecting humans and recognizing their interactions[END_REF][START_REF] Patron-Perez | Structured learning of human interactions in tv shows[END_REF][START_REF] Vu | Contextaware CNNs for person head detection[END_REF][START_REF] Ren | Finding people in archive films through tracking[END_REF], which provide ground truth head annotations of subjects in Hollywood movies. In the recent past, SCUT-Head [START_REF] Peng | Detecting heads using feature refine net and cascaded multi-scale architecture[END_REF] and CrowdHuman dataset [START_REF] Shao | Crowdhuman: A benchmark for detecting human in a crowd[END_REF] provide head annotations of humans in crowded scenes. Head detection is also of significant interest in the crowd counting and analysis literature [START_REF] Idrees | Detecting humans in dense crowds using locally-consistent scale prior and global occlusion reasoning[END_REF]. Rodriguez et al. [START_REF] Rodriguez | Datadriven crowd analysis in videos[END_REF] introduced the idea of tracking by head detection with their dataset consisting of roughly 2200 head annotations. In the recent years, there has been a surge in research works attempting to narrow the gap between detection and crowd counting [START_REF] Liu | Decidenet: Counting varying density crowds through attention guided detection and density estimation[END_REF][START_REF] Deepak | Locate, size and count: Accurately resolving people in dense crowds via detection[END_REF][START_REF] Liu | Point in, box out: Beyond counting persons in crowds[END_REF][START_REF] Zhang | Singleimage crowd counting via multi-column convolutional neural network[END_REF] which attempts to hallucinate pseudo head ground truth bounding boxes in crowded scenes.

Head Detection Methods: Fundamentally, the task of head detection is a combination of multi-scale and contextual object detection problem. Objects at multiple scales are detected based on image pyramids [START_REF] Hu | Finding tiny faces[END_REF][START_REF] Peng | Detecting heads using feature refine net and cascaded multi-scale architecture[END_REF][START_REF] Singh | An analysis of scale invariance in object detection-snip[END_REF][START_REF] Singh | SNIPER: Efficient multi-scale training[END_REF][START_REF] Yoo | EXTD: extremely tiny face detector via iterative filter reuse[END_REF] or feature pyramids [START_REF] He | Mask R-CNN[END_REF][START_REF] Lin | Feature pyramid networks for object detection[END_REF][START_REF] Zhu | Seeing small faces from robust anchor's perspective[END_REF]. The former is computationally intensive task requiring multiple forward passes of images while the latter generates multiple pyramids in a single forward pass. Contextual object detection has been widely addressed in the literature of face detection, such as [START_REF] Deng | Single-stage joint face detection and alignment[END_REF][START_REF] Najibi | SSH: single stage headless face detector[END_REF][START_REF] Tang | Pyramidbox: A context-assisted single shot face detector[END_REF] who show improved detection accuracy by using convolutional filters of larger receptive size to model context. Sun et al. [START_REF] Sun | Scale mapping and dynamic re-detecting in dense head detection[END_REF] employ such a contextual and scaleinvariant applied to head detection.

Tracking Benchmarks and Metrics: The task of Multiple Object Tracking (MOT) is to track an initially unknown number of targets in a video sequence. The first MOT dataset for tracking humans were the PETS dataset [START_REF] Ferryman | Pets2009: Dataset and challenge[END_REF], soon followed by [START_REF] Andriluka | People-tracking-bydetection and people-detection-by-tracking[END_REF][START_REF] Ess | Moving obstacle detection in highly dynamic scenes[END_REF][START_REF] Geiger | Vision meets robotics: The kitti dataset[END_REF][START_REF] Geiger | Are we ready for autonomous driving? the kitti vision benchmark suite[END_REF]. Standardization of MOT benchmarks were later proposed in [START_REF] Leal-Taixé | MOTChallenge 2015: Towards a benchmark for multitarget tracking[END_REF] and since then, it has been updated with yearly challenges with increasingly crowded scenes [START_REF] Dendorfer | CVPR19 tracking and detection challenge: How crowded can it get?[END_REF][START_REF] Milan | MOT16: A benchmark for multi-object tracking[END_REF]. Recently, the TAO dataset [START_REF] Dave | Tao: A large-scale benchmark for tracking any object[END_REF] was introduced for Multi-object tracking, which focuses on tracking 833 object categories across 2907 short sequences. Our dataset pushes the challenge of tracking in crowded environments with pedestrian density reaching 346 humans per frame. Other relevant pedestrian tracking dataset include [START_REF] Cai | Exploring context information for inter-camera multiple target tracking[END_REF][START_REF] Chavdarova | The WILDTRACK multicamera person dataset[END_REF][START_REF] Voigtlaender | MOTS: multi-object tracking and segmentation[END_REF]. To evaluate algorithms on MOTChallenge dataset, classical MOT metrics [START_REF] Wu | Tracking of multiple, partially occluded humans based on static body part detection[END_REF] and CLEAR MOT metrics [START_REF] Bernardin | Evaluating multiple object tracking performance: The CLEAR MOT metrics[END_REF] have been de facto established as standardised way of quantifying performances. The CLEAR Metric proposes two important scores MOTA and MOTP which concisely summarise the classical metrics based on cumulative per frame accuracy and precision of bounding boxes respectively. Recently, Ristani et al. [START_REF] Ristani | Performance measures and a data set for multi-target, multi-camera tracking[END_REF] propose the ID metric, which rewards a tracker based on its efficiency in preserving an identity for the longest duration of the Ground Truth trajectory.

Tracking Algorithms: Online Multi-object tracking algorithms can be summarised into: (i) Detection, (ii) Motion Prediction, (iii) Affinity Computation and, (iv) Association steps. R-CNN based networks have been common choice for the detection stage due to the innate advantage of proposal based detectors over Single-Stage detection methods [START_REF] Huang | Speed/accuracy trade-offs for modern convolutional object detectors[END_REF]. Amongst online Multiple Object Tracking algorithm, Chen et al. [START_REF] Chen | Online multi-object tracking with convolutional neural networks[END_REF] use Particle Filter framework and weigh the importance of each particle by their appearance classification score, computed by a separate network, trained independently. Earlier works such as [START_REF] Breitenstein | Robust tracking-by-detection using a detector confidence particle filter[END_REF][START_REF] Junliang Xing | Multi-object tracking through occlusions by local tracklets filtering and global tracklets association with detection responses[END_REF] use Sequential Importance Sampling (SIS) with Constant Velocity Assumption for assigning importance weight to particles. Henschel et al. [START_REF] Henschel | Improvements to frank-wolfe optimization for multi-detector multi-object tracking[END_REF] demonstrated the the limitation of single object detector for tracking and used a head detector [START_REF] Stewart | End-to-end people detection in crowded scenes[END_REF] in tandem with the pedestrian detector [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF]. However, in the recent past, research works in MOT have attempted to bridge the gap between tracking and detection through a unified framework [START_REF] Bergmann | Tracking without bells and whistles[END_REF][START_REF] Fang | Recurrent autoregressive networks for online multi-object tracking[END_REF][START_REF] Feichtenhofer | Detect to track and track to detect[END_REF][START_REF] Kieritz | Joint detection and online multi-object tracking[END_REF][START_REF] Liu | Looking fast and slow: Memory-guided mobile video object detection[END_REF][START_REF] Voigtlaender | MOTS: multi-object tracking and segmentation[END_REF]. Most notable amongst them is Tracktor [START_REF] Bergmann | Tracking without bells and whistles[END_REF], who demonstrated that an object detector alone is sufficient to predict locations of targets in subsequent frames, benefiting from the high-frame rates in video.

CroHD Dataset

Description: The objective of CroHD is to provide tracking annotation of pedestrian heads in densely populated video sequences. To the best of our knowledge, no such benchmark exists in the community and hence we annotated 2,276,838 human heads in 11,463 frames across 9 sequences of Full-HD resolution. We built CroHD upon 5 sequences from the publicly available MOTChallenge CVPR19 benchmark [START_REF] Dendorfer | CVPR19 tracking and detection challenge: How crowded can it get?[END_REF] to enable performance comparison of trackers in the same scene between two paradigms -head tracking and pedestrian tracking. We maintain the training set and test set classification of the aforementioned sequences to be the same in CroHD as the MOTChallenge CVPR19 benchmark. We further annotated 4 new sequences of higher crowd densities in two new scenarios. The new scenario centers on the Shibuya Train station and Shibuya Crossing, one of the busiest pedestrian crossings in the world. All sequences in CroHD have a framerate of 25f ps and are captured from an elevated viewpoint. The sequences involve crowded indoor and outdoor scenes, recorded across different lighting and environmental conditions. This ensures sufficient diversity in the dataset in order to make it viable for training and evaluating the comprehensiveness of modern Deep Learning based techniques. The maximum pedestrian density reaches approximately 346 persons per frame while the average pedestrian density across the dataset is 178. A detailed sequence-wise summary of CroHD is given in Table 1. We split CroHD into 4 sequences of 5740 frames for training and 5 sequences of 5723 frames for testing. They share three scenes in common, while the fourth scene is disparate to ensure generalization of trackers on this dataset. A representative frame from each sequence of CroHD and their respective training, testing splits are depicted in Figure 2. We will make our sequences and training set annotations publicly available. To preserve the fairness of the MOTChallenge CVPR19 benchmark, we will not release the annotations corresponding to the test set. Annotation: The annotation and data format of CroHD follows the standard guidelines outlined by MOTChallenge benchmark [START_REF] Dendorfer | CVPR19 tracking and detection challenge: How crowded can it get?[END_REF][START_REF] Milan | MOT16: A benchmark for multi-object tracking[END_REF]. We annotated all visible heads of humans in a scene with the visibility left to the best of discretion of annotators. Heads of all humans, whose shoulder is visible were annotated, including the heads occluded by head coverings such as hood, caps etc. For sequences inherited from MOTChallenge CVPR19 benchmark, the annotations were performed independent of pedestrian tracking ground truth in order to have no dependencies between the two modalities. Due to the high frame rate in our video sequences, we interpolate annotations in between keyframes and adjust a track only when necessary. CroHD constitutes four classes -Pedestrian, Person on Vehicle, Static and Ignore. Heads of statues or human faces on clothing have been annotated with an ignore label. Heads of pedestrians on vehicles, wheelchairs or baby transport have been annotated as Person on Vehicle. Pedestrians who do not move throughout the sequence are classified as static persons. Unlike the case of standard MOTChallenge benchmarks, we observe that overlap between bounding boxes are minimal since head bounding boxes from an elevated viewpoint are almost distinct. Hence, we limit our visibility flag to be binary -either visible (1.0) or occluded (0.0). We consider a proposal to be a match if the Intersection Over Union (IoU) with the ground truth is larger than 0.4.

Evaluation Metrics

For evaluation of head detection on CroHD, we follow the standard Multiple Object detection metricsmean Average Precision (mAP), Multiple Object Detection Accuracy (MODA), Multiple Object Detection Precision (MODP) [START_REF] Garofolo | Performance evaluation protocol for face, person and vehicle detection & tracking in video analysis and content extraction (vace-ii) clear -classification of events, activities and relationships[END_REF] and mAP COCO respectively. mAP COCO refers to a stricter metric which computes the mean of AP across IoU thresholds of {50%, 55%, 60%, . . . , 95%}. For evaluation of trackers, we adapt the well established Multiple Object Tracking metrics [START_REF] Bernardin | Evaluating multiple object tracking performance: The CLEAR MOT metrics[END_REF][START_REF] Ristani | Performance measures and a data set for multi-target, multi-camera tracking[END_REF], and extend with the proposed "IDEucl" metric.

IDEucl: While the event based metrics [START_REF] Bernardin | Evaluating multiple object tracking performance: The CLEAR MOT metrics[END_REF] and identity based metric (IDF1) [START_REF] Ristani | Performance measures and a data set for multi-target, multi-camera tracking[END_REF] are persuasive performance indicators of a tracking algorithm from a local and global perspective, they do not quantify the proportion of the ground truth trajectory a tracker in capable of covering. Specifically, existing metrics do not measure the proportion of ground truth trajectory in the image coordinate space a tracker is able to preserve an identity. It is important to quantitatively distinguish between trackers which are more effective in tracking a larger portion of ground truth pedestrian trajectories. This is particularly useful in dense crowds, for better understanding of global crowd motion pattern [START_REF] Dupont | Crowd-11: A dataset for fine grained crowd behaviour analysis[END_REF]. To that end, we propose a new evaluation metric, "IDEucl", which gauges a tracker based on its efficiency in maintaining consistent identity over the length of ground truth trajectory in image coordinate space. Albeit, IDEucl might seem related to the existing IDF1 metric which measures the fraction of frames of a ground truth trajectory in which consistent ID is maintained. In contrast, IDEucl measures the fraction of the distance travelled for which the correct ID is assigned.

To elucidate this difference, consider the example shown in Figure 3. Two trackers A and B compute four different identities for a ground truth trajectory G. Tracker A commits three identity switches in the first 150 frames while maintaining consistent identity for the remaining 150 frames. Tracker B, on the other hand, maintains consistent identity for the first 150 frames but commits 3 identity switches in the latter 150 frames. Our metric reports a score of 0.3 for Tracker A (Figure 3a) and a score of 0.67 for Tracker B (Figure 3b). Meanwhile, IDF1 and the classical metric reports a score of "0.5" and "3 identity switches" respectively for both the trackers. Following existing metrics, Tracker A and Tracker B are considered equally efficient. They neither highlight the ineffectiveness of Tracker A nor the ability of Tracker B in covering an adequate portion of ground truth trajectory with consistent identity. Therefore, IDEucl is more appropriate for judging the quality of the estimated pedestrian motion.

Thus, to formulate this metric, we perform a global hypothesis to ground truth matching by constructing a Bipartite Graph G = (U, V, E), similar to [START_REF] Ristani | Performance measures and a data set for multi-target, multi-camera tracking[END_REF]. Two "regular" nodes are connected by an edge e if they overlap in time, with the overlap defined by ∆,

∆ t,t-1 = 1, if δ > 0.5 0, otherwise (1)
Considering τ t , h t to be an arbitrary ground truth and hypothesis track at time t, δ is defined as,

δ = IoU(τ t , h t ) (2)
The cost on each edge E ∈ R N of this graph, M ∈ R N -1 is represented as the distance in image space between two successive temporal associations of "regular" node. In particular, cost of an edge is defined as ,

M = N t=1 m t = d(τ t , τ t-1 ), if ∆ t,t-1 = 1. 0, otherwise (3)
where d denotes the Euclidean distance in image coordinate space. A ground truth trajectory is assigned a unique hypothesis which maintains a consistent identity for the predominant distance of ground truth in image coordinate space. We employ the Hungarian algorithm to solve this maximum weight matching problem to obtain the best (longest) hypothesis. Once we obtain an optimal hypothesis, we formulate the metric C as the ratio of length of ground truth in image coordinates covered by the best hypothesis,

C = K i=1 M i K i=1 T i (4)
Note that this formulation of cost function naturally weighs the significance of each ground truth track based on its distance in image coordinate space.

Methods : Head Detection and Tracking

In this section, we elucidate the design and working of HeadHunter and HeadHunter-T.

HeadHunter

As detection is the pivotal step in object tracking, we designed HeadHunter differently from traditional object detectors [START_REF] Felzenszwalb | Object detection with discriminatively trained part-based models[END_REF][START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF][START_REF] Yang | Exploit all the layers: Fast and accurate cnn object detector with scale dependent pooling and cascaded rejection classifiers[END_REF] by taking into account the nature and size of objects we detect. HeadHunter is an end-to-end two stage detector, with three functional characteristics. First, it extracts feature at multiple scales using Feature Pyramid Network (FPN) [START_REF] Lin | Feature pyramid networks for object detection[END_REF] using a Resnet-50 [START_REF] He | Deep residual learning for image recognition[END_REF] backbone. Images of heads are homogeneous in appearance and often, in crowded scenes, resemble extraneous objects (typically background). For that reason, inspired by the head detection literature, we augmented on top of each individual FPNs, a Context-sensitive Prediction Module (CPM) [START_REF] Tang | Pyramidbox: A context-assisted single shot face detector[END_REF]. This contextual module consists of 4 Inception-ResNet-A blocks [START_REF] Szegedy | Inception-v4, inception-resnet and the impact of residual connections on learning[END_REF] with 128 and 256 filters for 3 × 3 convolution and 1024 filters for 1 × 1 convolution. As detecting pedestrian heads in crowded scenes is a problem of detecting many small-sized adjacently placed objects, we used Transpose Convolution on features across all pyramid levels to upscale the spatial resolution of each feature map. Finally, we used a Faster-RCNN head with Region Proposal Network (RPN) generating object proposals while the regression and classification head, each providing location offsets and confidence scores respectively. The architecture of our proposed network is summarised in Figure 4.

HeadHunter-T

We extended HeadHunter with two motion models and a color histogram based re-identification module for headtracking. Our motion models consist of Particle Filter to predict motion of targets and Enhanced Correlation Coefficient Maximization [START_REF] Evangelidis | Parametric image alignment using enhanced correlation coefficient maximization[END_REF] to compensate the Camera motion in the sequence. A Particle Filter is a Sequential Monte Carlo (SMC) process, which recursively estimates the state of dynamic systems. In our implementation, we represent the posterior density function by a set of bounding box proposals for each target, referred to as particles. The use of Particle Filter enables us to simultaneously model nonlinearity in motion occurring due to rapid movements of heads and pedestrian displacement across frames.

Notation: Given a video sequence I, we denote the ordered set of frames in it as {I 0 , • • • , I T -1 }, where T is the total number of frames in the sequence. Throughout the paper, we use subscript notation to represent time instance in a video sequence. In a frame I t at time t, the active tracks are denoted by T t = {b Particle Initialization: New tracks are initialized at the start of the sequence, I 0 from the detection provided by HeadHunter and at frame I t for detection(s) which cannot be associated with an existing track. A plausible association of new detection with existing track is resolved by Non-Maximal-Suppression (NMS). The importance weights of each particle are set to be equal at the time of initialisation. Each particles represent 4 dimensional state space, with the state of each targets modelled as (x c , y c , w, h, ẋc , ẏc , ẇ, ḣ), where, (x c , y c , w, h) denote the centroids, width and the height of bounding boxes.

Prediction and Update: At time t > 0, we perform RoI pooling on the current frame's feature map, F t , with the bounding box of particles corresponding to active tracks. Each particles' location in the current frame is then adjusted using the regression head of HeadHunter, given their location in the previous frame. The importance weights of each particle are set to their respective foreground classification score from the classification head of HeadHunter. Our prediction step is similar to the Tracktor [START_REF] Bergmann | Tracking without bells and whistles[END_REF], applied to particles instead of tracks. Given the new location and importance weight of each particle, estimated position of k th track is computed as weighted mean of the particles,

S k t = 1 M M i=1 w k,i t p k,i t (5)
Resampling: Particle Filtering frameworks are known to suffer from degeneracy problems [START_REF] Arulampalam | A tutorial on particle filters for online nonlinear/non-gaussian bayesian tracking[END_REF] and as a result we resample to replace particles of low importance weight. M particles corresponding to k th track are re-sampled when the number of particles which meaningfully contributes to probability distribution of location of each head, Nk eff exceeds a threshold, where, Cost Matching: Tracks are set to inactive when scores of their estimated state S a t falls below a threshold, λ reg nms . Positions of such tracks are predicted following Constant Velocity Assumption (CVA) and their tracking is resumed if it has a convincing similarity with a newly detected track. The similarity, C is defined as

Nk eff = 1 M i=1 (w k,i ) 2 (6)
C = α • IoU (L i t , N j t ) + β • d 1 (L i t , N j t ) (7)
where L i t and N j t are the i th lost track and j th new track respectively. And, d 1 denotes the Bhattacharyya distance between the respective color histograms in the HSV space [START_REF] Nummiaro | Object tracking with an adaptive color-based particle filter[END_REF]. Once tracks are re-identified, we re-initialize particles around its new position.

Experiments

HeadHunter

We first detail the experimental setup and analyse the performance of HeadHunter on two datasets -SCUT-HEAD [START_REF] Peng | Detecting heads using feature refine net and cascaded multi-scale architecture[END_REF] and CroHD respectively. For the Faster-RCNN head of HeadHunter, we used 8 anchors, whose sizes were obtained by performing K-means over ground truth bounding boxes from the training set. To avoid overlapping anchors, they were split equally across the four pyramid levels, with the stride of anchors given by max(16, s/d) where s is square-root of the area of an anchor-box and d is the scaling factor [START_REF] Najibi | Fa-rpn: Floating region proposals for face detection[END_REF]. For all experiments, we used Online Hard Example Mining [START_REF] Shrivastava | Training region-based object detectors with online hard example mining[END_REF] with 1000 proposals and a batch size of 512.

SCUT-Head is a large-scale head detection dataset consisting of 4405 images and 111,251 annotated heads split across Part A and Part B. We trained HeadHunter for 20 epochs with the input resolution to be the median image resolution of the training set (1000x600 pixels) and an initial learning rate of 0.01 halved at 5th, 10th and 15th epochs respectively. For a fair comparison, we trained HeadHunter only on the training set of this dataset and do not use any pre-trained models. We summarize the quantitative comparisons with other head detectors on this dataset in Table 2. HeadHunter outperforms other state-of-the-art head detectors based on Precision, Recall and F1 scores.

Methods

Precision Recall F1 Faster-RCNN [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF] 0.87 0.80 0.83 R-FCN+FRN [START_REF] Peng | Detecting heads using feature refine net and cascaded multi-scale architecture[END_REF] 0.91 0.84 0.87 SMD [START_REF] Sun | Scale mapping and dynamic re-detecting in dense head detection[END_REF] 0.93 0.90 0.91 HSFA2Net [START_REF] Shen | An indoor crowd detection network framework based on feature aggregation module and hybrid attention selection module[END_REF] 0.94 0.92 0.93 HeadHunter (Ours) 0.95 0.93 0.94 Table 2. Comparison between HeadHunter's and other state-ofthe-art head detectors on the SCUT-Head dataset.

CroHD: We first trained HeadHunter on the combination of training set images from SCUT-HEAD dataset and CrowdHuman dataset [START_REF] Shao | Crowdhuman: A benchmark for detecting human in a crowd[END_REF] for 20 epochs at a learning rate of 0.001. With variations well characterized, pre-training on large-scale image dataset improves the robustness of head detection. We then fine-tuned HeadHunter on the training set of CroHD, for a total of 25 epochs with an initial learning rate of 0.0001 using the ADAM optimizer [START_REF] Diederik | Adam: A method for stochastic optimization[END_REF]. The learning rate is then decreased by a factor of 0.1 at 10th and 20th epochs respectively.

Ablation: We examined our design choices for Head-Hunter, namely the use of context module and the anchor selection strategy by removing them. The head detection performance of HeadHunter and its variants on CroHD are summarised in Table 3. We threshold the minimum confidence of detection to 0.5 for evaluation. W/O Cont refers to the HeadHunter without Context Module. We further removed the median anchor sampling strategy and refer to as W/O Cont, mAn. We also provide baseline performance of Faster-RCNN with Resnet-50 backbone on CroHD, the object detector upon which we built HeadHunter. We followed the same training strategy for Faster-RCNN as Head-Hunter. All variants of HeadHunter significantly outperformed Faster-RCNN. Inclusion of the context module and the anchor initialisation strategy also has a noteworthy impact on head detection.

Method

Precision Recall F1 MODA MODP mAP COCO Faster-RCNN [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF]

HeadHunter-T

For the Particle Filtering framework, we used a maximum of N=100 particles for each object. The N particles were uniformly placed around the initial bounding box. To ensure that particles were not spread immoderately and were distinct enough, we sampled particles from a Uniform distribution whose lower and upper limit were ((x -1.5w, y -1.5h), (x + 1.5w, y + -1.5h)) respectively. Where, x, y, w, h denote the centroid, width and height of the initial bounding box. For the color based re-identification, we used 16, 16 and 8 bins for the H, S and V channels respectively, where the brightness invariant Hue [START_REF] Finlayson | Hue that is invariant to brightness and gamma[END_REF] was used instead of the standard Hue. α, β, which denotes the importance of IoU and color histogram matching, corresponding to Equation 7were set to 0.8 and 0.2 respectively. We deactivated a track if it remained inactive for λ age = 25 frames or if its motion prediction falls outside the image coordinates.

We evaluated three state-of-the-art trackers on CroHD, namely, SORT [START_REF] Bewley | Simple online and realtime tracking[END_REF], V-IOU [START_REF] Bochinski | Extending iou based multi-object tracking by visual information[END_REF] and Tracktor [START_REF] Bergmann | Tracking without bells and whistles[END_REF] to compare with HeadHunter-T. We chose methods which do not require any tracking specific training, whose implementations have been made publicly available and are top-performing on the crowded MOTChallenge CVPR19 benchmark [START_REF] Dendorfer | CVPR19 tracking and detection challenge: How crowded can it get?[END_REF]. For a fair comparison, we performed all experiments with head detection provided by HeadHunter, thresholded to a minimum confidence of 0.6. SORT is an online tracker, which uses a Kalman Filter motion model and temporally associates detection based on IoU matching and Hungarian Algorithm. V IOU associates detection based on IoU matching and employs visual information to reduce tracking inconsistencies due to missing detection. Parameters for V IOU and SORT were set based on fine-tuning on the training set of CroHD, as discussed in the supplementary material. We evaluated two variants of Tracktor, with and without motion model. Tracktor+MM denotes the Tracktor extended with Camera Motion Compensation [START_REF] Evangelidis | Parametric image alignment using enhanced correlation coefficient maximization[END_REF] and CVA for inactive tracks. For the two versions of Tracktor, we set tracking parameters similar to HeadHunter. tity preserved tracking. Although Tracktor [START_REF] Bergmann | Tracking without bells and whistles[END_REF] is similar to HeadHunter-T, there is a noticeable difference in its head tracking performance. We hypothesize the use of Particle Filter framework, which can handle arbitrary posteriors, as the reason for improvement. This claim is justified in the forthcoming section.

Ablation Experiments

HeadHunter-T: In this section, we analyse the design choices, in particular, the utility of re-identification module and Particle Filter of HeadHunter-T on the training set of CroHD. The results are summarised in Table 5. For variations in motion model, we removed the Particle Filter and used simple Camera Motion Compensation, denoted as HT w/o PF. We also experimented with a reduced number of particles initialized around the head, with n=10, denoted as HT + 10F. Introducing Particle Filter noticeably improved identity preserving scores (IDF1 and IDEucl) for HT + 10F. Further increasing the number of filters to 100 demonstrated the best performance. However, using more than 100 filters resulted in either duplicates or immoderate spreading, which are undesirable. We removed the re-identification module, to understand its influence, denoted as w/o ReID. Although color histogram is a modest image descriptor, yet it drastically reduced the number of identity switches and showed superior performance in identity preserving metrics -IDEucl, IDF1. We also experimented with α and β values corresponding to the importance of IoU and histogram matching (Equation 7). We set β to 0.8 and α to 0.2 and this configuration is denoted as HT + sReID. Surprisingly, we observed more identity switches and a slight decrease in performance across other tracking metrics. HeadHunter-T, our final model, outperformed all the other variants.

Choice of Filter: To further substantiate our choice of a multi-modal filter, we replaced the Particle Filter of HeadHunter-T with a Kalman Filter motion model [START_REF] Kalman | A new approach to linear filtering and prediction problems" transaction of the asme journal of basic[END_REF]. While both Kalman Filter and Particle Filter are recursive state estimation algorithms, Kalman Filter assumes the system to be linear with Gaussian noise [START_REF] Arulampalam | A tutorial on particle filters for online nonlinear/non-gaussian bayesian tracking[END_REF] inter-frame displacement of bounding boxes with CVA. The four states are x, y centroid coordinates, the height and aspect ratio of bounding boxes respectively, similar to the SORT [START_REF] Bewley | Simple online and realtime tracking[END_REF]. The performance of this tracker, denoted as HT + KF is summarised in Table 5. HeadHunter-T with Particle Filter demonstrates superior performance than its Kalman Filter variant with respect to all the tracking metrics reported and in particular, we observe major improvement in-terms of IDEucl metric. Motion of heads along with the pedestrian displacement induces non-linearities in the position of bounding boxes. Although pedestrian motion in general is non-linear, this issue is exacerbated with the small size of head bounding boxes. Hence, using a multimodal posterior state estimation is necessary to address the perceptible impact of non-linear motion. We remark this to be the reason behind improvement in performance while using a Particle Filter in comparison to the Kalman Filter. Comparison across paradigm : We compare pedestrian and head tracking performances on the common sequences between CroHD and MOTChallenge CVPR19 dataset. The sequence being the same ensures that trackers are evaluated on full body and head bounding boxes of the same pedestrians in the scene. For this comparison, we chose published state-of-the-art methods on the aforementioned dataset, namely, Tracktor++ [START_REF] Bergmann | Tracking without bells and whistles[END_REF], V IOU, DD TMA [START_REF] Yoon | Online multiple pedestrians tracking using deep temporal appearance matching association[END_REF] and HAM HI [START_REF] Yoon | Online multi-object tracking with historical appearance matching and scene adaptive detection filtering[END_REF]. We performed comparison in-terms of MOTA, IDF1, MT (Mostly Tracked in percentage) metrics. Since we used a different object detector than the rest, a straightforward comparison between performance metrics would not be fair. Hence, for each sequence, we measure the ratio of aforementioned performance metrics with their object detector's MODA score to obtain the scaled scores -s-MOTA, s-IDF1 and s-MT. The scaled scores, averaged across five common sequences are illustrated in Figure 5. Our approach substantially outperforms other methods indicating that tracking by head detection is more suited for tracking in environments involving high pedestrian density where preserving identity is important. It is also worthy to note that HeadHunter uses a ResNet-50 backbone in contrast to a Resnet-101 backbone used by other methods. Furthermore, Tracktor++, HAM HI and DD TMA all use Deep Networks for extracting appearance features, while HeadHunter-T uses a color histogram based appearance feature. By compromising our tracking space (size of bounding box) to avoid mutual occlusion, we observe notable performance gain and significantly reduce the computation cost. This suggests that tracking by head detection paradigm is more desirable for real-time tracking applications focused on identity preservation.

Conclusion

To advance algorithms to track pedestrians in dense crowds, we introduced a new dataset, CroHD, for tracking by head detection. To further quantify the efficacy of a tracker in describing pedestrian motion, we introduced a new metric, IDEucl. We developed two new baseline methods, HeadHunter, HeadHunter-T for head detection and head tracking on CroHD respectively. We demonstrated HeadHunter-T to be consistently more reliable for identity preserving tracking applications than existing state-of-theart trackers adapted for head tracking. Additionally, the adequacy of HeadHunter-T's performance with a modest computational complexity, opens up opportunities for future research focused on tracking methods adapted for low computational complexity and real-time applications. We also hope that CroHD will serve useful in contiguous fields, such as Crowd Counting and Crowd Motion Analysis.

Figure 2 .

 2 Figure 2. Depiction of a frame per each scene from our Crowd of Heads Dataset, CroHD. The top row shows frames from the training set while the bottom row illustrates frames from the test set. Scene 5 is kept exclusive for the test set.

Figure 3 .

 3 Figure 3. Identity prediction of two trackers -Tracker A (3a) and Tracker B (3b) for the same ground truth. A change of color implies an identity switch with both trackers registering 3 switches.

Figure 4 .

 4 Figure 4. An overview of the architecture of our proposed head detector, Headhunter. We augment the features extracted using FPN (C4. . . P4) with Context Sensitive feature extractor followed by series of transpose convolutions to enhance spatial resolution of feature maps. Cls and Reg denote the Classification and Regression branches of Faster-RCNN [48] respectively.

Figure 5 .

 5 Figure 5. Comparison between HeadHunter-T and state-of-theart trackers on common sequences of CroHD and MOTChallenge benchmark [13]. s-MOTA, s-IDF1, s-MT are scaled version of MOTA, IDF1 and Most Tracked (MT) metrics respectively.

Table 1 .

 1 Sequence-wise statistics CroHD. Sequences are named CroHD-XY, with X being either 0 or 1 depending on training set or testing set respectively. Y denotes the serial number of videos.

	Name	Frames	Scenario	Tracks Boxes Density
	CroHD-01 429	Indoor	85	21,456	50.0
	CroHD-02 3,315 Outdoor, night 1,276 733,622 222.0
	CroHD-03 1,000 Outdoor, day 811	258,012 258.0
	CroHD-04 997	Indoor	580	175,703 176.2
	CroHD-11 584	Indoor	133	38,492	65.8
	CroHD-12 2,080 Outdoor, night 737	383,677 185.0
	CroHD-13 1,000 Outdoor, day 725	257,828 258.0
	CroHD-14 1,050 Outdoor, day 562	258,227 246.0
	CroHD-15 1,008 Outdoor, day 321	149,821 149.0
	Total	11,463		5,230 2,276,838 178

  At time t, the i th particle corresponding to k th track is denoted by p k,i t and its respective importance weight by w k,i t . L t and N t denote the set of inactive tracks and newly initialized tracks respectively.

	1 t , b 2 t , . . . , b N t }, where b k t refers to bounding box of the k th active track, denoted as b k t = x k t , y k t , w k t , h k t .

Table 3 .

 3 Summary of various head detector's performances on the test set of CroHD.

		34.4	42.2 50.1 40.3	30.8	11.2
	W/O Cont, mAn	40.9	50.8 57.8 38.1	37.8	14.4
	W/O Cont	44.3	57.8 64.5 40.0	42.7	15.0
	HeadHunter	52.8	63.4 68.3 50.0	47.0	19.7

Table 6 .

 6 2 summarises the performance of aforementioned methods on the test set of CroHD. HeadHunter-T outperforms all the other methods, and furthermore demonstrates superiority in iden-

	Method	MOTA ↑ IDEucl ↑ IDF1 ↑ MT ↑ ML ↓ ID Sw. ↓
	SORT [5]	46.4	58.0	48.4 49 216	649
	V IOU [6]	53.4	34.3	35.4 80 182 1890
	Tracktor [3]	58.9	31.8	38.5 125 117 3474
	Tracktor+MM [3] 61.7	44.2	45.0 141 104 2186
	HeadHunter-T	63.6	60.3	57.1 146 93	892

Table 4 .

 4 Main tracking result comparing the performances of various state-of-the-art trackers and HeadHunter-T on the test set of CroHD. The direction of arrows indicate smaller or larger desired value for the metric.

Table 5 .

 5 while Particle Filter's multimodal posterior distribution enables it to model states of nonlinear systems. We replaced the Particle Filter with a four state Kalman Filter to model the Method MOTA ↑ IDEucl ↑ IDF1 ↑ MT ↑ ML ↓ ID Sw. ↓ Illustration of ablation studies of HeadHunter-T (denoted as HT) on the training set of CroHD. The direction of arrows indicate small or large desired metric values.

	HT w/o PF	60.6	40.1	43.9 200 102 3652
	HT + 10F	63.3	58.2	56.3 214 98	1534
	HT w/o ReID	59.5	57.7	57.5 225 91	1411
	HT + sReID	59.1	57.8	58.3 225 91	1280
	HT + KF	63.4	53.8	55.9 214 93	2451
	HeadHunter-T 64.0	61.5	58.5 225 91	1247

https://motchallenge.net/results/MOT20/

Acknowledgements: We are thankful to Dr. Vicky Kalogeiton and Prof. Dr. Bastian Leibe for their insightful feedback. We are also thankful to our annotators for their hard work. This work has received funding from the European Union's Horizon 2020 research and innovation programme under the grant agreement No. 899739. Experiments were performed on the IGRIDA computing cluster.


Patrice Meunier

Nicolas Boulanger

Xavier Riedinger

Stéphane Le Dizès

Instabilities of a columnar vortex in a stratified fluid

In this paper, we show that the presence of a vertical density gradient can destabilise a Gaussian vortex in two different ways. First, when the vortex axis is slightly tilted with respect to the vertical and if its Froude number is larger than one, the vortex exhibits strong axial flows inside a critical layer, which lead to a violent three-dimensional instability. Second, when the vortex axis is vertical, the vortex is weakly unstable by emission of internal gravity waves. Theoretical predictions for both instabilities are compared with experimental evidences.

Introduction

It is well known that rotating and stratified turbulence contains many vortices whose axes tend to align with the background rotation axis and the density gradient direction. It is thus of interest to study the stability of such vortices to gain comprehension on the mechanisms which drive rotating and stratified turbulence. As a first step, we will focus here on the effect of stratification only.

There exist many inviscid instabilities which can destabilise a vortex in a homogeneous fluid. When the velocity profile of the vortex satisfies the Rayleigh criterion, the vortex is subject to a centrifugal instability. Otherwise, all the eigenmodes of the vortex (called vortex Kelvin modes) are neutral and can only grow due to an external forcing. For example, the presence of a second vortex can induce two different instabilities: the large wavelength Crow instability [START_REF] Crow | Stability theory for a pair of trailing vortices[END_REF] and the short wavelength elliptical instability [START_REF] Kerswell | Elliptical instability[END_REF].

The presence of a stable stratification modifies the characteristics of these three instabilities. The Rayleigh criterion is not modified by stratification but the instability is weaker as soon as viscous effects are also considered [START_REF] Hopfinger | Vortices in rotating fluids[END_REF]. The Crow instability Patrice Meunier Institut de Recherche sur les Phenomenes Hors Equilibre, CNRS / Aix-Marseille Universites, 49 rue Joliot-Curie, 13384 Marseille, France, e-mail: meunier@irphe.univ-mrs.fr is transformed into a zig-zag instability [1; 6]. And the elliptical instability is weakened due to the appearance of critical points in the vortex Kelvin modes [START_REF] Dizès | Inviscid waves on a lamb-oseen vortex in a rotating stratified fluid: consequences for the elliptic instability[END_REF].

We show in this paper that the presence of stratification can also destabilise a centrifugally stable vortex. After describing the vortex and the density profile in section 2, we first focus on the effect of a tilt of the vortex in section 3. We then study the instabilities of a vertical vortex in section 4 and finally conclude in section 5.

A Gaussian vortex in a stratified fluid

We wish to study the stability of a Gaussian vortex, whose angular velocity profile is defined as:

Ω (r) = v θ r = Γ 2πr 2 (1 -e -r 2 /a 2 ) (1)
for a vortex of circulation Γ and core size a. Such a vortex is stable with respect to the centrifugal instability in a non-stratified fluid. It is also interesting because it is the attractive solution for any localised patch of vorticity evolving due to viscous diffusion. This vortex is placed in a vertical density stratification ρ(z), which is assumed to be linear and characterized by a constant buoyancy frequency N = -(g/ρ)(∂ ρ/∂ z). We consider that the vortex axis may be tilted with respect to the vertical with an angle α, as shown on Fig. 1(a).

This system is characterised by two non-dimensional parameters: the Froude number F = Γ /2πNa 2 and the Reynolds number Re = Γ /ν. Experimentally, it is possible to create a Gaussian vortex by rotating a flap impulsively in a fluid initially at rest. Some PIV measurements (shown on Fig. 1b) allow to check that the vortex is circular and close to Gaussian. The linear stratification of the fluid is obtained using the two-tank method with salted water.

The linear temporal stability characteristics of the vertical vortex have also been obtained by a numerical spectral method.

Instabilities of a tilted vortex

We show in this section that the presence of a small tilt angle α can drastically change the structure and the stability of the vortex when its Froude number is larger than one.

Spatial structure of a tilted vortex

It has been shown [START_REF] Boulanger | Structure of a stratified tilted vortex[END_REF] that tilting the vortex induces an axial velocity w proportional to the tilt angle α for small angles:

w = α rΩ (r) 3 Ω (r) -N 2 sin θ ( 2
)
where θ is the angle with respect to the direction of the tilt. At low Froude number, the vertical velocity decreases as 1/N 2 and thus vanishes, leading to the same model of vortex described by [START_REF] Canals | Tilted baroclinic tidal vortices[END_REF], where the circular streamlines are horizontal but with their centers aligned on a tilted axis. However, when the Froude number is larger than 1, this solution diverges at the critical radius r c where the angular velocity Ω (r c ) equals the buoyancy frequency N. This creates a critical layer with strong axial velocities.

Theoretically, the addition of the viscous terms smoothes the singularity in a viscous critical layer of width Re -1/3 . In this layer, the axial velocity reaches values of order αRe 1/3 . In Fig. 2, the analytic solution is compared to the experimental PIV measurement of the axial velocity. A good qualitative and quantitative agreement is observed. A cut along the direction of the tilt (x-axis) indicates that the velocity profile contains a strong positive jet on one side and a strong negative jet on the other side. On the contrary, the velocity profile in the direction perpendicular to the tilt (y-axis) exhibits two shear layers on each side of the vortex center.

Tilt-induced instabilities

The presence of strong axial flows inside the critical layers can destabilise the vortex. Indeed, when the Reynolds number or the tilt angle is increased, the vortex suddenly exhibits strong three-dimensional perturbations which lead to a turbulent flow inside the vortex. Such perturbations are shown on the photographs of Fig. 3 in two different sideviews.

Fig. 3(a) shows a shadowgraph visualisation in the vertical longitudinal plane (i.e. the plane of the figure in Fig. 1), where the vortex thus appears tilted. The visualisation highlights a row of clockwise vortices on each side of the vortex. These vortices correspond to Kelvin-Helmholtz billows issued from a Kelvin-Helmholtz instability of the shear present in thes critical layer. In the tilted longitudinal plane (i.e. normal to the plane of the figure in Fig. 1), on the other hand, the perturbation seems to contain a row of alternate vortices on each side (see Fig. 3b), which is due to the sinuous mode of a jet instability.

These mechanisms are further confirmed by PIV measurements in these two planes shown in Fig. 4. In the vertical longitudinal plane (Fig. 4a), the jet instability made of alternate vortices is clearly visible on each side of the vortex axis. In the tilted longitudinal plane (Fig. 4b), each mixing layer destabilizes in a row of vortices. It can be noted that the different instabilities do not correspond between visualisations and PIV measurements: they are in phase quadrature. This is due to the fact that density structures are advected by the vortex and they are thus maximal at the end of the forcing by the flow, and not in the plane where the forcing is maximal. A local stability analysis has been performed [START_REF] Boulanger | Tilt-induced instability of a stratified vortex[END_REF] and proves that the threedimensional instability of the vortex is indeed linked to a local shear or jet instability inside the critical layer. The growth rate (of the order of 0.1αRe 2/3 N) and the wavelength are in good quantitative agreement with the observations. Moreover, a qualitative criterion for the stability of a tilted stratified vortex can be obtained by saying that the growth rate must be larger than the angular frequency of the vortex. The vortex is thus unstable when αRe 2/3 1 and for F > 1.

Consequences

These local instabilities create a strong turbulent flow inside the vortex. However, the vortex is not broken and its circulation remains constant even in the presence of instabilities. The only effect of the instabilities on the global properties of the vortex is to increase its core size a if the Froude number F is slightly larger than 1. This can be understood by the fact that for Froude of order 1, the critical layer is close to the vortex center and the mixing created by the instabilities can thus diffuse the vorticity more efficiently than for very large Froude number (for which the critical layer is far from the vortex center). In the presence of instabilities, the vortex core size increases, leading to a decrease of the Froude number until it becomes smaller than one at which the critical layer suddenly disappears.

Since the control parameter of the instability is αRe 2/3 , for the large Reynolds numbers of geophysical flow, vortices are expected to be unstable for very small tilt angle. This instability should create a deficit of vortices with a Froude number slightly above one, which could be interesting to analyse in oceanic and atmospheric data.. The strong axial velocity created in the critical layer may also be a visible effect. It could have an important impact on the diapycnal mixing in geophysical flows. It would be interesting to quantify its real contribution in oceans.

Radiative instability

We now focus on the stability of a Gaussian vortex whose axis is aligned with the density gradient (α = 0).

Linear stability analysis

A linear stability analysis has been performed using a pseudo-spectral Chebyshev code. The boundary conditions, especially the condition of radiation at infinity, have been implemented by considering an integration in the complex plane. Further explanations on the numerical technique can be obtained in [START_REF] Riedinger | Viscous stability properties of a Lamb-Oseen vortex in a stratified fluid[END_REF].

Figure 5(a) shows the growth rate of the most unstable eigenmode as functions of Froude and Reynolds numbers. Its spatial structure exhibits two opposite sign vortices at the center surrounded by weak edges spiralling around the vortex, associated with internal gravity waves. This instablity is due to the gravity wave emission, as described in [START_REF] Ford | The response of a rotating ellipse of uniform potential vorticity to gravity wave radiation[END_REF] and [START_REF] Plougonven | Internal gravity wave emission from a pancake vortex: an example of wave-vortex interaction in strongly stratified flows[END_REF]. It can be understood as an interaction between the vortex Kelvin mode (creating the two central vortices) and internal gravity waves which radiate to infinity. Such a description has also been given in [START_REF] Schecter | A theory for the vertical alignment of a quasigeostrophic vortex[END_REF].

It can be noted that the vortex is unstable for any Reynolds number and is maximal for Froude numbers around 1. The damping of the instability at large Froude number is due to the appearance of a critical point in the eigenmode. The damping of the instability at small Froude numbers is due to viscous effects: the wavelength scales as F for small Froude numbers, and the viscous damping rate νk 2 thus increases as F -2 when F tends to zero. These two effects explain why the instability is only present on a band of Froude numbers around one.

Experimental evidence?

Experiments have been performed to try to observe this instability. However, the growth rate of this instability is very small (close to 0.1% of the maximal angular velocity of the vortex) and the experimental results are not clear. The vortex center, visualised by a dye, usually exhibits an undulation as shown on Fig. 6 in the unstable region of the parameter space. The wavelength of the undulation is of the order of the most unstable wavelength found theoretically although it is usually larger. For the parameters of Fig. 6, the wavelength was measured to be twice larger than the predicted wavelength.

Conclusion

We have shown that the presence of stratification can drastically change the stability properties of a Gaussian vortex. Two cases have been studied. We have demonstrated that the presence of small tilt between the vortex axis and the density gradients can induce strong axial flows inside a critical layer, which lead to three-dimensional perturbations due to a Kelvin-Helmholtz or a jet instability. These instabilities are extremely violent and appear for very small tilt angles at large Reynolds numbers. Theoretical predictions are confirmed by experimental visualisations and measurements.

In a second study, we have shown theoretically that a vertical Gaussian vortex is unstable in a stably-stratified fluid whatever the Reynolds number in an interval of Froude number around one. The unstable mode creates an undulation of the vortex centerline and is associated with an emission of internal gravity waves far from the vortex center.

We believe that these results are interesting for geophysical applications because they show that the coherent structures found in stratified turbulence may not be stable, especially at high Reynolds numbers. Indeed, it is natural to think that vortices in a turbulent flow might experience a small tilt angle due to a small vertical shear or to the presence of weak internal waves. These vortices would thus be subject to the tilt-induced instability which is extremely powerful, even for small tilt angles. This would mean that stratified turbulence is not similar to a two-dimensional turbulence composed of horizontal vortices but that it is composed of intense three-dimensional motions at least for Froude numbers around 1. Another interesting consequence is that there is a coupling between vortices and internal gravity waves since vortices can emit internal gravity waves.

Fig. 1

 1 Fig. 1 (a) Sketch of a vortex inclined of an angle α with respect to a stable linear density stratfication ρ(z). (b) Experimental velocity field (and vorticity contours) obtained for a vortex generated by a flat plate impulsively rotated in a stratified fluid measured by PIV.

Fig. 2

 2 Fig. 2 (a) Theoretical and (b) experimental contours of axial velocity in a horizontal plane, for a vortex tilted with an angle α = 0.07 rad. Re = 450 and F = 3.2.

Fig. 3

 3 Fig. 3 Shadowgraph visualisation of the tilt-induced instabilities in (a) the vertical longitudinal plane and (b) the tilted longitudinal plane. F = 3, Re = 720, α = 0.07rad.

Fig. 4

 4 Fig. 4 PIV measurements of the instantaneous azimuthal vorticity showing the tilt-induced instabilities in (a) the vertical longitudinal plane and (b) the tilted longitudinal plane.

Fig. 5

 5 Fig. 5 (a) Growth rate of the most unstable mode for a vertical vortex in a stratified fluid. (b) Structure of the most unstable mode for Re = 10 5 and F = 1.5.

Fig. 6

 6 Fig. 6 Dye visualisation of the instability of a vertical vortex in a stratified fluid (sideview). F = 1.15, Re = 610.

Acknowledgements

The support of ACI grant 'Prévention des catastrophes naturelles' and ANR grant 'Fronts, secondary cyclones and gravity waves: instability, mixing and transport' by the french ministry of research is most gratefully acknowledged.