Quarkonium production and polarization: where do we stand?

We review the current status of heavy quarkonium production phenomenology based on nonrelativistic effective field theories, focusing on spin-triplet $S$-wave states such as $J/\psi$, $\psi(2S)$, and $\Upsilon$. We present some representative examples for heavy quarkonium production mechanisms proposed in the literature, which vary significantly depending on the choice of data employed in analyses. We then discuss the r\^{o}le of polarization in discriminating between the different possible scenarios for quarkonium production. Other observables that may be useful in pinpointing the production mechanism are also introduced, such as the $\eta_c$ production, associated production of $J/\psi$ plus a gauge boson, and $J/\psi$ production at the Electron-Ion Collider.


Introduction
Heavy quarkonia are useful laboratories to study perturbative and nonperturbative aspects of QCD [1][2][3][4].An important class of observables involve inclusive production of heavy quarkonia, which are considered to be promising contenders for tools to study QCD in colliders.This requires a robust understanding of the quarkonium production mechanism based on first principles, which remains a formidable challenge.
Most theoretical studies of quarkonium production phenomenology rely on nonrelativistic effective field theories, which are based on the fact that the mass m of the heavy quark Q and antiquark Q that constitute a quarkonium is much larger than Λ QCD .This allows an interpretation of heavy quarkonium states as nonrelativistic bound states; the scales that appear in a Q Q bound state are the momentum mv and the binding energy mv 2 , where m is the heavy quark mass and v is the velocity of the Q or Q inside the quarkonium.Typical values of v are v 2 ≈ 0.3 for charmonia, and v 2 ≈ 0.1 for bottomonia.The nonrelativistic QCD (NRQCD) effective field theory [5,6] provides a factorization formalism that separates the perturbative short-distance physics of scales of order m and higher from the nonperturbative long-distance physics which is encoded in the NRQCD matrix elements.This formalism has been widely adopted in phenomenological studies of heavy quarkonium production.
A difficulty in the NRQCD factorization approach is that, while the short-distance part can be computed in perturbative QCD, it is generally not known how to compute the longdistance quantities from first principles, and they are usually determined phenomenologically.This approach has not lead to a satisfactory description of the heavy quarkonium production mechanism [7].Determinations of NRQCD matrix elements from different choices of data can disagree with one another, and none of the determinations have been able to give a comprehensive description of the important observables.Hence, it is fair to say that a QCD-based understanding of the production mechanism of heavy quarkonium still remains elusive.
In this paper, we review recent efforts made towards understanding the heavy quarkonium production mechanism based on nonrelativistic effective field theories.In Sect.2, we briefly introduce the nonrelativistic effective field theory formalisms that are used in heavy quarkonium phenomenology.We review the various NRQCD matrix element determinations for spin-triplet S -wave quarkonia, which include J/ψ, ψ(2S ), and Υ in Sect.3. We discuss the rôle of quarkonium polarization in understanding the production mechanism in Sect.4, and introduce other observables related to quarkonium production in Sect. 5. We conclude in Sect.6.
2 Nonrelativistic effective theories for quarkonium production 2.1 Nonrelativistic QCD factorization NRQCD provides a factorization formalism for inclusive production cross sections of a heavy quarkonium Q in the form [6] where the sum is over color, spin, and orbital angular momentum states are perturbatively calculable short-distance cross sections for inclusive production of Q Q in the n state, and ⟨O Q (n)⟩ are NRQCD matrix elements that correspond to the nonperturbative probabilities for the Q Q(n) to evolve into a quarkonium Q plus anything.The matrix elements have known scalings in v, so that the sum over n in Eq. ( 1) is organized in powers of v, and is in practice truncated at a desired order in v.The NRQCD matrix elements are universal quantities that depend on the nonperturbative nature of the heavy quarkonium state.Hence, Eq. ( 1) provides descriptions of inclusive quarkonium production rates with a few universal, process-independent nonperturbative parameters.A novel feature of NRQCD factorization is that the quarkonium Q can be produced from Q Q in color-octet states.In the case of χ cJ and χ bJ production, contributions at leading order in v come from both color-singlet (n = 3 P [1]  J , J = 0, 1, 2) and color-octet (n = 3 S [8]  1 ) channels.For production of spin-triplet S -wave quarkonia such as J/ψ, ψ(2S ), and Υ, the color-octet channel (n = 3 S [8]  1 , 1 S [8]  0 , 3 P [8]  J , J = 0, 1, 2) contributions 1 are suppressed by several powers of v compared to the contribution from the color-singlet channel (n = 3 S [1]  1 ).However, when the transverse momentum p T of the quarkonium is much larger than the heavy quarkonium mass, the short-distance cross sections for the color-octet channels are strongly enhanced compared to the color-singlet channel; moreover, the color-singlet channel contribution severely underestimates the large-p T cross section measured at hadron colliders, so that the cross section is dominated by color-octet contributions [8][9][10].Hence, precise determination of the coloroctet matrix elements is crucial for understanding heavy quarkonium production based on first principles.
It is worth mentioning that in general, NRQCD matrix elements are ultraviolet divergent quantities that require renormalization.In particular, the color-octet matrix element ⟨O Q ( 3 S [8]  1 )⟩ acquires dependence on the renormalization scheme and scale from one loop, in a way that contributions to the cross section from different channels mix under renormalization.In χ cJ and χ bJ production, 3 P [1]  J and 3 S [8]  1 channels mix under changes of the NRQCD scale, and in spin-triplet S -wave quarkonium production, 3 P [8]  J and 3 S [8]  1 channels mix under renormalization [6,11].In calculations of the short-distance cross sections at one-loop level, the matrix elements are usually renormalized in the MS scheme at the scale of the heavy quark mass m.In this case, the short-distance cross sections for the 3 P [1]  J and 3 P [8]  J channels are negative at values of p T much larger than the heavy quarkonium mass.Because of the mixing, only the sum of the contributions from all channels is physically meaningful, while the contribution from a single channel can in principle become negative.
It is generally not known how to compute the NRQCD matrix elements from first principles, except for the color-singlet matrix elements at leading order in v, which can be related to decay matrix elements or quarkonium wavefunctions at the origin [6].Because of this, the color-octet matrix elements are usually determined phenomenologically by comparing Eq. (1) with measured cross section data.As will be explained in a following section, the values of color-octet matrix elements extracted from data depend strongly on the choice of measurements employed in the determination.This completely phenomenological approach to NRQCD matrix element determination has not lead to a satisfactory description of the quarkonium production mechanism, as none of the determinations have been able to give a comprehensive description of important observables associated with inclusive quarkonium production [7].Some representative examples of matrix element determinations will be shown in Sect.3.

Potential NRQCD
Recently, attempts have been made towards computing the NRQCD matrix elements in the potential NRQCD (pNRQCD) effective field theory [12][13][14].For strongly coupled quarkonia, pNRQCD provides expressions for NRQCD matrix elements in terms of quarkonium wavefunctions at the origin and universal gluonic correlators [15][16][17].The gluonic correlators are defined by vacuum expectation values of products of gluon field strengths and Wilson lines.A similar formalism has previously been used to compute NRQCD matrix elements for quarkonium decays into light particles; in this case, the gluonic correlators have different definitions from the ones that appear in the production matrix elements [18,19].In the case of color-singlet matrix elements, the pNRQCD expressions reproduce at leading order in v the known results in terms of the wavefunctions at the origin.For the color-octet case, the gluonic correlators appear from leading order in v in the pNRQCD expressions for the matrix elements.
Because the gluonic correlators do not involve heavy quark fields or projection operators, they are universal quantities that do not depend on the specific heavy quarkonium state.In particular, the same gluonic correlators appear in expressions for color-octet matrix elements for production of heavy quarkonia with different radial excitation or heavy quark flavor.Based on this point, it has been argued that the gluonic correlators are more basic quantities that are better suited for lattice QCD evaluations than the original definitions for NRQCD matrix elements, although a lattice calculation of the correlators is yet to be done [15][16][17].
Even though first-principles determinations of NRQCD matrix elements through lattice calculations of gluonic correlators have not been made possible yet, the pNRQCD expressions imply universal relations between color-octet matrix elements for heavy quarkonium states with different radial excitation or heavy quark flavor.Hence, the pNRQCD formalism allows simultaneous inclusion of charmonium and bottomonium data in a single analysis of color-octet matrix elements.This provides a strong constraint on color-octet matrix elements.
Table 1.Phenomenological determinations of J/ψ color-octet matrix elements in units of 10 −2 GeV 3 from refs.[20][21][22][23][24][25][26].In ref. [21], only two linear combinations of the three matrix elements are determined, and the maximum ranges of the matrix elements are obtained from maximizing and minimizing ⟨O J/ψ ( 1 S [8]  0 )⟩ under the assumption made in ref. [21] that the matrix elements are positive definite.Refs.[24,26] provide covariance matrices from which the correlations in the uncertainties can be obtained.The pNRQCD result from ref. [26] is from the analysis with p T larger than 5 times the heavy quarkonium mass.

NRQCD matrix elements
We now list some representative examples of NRQCD matrix element determinations based on calculations of the short-distance cross sections at next-to-leading order in the strong coupling.The color-octet matrix elements for J/ψ production from refs.[20][21][22][23][24][25][26] are shown in table 1.In table 2 we show ψ(2S ) matrix elements from refs.[21,[24][25][26][27].For bottomonium, we show color-octet matrix elements for production of Υ(1S ), Υ(2S ), and Υ(3S ) states from refs.[26,28,29] in table 3.In all cases listed here, the color-singlet matrix elements employed in the analyses are obtained from potential models or quarkonium decay rates, and are consistent within uncertainties.The color-octet matrix elements are determined by comparing Eq. ( 1) with cross Table 3. Phenomenological determinations of Υ color-octet matrix elements in units of 10 −2 GeV 3 from refs.[26,28,29].The analysis in ref. [28] employed a smaller NRQCD scale of 1.5 GeV whereas in refs.[26,29] the scale was chosen to be the bottom quark mass m b .Due to the running of the 3 S [8]  1 matrix element, at the scale m b the values of ⟨O Υ ( 3 S [8]  1 )⟩ from ref. [28] will be more negative for Υ(1S ) and Υ(2S ), and more positive for Υ(3S ) compared to what are listed in this table.It is also worth noting that the treatment of P-wave feeddowns in ref. [28] are inconsistent with measurements in ref. [30], which became available after ref. [28] was published.In ref. [29], only two linear combinations of the three matrix elements are determined, and the maximum ranges of the matrix elements are obtained from maximizing and minimizing ⟨O Υ ( 1 S [8]  0 )⟩ under the assumption that the matrix elements are positive definite similarly to what has been done in ref. [21].Ref. [26] provides covariance matrices from which the correlations in the uncertainties can be obtained.The pNRQCD results from ref. [26] are from the analysis with p T larger than 5 times the heavy quarkonium mass; the pNRQCD results for Υ(1S ) matrix elements are obtained by assuming the results in the strongly coupled pNRQCD formalism developed in refs.[15][16][17] also applies to the 1S state.1.17 ± 0.37 −0.16 ± 0.81 0.84 ± 0.27 section data, taking into account the effect of feeddowns.We can see that the resulting values of the color-octet matrix elements differ wildly, and even the signs of the matrix elements can be different, although none of the color-octet matrix elements exceed the typical sizes expected from the nonrelativistic power counting: they are usually more than an order of magnitude smaller than the color-singlet matrix element.With the exception of ref. [20], the matrix element extractions are solely based on p T -differential cross sections from hadron colliders with various choices of lower p T cuts.In these cases, an approximate degeneracy in the p T shapes of the short-distance cross sections can prevent strongly constraining all three color-octet matrix elements (see, e.g., ref. [31]).In refs.[21,29], only two linear combinations of the color-octet matrix elements were extracted, and the ranges of matrix elements were determined by assuming positivity of the matrix elements; the results shown in the tables correspond to two extreme cases where ⟨O Q ( 1 S [8]  0 )⟩ is maximized or minimized.Other hadroproduction-based determinations from refs.[24,26,27] employed covariancematrix analyses to obtain linear combinations of matrix elements that are more suited for phenomenological determinations; in many of these cases, one of the three linear combinations is poorly determined compared to others, which corresponds to the undetermined linear combination of matrix elements in refs.[21,29].The hadroproduction-based approaches lead to predictions of the spin-triplet S -wave quarkonium production mechanism that lie somewhere between two extreme scenarios: in the 1 S [8]  0 dominance scenario, the cross section is dominated by the 1 S [8]  0 channel contribution, while the sum of the 3 S [8]  1 and 3 P [8]  0 channel , 01011 (2022) https://doi.org/10.1051/epjconf/202227401011t h Quark Confinement and the Hadron Spectrum EPJ Web of Conferences 274 XV contributions are small; in the opposite scenario, the bulk of the cross section comes from the sum of the 3 S [8]  1 and 3 P [8]  0 channel contributions, while the 1 S [8]  0 contribution is small.In the hadroproduction-based approaches, the 3 S [8]  1 and 3 P [8]  0 matrix elements have same signs, so that the contributions from the two channels tend to cancel at large p T , because there the short-distance cross sections for the two channels have opposite signs.
The J/ψ matrix elements in ref. [20] were obtained from a global fit of cross section data including hadroproduction, photoproduction, and p T -integrated cross section at B factories.This helps lift the approximate degeneracy in the p T shapes of the short-distance cross sections, which allows all three color-octet matrix elements to be well determined.However, the values of the matrix elements obtained in the global fit are very different from hadroproduction-based approaches.In refs.[21][22][23][24][25][26], the signs of the 3 S [8]  1 and 3 P [8]  0 matrix elements are same, leading to cancellations between the two channels at large p T .In contrast, in the global fit the 3 P [8]  0 matrix element is negative, while the 3 S [8]  1 matrix element is positive, so that the contributions from the two channels add at large p T .Because of this, the J/ψ hadroproduction cross sections from the global fit tend to be in tension with measurements at very large p T .
The ψ(2S ) matrix elements in ref. [27] were also obtained from a global fit of available cross section data; however, unlike the J/ψ case, availability of ψ(2S ) production data is mostly limited to hadron collider experiments.The analysis with the cut p T > 1 GeV shows a pattern of color-octet matrix elements that is similar to the J/ψ global fit, yielding a negative 3 P [8]  0 matrix element, while the other two remain positive; this leads to predictions that are in tension with measurements of the p T shape of the ψ(2S ) cross section, as well as the polarization, as we will see in the next section.In contrast, an alternative analysis with the cut p T > 7 GeV presented in the same work results in color-octet matrix elements that are similar to the 3 S [8]  1 plus 3 P [8]  0 dominance scenario from other hadroproduction-based approaches.The J/ψ matrix element extractions in refs.[22,23] are based on large-p T hadroproduction data of J/ψ and η c at the LHC.As will be explained in Sect.5, inclusion of the η c data gives additional constraints to J/ψ matrix elements through approximate heavy quark spin symmetry, and results in configurations where the 1 S [8]  0 channel contribution to the J/ψ production rate is small.That is, the analyses based on J/ψ and η c hadroproduction data prefer the scenario where the sum of 3 S [8]  1 and 3 P [8]  J channel contributions dominate the J/ψ cross section.
The pNRQCD-based analysis in ref. [26] employed hadroproduction data of J/ψ, ψ(2S ), Υ(2S ), and Υ(3S ) at the LHC, by using the universal relations between color-octet matrix elements for spin-triplet S -wave quarkonia.This results in values of color-octet matrix elements that are better constrained than some conventional NRQCD approaches such as refs.[21,29].This happens because the pNRQCD analysis includes the charmonium and bottomonium data simultaneously in the extraction of matrix elements, and acquires sensitivity to the running of the 3 S [8]  1 matrix element at scales ranging from the charm to the bottom quark masses.Because the one-loop anomalous dimension of the 3 S [8]  1 matrix element is proportional to the 3 P [8]  0 matrix element, this constrains the 3 P [8]  0 matrix element to a positive definite value.As a result, the pNRQCD analysis yields a configuration of color-octet matrix elements where the bulk of the cross sections come from the sum of 3 S [8]  1 and 3 P [8]  J channel contributions for all spin-triplet S -wave quarkonia.
As have been shown in this section, the phenomenological determinations of color-octet matrix elements result in values that vary wildly depending on the choice of data.Notably, the large-p T analyses based on J/ψ and η c hadroproduction data [22,23], as well as the pNRQCD analysis based on charmonium and bottomonium hadroproduction data [26]  [8]  1 plus 3 P [8]  J dominance scenario, in contrast with global fits including low-p T data and other hadroproduction-based approaches favoring 1 S [8]  0 dominance.

Polarization of J/ψ, ψ(2S), and Υ in hadron colliders
The polarization of spin-triplet S -wave heavy quarkonia has long been considered an important test of the color-octet matrix elements.Early analyses based on tree-level calculations of the short-distance cross sections predicted that the J/ψ will be strongly transverse at large p T [32][33][34].This has not been supported by experiment: measurements at the LHC show little or no evidence of any strong polarization of spin-triplet S -wave quarkonia (see for example refs.[35,36]).
The tree-level prediction of transversely polarized J/ψ was based on the observation that only the 3 S [8]  1 channel can contribute appreciably at large p T .This no longer holds at one loop: all three color-octet channels can contribute at large p T through gluon fragmentation [31,[37][38][39][40][41][42].The polarization can still discriminate between different color-octet channels, because the polarization of the quarkonium is affected by the spin and orbital angular momentum of the Q Q produced in gluon fragmentation.For both 3 S [8]  1 and 3 P [8]  J channels, the transverse polarization of the fragmenting gluon is mostly transferred to the Q Q, because the fragmentation can occur by emitting soft gluons: as a result, the Q Q produced in 3 S [8]   1 and 3 P [8]  J channels is mostly transverse, while the longitudinal production rate is small; note that, due to the subtraction of the infrared divergence, the large-p T transverse production rate of Q Q( 3 P [8]  J ) is negative, while the longitudinal production rates are positive.On the other hand, the 1 S [8]  0 channel is isotropic, so it cannot produce polarized final states.Unpolarized spin-triplet S -wave quarkonia can be produced in two ways: if the production rate is dominated by the 1 S [8]  0 channel, then the quarkonium cannot be strongly polarized, because Q Q( 1 S [8]  0 ) is isotropic.In the 3 S [8]  1 plus 3 P [8]  J dominance scenario, the color-octet matrix elements for the two channels have same signs, so the transverse production rate largely cancels between the two channels, while the longitudinal cross sections add; this way, unpolarized final states can be obtained even when the 1 S [8]  0 channel contribution is small.In contrast, the global fit analyses with small lower p T cuts that give negative values for the 3 P [8]  0 matrix elements yield transversely polarized quarkonia at large p T , because the transverse production rates from the 3 S [8]  1 and 3 P [8]  J channels add.As a result, the polarization measurements at the LHC are in tensions with the predictions based on the global fit analyses with small lower p T cuts [20,27], while the polarization results based on large-p T hadroproduction measurements agree with experiments [21][22][23][24][25][26]29].
A shortcoming of the use of polarization for discriminating color-octet matrix elements is that it can hardly distinguish between the 1 S [8]  0 dominance and the 3 S [8]  1 plus 3 P [8]  J dominance scenarios, because both cases lead to unpolarized quarkonia.In the case of the pNRQCD analysis, which favors the 3 S [8]  1 plus 3 P [8]  J dominance scenario, Υ is predicted to be more transverse than J/ψ or ψ(2S ) due to the running of the 3 S [8]  1 matrix element coming from the large and positive 3 P [8]  J matrix element; this running would not have a prominent effect to polarization in the 1 S [8]  0 dominance scenario.Even though this prediction agrees with measurements of Υ polarization at the LHC [35], which show slightly more transverse polarization than J/ψ, this effect is numerically small and diluted by feeddown effects, especially for 1S and 2S bottomonia.This makes it desirable to have more observables that may help distinguish between the two competing scenarios.

Comparison with other observables
There have been quite a few observables related to inclusive quarkonium production measured in collider experiments, but many have not been able to strongly scrutinize the heavy quarkonium production mechanism.For example, the Belle measurement for the total inclusive J/ψ production rate [43] involves an unknown branching fraction into four or more charged tracks; furthermore, it is unclear whether the form of NRQCD factorization given in Eq. ( 1) would hold for total inclusive production rates 2 .In the case of photoproduction, measurements at the DESY HERA [44,45] were made with kinematical cuts on the elasticity, which can make it difficult for NRQCD to make reliable predictions [46].Studies of J/ψ momentum distribution in jet [47,48] showed that the measured distribution from LHCb [49] is incompatible with the global fit results for J/ψ matrix elements, while the matrix elements in the 1 S [8]  0 dominance scenario lead to results that are in fair agreement with measurements.Although a calculation based on the 3 S [8]  1 plus 3 P [8]  J dominance scenario has not been done in ref. [48], we can expect that this will yield results that are qualitatively similar to the 1 S [8]  0 dominance scenario, based on the general behavior of the shapes of the 3 S [8]  1 and 3 P [8]   J channel contributions to the distribution.There are still several observables proposed in the literature that can help distinguish the different scenarios for quarkonium production mechanism.The η c production rate measured by LHCb [50,51] has been considered a good observable, as it gives additional constraints for J/ψ matrix elements based on heavy quark spin symmetry.Heavy quark spin symmetry implies that the J/ψ and η c matrix elements that differ by one unit of the Q Q spin are same at leading order in v, up to calculable spin multiplicity factors.That is, the 1 S [8]  0 matrix element for J/ψ determines the 3 S [8]  1 matrix element for η c .In the case of η c , the cross section is dominated by 1 S [1]  0 and 3 S [8]  1 channels, so that the measured cross section gives a strong constraint on the 1 S [8]  0 matrix element for J/ψ [22,23,52].The measurements imply that the 1 S [8]  0 contribution to the J/ψ cross section must be small, because the color-singlet contribution makes up for the bulk of the measured η c production rate.As a result, the analyses based on J/ψ and η c production data favor the 3 S [8]  1 plus 3 P [8]  J dominance scenario, as have been presented in the previous section.Similarly, the pNRQCD analysis leads to η c production rates that are compatible with measurements, albeit with large uncertainties due to the limited precision for ⟨O J/ψ ( 1 S [8]  0 )⟩.A shortcoming of the NRQCD description of the η c production rate currently adopted in the literature is that, unlike the J/ψ case, the contribution from the color-singlet channel at leading order in v is significant; recall that, in the case of S -wave quarkonia, the color-octet matrix elements are suppressed by several powers of v compared to the color-singlet one.This means that in the η c case, it may be necessary that the relativistic corrections to the color-singlet channel must be included up to relative order v 4 , because they can be the same order as the color-octet contributions.The tension between measurement and NRQCD calculations of exclusive production rates of η c at B factories, which only involve color-singlet contributions, may imply that the relativistic corrections to the color-singlet channel can be significant [53][54][55][56][57][58][59].This effect has so far not been taken into account in existing analyses of inclusive η c production.
Another observable that may help discriminate the quarkonium production mechanism is the associated production of a heavy quarkonium plus a gauge boson.Calculations of short-distance cross sections at one-loop level have been done for the production of J/ψ plus a photon [60], and the weak gauge bosons W and Z [61].Measurements have been made available by ATLAS for production of J/ψ + W [62,63] and J/ψ + Z [64].The data are 2 One reason would be that in such case, the Q and Q do not necessarily need to be produced within a distance of 1/m in order to produce a quarkonium, which would not allow for the usual form of NRQCD factorization to hold.
, 01011 (2022) https://doi.org/10.1051/epjconf/202227401011t h Quark Confinement and the Hadron Spectrum EPJ Web of Conferences 274 XV mostly available for J/ψ transverse momentum larger than the J/ψ mass.A recent analysis from ref. [61] shows that only the3 S [8]  1 plus 3 P [8]  J dominance scenario 3 results in predictions for the J/ψ + W and J/ψ + Z production rates that are compatible with measurements; the 1 S [8]  0 dominance scenario can even lead to negative direct cross sections, and the global fit that mainly comes from low-p T data gives cross sections that underestimate data.
Finally, predictions for the J/ψ production cross sections from electron-proton collisions at the Electron-Ion Collider have recently been made available [65].It has been shown that the large-p T hadroproduction-based analyses lead to predictions that are distinct from what is obtained from the global fit [26,65].While the 3 S [8]  1 plus 3 P [8]  J dominance scenario yields slightly larger p T -differential production rates than the 1 S [8]  0 dominance scenario, precise measurements of the cross sections at large p T will be needed to distinguish between the two scenarios.

Summary and outlook
In this paper we have presented a concise review of the current status of phenomenology of inclusive heavy quarkonium production and polarization based on nonrelativistic effective field theories.Theoretical calculations of heavy quarkonium production rates in the nonrelativistic QCD (NRQCD) factorization formalism require perturbative calculations of the short-distance cross sections as well as nonperturbative determinations of NRQCD matrix elements.While color-singlet matrix elements have been computed in potential models and lattice QCD or determined from decay rates, color-octet matrix elements have not been computed from first principles.In the case of the production of spin-triplet S -wave quarkonia, color-octet matrix elements for the 3 S [8]  1 , 1 S [8]  0 , and the 3 P [8]  J channels have significant contributions to the cross section.
While perturbative QCD calculations of short-distance cross sections have been carried out at one-loop accuracy for many important processes including hadroproduction and polarization at the LHC, results for phenomenological determinations of NRQCD matrix elements depend strongly on the choice of data.For J/ψ and ψ(2S ), large-p T hadroproduction-based determinations lead to scenarios where the cross section is dominated by either the 1 S [8]   0 channel or the remnant of the cancellation between 3 S [8]  1 and 3 P [8]  J channels that mix under renormalization.On the other hand, global fits that include data with p T similar or smaller than the heavy quarkonium mass result in values of color-octet matrix elements that make 3 S [8]  1 and 3 P [8]  J channel contributions add at large p T .The matrix element determinations from global fits including low-p T data [20,27] and large-p T hadroproduction-based analyses [21-26, 28, 29] lead to contrasting predictions for polarization.Because both the 3 S [8]  1 and 3 P [8]  J channels are strongly transversely polarized, global fits including low p T data predict transversely polarized charmonia at large p T , while hadroproduction-based approaches predict almost no polarization.Polarization measurements at the LHC disfavor the low-p T global fit predictions, showing no strong evidence of polarization.While LHC polarization measurements seem to agree with predictions from the hadroproduction-based analyses, polarization cannot strongly discriminate between the 1 S [8]  0 dominance and 3 S [8]  1 plus 3 P [8]  J dominance scenarios, because they both lead to similar near-zero polarization predictions.
This unfavorable situation could be improved by efforts from both theory and experiment.On the theory side, the potential NRQCD (pNRQCD) effective field theory has been employed to further factorize the NRQCD matrix elements into quarkonium wavefunctions at the origin and universal gluonic correlators [15][16][17]26].While first-principles determinations of NRQCD matrix elements through lattice calculations of the gluonic correlators are yet to be done, the universality of the gluonic correlators give rise to relations between color-octet matrix elements for different heavy quarkonium states, which provide additional constraints in phenomenological extractions of matrix elements.Analyses based on pNRQCD calculations of the color-octet matrix elements and large-p T hadroproduction data favor the 3 S [8]   1 plus 3 P [8]  J dominance scenario for all spin-triplet S -wave quarkonium states including J/ψ, ψ(2S ), and Υ.On the experimental side, measurements of additional observables such as η c production [50,51] and the associated production of J/ψ + W [62,63] and J/ψ + Z [64] at the LHC have also been shown to prefer the 3 S [8]  1 plus 3 P [8]  J dominance scenario [22,23,61].While it looks promising that analyses based on large-p T production seem to be converging to the 3 S [8]  1 plus 3 P [8]  J dominance scenario, it is well known that these approaches lead to bad descriptions of low-p T observables, including total inclusive production rates in lepton colliders and photoproduction cross sections at HERA [66].Even in hadroproduction, analyses based on large p T production have trouble describing low-p T data, as has been demonstrated in ref. [27].The fact that the heavy quarkonium production mechanism that correctly describes both high and low p T regions still remains out of reach suggests that there is much more to be understood in QCD and factorization formalisms.