Study of $D_J$ meson decays to $D^+\pi^-$, $D^0 \pi^+$ and $D^{*+}\pi^-$ final states in pp collision

A study of $D^+\pi^-$, $D^0 \pi^+$ and $D^{*+}\pi^-$ final states is performed using pp collision data, corresponding to an integrated luminosity of 1.0 $fb^{-1}$, collected at a centre-of-mass energy of 7 TeV with the LHCb detector. The $D_1(2420)^0$ resonance is observed in the $D^{*+}\pi^-$ final state and the $D^*_2(2460)$ resonance is observed in the $D^+\pi^-$, $D^0 \pi^+$ and $D^{*+}\pi^-$ final states. For both resonances, their properties and spin-parity assignments are obtained. In addition, two natural parity and two unnatural parity resonances are observed in the mass region between 2500 and 2800 MeV. Further structures in the region around 3000 MeV are observed in all the $D^{*+}\pi^-$, $D^+\pi^-$ and $D^0 \pi^+$ final states.


Introduction
Charm meson spectroscopy provides a powerful test of the quark model predictions of the Standard Model. Many charm meson states, predicted in the 1980s [1], have not yet been observed experimentally. The expected spectrum for the cū system is shown in Fig. 1 (the spectrum of the cd system is almost identical). The J P states having P = (−1) J and therefore J P = 0 + , 1 − , 2 + , ... are called natural parity states and are labelled as D * , while unnatural parity indicates the series J P = 0 − , 1 + , 2 − , .... The low-mass spectrum of the cū system is comprised of the ground states (1S), the orbital excitations with angular momentum L=1, 2 (1P, 1D), and the first radial excitations (2S). Apart from the ground states (D, D * ), only two of the 1P states, D 1 (2420) and D * 2 (2460) [2], are experimentally well established since they have relatively narrow widths (∼30 MeV). 1 In contrast, the broad L = 1 states, D * 0 (2400) and D 1 (2430), have been established by the Belle and BaBar experiments in exclusive B decays [3,4].
The theoretical predictions are in agreement (within 20-30 MeV) with observations for the 1S states and the J P = 2 + and J P = 1 + 1P states. In the cs system, the J P = 0 + and J P = 1 + states (both L = 1) have predicted masses about 100 MeV higher than the measured masses of the D sJ mesons. To quantitatively assess the accuracy of the quark model predictions, assumptions are needed to formulate a wave equation for quark-antiquark bound states starting from the QCD Lagrangian [5]. Nevertheless, the 1 We work in units where c = 1.  Figure 1: Modified Godfrey-Isgur mass predictions [1]. The figure shows the cū spectrum in which the masses have been scaled such that the ground state coincides with the D 0 mass. The 2 − states, not shown in the original publication, have been inserted following the splitting structure of the 1P states.
discrepancy between the predictions of various models and the mass measurements [6][7][8][9] suggests that some observed states in the cs case are not simple quark-antiquark configurations. Possible interpretations include more complex structures, such as bound states ("molecules") of other mesons [10], or mixtures of conventional quark-antiquark with four-quark components [11]. The properties of hadrons can be computed from the QCD Lagrangian using lattice calculations and the resulting cū, cd and cs mass spectra can be compared to measurements. However, the calculation of hadronic quantities for dynamical light quarks is still a challenging task and different results are obtained [12][13][14][15][16][17].
To search for excited charmed mesons, labelled D J , BaBar analyzed the inclusive production of the D + π − , D 0 π + and D * + π − final states in the inclusive reaction e + e − → cc → D ( * ) πX, where X is any additional system [18]. 2 They observe four signals, labelled D(2550) 0 , D * (2600) 0 , D(2750) 0 and D * (2760) 0 , and the isospin partners D * (2600) + and D * (2760) + . This paper reports a search for D J mesons in a data sample, corresponding to an integrated luminosity of 1.0 fb −1 , of pp collisions collected at a centre-of-mass energy of 7 TeV with the LHCb detector.

Detector
The LHCb detector [19] is a single-arm forward spectrometer covering the pseudorapidity range 2 < η < 5, designed for the study of particles containing b or c quarks. The detector includes a high precision tracking system consisting of a silicon-strip vertex detector surrounding the pp interaction region, a large-area silicon-strip detector located upstream of a dipole magnet with a bending power of about 4 Tm, and three stations of silicon-strip detectors and straw drift tubes placed downstream. The combined tracking system has momentum resolution that varies from 0.4% at 5 GeV to 0.6% at 100 GeV, and impact parameter resolution of 20 µm for tracks with high transverse momentum p T with respect to the beam direction. The impact parameter is defined as the perpendicular distance between the track path and the position of a pp collision. Charged hadrons are identified using two ring-imaging Cherenkov (RICH) detectors. Photon, electron and hadron candidates are identified by a calorimeter system consisting of scintillating-pad and preshower detectors, an electromagnetic calorimeter and a hadronic calorimeter. Muons are identified by a system composed of alternating layers of iron and multiwire proportional chambers. The trigger [20] consists of a hardware stage, based on information from the calorimeter and muon systems, followed by a software stage which applies a full event reconstruction.

Event selection
The search for D J mesons is performed using the inclusive reactions where X represents a system composed of any collection of charged and neutral particles. The charmed mesons in the final state are reconstructed in the decay modes D + → K − π + π + , D 0 → K − π + and D * + → D 0 π + . Charged tracks are required to have good track fit quality, momentum p > 3 GeV and p T > 250 MeV. These conditions are relaxed to lower limits for the pion originating directly from the D * + decay. Tracks pointing to a pp collision vertex (primary vertex) are rejected by means of an impact parameter requirement in the reconstruction of the D + and D 0 candidates. All tracks used to reconstruct the mesons must have a distance of closest approach to each other smaller than 0.5 mm. The cosine of the angle between the momentum of the D meson candidate and its direction, defined by the positions of the primary vertex and the meson decay vertex, is required to be larger than 0.99999. This ensures that the D meson candidates are produced at the primary vertex and reduces the contribution from particles originating from b-hadron decays. The reconstructed D + , D 0 and D * + candidates are combined with all the right-sign charged pions in the event. Each of the D + π − , the D 0 π + , and the D * + π − candidates are fitted to a common vertex with χ 2 /ndf < 8, where ndf is the number of degrees of freedom. The purity of the charmed meson candidates is enhanced by requiring the decay products to be identified by the RICH detectors, using the difference in the log-likelihood between the kaon and pion hypotheses ∆ ln L Kπ [21]. We require ∆ ln L Kπ > 3 for kaon tracks and a loose requirement of ∆ ln L Kπ < 10 for pions. In the reconstruction of D + → K − π + π + decays, a small D * + signal in the D 0 π + mass spectrum is removed by demanding ∆m ≡ m(K − π + π + ) − m(K − π + ) > 152 MeV.
In order to reduce combinatorial background, the cosine of the angle between the momentum direction of the charged pion in the D ( * ) π ± rest frame and the momentum direction of the D ( * ) π ± system in the laboratory frame is required to be greater than zero. Due to the possible presence of multiple primary vertex candidates, it is required that the D ( * ) and the π ± point to the same primary vertex.
To reduce any dependence on the mass scale, the invariant mass of the D ( * ) π ± system is calculated from the measured mass difference. For example, the D 0 π + invariant mass is given by where m D 0 is the known value of the D 0 mass [2]. Figure 2 shows the K − π + π + , K − π + and D 0 π + invariant mass spectra after the selection criteria are applied. The distributions are fitted by the sum of two Gaussian functions, with a common mean to describe the signal shape and a linear term to describe the combinatorial background. The mean mass resolutions for the three distributions are 8.  and contain 15.1×10 6 , 20.4×10 6 and 6.4×10 6 candidates for the D + π − , D 0 π + and D * + π − modes, respectively.

Mass spectra
The D + π − , D 0 π + and D * + π − mass spectra are shown in Fig. 3. The D + π − and D 0 π + mass spectra evidence strong D * 2 (2460) signals, while in the D * + π − mass spectrum clear D 1 (2420) 0 and D * 2 (2460) 0 signals are visible. A further reduction of the combinatorial background is achieved by performing an optimization of the signal significance and purity as a function of p T of the D ( * ) π ± system using the well known D 1 (2420) and D * 2 (2460) resonances. 3 For this purpose, we fit the three mass spectra as explained in Sect. 7 and Sect. 9 and obtain, for each resonance, the signal yield N S and background yield N B events. We compute the signal significance S = N S / √ N S + N B and signal purity P = N S /(N S + N B ) and find that the requirement p T (D ( * ) π) > 7.5 GeV provides a good compromise between significance and purity. After the optimization 7.9×10 6 , 7.5×10 6 and 2.1×10 6 D + π − , D 0 π + and D * + π − candidates are obtained. We also study the dependence of the signal to background ratio in the three mass spectra on the pseudorapidity of the D ( * ) π system and find a very weak correlation. We analyze, for comparison and using the same selections, the wrong-sign D + π + , D 0 π − and D * + π + combinations which are also shown in Fig. 3.
The D + π − mass spectrum, Fig. 3(a), shows a double peak structure around 2300 MeV due to cross-feed from the decay where the π 0 /γ is not reconstructed; the last number, in parentheses, indicates the branching fraction of D * + → D + π 0 /γ decays [2]. We observe a strong D * 2 (2460) 0 signal and weak structures around 2600 and 2750 MeV. The wrong-sign D + π + mass spectrum does not show any structure.
The D 0 π + mass spectrum, Fig. 3(b), shows an enhanced double peak structure around 2300 MeV due to cross-feed from the decays The D * 2 (2460) + signal and weak structures around 2600 and 2750 MeV are observed. In comparison, the wrong-sign D 0 π − mass spectrum does show the presence of structures in the 2300 MeV mass region, similar to those observed in the D 0 π + mass spectrum. These structures are due to cross-feed from the decay  The D * + π − mass spectrum, Fig. 3(c), is dominated by the presence of the D 1 (2420) 0 and D * 2 (2460) 0 signals. At higher mass, complex broad structures are evident in the mass region between 2500 and 2800 MeV.

Simulation
Simulated events are used to study the effects of the detector on the observed mass resolution. The pp collisions are generated using Pythia 6.4 [22] with a specific LHCb configuration [23]. Decays of hadronic particles are described by EvtGen [24] and the interaction of the generated particles with the detector and its response are implemented using the Geant4 toolkit [25,26] as described in Ref. [27].
Simulated events are reconstructed in the same manner as data. We analyze samples of full detector simulations to estimate the reconstruction efficiency, mass resolution and possible bias in the reconstruction chain. We also make use of simple event generator level simulations [28] to study kinematic effects. The tight trigger conditions and selection criteria have the effect of producing very low integrated efficiencies, which we calculate to be (0.149±0.004)%, (0.056±0.005)% and (0.064±0.003)% for D + π − , D 0 π + and D * + π − candidates, respectively.
To estimate the detector resolution we compare generated and reconstructed invariant masses and obtain experimental resolutions as functions of the reconstructed mass. The analysis of these simulated samples shows no bias in the reconstructed invariant masses. We estimate resolutions which, in the mass region between 2000 and 2900 MeV, are similar for the three mass spectra and range from 1.0 to 4.5 MeV as a function of the mass. Since the widths of the resonances appearing in the three mass spectra are much larger than the experimental resolutions, resolution effects are neglected.

Mass fit model
Binned χ 2 fits to the three mass spectra are performed. The D * 2 (2460) and D * 0 (2400) signal shapes in two-body decays are parameterized with a relativistic Breit-Wigner that includes the mass-dependent factors for a D-wave and S-wave decay, respectively. The radius entering in the Blatt-Weisskopf [29] form factor is fixed to 4 GeV −1 . Other resonances appearing in the mass spectra are described by Breit-Wigner lineshapes. All Breit-Wigner expressions are multiplied by two-body phase space. The cross-feed lineshapes from D 1 (2420) and D * 2 (2460) appearing in the D + π − and D 0 π + mass spectra are described by a Breit-Wigner function fitted to the data with the parameters given in Table 1. Resonances are included sequentially in order to test the χ 2 improvement when a new contribution is included. The background B(m) is described by an empirical shape [18] B(m) = P (m)e a 1 m+a 2 m 2 for m < m 0 , B(m) = P (m)e b 0 +b 1 m+b 2 m 2 for m > m 0 , where P (m) is the two-body phase space and m 0 is a free parameter.
The two functions and their first derivatives are required to be continuous at m 0 such that Therefore the background model has four free parameters: m 0 , a 1 , a 2 and b 2 .
7 Fit to the D * + π − mass spectrum Due to the three-body decay and the availability of the helicity angle information, the fit to the D * + π − mass spectrum allows a spin analysis of the produced resonances and a separation of the different spin-parity components. We define the helicity angle θ H as the angle between the π − and the π + from the D * + decay, in the rest frame of the D * + π − system. Full detector simulations are used to measure the efficiency as a function of θ H , which is found to be uniform.
It is expected that the angular distributions are proportional to sin 2 θ H for natural parity resonances and proportional to 1 + h cos 2 θ H for unnatural parity resonances, where h > 0 is a free parameter. The D * π decay of a J P = 0 + resonance is forbidden. Therefore candidates selected in different ranges of cos θ H can enhance or suppress the different spin-parity contributions. We separate the D * + π − data into three different categories, summarized in Table 2. The candidate yields for these categories are given in Table 3, which also reports the mass intervals, the number of bins, and the resulting χ 2 /ndf in the fits to the different mass spectra. The data and fit for the D * + π − enhanced unnatural parity sample are shown in Fig. 4 and the resulting fit parameters are summarized in Table 4. The mass spectrum is dominated by the presence of the unnatural parity D 1 (2420) 0 resonance. The fitted natural parity D * 2 (2460) 0 contribution is consistent with zero, as expected. To obtain a good fit to the mass spectrum, three further resonances are needed. We label them D J (2580) 0 , D J (2740) 0 , and D J (3000) 0 . The presence of these states in this sample indicates unnatural parity assignments.
The masses and widths of the unnatural parity resonances are fixed in the fit to the natural parity sample. The fit is shown in Fig. 5 and the obtained resonance parameters are summarized in Table 4. The mass spectrum shows that the unnatural parity resonance   Fig. 6. We observe, as expected, small contributions from the natural parity resonances. We also fit the total D * + π − mass spectrum, again with all the resonance parameters fixed. The data and fit are shown in Fig. 7. Table 4 summarizes the measured resonance parameters and yields. The resonance parameters are obtained from the fits to the enhanced unnatural parity sample and natural parity sample, apart for the parameters of the D 1 (2420) 0 resonance, which are extracted from the fit to the total sample. The significances are computed as ∆χ 2 where ∆χ 2 is the difference between the χ 2 values when a resonance is included or excluded from the fit while all the other resonances parameters are allowed to vary. All the statistical significances are well above 5σ.
8 Spin-parity analysis of the D * + π − system In order to obtain information on the spin-parity assignment of the states observed in the D * + π − mass spectrum, the data are subdivided into ten equally spaced bins in cos θ H . The ten mass spectra are then fitted with the model described above with fixed resonance parameters to obtain the yields as functions of cos θ H for each resonance.
The resulting distributions for D 1 (2420) 0 and D * 2 (2460) 0 are shown in Fig. 8. They have been fitted using the functions described in Table 5. A good description of the data is obtained in terms of the expected angular distributions for J P = 1 + and J P = 2 + resonances. Figure 9 shows this case we compare the distributions with expectations from natural parity, unnatural parity and J P = 0 − . In the case of unnatural parity, the h parameter, in 1 + h cos 2 θ H , is constrained to be positive and therefore the fit gives h = 0. In both cases, the distributions are best fitted by the natural parity hypothesis. Figure 10 shows the angular distributions for the D J (2580) 0 , D J (2740) 0 and D J (3000) 0 states. The distributions are fitted with natural parity and unnatural parity. The J P = 0 − hypothesis is also considered for D J (2580) 0 . The results from the fits are given in Table 5. In all cases unnatural parity is preferred over a natural parity assignment.
9 Fit to the D + π − and D 0 π + mass spectra The D + π − and D 0 π + mass spectra consist of natural parity resonances. However these final states are affected by cross-feed from all the resonances that decay to the D * π final state. Figures 3(a) and (b) show (in the mass region around 2300 MeV) cross-  feed contributions from D 1 (2420) and D * 2 (2460) decays. However we also expect (in the mass region between 2400 and 2600 MeV) the presence of structures originating from the complex resonance structure present in the D * π mass spectrum in the mass region between 2500 and 2800 MeV.
To obtain an estimate of the lineshape and size of the cross-feed, we normalize the D * + π − mass spectrum to the D + π − mass spectrum using the sum of the D 1 (2420) 0 and D * 2 (2460) 0 yields in the D * + π − mass spectrum (N sig ) and the sum of the cross-feed in the D + π − mass spectrum (N feed D + π − ). We estimate that each resonance appearing in the Table 4: Resonance parameters, yields and statistical significances. The first uncertainty is statistical, the second systematic.
Significances are evaluated using the method described in the text. D * + π − should also appear in the D + π − mass spectrum with a yield given by where R D + π − = N feed D + π − /N sig . Here N (D * + π − ) is the yield measured in the D * + π − final state, N (D + π − ) is the expected yield in the D + π − mass spectrum and R D + π − = 1.41±0.02 where the uncertainty is statistical only. Assuming similar yields for the D 1 (2420) + and D * 2 (2460) + resonances, we estimate for the D 0 π + channel, where R D 0 π + = N feed D 0 π + /N sig = 1.87±0.02 is the corresponding value for the D 0 π + channel. To obtain the expected lineshape of the cross-feed in the D + π − final state, we perform a study based on a generator level simulation. We generate D * J (2650) 0 , D * J (2760) 0 , D J (2580) 0 and D J (2740) 0 decays according to the chain described in Eq. (3). Given the small branching fraction of the D * + → D + γ decay, (1.6 ± 0.4)%, we only generate the D * + → D + π 0 decay. The parameters of the resonances are as reported in Table 4 and the decays to D * + π − are uniform over phase space. We then compute the resulting D + π − mass spectra and normalize each contribution to the measured yields. The overall resulting structures are then scaled by the factor R D + π − and superimposed on the D + π − mass spectrum shown in Fig. 11.
Similarly, to obtain the expected lineshape of the cross-feed in the D 0 π + final state, we generate the four resonances according to the decays shown in Eq. (4). We assume, for the charged modes, rates for the four states similar to that for the neutral modes. The overall resulting structures obtained for the D * 0 → D 0 π 0 and D * 0 → D 0 γ decays are scaled according to their branching fractions and the distribution is scaled by the factor R D 0 π + discussed above. The resulting contribution is superimposed on the D 0 π + mass spectrum shown in Fig. 12.
The cross-feed lineshapes obtained by the generator level simulation are not precise enough to be included in the fits to the D + π − and D 0 π + mass spectra. We therefore follow an empirical procedure to obtain good fits in this mass region. We first notice that these contributions produce a distortion of the D * 2 (2460) 0 and D * 2 (2460) + lineshapes. These are accommodated in the fit by means of a Breit-Wigner function, which we include to obtain a good description of the data. The parameters of the Breit-Wigner function are M = 2414.3 ± 1.4 MeV and Γ = 103.2 ± 2.7 MeV for the D + π − final state and M = 2435.1 ± 5.2 MeV and Γ = 106.9 ± 6.2 MeV for the D 0 π + final state. We consider these contributions as methods to improve the description of the cross-feeds.
We expect, in both the D + π − and D 0 π + mass spectra, the presence of D * J (2650) and D * J (2760) states. Enhancements in these mass regions can be seen in the two mass spectra shown in Fig. 11 and Fig. 12. However the D * J (2650) region is strongly affected by cross-feed. We include a simple Breit-Wigner function to describe these mass regions and obtain M = 2621.7 ± 1.4 MeV and Γ = 119.7 ± 6.5 MeV for the D + π − final state and M = 2599.9 ± 0.9 MeV and Γ = 72.3 ± 4.0 MeV for the D 0 π + final state. However the parameters so far obtained are strongly biased by the presence of the cross-feed and we therefore report, for the D * J (2650) resonance, only the results obtained from the fit to the D * + π − mass spectrum.
To obtain good quality fits we add broad structures around 3000 MeV, which we label D * J (3000) 0 and D * J (3000) + . Their parameters are derived from the fit to the D + π − mass spectrum and then fixed in the fit to the D 0 π + mass spectrum, where the effect is weaker.
The sensitivity of the fits to the presence of the broad D * 0 (2400) resonance is tested by performing simulations that include a D * 0 (2400) resonance with parameters fixed to their known values and vary the background lineshape within a wide range of values. We find a high correlation between the D * 0 (2400) parameters and the background lineshape and a failure of the fit to obtain correct estimates of its parameters and yields. Therefore this contribution is not included in the fit.
The fits to the D + π − and D 0 π + mass spectra are shown in Fig. 11 and Fig. 12, respectively. Masses, widths, yields, and significances for the all the fitted resonances are displayed in Table 4. All the statistical significances are well above 5σ.

Cross-checks and systematic uncertainties
Several cross-checks are performed to test the stability of the fits and their correct statistical behaviour. We first repeat all the fits, including the spin-parity analysis, lowering the p T requirement from 7.5 to 7.0 GeV. We find that all the resonance parameters vary within their statistical uncertainties and that the spin-parity assignments are not affected by this selection.
The fits stability and the uncertainties on the resonance parameters are tested using random variations of the histogram contents. For each histogram, we obtain and fit 500 new histograms by random Poisson variation of each bin content. We find in all cases a Gaussian behaviour of all the fit components with r.m.s. values that agree well with the statistical uncertainties given by the fits.
The systematic uncertainties on the resonance parameters and yields reported in Table 4 are estimated as follows. The background lineshape uncertainty is estimated using an alternative function B(m) = (m − m th ) a e −b 1 m−b 2 m 2 −b 3 m 3 , where m th is the threshold mass. This function gives acceptable fits for the Dπ mass spectra but generally a worse description of the threshold region.
The background lineshapes are additionaly tested by random variation of their parameters. For each mass spectrum, we generate and fit 500 new histograms where the resonance parameters and yields are fixed to the values obtained from the data, while the background yield is fixed but has parameters varying within ±3σ from the values obtained from the data. The distributions obtained from these fits are used to obtain systematic uncertainties due to the background lineshape. For the uncertainty due to the background lineshape the largest value between the estimates from the two methods described above is taken.
In the fits to the D * + π − mass spectra, where resonances have in some cases fixed parameters, we let the resonance parameters float sequentially. The procedure is repeated for each helicity sample and for the fit to the total mass spectrum. The Breit-Wigner shapes used to describe the D * J (2760) resonance in the D + π − and D 0 π + mass spectra are replaced by a relativistic Breit-Wigner functions with different spin assignments. We also include the D * 0 (2400) resonance with parameters fixed to the known values and obtain a small improvement in the fit to the D + π − mass spectrum but a yield consistent with zero in the fit to the D 0 π + mass spectrum.
The various estimated systematic uncertainties are added in quadrature. We do not report systematic uncertainties on the structures labelled as D J (3000) 0 and D * J (3000) + because, being at the limit of the mass spectra, they are strongly correlated with the background parameters.

Discussion and conclusions
Due to the correlation between the resonances parameters, this procedure also affects the properties of the other states appearing in the D * + π − mass spectrum. In the present analysis, as stated above, we measure important cross-feeds in the 2500-2600 MeV region of the D + π − and D 0 π + final states and therefore we obtain the D * J (2650) 0 parameters from the D * + π − final state only.
We compare the quark-model predictions given in Fig. 1 with our mass measurements and spin-parity analysis reported in Table 4 and Table 5, respectively. The BaBar analysis suggests a J P = 0 − assignment for the D J (2580) 0 state (labelled 2S D 0 (2558) in Fig. 1). Our results are consistent with the BaBar measurement, but cannot confirm it, due to the superposition of many relatively broad resonances in a limited mass region which complicates the extraction of the resonances parameters.
The D * J (2650) 0 resonance is observed to decay to D * + π − , has natural parity and therefore is expected to decay to Dπ. However the presence of this state in the Dπ mass spectra is obscured by the presence of cross-feeds from the D * π channels. We tentatively identify the D * J (2650) 0 resonance as a J P = 1 − state (2S D * 1 (2618)). The D * J (2760) 0 is observed in the D * + π − and D + π − decay modes with consistent parameters. We also observe the D * J (2760) + in the D 0 π + final state which can be identified as a J P = 1 − state (1D D * 1 (2796)). The D J (2740) 0 could be identified as the J P = 2 − (1D D 2 (2801)) resonance, although in this case the measured and predicted mass do not agree well. Definitive spin-parity assignments will be possible if these states are observed in B decays.