Derivation Notes on SPheRIO Implementation of Back to Back Correlation

Derivation Notes on SPheRIO implementation of Back to Back Correlation

Two-particle correlation
The definition of two-particle correlation reads


 * $$\begin{align}

R(p_1,p_2) \equiv \frac{\langle N\rangle^2}{\langle N(N-1)\rangle}\frac{P_2(p_1,p_2)}{P_1(p_1)P_1(p_2)} =\frac{\langle N\rangle^2\sigma_T}{\langle N(N-1)\rangle}\frac{d^6\sigma/d^3p_1d^3p_2}{(d^3\sigma/d^3 p_1)(d^3\sigma/d^3 p_2)} \end{align}$$

where $$P_n$$ are the inclusive n-particle distribution function. The expression can be written down in terms of differential cross section and $$\sigma$$, the latter is normalized to the total cross section $$\sigma_T$$ by the following relations


 * $$\begin{align}

&\int d^3p_1d^3p_2 (d^6\sigma/d^3p_1d^3p_2) = \langle N(N-1)\rangle\sigma_T \\ &\int d^3p (d^3\sigma/d^3 p) = \langle N\rangle\sigma_T \end{align}$$

Note that


 * $$\begin{align}

\int d^3p (d^3\sigma/d^3 p) = \int d\Omega p^2dp \frac{d^3\sigma}{p^2dpd\Omega}= \int d\Omega \frac{d\sigma}{d\Omega} \end{align}$$

The factor $$\langle N\rangle$$ and $$\langle N(N-1)\rangle$$ to r.h.s of the above expressions are due to the fact that the collision is not elastic, therefore the out-coming count is normalized to $$\langle N\rangle$$ and $$\langle N(N-1)\rangle$$ instead of $$1$$. The difference between $$\langle N\rangle^2$$ and $$\langle N(N-1)\rangle$$ can be ignored in big multiplicity, as in the case of nuclear collisions.

For simplicity, we first consider a state $$|\phi\rangle$$ at zero temperature. The number of particle with momentum $$p$$ can be obtained by calculating the expectation value of the number operator


 * $$\begin{align}

&P_1(p) = \langle\phi|\hat{n}_p|\phi\rangle \\ &\hat{n}_p=\hat{a_p}^+\hat{a_p} \end{align}$$

in quantum mechanics. This average number is also obtained by explicitly enumerating all the possible states where at least one particle of given momentum is observed and the probability to find such a state. To calculate this we evaluate the squared modulus of the probability amplitude to transit to any state $$X$$ after annihilation of one particle (If the corresponding transition amplitude to $$X$$ is zero, it gives no contribution), and the corresponding number of particles is taken into account by the property of annihilation operator.


 * $$P_1(p) = \sum_X |\langle X|\hat{a_p}|\phi\rangle|^2 =\langle\phi |\hat{a_p}^+\hat{a_p}|\phi\rangle$$

The above expression is consistent with our previous result.

One needs to take more careful when considering two-particle distribution function. If the two momenta $$p_1, p_2$$ are different, one uses $$\hat{n}_{p_1}\hat{n}_{p_2}$$.


 * $$\begin{align}

P_2(p_1,p_2)=\langle\phi|\hat{n}_{p_1}\hat{n}_{p_2}|\phi\rangle \end{align}$$

But if we are to detect two particles of the same momentum, one should use $$\hat{n}_{p_1}(\hat{n}_{p_1}-1)$$.


 * $$\begin{align}

P_2(p_1,p_2)=\langle\phi|\hat{n}_{p_1}(\hat{n}_{p_1}-1)|\phi\rangle \end{align}$$

This is because the number of pairs changes in the second case. The issue is automatically taken care of when one uses the second approach for Boson, for Fermi, one will never detect two particles of the same state


 * $$\begin{align}

P_2(p_1,p_2) =\sum_X |\langle X|a_{p_1}a_{p_2}|\phi\rangle|^2 =\langle\phi|a_{p_2}^+a_{p_1}^+a_{p_1}a_{p_2}|\phi\rangle \end{align}$$

Let us consider the following state of a non-interacting system
 * $$\begin{align}

\phi\rangle = \mathcal{C}\prod_p (\hat{a_p}^+)^{N_p} \end{align}$$

where $$\mathcal{C}$$ is the normalization constant. Using the explicit form of $$|\phi\rangle$$, one obtains trivially


 * $$\begin{align}

R(p_1,p_2)=1 \end{align}$$

Harons are emitted in collisions through interactions. To take into account the physics, these interactions have to be treated, to say the least, as external classical sources. Analytic solutions are found in a few of these models, therefore in such cases one may move on straightforward and discuss the two-particle correlations. One typical example was studied by Gyulassy et. al., surprisingly enough, the coherent state of spinless Bosons produced by a classical souce does not exibit any non-trivial two particle correlation. And non-trivial corerlatin only appreas when one imposes explicit charge conservation, or introduces some source variations (for instance, different hadrons may be emitted from independent sources). It is therefore understood that Boson statistics alone does not describe all the physics behind HBT correlation (see PRC20, 2267, (1979) for details). In general, one usually has to adopt the perturbative approach by using the Wick theorem discussed below.

Wick theorem
As shown above, two particle correlation involves terms like


 * $$\begin{align}

P_2(p_1,p_2)=\langle | a^+_{p_1}a^+_{p_2}a_{p_2}a_{p_1} |\rangle=\text{Tr}(\hat{\rho} a^+_{p_1}a^+_{p_2}a_{p_2}a_{p_1}) \end{align}$$

In this section, we will discuss the Wick theorem and show that these term can be written in terms of two-point Green functions, then later on we will further express the resulting expression in terms of hydrodynamical variables.

For two operators $$A, B$$, one has


 * $$\begin{align}

&A=A^++A^- \\ &B=B^++B^- \\ &A^\bullet B^\bullet \equiv AB-\text{N}(AB)=[A^+,B^-]_{\pm} \end{align}$$

where for instance $$A^+$$ is proportional to annihilation operator $$a_p$$ and it contains an integral in momentum $$p$$ with time $$t$$ as a parameter.

As a result, the term $$[A^+,B^-]_{\pm}$$ on the r.h.s. is a number, so it can be replaced by its vacuum expectation value


 * $$\begin{align}

\langle 0|[A^+,B^-]_{\pm}|0\rangle=\langle 0|AB|0\rangle \end{align}$$

Therefore


 * $$\begin{align}

AB=\text{N}(AB)+\langle 0|AB|0\rangle \end{align}$$

One then considers different time orders to obtain


 * $$\begin{align}

\text{T}\{AB\}=\text{N}(AB)+\langle 0|\text{T}\{AB\}|0\rangle \end{align}$$

This is simply because


 * $$\begin{align}

&\text{T}\{AB\}=\pm \text{T}\{AB\}\\ &\text{N}(AB)=\pm \text{N}(AB) \end{align}$$

This result can be generalize to N-point case to obtain the Wick theorem at zero temperature. See F.Mandl (6.30) for more details.

For finite temperature, the situation is more complicated. Although, instead of replacing the term $$[A^+,B^-]_{\pm}$$ with its expectation value on vacuum state $$|0\rangle$$, one may consider the vacuum state in finite temperature, namely, $$|\rangle$$.


 * $$\begin{align}

A^\bullet B^\bullet = \langle |[A^+,B^-]_{\pm}|\rangle\ne\langle |AB|\rangle \end{align}$$

But this does not help further with the evaluation, since generally


 * $$\begin{align}

\langle |\text{N}(AB)|\rangle \ne 0 \end{align}$$

A careful study takes a different approach and gives the Wick theorem at finite temperature, for 4-point Green function it reads


 * $$\begin{align}

\langle |a^+_{p_1}a^+_{p_2}a_{p_2}a_{p_1}|\rangle = \langle |a^+_{p_1}a_{p_1}|\rangle\langle |a^+_{p_2}a_{p_2}|\rangle +\langle |a^+_{p_1}a_{p_2}|\rangle\langle |a^+_{p_2}a_{p_1}|\rangle +\langle |a^+_{p_1}a^+_{p_2}|\rangle\langle |a_{p_2}a_{p_1}|\rangle \end{align}$$

which does not explicitly involves normal order neither time order. In fact, the normal order (not the time order, which is defined as at zero temperature) can be employed to facilitate the proof. In particular, one can conveniently split an operator $$A$$ in such a way, that the expectation value of normal order on thermal vacuum states is always zero (therefore $$A^\bullet B^\bullet$$ changes accordingly). One then employs the Wick theorem at zero temperature and evaluate the expectation value of both side of the equality in thermal vacuum state, the l.h.s is defined as N-point Green function at finite temperature, while the r.h.s. is expressed in terms of two point Green functions. Therefore the above expression is obtained. See arXiv: hep-ph/9601268 for details.

HBT correlation
In the case of hydrodynamics, it is more natural to try to express physical quantities in terms of particle distribution functions. The one-particle and two-particle distributions then can be expressed as integrals on the freeze-out surface.

We will first discuss the simplest case. In free space, one-particle distribution function can be intuitively obtained as the following


 * $$\begin{align}

I(0,p_1)\equiv {p_1^0}\frac{d^3N}{d^3{p}_1}={p_1^0}\langle |a^+_{p_1}a_{p_1} |\rangle = {p_1^0} \text{Tr}(\hat{\rho} a^+_{p_1}a_{p_1}) \end{align}$$

Note that the appearance of the factor $$p_1^0$$ is due to dimensional arguments. Owing to the commutation relation of $${a_p}, {a_{p'}}$$, $$\langle | a^+_{p_1}a_{p_1} |\rangle$$ has the dimension of $$[\delta(\vec p)]=[p^0]=[\text{Energy}]$$, where one makes use of the fact that $$\int\delta(\vec p)d^3p=1$$ and $$d^3p/p^0$$ is invariant and dimensionless due to $$\int\delta(p^2-m^2)d^4p=1$$.

Hydrodynamics assumes local equilibrium, which implies that we know how to obtain the expectation value such as $$\langle |a^+_{p_1}a_{p_1} |\rangle$$ once the system evolution is determined. For a static source, once the particle density $$n$$ is obtained, one can calculate the particle number by a spatial integral.


 * $$\begin{align}

N=\int n d^3x \end{align}$$

The above integral can be seen as carried out at constant time, which is not the case when hadrons freeze-out in the final stage of hydrodynamic evolution. From an hydrodynamical point of view, the freeze-out surface is not necessarily time-like, it can be any surface in 4-space. The hadrons measured by the detector only consist of those whose fluxes pass the surface element. Therefore, one has to generalize the above expression to a form which is covariant, namely, a contraction. Following this line of thought, the term $$d^3x$$ should be understand as an infinitesimal dual vector, in other words, the surface element is a four dual vector which only has time component. Therefore, the particle number is given by the contraction of surface element and particle four momentum.


 * $$\begin{align}

p^0N=\int_{\sigma} d\sigma_{\mu}p^{\mu}f(x,p) \end{align}$$

Using Wick theorem, one shows that the two-particle correlation can be simplifies to the summation of three terms, all of which only involves two point Green functions.


 * $$\begin{align}

\langle |a^+_{p_1}a^+_{p_2}a_{p_2}a_{p_1}|\rangle = \langle |a^+_{p_1}a_{p_1}|\rangle\langle |a^+_{p_2}a_{p_2}|\rangle +\langle |a^+_{p_1}a_{p_2}|\rangle\langle |a^+_{p_2}a_{p_1}|\rangle +\langle |a^+_{p_1}a^+_{p_2}|\rangle\langle |a_{p_2}a_{p_1}|\rangle \end{align}$$

The first term is, in fact, the product of one-particle distribution. We have discussed above how to evaluate it in hydrodynamics, however, it is not even necessary in this context, since it is exactly the same form as the denominator, therefore gives a constant of 1. The second term corresponds to HBT effect, and the third term corresponding to back to back correlation. Consider the following pure state


 * $$|\phi\rangle \equiv \mathcal{C}|(a^+_k)^2\rangle+\mathcal{D}|(a^+_k)^3\rangle+\mathcal{E}|(a^+_{k'})^2\rangle+\mathcal{F}|(a^+_{k'})^3\rangle$$

where $$\mathcal{C} \mathcal{D}\mathcal{E} \mathcal{F}$$ are normalization constants. It is straightforward to verify that none of the following expectation value vanishes,


 * $$\begin{align}

&\langle\phi|a_k|\phi\rangle \\ &\langle\phi|a^+_{k'}a_k|\phi\rangle \\ &\langle\phi|a_{k'}a_k|\phi\rangle \end{align}$$

In particular, they are non-vanishing for coherent state $$|\alpha\rangle$$, which exibits maximal coherence


 * $$\begin{align}

&|\alpha\rangle =e^{-{|\alpha|^2\over2}}\sum_{n=0}^{\infty}{\alpha^n\over\sqrt{n!}}|n\rangle =e^{-{|\alpha|^2\over2}}e^{\alpha a^+}|0\rangle \\ &a|\alpha\rangle = \alpha |\alpha\rangle \end{align}$$

The only problem with such coherent state is that they will tend to collapse to an ensemble of energy eigen-states in a thermal environment. Due to the second law of thermodynamics, entropy obtains its maximum when the system is in equilibrium, leading to the canonical distribution, which is known as thermal decoherence (see Sakurai Modern Quantum Mechanics(3.4.48)). The effect of decoherence, in terms of density matrix, is essentially the decay or rapid vanishing of the off-diagonal elements of the partial trace of the joint system's density matrix, i.e. the trace, with respect to any environmental basis, of the density matrix of the combined system and its environment. The decoherence irreversibly converts the "environmentally traced over" density matrix from a pure state to a reduced mixture (see arXiv:quant-ph/0105127). Due to this very reason, the third therm is usually understood to be zero, but we will come back with more details in the below.


 * $$\begin{align}

\langle |a^+_{p_1}a^+_{p_2}|\rangle = 0 \end{align}$$

Now let us focus on HBT effect. Sinyukov has shown (see arXiv: nucl-th/9909018) that such terms can be written using Wigner function as following


 * $$\begin{align}

&\langle a^+_{p_1}a_{p_2}\rangle =\int d\sigma_{\mu} p^{\mu}e^{iq\cdot x} f_W(x,p) \\ &f_W(x,p) \equiv \int d^4u (2\pi)^{-3}e^{-iu\cdot x}\delta(p\cdot u)\langle a^+(p-\frac{u}{2})a(p+\frac{u}{2})\rangle \end{align}$$

where $$f_W(x,p)$$ is the single particle Wigner function. The above expression is general, except it is written in a covariant form as we did for one-particle distribution function. As a result, the l.h.s. which is equal time two-point Green function into the r.h.s. where the integral is carried out on the freeze-out surface (the two-point Green function is still at equal time).

Then it is argued that at the limit when the density matrix is diagonal enough, one can replace the Wigner function with Bose-Einstein distribution for bosons.


 * $$\begin{align}

&f_W(x,p)\to f_{B}(p\beta u) = \frac{(2\pi)^{-3}}{exp(\beta p\cdot u) - 1}\\ &\beta = \frac{1}{T} \end{align}$$

In hydrodynamics, the most simple scenario for the freeze-out surface is the Cooper-Frye sudden freeze-out at $$T=T_f=\text{const.}$$. This freeze-out temperature should be the same as that one uses to study the $$p_T$$ spectra. We neglect the resonance decays in this work. Since it was argued that, since resonance decays contribute to the correlations with very small $$q$$ values, the experimentally determined HBT radii are essentially due to the direct pions.

Putting all the pieces together, we obtains the expression of HBT correlation in hydrodynamics


 * $$\begin{align}

R(p_1,p_2) \equiv C_2^{HBT}(q,P)=1+\frac{|I(q,P)|^2}{I(0,p_1)I(0,p_2)} \end{align}$$

where


 * $$\begin{align}

&I(q,P)\equiv \int_{\sigma=T_f}d\sigma_\mu P^\mu f(x,P)e^{iqx}\\ &P=(p_1+p_2)/2 \\ &q=(p_1-p_2) \end{align}$$

In SPH representation, we write $$I(q,P)$$ as


 * $$\begin{align}

I(q,P)=\sum_j \frac{\nu_j n_{j\mu} P^\mu}{s_j |n_{j\mu} u_j^\mu|} {\mathrm{e}}^{iq_\mu x_j^\mu}f(u_{j\mu} P^\mu) \end{align}$$

where the summation is over all the SPH particles. In the Cooper-Frye freeze-out, these particles should be taken where they cross the hyper-surface $$T=T_f$$ and $$n_{j\mu}$$ is the normal to this hyper-surface. Notice that, if we put $$p_1=p_2=P$$, so that $$q=0$$, $$I(0,P)$$ is reduced to the inclusive one-particle distribution.

Back to back correlation
As mentioned above, the contribution from the third on the r.h.s. of the Green function gives rise to the so called "Back to back correlation".


 * $$\begin{align}

\langle |a^+_{p_1}a^+_{p_2}a_{p_2}a_{p_1}|\rangle = \langle |a^+_{p_1}a_{p_1}|\rangle\langle |a^+_{p_2}a_{p_2}|\rangle +\langle |a^+_{p_1}a_{p_2}|\rangle\langle |a^+_{p_2}a_{p_1}|\rangle +\langle |a^+_{p_1}a^+_{p_2}|\rangle\langle |a_{p_2}a_{p_1}|\rangle \end{align}$$

The basic idea was introduced by Asakawa and Csorgo (see arXiv:hep-ph/9612331) that in the ultra-relativistic collisions the hot and dense media may have the hadrons dressed up so that their effective masses are modified, on the other hand, when the hadrons finally go cross the freeze-out surface, they become free particles and are subsequently detected. If one desires to evaluate the ensemble average value of terms such as
 * $$\langle |a_{p_2}a_{p_1}|\rangle = \text{Tr}(\hat{\rho} a_{p_2}a_{p_1})$$

one should use the density matrix $$\hat{\rho}$$ of the dressed up hadrons, while $$a_{p_2}, a_{p_1}$$ are particle creation and annihilation operators in the vacuum. It turns out naturally these terms are non-vanishing in such case. In fact, one can establish two sets of particles creation and annihilation operators, those for the vacuum ($$a^+_p,a_{p'}$$) and those for the hot and dense media ($$b^+_p,b_{p'}$$), they are connected by the Bogoliubov transformation.

However, their calculation only treats homogeneous system, the model Hamiltonian does not contain any coordinate dependent term souce term which describes the spatial distribution of the system in question. As a results, the correlation is non-trivial only at parallel and anti-parallel directions, namely


 * $$\begin{align}

&C_2^{HBT}(p,p')=2\delta(p-p')\\ &C_2^{B2B}(p,p')=\delta(p+p')\left(1+\frac{|s_pc_p^*n_p+s_{-p}c_{-p}^*(n_{-p}+1)|^2}{n_1(p)n_1(-p)}\right) \end{align}$$

where


 * $$\begin{align}

&c_p = \cosh(r_p)\\ &s_p = \sinh(r_p)\\ &r_p =\frac{1}{2}\log\left[\frac{\omega_p}{\Omega_p}\right]\\ &\omega_p = \sqrt{m^2+{\vec p}^2}\\ &\Omega_p= \sqrt{m_*^2+{\vec p}^2}\\ &n_p=[\exp(\Omega_p/T -1)]^{-1}\\ &n_1(p)=|c_p|^2n_p+|s_k|^2(n_{-p}+1) \end{align}$$

and $$m_*$$ is the modified mass.

However, as it is well known, in real collisions, hadrons are emitted from freeze-out surface due to the hydrodynamical nature of the system, and the size of the fireball is finite and the system is far from homogeneous. As pointed out by Gastão and Sandra, such finite volume effect can be effectively taken care of by the Wigner function procedure introduced by Sinyukov. The trick to do this is to try to express the correlation in term of creation and annihilation operator pairs of the hot media in a form similar to that in HBT correlation, namely, $$\langle|b^+_pb_{p'}|\rangle$$, then applies Sinyukov's procedure to rewrite it as an integral on the freeze-out surface. Meanwhile one also need to rewrite the remaining coefficients in a covariant fashion. The resulting expressions in terms of $$\langle|b^+_pb_{p'}|\rangle$$ are


 * $$\begin{align}

&R(p_1,p_2)\equiv C_2(p_1,p_2) =1+\frac{|G_c(p_1,p_2)|^2}{G_c(p_1,p_1)G_c(p_2,p_2)}+\frac{|G_s(p_1,p_2)|^2}{G_c(p_1,p_1)G_c(p_2,p_2)}\\ &G_c(p_1,p_2)=\sqrt{\omega_{p_1}\omega_{p_2}}\left[c_{p_1}^*c_{p_2}\langle |b_{p_1}^+b_{p_2}|\rangle+s_{-p_1}s_{-p_2}^*(\langle |b_{-p_2}^+b_{-p_1}|\rangle+1)\right] \\ &G_s(p_1,p_2)=\sqrt{\omega_{p_1}\omega_{p_2}}\left[s_{-p_1}^*c_{p_2}\langle |b_{-p_1}^+b_{p_2}|\rangle+c_{p_1}s_{-p_2}^*(\langle |b_{-p_2}^+b_{p_1}|\rangle+1)\right] \end{align}$$

From which, one obtains the expressions which can be used in hydrodynamics (see arXiv:nucl-th/9810034 for details.)


 * $$\begin{align}

&G_c(p_1,p_2)=\frac{1}{(2\pi)^3}\int_{\sigma}d\sigma_{\mu}P^{\mu}e^{iq\cdot x} \left[|c_{(p_1,p_2)}|^2n_{(p_1,p_2)}+|s_{(-p_1,-p_2)}|^2(n_{(-p_1,-p_2)}+1)\right] \\ &G_s(p_1,p_2)=\frac{1}{(2\pi)^3}\int_{\sigma}d\sigma_{\mu}P^{\mu}e^{iP\cdot x} \left[s_{(-p_1,p_2)}^*c_{(p_2,-p_1)}n_{(-p_1,p_2)}+c_{(p_1,-p_2)}s_{(-p_2,p_1)}^*(n_{(p_1,-p_2)}+1)\right] \end{align}$$

where


 * $$\begin{align}

&P=(p_1+p_2)/2 \\ &q=(p_1-p_2)\\ &\omega_i=p_i \cdot u\\ &\tilde{p}_i= p_i - \omega_i u \\ &\Omega_i=\sqrt{m_*^2-\tilde{p}_i \cdot \tilde{p}_i}\\ &p_i^*=\Omega_i u+ \tilde{p_i}\\ &P^*=(p_1^*+p_2^*)/2\\ &c_{(\pm p_1,\pm p_2)} \equiv c_{(p_1, p_2)}= \cosh(r_{(p_1, p_2)})\\ &s_{(\pm p_1,\pm p_2)} \equiv s_{(p_1, p_2)} = \sinh(r_{( p_1, p_2)})\\ &r_{(\pm p_1,\pm p_2)} \equiv r_{(p_1, p_2)} =\frac{1}{2}\log\left[\frac{P \cdot u}{P^*\cdot u}\right]\\ &n_{(\pm p_1,\pm p_2)} \equiv n_{(p_1,p_2)}=[\exp((P^*\cdot u)/T -1)]^{-1} \end{align}$$

for further discussion, one may also define


 * $$\begin{align}

&p_{\pm i}=\omega u\pm \tilde{p_i}\\ &p_{\pm i}^*=\Omega u\pm \tilde{p_i} \end{align}$$

where $$i=1,2$$, $$u$$ is the four-velocity of the fluid on the freeze-out surface, $$\tilde{p}$$ measures the part of the four-momentum orthogonal to the four-velocity, $$p_{-i}$$ is the momentum which inverts the part of $$p_i$$ orthogonal to the four-velocity. Here I am not sure what was the purpose to introduce symbols such like $$c_{(p_1,-p_2)}, n_{(p_1,-p_2)}$$ where they only depend on the average momentum $$P$$ as a natural generalization of Sinyukov's formulae. The last step is to write everything in terms of SPheRIO variables, which reads


 * $$\begin{align}

&G_c(p_1,p_2)=\frac{1}{(2\pi)^3}\sum_j \frac{\nu_j n_{j\mu} P^\mu}{s_j |n_{j\mu} u_j^\mu|}e^{iq\cdot x} \left[|c_{(p_1,p_2)}|^2n_{(p_1,p_2)}+|s_{(-p_1,-p_2)}|^2(n_{(-p_1,-p_2)}+1)\right] \\ &G_s(p_1,p_2)=\frac{1}{(2\pi)^3}\sum_j \frac{\nu_j n_{j\mu} P^\mu}{s_j |n_{j\mu} u_j^\mu|}e^{iP\cdot x} \left[s_{(-p_1,p_2)}^*c_{(p_2,-p_1)}n_{(-p_1,p_2)}+c_{(p_1,-p_2)}s_{(-p_2,p_1)}^*(n_{(p_1,-p_2)}+1)\right] \end{align}$$

Particle and anti-particle
The above results hold only for particles that are their own antiparticles, e.g. the $$\phi$$ meson and $$\pi^0$$. In a general case, the Bogoliubov transformation will mix particles with antiparticles, but the extension to particle-antiparticle correlations is quite straightforward. For instance, the particle-antiparticle involves the evaluation of the following term


 * $$\begin{align}

&\langle |a^+_{p_1}\bar{a}_{p_2}|\rangle \\ &=\langle |(c_{p_1}^*b_{p_1}^++s_{-p_1}\bar{b}_{-p_1})(c_{p_2}\bar{b}_{p_2}+s_{-p_2}^*b_{-p_2}^+)|\rangle\\ &=\langle |c_{p_1}^*c_{p_2}b_{p_1}^+\bar{b}_{p_2}|\rangle +\langle |c_{p_1}^*s_{-p_2}^*b_{p_1}^+b_{-p_2}^+|\rangle +\langle |s_{-p_1}c_{p_2}\bar{b}_{-p_1}\bar{b}_{p_2}|\rangle +\langle |s_{-p_1}s_{-p_2}^*\bar{b}_{-p_1}b_{-p_2}^+|\rangle\\ &=\langle |c_{p_1}^*c_{p_2}b_{p_1}^+\bar{b}_{p_2}|\rangle +\langle |s_{-p_1}s_{-p_2}^*\bar{b}_{-p_1}b_{-p_2}^+|\rangle\\ &=\langle |c_{p_1}^*c_{p_2}b_{p_1}^+\bar{b}_{p_2}|\rangle +\langle |s_{-p_1}s_{-p_2}^*(b_{-p_2}^+\bar{b}_{-p_1}+1)|\rangle\\ &=c_{p_1}^*c_{p_2}\langle |b_{p_1}^+\bar{b}_{p_2}|\rangle +s_{-p_1}s_{-p_2}^*\left(\langle |b_{-p_2}^+\bar{b}_{-p_1}|\rangle+1\right)\\ &=0 \end{align}$$

which was non-vanishing in the previous case. As a result, the identical particle correlation does not have the contribution from $$G_s$$, while the particle-antiparticle correlation does not receive contribution from $$G_c$$. Besides, one need to take bit more care to identify the particle and substitute correct chemical potential when applying integral on the freeze-out surface. Other than that, things remain unchanged, since they have the same mass, and the formulae do not involve charges.

Numerical results
To be updated.