text:
# Janossy Densities I. Determinantal Ensembles. 11footnote 1 AMS 2000 subject classification: 15A52; 60G55 keywords and phrases: random matrices, orthogonal polynomials, Janossy densities, Riemann-Hilbert problem

## 1 Introduction

We consider an ensemble of $`n`$ particles on a measure space $`(X,\mu )`$ with the joint distribution density (with respect to the product measure $`\mu ^n`$) given by the formula

$$p(x_1,\mathrm{},x_n)=const_ndet(\varphi _j(x_k))_{j,k=1,\mathrm{},n}det(\psi _j(x_k))_{j,k=1,\mathrm{},n}.$$
(1)
Here $`\varphi _k(x),\psi _k(x),k=1,\mathrm{},n`$, are some functions on $`X`$ and $`const_n`$ is the normalization constant

$`const_n^1`$ $`=`$ $`{\displaystyle _{X^n}}det(\varphi _j(x_k))_{j,k=1,\mathrm{},n}det(\psi _j(x_k))_{j,k=1,\mathrm{},n}{\displaystyle \underset{j=1,\mathrm{},n}{}}\mu (dx_j)`$ (2)
$`=`$ $`n!det\left({\displaystyle _X}\varphi _i(x)\psi _j(x)\mu (dx)\right)_{i,j=1,\mathrm{},n}`$
where $`X^n=X\times \mathrm{}\times X`$ ($`n`$ times). Ensembles of this form were introduced in , . In the special cases when $`X=`$, $`\varphi _i=\psi _i=x^{i1}`$, and $`X=\{z||z|=1\}`$, $`\varphi _i=\overline{\psi }_i=z^{i1}`$, such ensembles were extensively studied in Random Matrix Theory much earlier under the general name unitary ensembles, see  for details. An example of the form (1) which is different from random matrix ensembles was considered in .

Let us assume that we can biorthogonalize $`\{\varphi _j\}_{j=1,\mathrm{},n}`$ and $`\{\psi _j\}_{j=1,\mathrm{},n}`$ with respect to the pairing

$$\varphi ,\psi =_X\varphi (x)\psi (x)\mu (dx).$$

In other words, suppose that we can find functions $`\xi _k(x),\eta _k(x),k=1,\mathrm{},n`$ such that

$$\xi _k\text{Span}(\varphi _j,j=1,\mathrm{},n),\eta _k\text{Span}(\psi _j,j=1,\mathrm{},n),\xi _k,\eta _m=\delta _{km}.$$

The families $`\{\xi _j\}`$ and $`\{\eta _j\}`$ are called biorthogonal bases in $`\text{Span}(\varphi _j,j=1,\mathrm{},n)`$ and $`\text{Span}(\psi _j,j=1,\mathrm{},n)`$ considered as subspaces in $`L^2(X,\mu ).`$ Then the distribution (1) can be rewritten as (, )

$$p_n(x_1,\mathrm{},x_n)=\frac{1}{n!}det(K(x_i,x_j))_{i,j=1,\mathrm{},n},$$
(3)
with

$$K(x,y)=\underset{j=1}{\overset{n}{}}\xi _j(x)\eta _j(y).$$
(4)
One of the particularly nice properties of the ensemble (1), (3) is that one can explicitly calculate the correlation functions

$$\rho _k(x_1,\mathrm{},x_k):=\frac{n!}{(nk)!}_{X^{nk}}p(x_1,\mathrm{},x_k,x_{k+1},\mathrm{},x_n)\mu (dx_{k+1})\mathrm{}\mu (dx_n)$$

which still have a determinantal form with the same kernel $`K(x,y)`$:

$$\rho _k(x_1,\mathrm{},x_k)=det(K(x_i,x_j))_{i,j=1,\mathrm{}k}.$$
(5)
If $`\mu `$ is supported by a discrete set of points then the probabilistic meaning of the $`k`$-point correlation function is that of the probability to find a particle at each of $`k`$ sites $`x_1,x_2,\mathrm{}x_k`$. In other words,

$$\rho _k(x_1,x_2,\mathrm{}x_k)\mu (x_1)\mathrm{}\mu (x_k)=\mathrm{Pr}\{\text{there is a particle at each of the points}x_i,i=1,\mathrm{},k\}.$$

Analogously, if $`X`$ and $`\mu `$ is absolutely continuous with respect to the Lebesgue measure then

$$\rho _k(x_1,x_2,\mathrm{}x_k)\mu (dx_1)\mathrm{}\mu (dx_k)=\mathrm{Pr}\{\text{there is a particle in each interval}(x_i,x_i+dx_i)\}.$$

In general, random point processes with the $`k`$-point correlation functions of the determinantal form (5) are called determinantal or fermion (see e.g. ).

So-called Janossy densities $`𝒥_{k,I}(x_1,\mathrm{},x_k)`$, $`k=0,1,2,\mathrm{}`$, describe the distribution of the particles in a subset $`I`$ of $`X`$. If $`X`$ and $`\mu `$ is absolutely continuous with respect to the Lebesgue measure then

$`𝒥_{k,I}(x_1,\mathrm{}x_k)\mu (dx_1)\mathrm{}\mu (dx_k)=\mathrm{Pr}`$ $`\{\text{there are exactly }k\text{ particles in }I\text{,}`$
$`\text{one in each of the }k\text{ infinitesimal intervals }(x_i,x_i+dx_i)\}.`$
If $`\mu `$ is discrete then

$`𝒥_{k,I}(x_1,\mathrm{}x_k)\mu (x_1)\mathrm{}\mu (x_k)=\mathrm{Pr}`$ $`\{\text{there are exactly }k\text{ particles in }I\text{,}`$
$`\text{one at each of the }k\text{ points }x_i\}.`$
See  for details.

For determinantal point processes Janossy densities also have a determinantal form (see , p. 140, or , section 2):

$$𝒥_{k,I}(x_1,\mathrm{},x_k)=const(I)det(L_I(x_i,x_j))_{i,j=1,\mathrm{}k},$$
(6)
where

$$L_I=K_I(IdK_I)^1.$$
(7)
Here the kernel of $`K_I`$ is the restriction of the kernel $`K(x,y)`$ to $`I`$: $`K_I(x,y)=\chi _I(x)K(x,y)\chi _I(y)`$, where $`\chi _I()`$ is the characteristic function of $`I`$, and $`const(I)`$ is the Fredholm determinant

$$const(I)=det(IdK_I)=det(Id+L_I)^1.$$

The main result of this paper is

###### Theorem 1.

Let $`\stackrel{~}{\xi }_j,j=1,\mathrm{},n`$ and $`\stackrel{~}{\eta }_j,j=1,\mathrm{},n`$ be biorthonormal bases in $`\text{Span}\{\varphi _j,j=1,\mathrm{},n\}`$ and $`\text{Span}\{\psi _j,j=1,\mathrm{},n\}`$ considered as subspaces of $`L^2(XI,\mu )`$:

$$\stackrel{~}{\xi }_k\text{Span}(\varphi _j,j=1,\mathrm{},n),\stackrel{~}{\eta }_k\text{Span}(\psi _j,j=1,\mathrm{},n),_{XI}\stackrel{~}{\xi }_k(x)\stackrel{~}{\eta }_m(x)\mu (dx)=\delta _{km}.$$

Then the kernel of $`L_I=K_I(IdK_I)^1`$ is equal to

$$L_I(x,y)=\underset{j=1}{\overset{n}{}}\stackrel{~}{\xi }_j(x)\stackrel{~}{\eta }_j(y).$$
(8)
The above result readily applies to the so-called $`\beta =2`$ polynomial ensembles. Such ensembles arise, in particular, in random matrix theory , , directed percolation and tiling models , and representation theory , , . The definition is as follows.

Assume that $`X`$ is a subset of $``$ and take

$$p(x_1,\mathrm{},x_n)=const_n\underset{1j<kn}{}(x_jx_k)^2.$$
(9)
(Recall that this formula gives the joint distribution density with respect to $`\mu (dx_1)\mathrm{}\mu (dx_n)`$.)

This is a special case of (1) with $`\varphi _j(x)=\psi _j(x)=x^{j1},j=1,\mathrm{},n`$. Then we have $`\xi _j=\eta _j=p_{j1}`$, where $`\{p_j(x)\}`$ are normalized orthogonal polynomials on $`(X,\mu (dx))`$, and $`\mathrm{deg}(p_j)=j`$. The kernel $`K(x,y)`$ is the $`n`$th Christoffel-Darboux kernel

$$K_n(x,y)=\underset{j=0}{\overset{n1}{}}p_j(x)p_j(y)=\frac{k_{n1}}{k_n}\frac{p_n(x)p_{n1}(y)p_n(y)p_{n1}(x)}{xy},$$

where $`k_j`$ is the coefficient of $`x^j`$ in $`p_j(x)`$. It should be noted that the kernel $`K`$ depends on $`n`$, but in what follows we will usually omit the subscript $`n`$ unless this may lead to a confusion.

Clearly, in the case of the polynomial ensemble (9), Theorem states that the kernel of $`L_I=K_I(1K_I)^1`$ is the $`n`$th Christoffel-Darboux kernel computed for the measure $`\mu `$ restricted to $`XI`$. That is,

$$L_I(x,y)=\underset{j=0}{\overset{n1}{}}\stackrel{~}{p}_j(x)\stackrel{~}{p}_j(y)=\frac{\stackrel{~}{k}_{n1}}{\stackrel{~}{k}_n}\frac{\stackrel{~}{p}_n(x)\stackrel{~}{p}_{n1}(y)\stackrel{~}{p}_n(y)\stackrel{~}{p}_{n1}(x)}{xy},$$

where

$$\stackrel{~}{p}_j(x)=\stackrel{~}{k}_jx^j+\{\text{lower order terms}\},_{XI}\stackrel{~}{p}_k(x)\stackrel{~}{p}_m(x)\mu (dx)=\delta _{km}.$$

One of the particulary nice properties of the Janossy densities is that for any interval (or, more genertally, a measurable set) $`I`$ and non-negative integer $`k`$ one has

$$\mathrm{Pr}(\text{ there are exactly }k\text{ particles in }I)=\frac{1}{k!}_{I^k}𝒥_{k,I}(x_1,\mathrm{},x_k)\mu (dx_1)\mathrm{}\mu (dx_k)$$

The Janossy densities can be particularly useful in calculating the distribution of the left-most (right-most) particles when the particle space $`X`$ is a subset of the real line. Indeed, let us denote by $`\lambda _1\lambda _2\mathrm{}\lambda _n`$ the locations of the particles in the increasing order. Then it is easy to see that

$`\mathrm{Pr}(\lambda _k(s,s+ds))=\left({\displaystyle \frac{1}{(k1)!}}{\displaystyle _{(\mathrm{},s)^{k1}}}𝒥_{k,(\mathrm{},s)}(x_1,\mathrm{},x_{k1},s)\mu (dx_1)\mathrm{}\mu (dx_{k1})\right)\mu (ds)`$
$`=`$ $`\mathrm{Pr}(\lambda _1s){\displaystyle \frac{1}{(k1)!}}{\displaystyle _{(\mathrm{},s)^{k1}}}det(L_{(\mathrm{},s)}(x_i,x_j))_{i,j=1,\mathrm{}k}\mu (dx_1)\mathrm{}\mu (dx_{k1})\mu (ds)`$
$`\mathrm{Pr}(\lambda _1s)=\left(det(Id+L_{(\mathrm{},s)})\right)^1=det(IdK_{(\mathrm{},s)}),`$
(where in (1) we put $`x_k=s`$).

This observation and the Theorem above allow us to compute explicitly the distribution functions of the left-most particles in the hard-edge scaling limit of random matrix models when the parameter (charge at the edge) is equal to zero. We refer to Section 4 below for the details.

The result of Theorem 1 was initially discovered in the case of polynomial ensembles using the techniques of Riemann-Hilbert problems. Later on, it was realized that Theorem 1 has a simpler linear algebraic proof. However, since the main idea of the “Riemann-Hilbert” computation is very useful in deriving Painlevé equations for the distribution of the left- or right-most particles in determinantal point processes, see , , , we decided to include the argument into this paper; it can be found in Section 2. Section 3 contains the simpler proof. Concluding remarks are given in Section 5.

To conclude the Introduction, let us note that Theorem 1 has a counterpart for the so-called pfaffian ensembles. (The $`\beta =1`$ and 4 (or “orthogonal” and “symplectic”) random matrix ensembles are the most known examples of the pfaffian ensembles.) See the companion paper  for details.

## 2 Riemann-Hilbert Problem

In this section we will briefly describe two applications of the Riemann-Hilbert problem (to computing orthogonal polynomials and to inverting integrable integral operators) and use them to derive Theorem 1 in the case of polynomials ensembles. Since we use the Riemann-Hilbert problem (RHP, for short) mainly for instructional purposes, we avoid the discussion of any technical issues involved.

Let $`\mathrm{\Sigma }`$ be an oriented contour in $``$. We agree that when we go along the contour in the direction of orientation, the positive side lies to the left and the negative side lies to the right. Let $`v`$ be a map from $`\mathrm{\Sigma }`$ to $`\text{GL}(l,)`$, where $`l=1,2,\mathrm{}`$ . We say that an $`l\times l`$ matrix function $`m=m(z)`$ is a solution of the RHP $`(\mathrm{\Sigma },v)`$ if (, )

$`(i)`$ $`m(z)\mathrm{is}\mathrm{analytic}\mathrm{in}\mathrm{\Sigma },`$ (11)
$`(ii)`$ $`m_+(z)=m_{}(z)v(z),z\mathrm{\Sigma }.`$ (12)
Here $`m_+(z),m_{}(z)`$ stand for the limiting values of $`m(z)`$ as $`z`$ approaches $`\mathrm{\Sigma }`$ from the positive (negative) side. If, in addition, $`m(z)Id`$ as $`z\mathrm{}`$ then we say that $`m(z)`$ solves the normalized RHP $`(\mathrm{\Sigma },v)`$. The matrix $`v`$ is usually called the jump matrix for the RHP.

First we describe the connection of RHP to orthogonal polynomials, see , . Let $`d\mu (x)=\omega (x)dx`$ be an absolutely continuous measure on the real line, such that the non-negative density $`\omega `$ decays at infinity sufficiently fast (in particular, all moments exist). Consider an RHP on $``$ oriented from left to right with the jump matrix

$$v(z)=\left(\begin{array}{cc}1& \omega (z)\\ 0& 1\end{array}\right).$$
(13)
Fix a non-negative integer $`n`$. We are looking for a solution of the RHP $`(,v)`$ satisfying

$$m(z)=\left(Id+O(z^1)\right)\left(\begin{array}{cc}z^n& 0\\ 0& z^n\end{array}\right),z\mathrm{}.$$

It appears that this RHP has a unique solution given by

$$m(z)=\left(\begin{array}{cc}\pi _n(z)& (C(\omega \pi _n))(z)\\ \gamma _{n1}\pi _{n1}(z)& \gamma _{n1}(C(\omega \pi _{n1}))(z)\end{array}\right),z,$$
(14)
where $`\pi _n(z)=z^n+\mathrm{}`$ is the $`n`$th monic orthogonal polynomial corresponding to the weight function $`\omega (x)`$,

$$(Ch)(z)=\frac{1}{2\pi i}_{}\frac{h(\xi )}{\xi z}𝑑\xi $$

is the Cauchy transform, $`\gamma _n=2\pi ik_n^2,`$ and $`k_n`$ is the leading coefficient of the $`n`$th orthonormal polynomial $`p_n`$, i.e. $`p_n(z)=k_n\pi _n(z)`$. Thus, computing the orthogonal polynomials with the weight $`\omega (z)`$ is equivalent to solving RHP of the form above.

Now let us explain the relation of RHP to integrable operators. Let $`I`$ be a subset of $``$ (typically, a disjoint union of finitely many intervals). We recall that an integral operator $`M`$ in $`L^2(I,dx)`$ with the kernel $`M(x,y)`$ is called integrable (, , ) if

$$M(x,y)=\frac{_{i=1}^lf_i(x)g_i(y)}{xy}$$
(15)
for some $`l=2,3,\mathrm{}`$ and some functions $`f_i`$, $`g_i`$ on $`I`$. We assume that $`_{i=1}^lf_i(x)g_i(x)=0`$ so that the kernel has no singularity on the diagonal.

In particular, the formula for the Christoffel-Darboux kernel in the case of polynomial ensembles discussed in §1 means that the operators $`K`$ and $`K_I`$ can be viewed as integrable operators in $`L^2(,dx)`$ and $`L^2(I,dx)`$ with $`l=2`$, and we may take

$$f_1(x)=\frac{1}{2\pi ik_n}p_n(x)\sqrt{\omega (x)},f_2(x)=k_{n1}p_{n1}(x)\sqrt{\omega (x)},$$
(16)
$$g_1(x)=2\pi ik_{n1}p_{n1}(x)\sqrt{\omega (x)},g_2(x)=\frac{1}{k_n}p_n(x)\sqrt{\omega (x)}.$$
(17)
The appearance of $`\sqrt{\omega (x)}`$ has to do with the fact that we consider the Lebesgue measure rather than $`\mu (dx)`$ as our reference measure for the $`L^2`$-space.

It turns out that if the operator $`IdM`$ is invertible then the resolvent $`R=M(IdM)^1`$ is also an integrable operator and

$$R(x,y)=\frac{_{i=1}^lF_i(x)G_i(y)}{xy},$$
(18)
$$F_i=(IdM)^1f_i,G_i=(IdM^t)^1g_i,i=1,2,\mathrm{},l,$$
(19)
see ,  and also  for a very nice exposition. Furthermore, the functions $`F_i`$ and $`G_i`$ can be obtained through solving a RHP as follows. Let $`v^{}`$ be an $`l\times l`$ matrix valued function on $`I`$ given by

$$v^{}=Id2\pi ifg^t,f=(f_1,\mathrm{},f_l)^t,g=(g_1,\mathrm{},g_l)^t.$$
(20)
One can prove (, , ) that the normalized RHP $`(I,v^{})`$ has a unique solution $`m^{}(z)`$, and

$`F`$ $`=`$ $`(F_1,\mathrm{},F_l)^t=(m^{})_\pm f,`$ (21)
$`G`$ $`=`$ $`(G_1,\mathrm{},G_l)^t=(m^{})_\pm ^tg.`$ (22)
The following observation is crucial.

###### Lemma 1.

Let $`m`$ be the solution (14) of the RHP $`(,v)`$ and let $`m^{}`$ be the solution of the normalized RHP $`(I,v^{})`$ with $`v^{}`$ given by (20) and (16), (17). Then $`M=m^{}m`$ solves the RHP $`(I,v)`$ with the asymptotics $`\text{diag}(z^n,z^n)`$ as $`z\mathrm{}`$, and hence

$$M(z)=\left(\begin{array}{cc}\stackrel{~}{\pi }_n(z)& (C(\omega \chi _I\stackrel{~}{\pi }_n))(z)\\ \stackrel{~}{\gamma }_{n1}\stackrel{~}{\pi }_{n1}(z)& \stackrel{~}{\gamma }_{n1}(C(\omega \chi _I\stackrel{~}{\pi }_{n1}))(z)\end{array}\right),zI,$$
(23)
where $`\stackrel{~}{}`$ signifies that the corresponding polynomials are orthogonal on $`I`$ with respect to the same weight function $`\omega `$.

The proof of this lemma is based on Lemma 4.3 of  (see also Lemma 2.4 of  for a discrete analog). The analog of Lemma 1 for weight functions with discrete support was one of the basic tools used in .

###### Proof.

A straightforward calculation shows that on $`I`$ we have $`v^{}=m_+v^1m_+^1=m_{}v^1m_{}^1`$. Thus on $`I`$

$$M_{}^1M_+=m_{}^1(m_{}^{})^1m_+^{}m_+=m_{}^1v^{}m_+=v^1m_{}^1m_+=Id.$$

On the other hand, since $`m^{}(z)`$ is holomorphic away from $`I`$ and tends to $`Id`$ as $`z\mathrm{}`$, it is clear that on $`I`$, $`M(z)`$ satisfies the same jump condition as $`m(z)`$, and that it also has the same asymptotics as $`m(z)`$ when $`z\mathrm{}`$. ∎

###### Proof of Theorem 1.

We apply the formalism described above to the Christoffel-Darboux kernel with $`f_i`$, $`g_i`$ specialized by (16), (17) above. As before, we will use the notation $`\stackrel{~}{p}_k`$, $`\stackrel{~}{\pi }_k`$ for the $`k`$th orthonormal and monic orthogonal polynomials corresponding to the weight $`\omega `$ on $`I`$, and we also denote

$$q_k=C(\omega \pi _k),\stackrel{~}{q}_k=C(\omega \chi _I\stackrel{~}{\pi }_k).$$

In the calculations below we use the identity $`detm(z)detM(z)1`$. Indeed, Liouville’s theorem readily implies that if the jump matrix of a RHP has determinant 1 and the determinant of the asymptotics of a solution at infinity is also equal to 1, then the determinant of any solution of this RHP (having the corresponding asymptotics at infinity) must equal 1 identically.

We have

$`F_1`$ $`=`$ $`(m^{}f)_1=(Mm^1)_{11}f_1+(Mm^1)_{12}f_2=M_{11}m_{22}f_1M_{12}m_{21}f_1M_{11}m_{12}f_2+M_{12}m_{11}f_2`$
$`=`$ $`{\displaystyle \frac{k_{n1}}{k_n\stackrel{~}{k}_n}}\left(k_{n1}\stackrel{~}{p}_nq_{n1}p_n+\stackrel{~}{k}_n\stackrel{~}{q}_np_{n1}p_n+k_{n1}\stackrel{~}{p}_nq_np_{n1}\stackrel{~}{k}_n\stackrel{~}{q}_np_np_{n1}\right)\sqrt{\omega }`$
$`=`$ $`{\displaystyle \frac{k_{n1}^2}{k_n\stackrel{~}{k}_n}}(q_np_{n1}q_{n1}p_n)\stackrel{~}{p}_n\sqrt{\omega }={\displaystyle \frac{1}{2\pi i\stackrel{~}{k}_n}}det(m)\stackrel{~}{p}_n\sqrt{\omega }={\displaystyle \frac{1}{2\pi i\stackrel{~}{k}_n}}\stackrel{~}{p}_n\sqrt{\omega }.`$
Similar calculations yield

$$G_1=2\pi i\stackrel{~}{k}_{n1}\stackrel{~}{p}_{n1}\sqrt{\omega },F_2=\stackrel{~}{k}_{n1}\stackrel{~}{p}_{n1}\sqrt{\omega },G_2=\frac{1}{\stackrel{~}{k}_n}\stackrel{~}{p}_n\sqrt{\omega }.$$

Hence, the kernel of $`L_I=K_I(IdK_I)^1`$ equals

$$L_I(x,y)=\frac{F_1(x)G_1(y)+F_2(x)G_2(y)}{xy}=\frac{\stackrel{~}{k}_{n1}}{\stackrel{~}{k}_n}\frac{\stackrel{~}{p}_n(x)\stackrel{~}{p}_{n1}(y)\stackrel{~}{p}_n(y)\stackrel{~}{p}_{n1}(x)}{xy}\sqrt{\omega (x)\omega (y)}.$$

Recall that the factor $`\sqrt{\omega (x)\omega (y)}`$ is due to the fact that we are working in $`L^2(I,dx)`$ rather than $`L^2(I,\omega (x)dx)`$. The proof of Theorem 1 for polynomial ensembles is complete.

## 3 Linear algebraic proof

We use the notation of §1. Consider the integral operators $`K_I`$ and $`L_I`$ in $`L^2(I,\mu )`$ with the kernels

$$K_n(x,y)=\underset{j=1}{\overset{n}{}}\xi _j(x)\eta _j(y)$$
(24)
and

$$L_I(x,y)=\underset{j=1}{\overset{n}{}}\stackrel{~}{\xi }_j(x)\stackrel{~}{\eta }_j(y).$$
(25)
Both operators are finite-dimensional:

$$\begin{array}{c}\text{Ran}(K_I)=\text{Ran}(L_I)=_1:=\text{Span}(\varphi _j)=\text{Span}(\xi _j)=\text{Span}(\stackrel{~}{\xi }_j),\\ \text{Ker}(K_I)=\text{Ker}(L_I)=_2^{},_2:=\text{Span}(\psi _j)=\text{Span}(\eta _j)=\text{Span}(\stackrel{~}{\eta }_j),\end{array}$$

where the index $`j`$ ranges over $`\{1,\mathrm{},n\}`$, and the subspaces are taken inside $`L^2(I,\mu )`$. Therefore, in order to prove that $`L_I=K_I(IdK_I)^1,`$ it is enough to prove this relation for the restrictions of $`K_I`$ and $`L_I`$ to the $`n`$-dimensional space $`_1.`$ To this end we compute the matrices of the restrictions of the operators $`K_I,L_I`$ to $`_1`$ in the basis $`\{\xi _j\}_{j=1,\mathrm{}n}`$.

Let us denote by $`G_I`$ and $`G_{XI}`$ the $`n\times n`$ matrices with entries

$$(G_I)_{jk}=_I\xi _j(x)\eta _k(x)\mu (dx),(G_{XI})_{jk}=_{XI}\xi _j(x)\eta _k(x)\mu (dx)j,k=1,\mathrm{},n.$$

Since $`\{\xi _j\},\{\eta _k\}`$ are biorthonormal on $`X`$, we have $`G_I+G_{XI}=Id.`$ The matrix of $`K_I`$ on $`_1`$ in the basis $`\{\xi _j\}`$ is given by $`G_I.`$ To calculate the matrix of the restriction of $`L_I`$ on $`_1`$ we biorthonormalize the functions $`\xi _j`$, $`\eta _j`$, $`j=1,\mathrm{}n`$, in $`L^2(XI,\mu )`$. This gives (cf. Proposition 2.2 in )

$$\underset{j=1}{\overset{n}{}}\stackrel{~}{\xi }_j(x)\stackrel{~}{\eta }_j(y)=\underset{j,k=1}{\overset{n}{}}\xi _j(x)\eta _k(y)(G_{XI})_{kj}^1=\underset{j,k=1}{\overset{n}{}}\xi _j(x)\eta _k(y)(IdG_I)_{kj}^1.$$
(26)
It immediately follows from (26) that the matrix of the restriction of $`L_I`$ to $`_1`$ in the basis $`\{\xi _j\}`$ is equal to $`G_I(IdG_I)^1`$. The proof is complete.

## 4 Hard edge with zero charge

In a few special cases, the polynomials orthogonal with respect to $`\omega `$ on $`XI`$ can be easily expressed in terms of the orthogonal polynomials on $`X`$. We consider the Laguerre ensemble of positive definite matrices as an example.

Every positive definite $`n\times n`$ matrix $`M`$ can be written (in a non-unique way) as $`M=AA^{}`$, where $`A`$ is an $`n\times n`$ matrix with complex entries and $`A^{}`$ is the adjoint matrix. The probability measure in the Laguerre ensemble (also called Wishart ensemble in statistics) is defined as ():

$$P(dM)=Z_n^1\mathrm{exp}(\text{Tr}(AA^{}))det(AA^{})^\alpha dA,$$
(27)
where $`dA`$ is the Lebesgue measure on the $`2n^2`$-dimensional space of $`n\times n`$ complex matrices, $`Z_n^1`$ is a normalization constant, and $`\alpha >1.`$ The joint probability density of the distribution of the eigenvalues of $`M=AA^{}`$ is equal to

$$p(x_1,\mathrm{},x_n)=const_n\underset{1i<jn}{}(x_ix_j)^2\underset{j=1}{\overset{n}{}}x_j^\alpha e^{x_j},x_j(0,+\mathrm{}),j=1,\mathrm{},n.$$
(28)
The polynomials orthogonal with the weight $`\omega (x)=x^\alpha e^x`$ on $`_+=(0,+\mathrm{})`$ are the classical Laguerre polynomials (see e.g. ). In the special case of $`\alpha =0`$ and $`I=(0,t)`$, the orthogonal polynomials on $`XI=_+(0,t)=[t,+\mathrm{})`$ are obtained from the Laguerre polynomials by the simple shift of variable $`xxt,`$ i.e. $`\stackrel{~}{p}_j(x)=p_j(xt)`$, $`j=0,1,\mathrm{}`$ . Thus, by Theorem 1, the kernel of $`L_I^{Lag(n)}=K_I^{Lag(n)}\left(IdK_I^{Lag(n)}\right)^1`$ is equal to

$$L_I^{Lag(n)}(x,y)=K^{Lag(n)}(xt,yt),$$
(29)
where $`K^{Lag(n)}(x,y)`$ is the order $`n`$ Christoffel-Darboux kernel for Laguerre polynomials with $`\alpha =0`$.

It is well known, see , , , that when $`n`$ becomes large, the smallest eigenvalues in the Laguerre ensemble are of order $`n^1`$. Moreover, if we rescale all the eigenvalues of the $`n`$th Laguerre ensemble by $`n^1`$ then there exists a scaling limit as $`n\mathrm{}`$ of all the correlation functions. (In the random matrix theory this procedure is usually referred to as “hard edge scaling limit”.) The limit correlation functions also have the determinantal form (5) with the so-called Bessel kernel:

$$\begin{array}{c}\underset{n\mathrm{}}{lim}n^k\rho _k^{Lag(\alpha ,n)}(\frac{x_1}{n},\mathrm{},\frac{x_k}{n})=det(K^{(\alpha )}(x_i,x_j))_{i,j=1,\mathrm{},k},\\ K^{(\alpha )}(x,y)=\frac{J_\alpha (2\sqrt{x})\sqrt{y}J_\alpha ^{}(2\sqrt{y})J_\alpha (2\sqrt{y})\sqrt{x}J_\alpha ^{}(2\sqrt{x})}{xy}\\ =_0^1J_{\alpha +1}(2\sqrt{\tau x})J_{\alpha +1}(2\sqrt{\tau y})𝑑\tau =\underset{k,l=0}{\overset{\mathrm{}}{}}\frac{(1)^kx^{k+\alpha /2}}{k!\mathrm{\Gamma }(\alpha +k+1)}\frac{(1)^ly^{l+\alpha /2}}{l!\mathrm{\Gamma }(\alpha +l+1)}\frac{1}{\alpha +k+l+1}.\end{array}$$

Here $`J_\nu ()`$ is the J-Bessel function, see e.g. . Note that if $`\alpha =0`$ then the above formula makes sense for any $`x,y`$.

###### Proposition 1.

For any $`s>0`$, let $`K_s^{(0)}`$ be the (bounded) integral operator in $`L^2((0,s),dx)`$ defined by the restriction of the Bessel kernel $`K^{(\alpha )}(x,y)`$ with $`\alpha =0`$ to $`(0,s)\times (0,s)`$. Then the operator $`K_s^{(0)}\left(1K_s^{(0)}\right)^1`$ is bounded and has a kernel which is equal to $`K^{(0)}(xs,ys)`$.

###### Proof.

The relation (29) implies

$$K^{Lag(n)}(x,y)=K^{Lag(n)}(xt,yt)_0^tK^{Lag(n)}(xt,ut)K^{Lag(n)}(u,y)𝑑u,x,y(0,t).$$
(30)
Since $`n^1K^{Lag(n)}(xn^1,yn^1)`$ tends to $`K^{(0)}(x,y)`$ as $`n\mathrm{}`$ uniformly on compact subsets of $``$ (this follows, e.g., from the proof of Theorem 4.5 in ), taking the scaling limit in (30) yields

$$K^{(0)}(x,y)=K^{(0)}(xs,ys)_0^sK^{(0)}(xs,us)K^{(0)}(u,y)𝑑u,x,y(0,s).$$
(31)
Denote by $`L_s^{(0)}`$ the operator in $`L^2((0,s),dx)`$ with the kernel $`L_s^{(0)}(x,y)=K^{(0)}(xs,ys)`$. Since the kernel of this operator is the uniform limit of the kernels of the nonnegative operators $`L_I^{Lag(n)}=K^{Lag(n)}(xt,yt)`$, we have $`L_s^{(0)}0`$, and hence $`1`$ does not belong to the spectrum of $`L_s^{(0)}`$. Thus, (31) can be rewritten in the form $`K_s^{(0)}=L_s^{(0)}\left(Id+L_s^{(0)}\right)^1`$, and this is equivalent to the statement of the proposition. ∎

###### Corollary 1.

Let $`\lambda _1^{(n)}\lambda _2^{(n)}\mathrm{}\lambda _n^{(n)}`$ be the ordered eigenvalues of the Laguerre ensemble (28) with $`\alpha =0`$. Then

$$\mathrm{Pr}\left(\lambda _1^{(n)}\frac{s}{n}\right)=e^s,$$
(32)
$$\underset{n\mathrm{}}{lim}\mathrm{Pr}\left(\lambda _{k+1}^{(n)}\frac{s}{n}\right)=e^s_{(s,0)^k}det(K^{(0)}(x_i,x_j))_{i,j=1,\mathrm{},k}dx_1\mathrm{}dx_k,k2.$$
(33)
In particular,

$$\underset{n\mathrm{}}{lim}\mathrm{Pr}\left(\lambda _2^{(n)}\frac{s}{n}\right)=\frac{e^s}{2}_0^{2\sqrt{s}}x(I_0^2(x)I_1^2(x))𝑑x,$$
(34)
where $`I_\nu ()`$ is the I-Bessel function.

The formula (32) was first observed in . The limiting distribution (34) of the second smallest eigenvalue was computed in  and . <sup>2</sup><sup>2</sup>2Note that the integral (34) can be evaluated in terms of Bessel functions, see e.g. (2.30) in . Further results in this direction, including formulas similar to (33) can be found in .

###### Proof.

The relation (32) is easy:

$$\begin{array}{c}\mathrm{Pr}\left(\lambda _1^{(n)}\frac{s}{n}\right)=\frac{_{(s/n,+\mathrm{})^n}\underset{i<j}{}(x_ix_j)^2\underset{j}{}e^{x_j}dx_j}{_{(0,+\mathrm{})^n}\underset{i<j}{}(x_ix_j)^2\underset{j}{}e^{x_j}dx_j}=\frac{_{(0,+\mathrm{})^n}\underset{i<j}{}(x_ix_j)^2\underset{j}{}e^{x_js/n}dx_j}{_{(0,+\mathrm{})^n}\underset{i<j}{}(x_ix_j)^2\underset{j}{}e^{x_j}dx_j}=e^s.\end{array}$$

The relation (33) follows from (1) applied to the Laguerre ensemble and the uniform convergence of kernels mentioned in the proof of the proposition above. Finally, using the L’Hôpital rule we obtain

$$_0^sK^{(0)}(x,x)𝑑x=_0^s\left(\left(J_0^{}(2i\sqrt{x})\right)^2\frac{J_0(2i\sqrt{x})J_0^{}(2i\sqrt{x})}{2i\sqrt{x}}J_0(2i\sqrt{x})J_0^{\prime \prime }(2i\sqrt{x})\right)𝑑x.$$

The formulas

$$J_0^{}(z)=J_1(z),J_0^{\prime \prime }(z)=z^1J_1(z)J_0(z),I_0(z)=J_0(iz),I_1(z)=iJ_1(iz),$$

see , and the change of variable $`xx^2/2`$ bring the last integral to the form (34). ∎

Of course, if $`\lambda _1<\lambda _2<\mathrm{}`$ are the ordered particles of the determinantal point process with the correlation functions given by the Bessel kernel with $`\alpha =0`$ then the right–hand sides of (32) and (33) are equal to $`\mathrm{Pr}\left(\lambda _1s\right)`$ and $`\mathrm{Pr}\left(\lambda _{k+1}s\right)`$, respectively.

The calculations similar to those above can be done for the Jacobi ensemble corresponding to $`\omega (x)=(1x)^\alpha (1+x)^\beta ,x(1,1)`$, in the special cases $`I=(t,1)`$, $`\alpha =0`$; $`I=(1,t)`$, $`\beta =0`$. After appropriate rescaling one again obtains the limit relations of the form (32), (33), (34).

## 5 Concluding remarks

Take two $`n`$-point ensembles with joint probability densities of the form $`constp_n(x_1,\mathrm{},x_n)`$ with the same $`p_n`$, but assume that these two ensembles are supported by different sets — the first one lives on $`(X,\mu )`$ while the second one lives on $`(XI,\mu )`$, where $`I`$ is a subset of $`X`$. Of course, the normalization constants for these two ensembles will be different. The $`k`$th Janossy density $`𝒥_{k,I}(x_1,\mathrm{}x_k)`$ of the first ensemble is given by the formula

$$𝒥_{k,I}(x_1,\mathrm{}x_k)=const^{}_{(XI)^{nk}}p_n(x_1,\mathrm{},x_n)\mu (dx_{k+1})\mathrm{}\mu (dx_n),x_1,\mathrm{},x_kI,$$

while the $`k`$th correlation function $`\stackrel{~}{\rho }_k(x_1,\mathrm{},x_k)`$ of the second ensemble equals

$$\stackrel{~}{\rho }_k(x_1,\mathrm{},x_k)=const^{\prime \prime }_{(XI)^{nk}}p_n(x_1,\mathrm{},x_n)\mu (dx_{k+1})\mathrm{}\mu (dx_n),x_1,\mathrm{},x_kXI.$$

The only difference between the two formulas above is in the constant prefactor and in the domain where $`x_1,\mathrm{},x_k`$ are allowed to vary. Clearly, this suggests that there should be a direct relation between $`𝒥_{k,I}`$ and $`\rho _k`$, and in the case of determinantal ensembles such relation is provided by Theorem 1.<sup>3</sup><sup>3</sup>3Note, however, that a simple comparison of the two formulas above does not give a proof of Theorem 1 because these formulas only provide the symmetric minors of the corresponding kernels.

Since the argument above does not depend on the specific form of the density $`p_n`$, one might expect that Theorem 1 should have an analog for the pfaffian ensembles (see e.g. , ,  for definitions). This is exactly the case, and the corresponding result is presented in the companion paper .


format: text

text:
# Angular Distribution of Particle Flows and Perturbative QCD Predictions

## 1 Introduction

We study the angular distribution of particles between the hadronic jets in the final state of a high energy collision. The angular flow of such particles, most of them with low momentum $`<1`$ GeV, is expected to depend characteristically on the colour connections of the primary partons  (for recent reviews, see, for example, Ref. ). An early example of such a phenomenon concerns the final state $`e^+e^{}q\overline{q}g`$. The soft radiation is not distributed symmetrically between the three jets, but rather is depleted in the angular region between the $`q`$ and $`\overline{q}`$ directions. This effect was first predicted within the string hadronization model . Subsequently, it has been derived within perturbative QCD  where the particle flow is derived directly from the soft gluon bremsstrahlung, which is emitted coherently from all primary partons but with different strength from gluon and quark emitters according to the QCD colour factors.

The first observation of the effect was due to the JADE collaboration  with many of the subsequent details having been studied by various experimental groups, for example . This “string/drag effect” is by now well established and reproduced by the popular Monte Carlo models. The analytic calculations have been verified mainly for angles midway between the jets. Until now, the full angular pattern predicted by perturbative QCD had not been compared systematically with data. At PEP energies the TPC collaboration  found that their data – though in qualitative agreement with QCD expectations – showed some quantitative differences with the asymptotic analytical results, which could have possibly been caused by the low purity of the quark/gluon identification of the jets. At LEP there has been considerable progress in jet identification, in particular by DELPHI ; their result for the “Mercedes configuration” has been compared successfully  with the analytic calculation . Comparison with the limited results by OPAL , averaged over a large range of jet angles, has been moderately successful.

In this paper we compare the experimental results for the various symmetric and asymmetric angular configurations of $`q\overline{q}g`$ presented by DELPHI and OPAL for the full angular region, along with the results for $`q\overline{q}\gamma `$ and two jet events with theoretical predictions. Furthermore, we discuss the energy and jet resolution dependence of the effect which has been ignored in the previous analyses.

## 2 Angular pattern of particle flows

### 2.1 Perturbative results and dependence on the jet selection

Let us begin by recalling the main ideas  in the perturbatively-based calculation of particle flows in $`e^+e^{}q\overline{q}g`$ at high energies in comparison with $`e^+e^{}q\overline{q}\gamma `$. We consider the configurations where all the angles $`\mathrm{\Theta }_{ij}`$ between jets are large ($`i=\{+1\}\{q\overline{q}g\}`$ or $`i=\{+\}\{q\overline{q}\}`$). Within the perturbative picture, the angular distribution of soft inter-jet hadrons is calculated from the distribution of soft gluons radiated coherently off the colour antenna formed by the primary emitters $`(q,\overline{q}`$ and $`g`$) or $`(q`$ and $`\overline{q})`$.

The angular distribution of a secondary soft gluon, $`g_2`$, is derived in lowest order perturbation theory from the corresponding Feynman diagrams of single gluon emission off the primary partons. For the $`q\overline{q}\gamma `$ process one finds neglecting the recoil

$$\frac{8\pi dN_{q\overline{q}\gamma }}{d\mathrm{\Omega }_{\stackrel{}{n_2}}E_2dE_2}=\frac{4C_F\alpha _s(k_t)}{\pi }\frac{(p_+p_{})}{(p_+k_2)(p_{}k_2)}$$
(1)
with 4-momenta $`p_+,p_{}`$ and $`k_2`$ and soft gluon energy and transverse momentum $`E_2`$ and $`k_t`$.

A modification is necessary to take into account the fact that the “detected” gluon belongs to a parton jet. The results of the calculation for $`q\overline{q}\gamma `$ can be written as

$$\frac{8\pi dN_{q\overline{q}\gamma }}{d\mathrm{\Omega }_{\stackrel{}{n_2}}}=\frac{1}{N_C}W_+(\stackrel{}{n}_2)N_g^{}(Y)=\frac{2C_F}{N_C}(\widehat{+})N_g^{}(Y),$$
(2)
and for $`q\overline{q}g`$ one finds

$$\frac{8\pi dN_{q\overline{q}g}}{d\mathrm{\Omega }_{\stackrel{}{n}_2}}=\frac{1}{N_C}W_{\pm 1}(\stackrel{}{n}_2)N_g^{}(Y)=\left[(\widehat{1+})+(\widehat{1})\frac{1}{N_C^2}(\widehat{+})\right]N_g^{}(Y).$$
(3)
Here the angular distribution of soft gluon bremsstrahlung from the “antenna”-dipole $`(\widehat{ij})`$ is obtained from (1) and reads

$$(\widehat{ij})=\frac{a_{ij}}{a_ia_j},a_{ij}=(1\stackrel{}{n}_i\stackrel{}{n}_j),a_i=(1\stackrel{}{n}_2\stackrel{}{n}_i).$$
(4)
This angular distribution is the same as that for photon bremsstrahlung in QED – as it occurs, for example, in $`e^+e^{}`$ pair creation. Unlike the QED case, however, there are also the QCD-specific colour factors and the “cascading factor” $`N_g^{}(Y)`$, which takes into account the fact that $`g_2`$ is part of a jet. It represents the derivative of particle multiplicity in a gluon jet with respect to $`Y=\mathrm{ln}(K_T/Q_0)`$ at the appropriate maximal transverse momentum scale $`K_T`$ and cut-off $`Q_0`$. The parameter $`Q_0`$ is determined from a fit to the multiplicity data together with the QCD scale $`\mathrm{\Lambda }`$ ($`Q_0>\mathrm{\Lambda }`$, see, for example ).

The scale $`K_T`$ depends on the way the jets are selected. To see this, we write down the perturbative expansion for the emission of gluon $`g_2`$ in (2) from the primary parton $`p`$ as follows (see also )

$$\begin{array}{cc}\hfill \frac{8\pi dN_{q\overline{q}\gamma }}{d\mathrm{\Omega }_{\stackrel{}{n_2}}}& =\frac{1}{N_C}_{E_{min}}^{E_{max}}\frac{dE_2}{E_2}W_+(\stackrel{}{n_2})\gamma _0^2(k_{2p}^t)\hfill \\ & +\frac{1}{N_C}_{E_{min}}^{E_{max}}\frac{dE_2}{E_2}_{E_2}^{E_{max}}\frac{dE_a}{E_a}\frac{d\mathrm{\Omega }_{ap}}{2\pi \mathrm{\Theta }_{ap}^2}\gamma _0^2(k_{ap}^t)W_+(\stackrel{}{n_a})\gamma _0^2(k_{2a}^t)+\mathrm{}\hfill \\ & \frac{1}{N_C}W_+(\stackrel{}{n_2})_{E_{min}}^{E_{max}}\frac{dE_a}{E_a}\gamma _0^2(k_{2p}^t)[1+_{E_{min}}^{E_a}\frac{dE_2}{E_2}\frac{d\mathrm{\Omega }_{a2}}{2\pi \mathrm{\Theta }_{a2}^2}\gamma _0^2(k_{a2}^t)+\mathrm{}].\hfill \end{array}$$
(5)
The first term corresponds to the Born result for direct gluon emission in (1), the second one to the emission through the intermediate parton $`a`$. We denote by $`k_{ab}^t`$ and $`\mathrm{\Theta }_{ab}`$ the transverse momentum and angle respectively of parton $`a`$ with respect to parton $`b`$, in addition $`\mathrm{\Theta }_{ap}\mathrm{\Theta }_a`$; the leading order multiplicity anomalous dimension is given by $`\gamma _0^2(k_t)=2N_C\alpha _s(k_t)/\pi `$. Eq. (5) is written in DLA, but it can also be generalized to MLLA after inserting the appropriate splitting functions in all intermediate branchings.

In the second equation in (5) we have replaced the integral over the angle $`\mathrm{\Theta }_{ap}`$ by the integral over $`\mathrm{\Theta }_{a2}`$, as the leading contribution to the $`a`$ integral comes from the region of quasi-collinear emission $`\stackrel{}{n_a}\stackrel{}{n_2}`$; then we also replace the corresponding angles involving particle $`a`$ in $`W_+`$. Furthermore, the order of integration between $`E_a`$ and $`E_2`$ is interchanged ($`E_2<E_a`$). Now, the angular integral in the second term has to satisfy the “angular ordering” requirement, i.e. the emission angle of the intermediate parton $`a`$ should respect

$$\mathrm{\Theta }_{2a}<\mathrm{\Theta }_{ap};\text{therefore, also}\mathrm{\Theta }_{ap}>\mathrm{\Theta }_{2p}/2.$$
(6)
Then the angular integral can be performed with $`\mathrm{\Theta }_{a2}<\mathrm{\Theta }_2`$ in pole approximation.111The discussion of boundaries in these integrals is the same as in the derivation of the inclusive spectra, see for example, Ref. .

Next we discuss the limits of the energy integrations in (5). The lower limit $`E_{min}`$ is determined by the transverse momentum cut-off. For gluon emission inside a well separated jet one requires $`k_t>Q_0`$, therefore $`E_{min}Q_0/\mathrm{\Theta }_2`$ in the small angle approximation. In case of emission from a boosted $`q\overline{q}`$ “dipole” as in (1) an azimuthal dependence of the effective cut-off $`Q_0`$ around the jet axis is generated. To estimate this effect we consider a generalization of the cut-off restriction following from (1)

$$\stackrel{~}{k}_t^2=\frac{2(p_+k)(p_{}k)}{p_+p_{}}=\frac{2E_2^2(1\mathrm{cos}\mathrm{\Theta }_{+2})(1\mathrm{cos}\mathrm{\Theta }_2)}{1\mathrm{cos}\mathrm{\Theta }_+}>Q_0^2,$$
(7)
which approaches the limit $`\stackrel{~}{k}_t=k_t`$ for the $`q\overline{q}`$ back-to-back configuration. For our applications the effective $`k_t`$ cut-off $`Q_0`$ will then change by up to a factor 2. In the asymptotic expansion of multiplicity  $`\mathrm{ln}N_g(Y)=c_1\sqrt{Y}+c_2\mathrm{ln}Y+c_3/\sqrt{Y}`$ such a rescaling $`Q_0\beta Q_0`$ would modify only the $`c_3`$ term, therefore the leading (DLA) and next-to-leading (MLLA) results remain unchanged. Numerically, for the present analysis we estimate the effect on $`N^{}`$ in (2),(3) from the variation of $`Q_0`$ to be at the 10-20% level where we take into account the data in Fig. 4a below. In our application of analytic high energy approximations (large $`K_T/Q_0`$) we continue therefore with the cut-off $`k_t>Q_0`$ in the following.

The upper limit $`E_{max}`$ of the energy integral depends on the jet selection procedure and we consider two cases
a) no momentum restriction in jet selection
In the simplest case, we require 2 or 3 respectively energetic particles in the given angular directions which define the jet directions. A more precise definition is possible using the energy-energy-multiplicity or energy-energy-energy-multiplicity correlations for 2 and 3 jets . For two jets one considers all pairs of particles (ij) weighted with their energies, and for a given relative angle $`\mathrm{\Theta }_{ij}`$, one studies the angular distribution of all soft particles in the event. In the same way one proceeds for 3-jet events with triples of energy weighted particles and two relative angles.

In these cases there is no specific restriction on the energy of the triggered soft particles whose angular distribution is investigated. The maximum energy $`E_{max}`$ of these particles is then of order of the respective jet energies $`E_{jet}`$, say $`\sqrt{s}/2`$ or $`\sqrt{s}/3`$ in 2 and 3 jet events at total cms energy $`\sqrt{s}`$ of the jet system, and the same applies to the intermediate partons $`a`$.

b) selection of jets for a given resolution parameter $`y_{cut}`$
We consider here the “Durham algorithm”  which selects the jets according to a predefined minimal relative transverse momentum, approximately $`K_T\sqrt{y_{cut}s}`$. In this case, the $`y_{cut}`$ selection restricts the transverse momenta of emitted particles $`g_2`$ as well as the transverse momenta of the intermediate partons $`a`$ in (5) in the same way ($`k_{ap}^t<K_T`$), whereas the longitudinal momenta of the latter partons $`a`$ are only limited by $`𝒪(\sqrt{s})`$. However, an additional restriction comes from the angular ordering (6), which is violated for intermediate partons $`a`$ with large energies and, thus, small angles $`\mathrm{\Theta }_{ap}`$. In our logarithmic approximation we take the same upper limit for the $`E_a`$ and $`E_2`$ integrals. Similarly, the angular integral is limited by $`\mathrm{\Theta }_{a2}<\mathrm{\Theta }_2`$ within the same accuracy. Then, the second equation in (5) can be resummed, resulting in the multiplicity $`N_g(K_T/Q_0)`$

$$\frac{8\pi dN_{q\overline{q}\gamma }}{d\mathrm{\Omega }_{\stackrel{}{n_2}}}\frac{1}{N_C}W_+(\stackrel{}{n_2})_{E_{min}}^{E_{max}}\frac{dE_2}{E_2}\gamma _0^2(k_{2p}^t)N_g(k_{2p}^t/Q_0).$$
(8)
This integral can be expressed in terms of $`N^{}(Y)`$ using the evolution equation for multiplicity, which then leads to Eq. (2). We note that this replacement is also possible in MLLA.222The MLLA evolution equation  is usually written in small angle approximation with $`K_T=K\mathrm{\Theta }`$, the continuation to larger angles is not unambiguous.

The upper limit of integration $`E_{max}`$ is either given by the maximal energy $`E_{jet}`$ of the jet or depends on the jet resolution parameter $`y_{cut}`$

case a (no restriction) $`Y=\mathrm{ln}(E_{jet}\mathrm{\Theta }_m/Q_0)`$ (9)
$`\text{case b (}y_{cut}\text{ restriction})`$ $`Y=\{\begin{array}{cc}\mathrm{ln}(E_{jet}\mathrm{\Theta }_m/Q_0)\hfill & \mathrm{\Theta }_m<\sqrt{y_{cut}s}/E_{jet}\hfill \\ \mathrm{ln}(\sqrt{y_{cut}s}/Q_0)\hfill & \text{otherwise}\hfill \end{array}`$ (10)
Here $`\mathrm{\Theta }_m`$ is the angle between the soft gluon $`g_2`$ and the jet closest in angle. As the multiplicity $`N_g`$ depends on the maximum transverse momentum scale, the Durham $`K_T`$ \- algorithm naturally provides simple results. Note that in case of $`y_{cut}`$ selection (case b), the $`N^{}(Y)`$ factor becomes independent on the emission angle $`\mathrm{\Theta }_m`$ away from the jet directions, and, therefore, the soft gluon angular distribution is given simply by the Born-term factors $`W(\stackrel{}{n_2})`$ in (2,3).

### 2.2 Projection onto the event plane

We are now interested in the projection of the soft gluon $`g_2`$ angular distribution onto the event plane defined by the momentum vectors of the $`q\overline{q}g(\gamma )`$ systems. The azimuthal angle $`\varphi `$ in the plane is defined to be zero in the direction of $`q`$, becoming positive in the direction of $`\overline{q}`$ for definiteness (in our applications $`q`$ and $`\overline{q}`$ are not distinguishable), see Fig. 1. The projection of the particle density onto the event plane is obtained in the case of $`q\overline{q}\gamma `$ from the integral over $`\mathrm{cos}\mathrm{\Theta }_2`$, where $`\mathrm{\Theta }_2`$ is the polar angle with respect to the normal to the event plane,

$`W_+(\varphi )`$ $`=`$ $`2C_F{\displaystyle \frac{d\mathrm{cos}\mathrm{\Theta }_2}{2}(\widehat{+})}=\mathrm{\hspace{0.33em}2}C_Fa_+V(\alpha ,\beta ),`$ (11)
$`V(\alpha ,\beta )`$ $`=`$ $`{\displaystyle \frac{2}{\mathrm{cos}\alpha \mathrm{cos}\beta }}\left({\displaystyle \frac{\pi \alpha }{\mathrm{sin}\alpha }}{\displaystyle \frac{\pi \beta }{\mathrm{sin}\beta }}\right),`$
the angles $`\alpha ,\beta `$ are given below in terms of $`\varphi `$.

Replacing $`\gamma `$ by the hard gluon $`g_1`$ one obtains an additional particle flow in the $`g_1`$ direction, and from the integration of Eq. (3) over cos$`\mathrm{\Theta }_2`$ one derives

$$W_{\pm 1}(\varphi )=N_C\left[a_{+1}V(\alpha ,\gamma )+a_1V(\beta ,\gamma )\frac{1}{N_C^2}a_+V(\alpha ,\beta )\right].$$
(12)
For an arbitrary choice of $`\varphi `$ the angles $`\alpha ,\beta ,\gamma `$ in the projection formulae (11),(12) have to be defined as minimal angles between the directions of the soft gluon ($`g_2`$) and the directions of $`q,\overline{q}`$ and $`g`$ respectively, with all angles within the interval $`[0,\pi ]`$. We find the following definition of the angles $`\alpha ,\beta ,\gamma `$ to satisfy this requirement and to apply in all angular sectors with $`\varphi `$ running from $`0`$ to $`2\pi `$ (see also Fig. 1)

$`\alpha `$ $`=`$ $`\mathrm{min}(\varphi ,2\pi \varphi );\beta =\mathrm{min}(|\mathrm{\Theta }_+\varphi |,2\pi \varphi +\mathrm{\Theta }_+);`$ (13)
$`\gamma `$ $`=`$ $`\mathrm{min}(\varphi +\mathrm{\Theta }_{1+},|2\pi \mathrm{\Theta }_{1+}\varphi |).`$ (14)
Eq. (12) corresponds to an incoherent superposition of two Lorentz-boosted dipoles, one between the gluon and the quark and one between the gluon and the anti-quark, modified by the negative correction of $`𝒪(1/N_c^2)`$. This expression shows a depletion opposite to the gluon direction.

### 2.3 Comparison with experiment

Let us now compare these analytical formulae with experimental data. The DELPHI collaboration  has presented results on three jet events with two identified $`b`$-quark jets allowing for a high gluon jet purity ($`94\%`$). Symmetric angular configurations have been selected – either “Mercedes” or “Y-symmetric” events with angles confined to an angular region of $`\pm 15^{}`$. Furthermore, the data on the final state $`q\overline{q}\gamma `$ in the Y-configuration have been presented. These well defined event classes are compared in Fig. 2 with the above analytical formulae.

The angle $`\varphi =0`$ is in the direction of the most energetic quark jet. For Mercedes events, we chose the other jet angles at the observed peak positions – slightly different ($`<8^{}`$) from the symmetric values. The events had been selected with fixed $`y_{cut}0.01`$ according to the Durham algorithm. As a result, we can neglect the variation of the cascading factor $`N_g^{}`$ in (2) and (3) (for angles $`>20^{}`$ away from the jet directions). The shape of the angular distribution is then predicted without a free parameter; in Fig. 2 only one overall normalization factor has been fitted in the description of the three distributions.

The general characteristics of the distributions are well reproduced, in particular the different heights of the inter-jet valleys. In addition, the larger width of the gluon jet in comparison with the central quark jet width becomes visible. The quark-jets in a) and b) are a bit broader than perturbatively predicted, as can be expected from the b-quark jet selection; the quark-jets at the larger angle $`\varphi `$ in a) and b) are even broader because of the angular fluctuation in the selection of the jets. This successfull description nicely demonstrates the similarity of particle flow in “strong interaction” processes to photon flow in electromagnetic interactions – both are caused by the underlying gauge boson bremsstrahlung in QED and QCD. In c) the formulae are equivalent, in a) and b) they are modified by the different colour charges for quarks and gluons – a degree of freedom not available in QED. An interesting QCD interference effect is the negative $`1/N_C^2`$ term in (3). Unfortunately, according to our analysis, the presented data are not sufficient for a reliable observation of this term.

A comparison of $`q\overline{q}g`$ and $`q\overline{q}\gamma `$ distributions between the quark and anti-quark jets has been presented by OPAL . In the analysis of $`q\overline{q}g`$ events, the jets are ordered according to their energies. The jet with the lowest energy in the $`q\overline{q}g`$ system is taken as the gluon jet and two intervals of this energy $`E_3`$ are selected around $`E_3=10`$ and $`E_3=20`$ GeV, furthermore, $`\mathrm{\Theta }_+165^{}`$ from which we can estimate $`\mathrm{\Theta }_1`$ and $`\mathrm{\Theta }_{1+}`$. In the calculation we take into account the estimated purities of gluon jets, $`p_g`$, by superimposing the distribution (2) with the one obtained by swapping the second and third jet, i.e. by swapping $`\mathrm{\Theta }_+`$ and $`\mathrm{\Theta }_{1+}`$ and changing $`\varphi `$ to $`2\pi \varphi `$. The original distribution gets the weight $`p_g=0.92`$ for the sample with $`E_3=10`$ and $`p_g=0.74`$ for $`E_3=20`$ GeV, according to the reported purities.

In Fig. 3 we show the distribution in the rescaled azimuthal angle $`X=\varphi /\mathrm{\Theta }_+`$ for the two processes as well as for their ratios in the two intervals of $`E_3`$. The different heights of both distributions near $`X0.5`$ is rather well reproduced, as is the variation with energy $`E_3`$ and jet angle $`\mathrm{\Theta }_1`$. In particular, there is a sizeable asymmetry in the distribution which would be even more pronounced if the gluon jets were identified with higher purity. Since the data in both intervals still average a sizable range of different angular configurations, some deviations from the predictions are to be expected. We therefore allowed an increase of relative normalization by 15% of the curves in Fig 3 b) over Fig 3 a) to improve the description. The predictions are seen to fall below the data for $`dn/dX`$ near the jet directions, which may be due to the averaging over the jets with different energies and angular widths. On the other hand, the difference in the splitting of the two distributions, in their asymmetry, and in their absolute height for the two angular configurations in a) and b) are well reproduced. The effects related to configuration asymmetry have been studied here for the first time; these effects could be investigated more precisely in the future if smaller intervals in the relative jet angles were chosen.

## 3 Energy dependence of particle angular flow

The energy ($`\sqrt{s}`$) and $`y_{cut}`$ dependences are determined by the cascading factor $`N^{}(Y)`$ according to (9),(10). Angular distributions for particles produced in the hemispheres of quark and gluon jets without $`y_{cut}`$ restriction are available at different energies from $`e^+e^{}q\overline{q}`$ and at LEP-1 energy from $`e^+e^{}b\overline{b}g`$. Instead of the distribution in angle $`\mathrm{\Theta }`$ we may also study the distribution in pseudo-rapidity $`y=\mathrm{ln}\mathrm{tan}(\mathrm{\Theta }/2)`$. For quark and gluon jets ($`A=q,g)`$ we have

$$\frac{dn_A}{dy}(E_{jet},y)=N_A^{}(Y),Y=\mathrm{ln}\left(\frac{E_{jet}\mathrm{\Theta }}{Q_0}\right).$$
(15)
This relation is not limited to the leading order result $`N_q/N_g=C_F/C_A`$ as in Eq. (2) but allows for higher order logarithmic effects. Relation (15) follows from the MLLA result that the multiplicity depends only on the maximum $`K_T`$ in the jet, i.e. on the product $`E_{jet}\mathrm{\Theta }`$ in the small angle approximation.333 Note that the remarkable MLLA prediction  of such $`E_{jet}\mathrm{\Theta }`$ scaling behaviour has been recently well confirmed experimentally by the CDF collaboration  in the measurements of inclusive charged particle production in restricted cones around the jet direction. Then, the variation of multiplicity with $`y\mathrm{ln}(\mathrm{\Theta }/2)`$ is the same as the one with $`\mathrm{ln}E_{jet}y_{max}=\mathrm{ln}(2E_{jet}/m)`$ or, equivalently, the particle density depends on rapidity and energy only through the scaling variable

$$x=y\mathrm{ln}(\sqrt{s}/\mu )=\mathrm{ln}(E_{jet}\mathrm{\Theta }/\mu )$$
(16)
where we set $`\mu =1`$ GeV and $`\sqrt{s}=2E_{jet}`$. Eq. (15) predicts both the energy dependence at fixed angle (rapidity) and the angular dependence at fixed energy. Recall that within the MLLA approach  the multiplicity $`N`$ depends only on the phenomenological parameter $`Q_0\mathrm{\Lambda }`$ and a normalization factor, then the angular distribution is uniquely predicted by (15).

These predictions have been tested by comparing the QCD results for $`N_A^{}`$ with data from TASSO and OPAL (here from $`uds`$-quark jets only to avoid the heavy-quark contribution at the $`Z`$ mass) in Fig. 4. The curves in Fig. 4 represent the high energy QCD 3NLO prediction (DG )444we selected the fit with parameters $`n_f`$=3, $`\mathrm{\Lambda }`$=0.67 GeV and K=0.288 for a gluon jet and the numerical solution (LO ) of the full MLLA evolution equations for quark and gluon jets  . In both calculations the parameters have been determined from a fit to $`e^+e^{}`$ data, and the two results agree for quark jets within $`3\%`$ in the energy region considered ($`10<\sqrt{s}<200`$ GeV). The data in Fig. 4a show the predicted increase with energy. In Fig. 4b we show the comparison for the rapidity distributions at different energies using the scaling variable $`x`$ in (16). At each energy, the lower bound of the rapidity interval is $`y=0.8`$. Above this value the scaling property is approximately satisfied, and, again, a good description by the theoretical calculation of $`N^{}(y)`$ is obtained in absolute terms. However, some discrepency at low $`y`$ is seen for the OPAL data.555OPAL defines rapidity with respect to the sphericity axis, while TASSO defines rapidity with respect to the thrust axis; without the $`uds`$ selection and using the thrust axis to define rapidity ALEPH  data would lie above the curves. Approximate agreement is also obtained with the distribution in gluon jets from OPAL  in the given rapidity range. For larger rapidities, the distibution would fall below the distribution for quarks as one may expect from energy conservation and recoil effects given the increased central production of particles. The curve shown corresponds to the common solution of the evolution equations for quark and gluon jets (LO ); the solution DG  with the quark jet parameters would result in a prediction larger by 20%; with readjusted parameters a good description is obtained again.

We conclude that the energy dependence of the angular distribution, given by $`N^{}(Y)`$ in absolute terms, is reasonably well reproduced. A further improvement would require the treatment of (i) large angle and recoil corrections, (ii) dependence on jet axis definition (iii) heavy quark contribution and hadron mass effects.

The energy dependent effects discussed here are also expected in multi-jet events as in Eqs. (2, 3). Different analysis methods and configurations selected by different groups do not allow a meaningful comparison yet. Another interesting test will be the dependence of angular distributions on the jet resolution parameter $`y_{cut}`$ in (10).

## 4 Conclusion

The angular distribution of particles in two and three jet events follows closely the expectations based on the lowest order bremsstrahlung formulae of QCD. These formulae can be represented as sums over dipole radiation terms (4), the same as in QED, modified by the proper colour factors for quarks and gluons. In addition, there is an effect of parton branching which modifies these formulae by the “cascading factor” $`N^{}(Y)`$, the logarithmic derivative of the multiplicity. For fixed resolution parameter $`y_{cut}`$, this factor does not modify the basic angular distribution in a wide range. On the other hand, it implies an additional dependence on energy and/or jet resolution. The predicted dependence on $`N^{}`$ has been clearly established for two-jet events. These results demonstrate that soft multi-hadron production follows the simple expectations of perturbative QCD bremsstrahlung and confirm the close connection between the colour flows and observable flows of hadrons.

## Acknowledgements

We thank Bill Gary for useful discussions. VAK thanks the Leverhulme Trust for a Fellowship and the theory group of the Max-Planck-Institute, Munich for their hospitality.