  Neutral meson cross-section and direct photon HBT measurements in high energy
pp, pA, and AA collisions rely on either subtracting or dividing out the
background underlying the $M_\mathrm{inv}$ or $Q_\mathrm{inv}$ distributions.
  In this paper we investigate the photon pair distributions for simulated
proton-proton collisions at $\sqrt{s}=13~$TeV, and show that the correlated
part of the background in $M_\mathrm{inv}$ gets poorly described by the
conventional event-mixing method implemented in experimental analyses. We show
that it fails to describe the background as it gives a qualitatively different
shape especially at invariant masses around the $\pi^{0}$ and $\eta$ meson
mass. We show that this is caused by the correlations originating from parton
fragmentation, and discuss a new method that attempts to describe the
correlated part of the combinatorial background by only pairing photons when
the parton ancestors are close in $\eta$ and $\varphi$.
