1 Introduction

In the framework of Regge theory of soft hadronic interactions, the energy dependence of total hadron–hadron scattering cross sections is described only after taking into account a specific type of effective exchange with vacuum quantum numbers [1]. Although it is used in various contexts, such an exchange is often referred to as a ‘pomeron’ (\({I\!\!P}\)) [2]. Pomeron exchange is a tool to describe diffractive processes, which are characterised by large gaps, devoid of activity, in the rapidity distribution of final state particles.

Diffractive processes in electron–protonFootnote 1 deep inelastic scattering were observed already in the very early part of the HERA experimental program [3, 4] and lead to revived interest in this class of soft peripheral hadronic interactions [5]. In reactions of the type \(ep \rightarrow eXY\) they are characterised by a large gap in rapidity between the systems X and Y. The system X can be considered as resulting from a diffractive dissociation of the virtual photon, while the system Y consists of the initial state proton or its low mass hadronic excitation, scattered at a small momentum transfer squared t relative to the initial state proton.

Perturbative quantum chromodynamics (pQCD) calculations are applicable in deep inelastic scattering even though the partonic structure of the proton is a priori unknown. In order to overcome this difficulty, the collinear factorisation theorem [6] is used, where the calculation of deep inelastic scattering (DIS) cross sections is described by a process-dependent partonic hard scattering part convoluted with a universal set of parton distribution functions of the proton (PDF). Collinear factorisation, therefore, opens the possibility to extract PDFs from one process and use them to predict cross sections for another process. For the PDF extraction the validity of the DGLAP evolution equations [7,8,9] is assumed.

A similar strategy can also be adapted to diffractive deep inelastic scattering (DDIS), where collinear factorisation is expected to be valid as well [10]. Assuming in addition the validity of proton vertex factorisation [11], diffractive processes are described by the exchange of collective colourless partonic states, such as the pomeron. Diffractive parton distribution functions (DPDFs) are extracted from diffractive data [12, 13]. Similarly to the normal PDFs, the DPDFs are expected to evolve as a function of the scale as predicted by the DGLAP equations. QCD analyses of diffractive data show that gluons constitute the main contribution to the DPDFs [12, 13]. To date, analyses of HERA data support the validity of the collinear factorisation theorem in DDIS as evidenced by experimental results on inclusive production [12, 13], dijet production [13,14,15,16,17,18,19] and \(D^{*}\) production [20,21,22,23,24].

Here, a new measurement of \(D^{*}(2010)\) meson production in DDIS is presented, where the \(D^{*}\) is reconstructed in the \(D^{*}\rightarrow K \pi \pi \) decay channel. The \(D^{*}\) meson originates from the fragmentation of a charm quark, which is produced at HERA energies mainly via the boson-gluon-fusion (\(\gamma ^{*}g \rightarrow c\bar{c}\)) process. Hence, the gluon content of the pomeron can be accessed directly, and allows the collinear factorisation to be tested. Compared to the previous H1 publication [23] the analysis presented corresponds to a sixfold increase in the integrated luminosity.

Fig. 1
figure 1

The leading order diagram for open charm production in diffractive DIS at HERA in the picture of collinear and proton vertex factorisation

2 Kinematics of diffractive deep inelastic scattering

The standard DIS kinematics is described in terms of the invariants

$$\begin{aligned}&s = (k+P)^{2},\quad Q^2 = -q^{2},\quad y = \frac{q\cdot P}{k\cdot P},\nonumber \\&W^2 = (q+P)^{2},\quad x = \frac{Q^{2}}{2\,q\cdot P}, \end{aligned}$$
(1)

where the four-vectors are indicated in Fig. 1. Here s is the square of the total centre-of-mass energy of the collision, \(Q^{2}\) the photon virtuality, y the scattered electron inelasticity, \(W^2\) the centre-of-mass energy squared of the \(\gamma ^{*}p\) system and x the Bjorken scaling variable.

Given the two hadronic systems X and Y, separated by a large rapidity gap, diffractive kinematic variables are defined as follows:

$$\begin{aligned}&M_{X}^2 = (P_{X})^2,\quad M_{Y}^2 = (P_{Y})^2,\quad t = (P-P_{Y})^{2},\nonumber \\&x_{{I\!\!P}} = \frac{q\cdot (P-P_{Y})}{q\cdot P}. \end{aligned}$$
(2)

In inclusive DDIS, where \(M_{X}\) and \(M_{Y}\) are the invariant masses of the systems X and Y, respectively, t is the squared four-momentum transfer at the proton vertex and \(x_{{I\!\!P}}\) the fraction of the proton’s longitudinal momentum transferred to the system X. In open charm production, the \(z_{{I\!\!P}}\) variable is defined as

$$\begin{aligned} z_{{I\!\!P}} = \frac{\hat{s} + Q^{2}}{M_{X}^{2} + Q^{2}}. \end{aligned}$$
(3)

It represents, in leading order, the pomeron’s momentum fraction participating in the \(\gamma ^{*}g \rightarrow c\bar{c}\) hard process. The variable \(\hat{s}\) denotes the centre-of-mass energy squared of the hard process, corresponding to the centre-of mass energy of the \(c\bar{c}\) quark pair in Fig. 1.

3 Monte Carlo models and fixed order QCD calculations

The diffractive and non-diffractive processes are modelled with the RAPGAP Monte Carlo event generator [25]. The generated Monte Carlo events are subjected to a detailed H1 detector response simulation based on GEANT-3 [26]. These simulated samples are passed through the same analysis chain as used for data and are used to correct the data for detector effects.

Diffractive events are simulated with leading (pomeron) and sub-leading (reggeon, \({I\!\!R}\)) exchanges based on the H1 2006 DPDF Fit B [12] diffractive parton density parameterisation obtained from a previous QCD analysis of inclusive diffractive data, convoluted with leading order matrix elements for open charm production via photon-gluon fusion. The contribution of non-diffractive processes to open charm production is simulated using RAPGAP in non-diffractive mode with the CTEQ6L PDF set [27]. Higher order QCD effects are modelled through parton showers in the leading-log approximation. QED radiation effects are simulated with the HERACLES program [28] interfaced to RAPGAP. To assess the effect of QED radiation a diffractive sample without QED radiation was also generated. Fragmentation is performed using the Lund string model [29] where all decay channels of the charm quark are included. The longitudinal part of the fragmentation function is reweighted according to the Kartvelishvili parameterisation \(D(z) \sim z^{\alpha }(1-z)\) [30] with an appropriate choice of \(\alpha \) [31]. The simulated events contain both elastic (\(ep \rightarrow eXp\)) and proton dissociative (\(ep \rightarrow eXY\)) processes. The two fractions are normalised relative to each other [32],

$$\begin{aligned} {\sigma (M_{Y} < 1.6~\mathrm{GeV})} / {\sigma (M_{Y} = m_{p})} = 1.20 \pm 0.11. \end{aligned}$$
(4)

Here \(\sigma (M_{Y} = m_{p})\) denotes the contribution in which the system Y contains only a leading proton, whereas \(\sigma (M_{Y} < 1.6~\mathrm{GeV})\) also includes contributions from proton dissociation processes integrated up to mass \(M_{Y} = 1.6~\mathrm{GeV}\). The simulated physics events are reweighted in the generated kinematics in order to reach good agreement with data at the reconstructed level as will be shown in Sect. 4.2.

Predictions for \(D^{*}\) cross sections in next-to-leading-order (NLO) QCD precision are obtained from the HVQDIS [33, 34] program adapted for diffraction. The calculation relies on collinear factorisation using H1 2006 DPDF Fit B NLO parton density functions involving gluons and light quarks in the quark singlet (fixed-flavour-number scheme). Massive charm quarks are produced via \(\gamma ^{*}\)-gluon fusion with the QCD scale parameter set to \(\Lambda _{5} = 0.228~\mathrm{GeV}\), which corresponds to a 2-loop \(\alpha _{s}(M_{Z})=0.118\), as was used in the DPDF extraction. The charm quarks are fragmented independently into \(D^{*}\) mesons with \(f(c\rightarrow D^{*}) = 0.235 \pm 0.007\) [35] in the \(\gamma ^{*}p\) rest frame using the Kartvelishvili parameterisation with parameters suited for use with HVQDIS [31]. The factorisation and renormalisation scales are set to \(\mu _{r} = \mu _{f} = \sqrt{Q^{2} + 4m_{c}^{2}}\) with the value \(m_{c} = 1.5\mathrm{~GeV}\) for the charm pole mass. The uncertainties arising from the choice of scales are estimated by simultaneously varying them by factors of 0.5 and 2. The uncertainty introduced in the calculation caused by the uncertainty of \(m_c\) is evaluated by varying \(m_{c}\) to 1.3 and 1.7 GeV. The Kartvelishvili parameters are varied within their uncertainties [31]. The DPDF uncertainties are estimated by propagating the eigenvector decomposition of the errors of the DPDF parameterisation. The individual sources of uncertainties are added in quadrature separated for up and down variations of the cross sections. The contribution of B-hadrons due to beauty fragmentation to the diffractive \(D^{*}\) cross section is neglected; it is expected to be less than 3% for non-diffractive DIS (see [36]) and even smaller for the diffractive production.

The HVQDIS calculation is performed also in the non-diffractive mode using the CT10F3 proton PDF set [37]. It is used for comparisons of predictions with measurements of the diffractive to inclusive cross section ratio (Sect. 5.2). The calculation is done following the one used for comparison with inclusive \(D^{*}\) data [38]. The uncertainties from the choice of scales and \(m_{c}\) as well as the fragmentation uncertainty are evaluated in the same manner as for the diffractive calculation. The uncertainty of the CT10F3 PDF set is not considered for this analysis but is expected to be small in comparison to the DPDF uncertainties.

4 Experimental technique

4.1 The H1 detector

A detailed description of the H1 detector can be found elsewhere [39]. Here, a brief account of the detector components most relevant to the present analysis is given. The H1 coordinate system is defined such that the origin is at the nominal ep interaction point and the polar angle \(\theta = 0\) and the positive z axis correspond to the direction of the outgoing proton beam. The region \(\theta < 90^\circ \), which has positive pseudorapidity \(\eta = - \ln \tan \theta /2\), is referred to as the ‘forward’ hemisphere.

The ep interaction point in H1 is surrounded by the central tracking system, which includes silicon strip detectors [40] as well as two large concentric drift chambers. These chambers cover a region in polar angle \(20^{\circ }< \theta < 160^{\circ }\) and provide a resolution of \(\sigma (P_{T})/P_{T} = 0.006P_{T}/\mathrm{~GeV} \oplus 0.02\). They also provide triggering information [41, 42]. The forward tracking detector, a set of drift chambers with sense wires oriented perpendicular to the z axis, extends the acceptance of the tracking system down to \(7^{\circ }\) in polar angle. The central tracking detectors are surrounded by a finely segmented liquid argon (LAr) sampling calorimeter covering \(-1.5< \eta < 3.4\). Its resolution is \(\sigma (E)/E = 0.11/\sqrt{E/\mathrm{~GeV}} \oplus 0.01\) in its electromagnetic part and \(\sigma (E)/E = 0.50/\sqrt{E/\mathrm{~GeV}} \oplus 0.02\) for hadrons, as measured in test beams [43, 44].

The central tracker and LAr calorimeter are placed inside a large superconducting solenoid, which provides a uniform magnetic field of 1.16 T. The backward region \(-4< \eta < -1.4\) is covered by a lead-scintillating fiber calorimeter (SpaCal [45]) with electromagnetic and hadronic sections. In the present analysis the energy and angle of the scattered electron is measured in the electromagnetic section of the SpaCal. It has an energy resolution for electrons \(\sigma (E)/E = 0.07 / \sqrt{E/\mathrm{~GeV}} \oplus 0.01 \) as measured in test beams [46].

Information from the central tracker and the LAr and SpaCal calorimeters is combined in an energy flow reconstruction algorithm which yields a list of hadronic final state objects [47, 48]. For these objects a calibration is applied ensuring the relative agreement of hadronic energy scale between the data and simulations at 1% accuracy [49].

In the forward region the H1 detector is equipped with drift chambers comprising the forward muon detector (FMD, \(1.9< \eta < 3.7\)). The forward tagger system (FTS) is a set of scintillators surrounding the beam pipe at several locations along the proton beamline, downstream of the H1 main detector. The FTS station at \(28 \mathrm{~m}\) covering the range \(6.0< \eta < 7.5\) is used in this analysis. Both FMD and FTS are sensitive to the very forward energy flow and improve the selection of large rapidity gap events.

The luminosity is measured via the Bethe–Heitler bremsstrahlung process \(ep \rightarrow ep\gamma \), with the final state photon detected by a photon detector located close to the beam pipe at position \(z=-103~\mathrm{m}\). The precision of the integrated luminosity determination is improved in a dedicated analysis of the QED Compton process [50].

4.2 Event selection

The analysis is based on data collected by H1 in the 2005–2006 electron and the 2006–2007 positron running periods with \(\sqrt{s}=319~\mathrm{GeV}\), where the proton and lepton beam energies are 920 and 27.6 GeV, respectively. The corresponding integrated luminosity is 287 pb\(^{-1}\). The events are triggered on the basis of a scattered electron signal in the SpaCal calorimeter together with at least one track above the 900 MeV transverse momentum threshold in the drift chambers of the central tracker.

4.2.1 Diffractive DIS selection and reconstruction of kinematics

The momentum transfer \(Q^{2}\) and inelasticity y are reconstructed using the electron-\(\Sigma \) method [51] which combines information on the scattered electron candidate and hadronic final state (HFS) kinematics. This choice optimises the resolution for these observables. The measurement phase space in \(Q^{2}\) and y is chosen to be

$$\begin{aligned} 5< Q^{2}< 100~\mathrm{GeV}^{2},\quad 0.02< y < 0.65. \end{aligned}$$
(5)

The selection of diffractive events is based on the presence of a forward large rapidity gap (LRG), which is primarily provided by a cut on the pseudorapidity of the most forward cluster in the LAr calorimeter above the 800 MeV energy threshold, \(\eta _\mathrm{{ max}} < 3.2\).

The variable \(x_{{I\!\!P}}\) is reconstructed as

$$\begin{aligned} x_{{I\!\!P}} = \frac{M_{X}^{2} + Q^{2}}{W^{2} + Q^{2}}, \end{aligned}$$
(6)

where W is calculated as \(W=\sqrt{ys - Q^{2}}\). The invariant mass of the hadronic final state, \(M_{X}\), is determined as follows

$$\begin{aligned} M_{X} = f_\mathrm{{corr}}(\eta _\mathrm{{ max}}) \sqrt{(P_{X})^{2}}, \end{aligned}$$
(7)

where \(P_{X}\) is the reconstructed four-momentum of the hadronic final state and \(f_\mathrm{{ corr}}\) is an \(\eta _\mathrm{{ max}}\) dependent factor introduced in order to correct for detector losses at large \(\eta \). It is determined from simulations yielding 16% enhancement factor on average. The range of the reconstructed \(x_{{I\!\!P}}\) values is limited to \(x_{{I\!\!P}} < 0.03\).

Fig. 2
figure 2

The \(\Delta m\) distributions in the data for the right charge sample with the combined signal and background fit indicated by the solid and dotted line, respectively, is shown in the left figure. The wrong charge sample with the background-only fit, performed simultaneously under the assumption of identical background shape in the right charge combinations, is shown in the right figure as the dotted line

The variable \(z_{{I\!\!P}}\) is reconstructed using in addition the \(D^{*}\) candidate four-momentum. This variable is denoted \(z_{{I\!\!P}}^{obs}\) and is defined as

$$\begin{aligned}&z_{{I\!\!P}}^{obs} = \frac{(M_{c\bar{c}}^{2})^{obs} + Q^{2} }{ {M_{X}^{2} + Q^{2}} },\quad \text {with}\; (M_{c\bar{c}}^{2})^{obs} = \frac{1.2\,p_{\perp D^{*}}^{* 2} + m_{c}^{2}}{z(1-z)}\quad \nonumber \\&\quad \text {and } z = \frac{(E-p_{z})_{D^{*}}^{\text {(lab)}}}{2yE_{e}} \end{aligned}$$
(8)

where \((M_{c\bar{c}}^{2})^{obs}\) is an estimate of \(\hat{s}\) in Eq. 3. \((M_{c\bar{c}}^{2})^{obs}\) is reconstructed from the \(D^{*}\) kinematics. This is done in close analogy to the \(x_g^{\text {obs}}\) measurement in inclusive \(D^{*}\) production [52]. The term \(1.2\,p_{\perp D^{*}}^{* 2}\) is an approximation to the value of the transverse momentum squared of the charm quark in the \(\gamma ^{*}p\) rest frame. The observable z denotes the inelasticity of the \(D^{*}\) meson which is calculated in the laboratory frame using the difference of the energy and the longitudinal momentum, \((E-p_{z})_{D^{*}}\), of the D* meson. The factor 1.2 is introduced to ensure \(z_{{I\!\!P}}^{\text {true}} \approx z_{{I\!\!P}}^{\text {obs}}\) on average, as deduced in studies of generated events using RAPGAP.

The activity in the FTS and the FMD is required not to exceed the noise levels monitored with an event sample triggered independently of detector activity. Noise effects are also propagated into the simulation of the detector response in a similar manner. The diffractive event selection requirements ensure that the analysed sample is dominated by \(ep \rightarrow eXp\) processes at small |t| with an intact proton in the final state, often called proton elastic processes. A small admixture of events is present with leading neutrons or low \(M_{Y}\) baryon excitations, referred to as proton dissociation contributions (PD). The values of \(M_Y\) and t are not reconstructed explicitly. However, as the diffractive selection rejects events at large \(M_{Y}\) or large \(\left| t \right| \), the measurement is corrected to the \(M_{Y} < 1.6 \mathrm{~GeV}\) and \(\left| t \right| < 1 \mathrm{~GeV}^{2}\) range.

4.2.2 \(D^{*}\) selection

The detection of \(D^{*}\) mesons is based on the full reconstruction of its decay products in the ‘golden channel’

$$\begin{aligned} D^{*+} \rightarrow D^{0} \pi ^{+}_{slow} \rightarrow (K^{-} \pi ^{+}) \pi ^{+}_{slow} \quad (+C.C.), \end{aligned}$$
(9)

with a branching ratio of \((2.66 \pm 0.03)\%\) [53]. Tracks reconstructed in the central tracker are used to identify the decay products. The kaon and pion candidate tracks from \(D^{0}\) decays are required to satisfy \(p_{t} > 0.3 \mathrm{~GeV}\) while the slow pion candidate track is required to have \(p_{t} > 0.12~\mathrm{GeV}\), where \(p_{t}\) is the transverse momentum of the reconstructed track in the laboratory frame. In order to suppress combinatorial background as well as contributions of other decay channels with at least three charged decay products, called reflections, the invariant mass of the \(K^{\mp } \pi ^{\pm }\) pair is required to be in agreement with the nominal \(D^{0}\) mass (\(1864.83 \pm 0.05 \mathrm{~MeV}\) [53]) within a window of \(\pm 80\) MeV. The kinematics of the \(D^{*}\) meson candidate reconstructed from the \(K^{\mp } \pi ^{\pm } \pi ^{\pm }_{slow}\) system is restricted to the range \(p_{t,D^{*}} > 1.5 \mathrm{~GeV}\) and \(| \eta _{D^{*}} | < 1.5\).

Fig. 3
figure 3

The differential \(N(D^{*})\) distributions obtained from \(\Delta m\) fits to the data and simulation as a function of \(Q^2\), y, \(\text {log}_{10}(x_{{I\!\!P}})\), \(z_{{I\!\!P}}^\mathrm{obs}\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\). The data are represented with dots. The contributions of individual processes in the simulation, reweighted RAPGAP, are indicated by filled histograms as follows; elastic proton pomeron exchange (light orange), proton dissociation (dark orange) and reggeon exchange (dark red)

The mass difference \(\Delta m = m(K^{\mp } \pi ^{\pm } \pi ^{\pm }_{slow})-m(K^{\mp } \pi ^{\pm })\) is used to determine the \(D^{*}\) signal. It is expected to peak near \(\Delta m=0.145~\mathrm{GeV}\) [53]. The wrong charge combinations \(K^{\pm } \pi ^{\pm } \pi ^{\mp }_{slow}\) selected with otherwise unchanged criteria do not contribute to the signal, they are, however, utilised to constrain the shape of the combinatorial background. The right and wrong charge \(\Delta m\) distributions are fitted simultaneously by means of an unbinned extended likelihood fit using RooFit [54] and ROOT [55]. The Crystal Ball [56] and Granet [57] probability distribution functions are used for modelling the signal and background, respectively. The fit to the total sample of selected \(D^{*}\) candidates is shown for the right and wrong charge combinations in Fig. 2. The fit to the total number of \(D^{*}\) mesons in the data yields \(N(D^{*}) = 1169 \pm 58\). The observed width is dominated by experimental effects.

The fits are repeated in bins of reconstructed kinematic quantities. Figure 3 shows the \(D^{*}\) yields determined as a function of the variables \(Q^2\), y, \(\text {log}_{10}(x_{{I\!\!P}})\), \(z_{{I\!\!P}}^{obs}\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\). The \(N(D^{*})\) distributions are well described by the reweighted simulation. The fraction of proton dissociation processes adjusted globally (as given by Eq. 4) is largely independent of the kinematics. The reggeon contribution is generally small and reaches 5% at large \(x_{I\!\!P}\). The non-diffractive background contribution is below 1% level and is not shown.

4.3 Cross section measurement

The number of fitted \(D^{*}\) mesons is corrected for trigger inefficiency, detector effects due to limited acceptance and resolution, the branching ratio of the golden channel, and the contribution of reflections and higher order QED processes at the lepton vertex. The bin averaged \(D^{*}\) cross section in bin i of a generic variable x in the phase space defined in Table 1 is measured as

$$\begin{aligned} \left( \frac{ \mathrm{d} \sigma }{ { \mathrm d} x } \right) _{i} = \frac{N_{i}^\mathrm{data} - N_{i}^\mathrm{sim, bgr} }{\mathcal {L}_\mathrm{int} \,\, \Delta ^{x}_{i} \,\, B_{r} \,\, \varepsilon _\mathrm{trigg} \,\, A_{i}} \, C^\mathrm{{ QED}}_\mathrm{{corr},i}, \end{aligned}$$
(10)

where

  • \(N_{i}^\mathrm{data}\) is the number of \(D^{*}\) mesons from the fit passing the experimental cuts in the data.

  • \(N_{i}^\mathrm{sim, bgr}\) is the number of \(D^{*}\) mesons from the fit to simulated events passing the experimental cuts while being generated outside the phase space (Table 1) of the measurement.

  • \(A_{i}\) is the acceptance correction factor accounting for effects related to the transition from the hadron level to the detector level determined from MC simulations.

  • \(\mathcal {L}_\mathrm{int}\) is the integrated luminosity of the data.

  • \(B_{r}\) is the branching ratio of the golden decay channel.

  • \(\varepsilon _\mathrm{trigg}\) is the trigger efficiency.

  • \(C^\mathrm{{QED}}_\mathrm{{corr},i}\) are corrections for QED radiation defined as \(\sigma ^\mathrm{{QED-off}}/\sigma ^\mathrm{{QED}{\text {-}}\mathrm{on}}\) as obtained from Monte Carlo generated events, where \(\sigma ^\mathrm{\mathrm{QED}{\text {-}}\mathrm{off}}\) (\(\sigma ^\mathrm{{QED}{\text {-}}\mathrm{on}}\)) is the bin-integrated cross section predicted by RAPGAP with QED radiation turned off (turned on) as described in Sect. 3.

  • \(\Delta _{i}^{x}\) is the bin width of the i-th bin of x.

The acceptance corrections, \(A_{i}\), are defined as

$$\begin{aligned} A_{i} = \frac{ N_{i}^\mathrm{sim} - N_{i}^\mathrm{sim, bgr} }{n_{i}^\mathrm{sim}}, \end{aligned}$$
(11)

where, for a given bin i, \(N_{i}^\mathrm{sim}\) is the fitted number of \(D^{*}\) mesons passing the experimental cuts in the simulation of MC generated events encompassing all charm quark decay channels as well as all \(D^{*}\) decay channels, \(N_{i}^\mathrm{sim, bgr}\) is defined above and \(n_{i}^\mathrm{sim}\) is the MC generated number of \(D^{*}\) mesons decaying solely via the golden channel with the event kinematics inside the phase space defined in Table 1.

Table 1 Definition of the phase space of the cross section measurement
Table 2 Bin averaged hadron level \(D^{*}\) production cross sections in diffractive DIS as a function of y, \(Q^{2}, \mathrm{log}_{10}(x_{{I\!\!P}})\), \(z_{{I\!\!P}}^{obs}\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\) together with statistical (\(\delta _\mathrm{stat}\)), systematic (\(\delta _\mathrm{syst}\)) and total (\(\delta _\mathrm{tot}\)) uncertainties. The total uncertainties are obtained as the statistical and systematic uncertainties added in quadrature

4.4 Systematic uncertainties

Experimental and model uncertainties are propagated to the differential and the integrated \(D^{*}\) cross section measurements. In the following only the effects on the integrated \(D^{*}\) cross sections are quantified.Footnote 2

Fig. 4
figure 4

Bin averaged single-differential \(D^{*}\) cross sections as a function of \(Q^2\) and y. Data are shown as dots, where the inner error bars indicate statistical uncertainties and the outer error bars represent the statistical and the full set of systematic uncertainties added in quadrature. The central NLO QCD prediction by HVQDIS is shown as a white line inside the coloured bands. The inner band represents the DPDF and fragmentation uncertainties added in quadrature. The outer band represents DPDF, fragmentation, charm mass, factorisation and renormalisation scale uncertainties added in quadrature

  • The energy scale (polar angle) of the scattered lepton is known to the 1% (1 mrad) level resulting in a 0.5% (1.5%) uncertainty.

  • The relative energy scale of the hadronic final state is known with a precision of 1% resulting in a 0.06% uncertainty.

  • Changing the function \(f_\mathrm{{corr}}\) (Eq. 7) to the constant 1.16 results in 2.7% uncertainty.

  • There is a certain ambiguity in describing the tails of the \(\Delta m\) signal distribution. Choosing a modified Crystal Ball function with an extra Gaussian component for the fits to the \(D^{*}\) signal has 3.8% effect.

  • The normalisation of the proton dissociative contribution (Eq. 4) introduces an uncertainty of 7.1%.

  • In a dedicated study [58], using forward proton tagging data, a 10% uncertainty on the large rapidity gap selection inefficiency is determined, which translates to a 2.4% uncertainty.

  • The shapes of the generated spectra of \(Q^{2}\), y, \(x_{{I\!\!P}}\), \(p_{t,D^{*}}\) and t are varied independently with the help of multiplicative weights of \(\mathrm{e}^{{}^{+0.007 }_{-0.013}\,Q^{2}/\mathrm{GeV}^{2}}\), \(y^{^{+0.9}_{-1.1}}\), \((x_{{I\!\!P}})^{^{+0.13}_{-0.16}}\), \(\mathrm{e}^{^{+0.06}_{-0.15}{\,}p_{t,D^{*}}/\mathrm{GeV}}\) and \(\mathrm{e}^{^{+0.8}_{-0.9}{\,}t/\mathrm{GeV}^{2}}\) resulting in variations of the fitted differential distributions compatible with the data control distributions (Fig. 3). The reweighting is an approach to assess the uncertainties on the data correction procedure stemming from the Monte Carlo model. The resulting uncertainties are 0.5, 0.9, 0.4, 3.7 and 1.1%, respectively.

Fig. 5
figure 5

Bin averaged single-differential \(D^{*}\) cross sections as a function of the diffractive variables \(\text {log}_{10}(x_{{I\!\!P}})\) and \(z_{{I\!\!P}}^{obs}\). Further details are indicated in the caption of Fig. 4

Fig. 6
figure 6

Bin averaged single-differential \(D^{*}\) cross sections as a function of \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\). Further details are indicated in the caption of Fig. 4

The following uncertainties affect only the normalisation of the measurement.

  • The integrated luminosity is known to 2.7% and the golden channel branching ratio to 1.1%.

  • The uncertainty on the trigger efficiency (98% on average) is covered by a 2% variation.

  • The impact of the restriction to the \(D^{0}\) mass window in terms of \(N(D^{*})\) yield loss caused by the choice of the 80 MeV value is evaluated. A systematic uncertainty of 2% covers the observed difference between data and simulation.

  • The reflections contribute about 3% to the fitted \(N(D^{*})\). The branching fractions of \(D^{*}\) decaying to reflections are not precisely reproduced in the simulation. The integrated cross section increases by about 1.2% if recent branching ratios of reflections are used [53].

  • The track reconstruction efficiency is known with 1% uncertainty resulting in 3% per \(D^{*}\).

  • The contribution of non-diffractive processes is suppressed by the diffractive selection to a level of less than 1%. A conservative uncertainty of 1% is assigned.

The contributions of the individual systematic uncertainties are added in quadrature both for the integrated and differential cross section measurements.

5 Results

In the first part of this section the measured integrated and differential cross sections for \(D^{*}\) production in diffractive deep inelastic scattering are presented. Theoretical predictions based on next-to-leading order QCD calculations are compared with the data. In the second part ratios of diffractive to non-diffractive \(D^{*}\) production cross sections are extracted and confronted with theoretical predictions as well as with previous results from HERA.

5.1 Diffractive \(D^{*}\) production cross sections

The integrated cross section of \(D^{*}\) production for the phase space region given in Table 1 is measured to be

$$\begin{aligned} \sigma _{ep \rightarrow \,e Y \!X(D^{*}) } = 314 \pm 23\,\text {(stat.)} \pm 35\,\text {(syst.)}\,\text {pb}. \end{aligned}$$
(12)

This can be compared with the theoretical value calculated in next-to-leading order QCD with the HVQDIS code [33, 34] adapted to diffraction using H1 2006 DPDF Fit B and Kartvelishvili fragmentation as described in Sect. 3.

$$\begin{aligned} \sigma _{ep \rightarrow \,e Y \!X(D^{*}) }^{\text {theory}} = 265 \text { }^{+54}_{-40} \text { (scale) } \text { }^{+68}_{-54} (m_{c}) \text {}^{+7.0}_{-8.2} \text { (frag.) } \text { }^{+31}_{-35} \text { (DPDF) } \text { pb.}\nonumber \\ \end{aligned}$$
(13)

Within its large uncertainties the prediction is compatible with the measured value, which supports collinear factorisation. However, the prediction depends substantially on the choice of the factorisation and renormalisation scale as well as on the value of the charm mass. Similar conclusions were reached in a previous H1 publication [23] albeit within larger uncertainties and in a slightly different kinematic domain.

Table 3 The values of diffractive fraction for \(D^{*}\) production cross sections together with statistical (\(\delta _\mathrm{stat}\)), systematic (\(\delta _\mathrm{syst}\)) and total uncertainties (\(\delta _\mathrm{tot}\)) as a function of y, \(Q^{2}\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\). The total uncertainties are obtained as the statistical and systematic uncertainties added in quadrature

The measured bin averaged single-differential cross sections as a function of y, \(Q^{2}\), \(\text {log}_{10}(x_{{I\!\!P}})\), \(z_{{I\!\!P}}^{obs}\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\) are given in Table 2 and are shown in Figs. 4, 5 and 6 together with the NLO predictions. In order to compare the shapes between data and theory the ratios of data to NLO calculations are also shown.

Figure 4 shows that the shape of the measured \(\mathrm{d}\sigma /\mathrm{d}y\) is well described by the theory. The measured \(\mathrm{d}\sigma /\mathrm{d}Q^{2}\) might indicate a slightly harder dependence in the data, however, within the large uncertainties the shape is in agreement with the theory. The shape of the \(\mathrm{d}\sigma /\mathrm{dlog}_{10}(x_{{I\!\!P}})\) shown in Fig. 5 is satisfactorily described by the prediction given the large relative uncertainties at low \(x_{{I\!\!P}}\) values. The shape of \(\mathrm{d}\sigma /\mathrm{d}z_{{I\!\!P}}^{obs}\) shown in Fig. 5 is not described as well by the prediction, however the experimental uncertainties at low \(z_{{I\!\!P}}^{obs}\) are sizeable. The shapes of \(\mathrm{d}\sigma /\mathrm{d} p_{t,D^{*}}\) and \(\mathrm{d}\sigma /\mathrm{d}\eta _{D^{*}}\) are well described by the theory (see Fig. 6). For \(\eta _{D^{*}} > 1\), however, the theory predicts a value which underestimates the data by about \(50\%\) with a large uncertainty. There is an indication of a similar effect in the corresponding non-diffractive \(D^{*}\) cross section measurement [38].

The differential comparison profits from the substantial increase of statistics in the present analysis as compared with previous measurements at HERA.

Fig. 7
figure 7

The diffractive fraction, \(R_{D}\), measured as a ratio of bin averaged diffractive to non-diffractive \(D^{*}\) production single differential cross sections in deep inelastic scattering as a function of y, \(Q^2\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\). The data ratios are represented with dots, where the inner error bars indicate statistical uncertainties and the outer error bars represent the statistical and systematic uncertainties added in quadrature. The central NLO QCD prediction of \(R_D\) by HVQDIS is shown as a white line inside the coloured bands. The inner band represents DPDF uncertainty. The outer band represents effect of the DPDF uncertainty and simultaneous variations of scales, charm mass and fragmentation settings in the diffractive and non-diffractive calculations added in quadrature

5.2 Diffractive fractions

The \(D^{*}\) and DIS selection criteria given in Table 1 are close to those used in the corresponding non-diffractive analysis [38]. The non-diffractive \(D^{*}\) differential cross sections thus can be used to calculate the diffractive fraction, \(R_{D} = {\sigma ^{\mathrm{diff}}_{{D^{*}}}} / {{\sigma ^{\text {non-diff}}_{{D^{*}}}}}\), in the phase space defined in Table 1.

The non-diffractive cross sections [38], originally given for \(0.02< y < 0.7\), are interpolated to \(0.02< y < 0.65\) using small corrections calculated with HVQDIS. The correction factors reduce the non-diffractive cross sections by about 1.5–3.5% differentially in \(Q^2\), \(p_{t,D^{*}}\) and \(\eta _{D^{*}}\) with typical uncertainties of \(0.2\%\). The uncertainties of both the diffractive and non-diffractive cross sections are accounted for in the \(R_{D}\) measurement. Integrated over the whole phase space the results are \(R_{D} = 6.6 \pm 0.5 \mathrm{(stat)} \text { }^{+0.9}_{-0.8} \mathrm{(syst)}\%\) for the data and \(R_{D}^\mathrm{\, theory} = 6.0 {}^{+1.0}_{-0.7} \mathrm{(scale)} {}^{+0.5}_{-0.4} \mathrm{(}m_c\mathrm{)} {}^{+0.7}_{-0.8} \mathrm{(DPDF)} {}^{+0.02}_{-0.04} \mathrm{(frag)}\%\) for the theoretical prediction. The uncertainties of the theoretical predictions are obtained from simultaneous variations of \(m_c\), fragmentation parameters and the factorisation and renormalisation scales. The DPDF uncertainty is also propagated to the prediction.

The differential fractions \(R_{D}(y)\), \(R_{D}(Q^2)\), \(R_{D}(p_{t,D^{*}})\) and \(R_{D}(\eta _{D^{*}})\) are listed in Table 3 and are shown in Fig. 7. Within uncertainties the data do not provide strong evidence for kinematic dependencies of \(R_D\) on \(Q^2\) or y, while at the same time they are also consistent with the kinematic dependencies predicted by theory. The diffractive fraction decreases from \(8\%\) to \(3\%\) with \(p_{t,D^{*}}\) increasing. The measured dependence of the diffractive fraction on \(\eta _{D^{*}}\) decreases from \(10\%\) to about \(5\%\) for the highest \(\eta _{D^{*}}\) values. These shapes are well reproduced by the NLO QCD predictions within the uncertainties. The shapes can be qualitatively understood as follows. Due to the high energy of the leading proton in diffraction (\(x_{{I\!\!P}} < 0.03\)) the system X is produced with low masses \(M_X\). Less energy is available from the proton side to produce the hard system containing the \(D^{*}\) meson as compared to the non-diffractive case. Similarly, the fraction is suppressed for small y, i.e. for small energy of the exchanged virtual photon. The \(Q^2\) dependence of \(R_D\) can be explained by the fact that at high \(Q^2\) (higher x) the diffractive cross section is suppressed due to a limited \(x_{I\!\!P}\) range. Likewise, due to the lack of energy, the events with higher \(p_{t,D^{*}}\) can be expected to be suppressed in diffraction. The diffractive fraction as a function of \(\eta _{D^{*}}\) indicates that in diffraction the hard system tends to be produced backwards, due to the kinematics constrained by the presence of a large rapidity gap, or equivalently the \(x_{{I\!\!P}} < 0.03\) condition.

Fig. 8
figure 8

Integrated diffractive fractions measured in \(D^{*}\) production in the deep-inelastic and the photoproduction (\(Q^{2} < 1 \mathrm{~GeV}^2\)) regime as measured at HERA previously [20,21,22, 24] and in the present analysis. The inner error bars represent statistical uncertainties, the outer ones the statistical and systematic uncertainties added in quadrature. The dashed line and the shaded band indicate the central value and the total experimental uncertainty of \(R_{D}\) of the measurement presented here, respectively

In Fig. 8 the diffractive fraction, integrated over the full phase space, is compared with previous measurements performed at HERA both in the DIS regime [20,21,22] and in photoproduction [24]. The average value of \(R_{D}\) measured in this article is in agreement with the previous results and within the sizeable experimental uncertainties is observed to be largely independent of the varying phase space constraints in \(x_{{I\!\!P}}\), \(Q^2\) and \(p_{t,D^{*}}\). In particular, the ratios observed in DIS and in photoproduction are compatible with each other.

6 Conclusions

Integrated and differential cross sections of \(D^{*}(2010)\) production in diffractive deep inelastic scattering are measured. The analysis is based on a data sample taken by the H1 experiment at the HERA collider corresponding to an integrated luminosity of 287 pb\(^{-1}\). The measured cross sections are compared with theoretical predictions in next to leading order QCD. Compared to the previous measurement in a similar kinematic domain the precision is improved by a factor of two. The new measurements are well described by the predictions within the large theoretical uncertainties which are dominated by variations of scales and the charm quark mass. This supports the validity of collinear factorisation in diffraction.

Measurements of diffractive fractions of \(D^{*}\) production cross section in deep inelastic scattering are also presented, using non-diffractive cross sections published earlier by H1. The fractions are in agreement with theoretical predictions in next-to-leading order QCD. Although the value of the diffractive fraction is found to decrease at high \(p_{t,D^{*}}\) and at high \(\eta _{D^{*}}\) due to limitations of the diffractive phase space, it is observed to be largely independent of other details of the phase space definition. This is confirmed by comparisons to previous measurements of the diffractive fraction.