Predicting event weights at next-to-leading order QCD for jet events defined by $2\rightarrow 1$ jet algorithms

As the differential cross section depicts the probability to observe an event in a specific phase space region, its calculation enables the prediction of event weights for events observed in the experiment. At leading-order, calculating event weights is straightforward thus allowing for applications such as the generation of unweighted events or the Matrix Element Method. Because of the non-trivial interplay between real and virtual corrections, predicting event weights including next-to-leading order corrections is more involved in general and dependent on the jet algorithm in particular. In this article, we present a method to calculate event weights including NLO QCD corrections for jet events defined using the (anti-)$k_t$-algorithm together with the conventional $2\rightarrow 1$ recombination. This is a major extension compared to approaches used for example in existing MC tools which only allow to produce jet events together with their associated weight but are unable to calculate the corresponding weight at NLO accuracy for a given event. We validate the approach by generating unweighted events for single top quarks produced in hadronic collisions. As a further sample application we apply the Matrix Element Method to single top-quark events generated with POWHEG in combination with Pythia.


Introduction
The steadily improving precision achieved in collider experiments like ATLAS and CMS requires an equal precision in the theoretical predictions to make optimal use of the experimental results. In recent years tremendous progress has been made concerning the calculation of nextto-and next-to-next-to-leading order QCD corrections (see e.g. refs. [1][2][3][4] for overviews of the theoretical progress). Meanwhile, next-to-leading order (NLO) corrections are considered a solved problem and are calculable for many processes using publicly available tools (see e.g. refs. [5][6][7][8][9][10]). As far as next-to-next-to-leading (NNLO) corrections are concerned not the same level of maturity has been achieved yet. However, many 2-to-2 processes have been calculated recently (see e.g. refs. [11,12], further examples can be found in refs. [4,13]). The forefront of current research is the application to multi-scale problems as they occur for example in 2-to-3 reactions or 2-to-2 processes involving particles with many different masses. For some quantities like for example the Higgs cross section even higher-order corrections have been calculated recently [14][15][16][17].
In collider experiments no smoking gun as a clear sign for physics beyond the Standard Model has been observed so far. This has led to an increasing interest in sophisticated analysis methods to compare theoretical and experimental results. Multivariate methods allow to utilise most of the information contained in the recorded events. Thus, they present promising tools to search for even smallest hints for New Physics. Among these methods the so called Matrix Element Method sticks out since it provides a very general and at the same time optimal approach to compare theory and experiment. Based on the principle of Maximum Likelihood, it allows for an unambiguous interpretation of the findings. Briefly worded, the cumulative Likelihood for a sample of recorded events is calculated by interpreting the differential cross section evaluated for each event as the theoretically predicted probability of having measured this particular event. However, until recently the use of the Matrix Element Method has been limited by the fact that only leading-order cross sections could be used to calculate weights for given events. Considering the progress in the calculation of higher-order corrections mentioned above this has been a major drawback of the otherwise promising method.
Even though packages like aMC@NLO [8] and POWHEG [23] already allow the calculation of weights for events beyond NLO accuracy, the 'direction' of these calculations is unfortunately not as needed for the Matrix Element Method: In MC@NLO and POWHEG one starts with partonic momenta which are mapped to jet momenta according to a jet algorithm which is applied to the partonic final state after a parton shower. The corresponding weight is calculated at NLO accuracy including parton shower effects along the way. After un-weighting, the jet momenta are distributed according to the NLO+parton shower predictions-as needed for most applications. In contrast, the starting point for the Matrix Element Method is a particular set of recorded events (e.g. certain values of energy or angular variables attributed to resolved final-state objects-typically jets) with the corresponding weights to be calculated a posteriori: The differential cross section is required as a function of the variables in the event definition. In perturbation theory the experimentally resolved jets are modelled by mapping partonic momenta to jet momenta according to the same jet algorithm that is used by the experiment. The calculation of an event weight in perturbative QCD therefore requires the differential cross section to be differential in variables of the jet momenta modelling the recorded events. At leading order partonic and jet momenta are uniquely identified allowing for a straightforward evaluation of the differential cross section calculated for partonic final states. However, when including higher-order corrections the jet algorithm dictates non-trivial mappings from partonic momenta to jet momenta depending on their kinematic configurations. This means that the simple identification of partonic momenta and jet momenta cannot be justified beyond the leading order. The calculation of weights for measured events in terms of cross sections which are differential in jet variables is thus non-trivial when higher-order corrections are taken into account.
In refs. [18,19] different approaches to include NLO QCD corrections have been investigated. However, both articles focus on special aspects and make no attempt to present a general solution. Ref. [18] concentrates on the effect of initial state radiation while ref. [19] excludes strongly interacting particles in the final state. In ref. [20] a general algorithm has been proposed relying on a modification of the recombination procedure used in jet algorithms. As a proof-of-concept the method has been applied in ref. [21] to hadronic single top-quark production-a process where the Matrix Element Method has been used from early days on to disentangle the signal from an overwhelming background (see e.g. ref. [22]).
The evaluation of fully differential event weights incorporating NLO QCD corrections also allows the generation of unweighted events exactly following the NLO differential cross sections. In fact, this possibility has been used already in ref. [20] and ref. [21] to simulate a toy experiment and to validate the approach. These events obviously depend on the jet algorithm used to cluster additional radiation. To distinguish them from partonic events which are ill-defined beyond leading order we call them jet events in the following. As outlined above, these jet events simulate-within the fixed-order NLO approximation-the events observed in a real experiment. It is well known that fixed-order predictions can be further improved by resumming certain logarithmically enhanced corrections through parton showers. As mentioned before, packages like aMC@NLO and POWHEG also allow the generation of unweighted events including these effects. Although unweighted events following the fixed-order NLO predictions may be considered as only half way on the way to unweighted events including parton shower corrections they are interesting on their own right, since they allow a detailed study of parton shower effects. In addition, as has been pointed out in ref. [24], the possibility to produce unweighted jet events may also help to improve the numerical integration over the real corrections-a major bottleneck in the evaluation of NLO QCD cross sections.
At NLO QCD real and virtual corrections are combined to calculate the NLO contribution to differential cross sections. Because of additional real radiation the dimensionality of the phase space is different for real and virtual corrections. As a consequence the definition of a "fully differential" event weight requires to integrate out all unobserved radiation and add this contribution to the virtual corrections for the specific event. While this approach is straightforward in theory, in practice two complications arise: 1.) The variables used to describe the event must not allow to distinguish between real and virtual contributions, since this would prevent the one-toone correspondence of the two contributions. For example, in case of standard jet algorithms, relying on the summation of 4-momenta to define the momentum of a jet obtained through recombination, the full 4-momenta cannot be used to describe the event: the jets obtained from recombination will in general acquire a jet mass and the point-wise correspondence of real and virtual corrections is lost. 2.) For a given event all phase space regions/additional radiation which after recombination contribute to the event need to be identified and integrated in an efficient way. While the combinatorial part (i.e. the different possibilities to cluster the additional radiation to obtain the given event) is easy to solve, the efficient numerical integration is nontrivial. Because of these two complications the standard approach is to avoid the definition of a fully differential event weight and combine finite phase space regions in terms of histogrammed results as approximations to differential distributions. However, drawbacks of this procedure are potential numerical instabilities encountered at bin boundaries and loss of correlations between different variables due to the bin-wise integration.
The algorithm presented in refs. [20,21] to calculate weights for jet events at NLO accuracy modifies the recombination procedure used to cluster two primary objects into a resulting jet. It is thus possible to keep the kinematics Born-like leading to a straightforward identification of real and virtual contributions. Using in addition the factorisation of the phase space measure in terms of a phase space for the recombined jets and a part due to the additional radiation, the second complication is also solved. Although the modification which clusters on-shell objects into on-shell jets is theoretically well motivated, this recombination procedure is not yet used in the experimental analysis since it would require a major effort in re-calibration and re-tuning of existing Monte Carlo tools. To circumvent this problem it has been shown in ref. [25] for the example of single top-quark production how the modification of the jet algorithm can be avoided provided the variables used to describe the event are carefully chosen. As mentioned above, the basic idea is, that the variables should not allow to reconstruct the invariant mass of the jets since outside soft and collinear regions this precludes a one-to-one correspondence of Born-like virtual corrections and contributions with additional real radiation-which is required to uniquely define an event weight incorporating NLO QCD corrections. In fact as we will show in the next section the method proposed in ref. [25] is rather general and can be applied to arbitrary processes.
We note that similar ideas have been presented in ref. [26], refs. [24,27] and ref. [28]. In ref. [26] the problem is analysed from a mathematical point of view and a formal solution is given. However, the method requires the numerical solution of a non-linear system of equations together with the numerical computation of the Jacobian for the transformation mapping of the phase space variables to the variables used to describe the event. No proof is given in ref. [26] that this can be done in a numerically stable and efficient way. The method presented here is very similar to the approach developed in parallel in ref. [24]. The major difference is that in ref. [24] an additional prescription to balance the transverse momentum is used. In ref. [28], the generation of unweighted events of 'resolved' pseudo-partons is described which can be used to calculate infrared observables at NLO accuracy.
In fact, factorising the real phase space in terms of a Born-like phase space times the integration over the additional radiation is not a novel idea. In different contexts it has been widely used to improve the efficiency of numerical phase space integrations (see e.g. refs. [9,29,30]).
The article is organised as follows. In the next section we describe the method to calculate event weights including NLO QCD corrections. In section 3 the approach is validated by generating unweighted t-channel single top-quark events which are distributed according to the NLO cross sections. In section 4 the ability to predict event weights at NLO accuracy for jet events defined by conventional jet algorithms is employed in the Matrix Element Method. Examplarily, events obtained from a state-of-the art NLO+parton shower event generator are analysed with the Matrix Element Method at NLO. A brief summary and the conclusions are given in the last section.

NLO event weights for jet events defined by a → 1 recombination procedure
In this section we describe the calculation of a fully differential event weight including NLO QCD corrections. Fully differential means in this context that the number of variables in which the cross section is differential is maximal: r = 3n − 4 variables for the production of n jets in e + e − annihilation and r = 3n − 2 variables for the hadronic production of n jets. 1 We assume that an event for the process under consideration is described by a set of r variables {x 1 , . . . , x r }. On the one hand, it is convenient although not necessary to think of these variables as functions of the momenta of the n jets: x i = x i (J 1 , . . . , J n ) or x = x(J 1 , . . . , J n ) in vector notation, where the J i are the 4-momenta of the n jets which are calculated according to the chosen jet algorithm. Note that we keep the mass as a free parameter. On the other hand, experiments usually record values for variables related to energy depositions and particle tracks in the detectors which are attributed to jets. It is therefore also natural to think of {x 1 , . . . , x r } as experimentally accessible variables used to describe the 4-momenta of these resolved jets by imposing certain kinematics: The hadronic differential cross section evaluated for r specific values {y 1 , . . . , y r } of the variables in x is then given by: i . The phase space measure is given by To allow an efficient numerical evaluation, the delta-functions in eq. (2.1) need to be integrated out analytically. In the Born and virtual corrections this corresponds to the solution of a system of equations implied by the delta-functions, since all integrations are fixed through the deltafunctions. Another way to phrase this is to rewrite the phase space dΦ n in terms of the variables denotes the Jacobian of the transformation) which allows for an integration over the deltafunctions. In case of the real corrections a similar factorisation of the form is required, where dΦ i denotes the phase space associated with the 'unresolved' radiation of parton i (andJ( x, p i ) is the respective Jacobian). In general, for arbitrary variables {x 1 , . . . , x r } such a factorisation does not exist, because the jets J (n+1) k obtained from the recombination of two partons do not necessarily satisfy the same kinematical constraints as the ones used in the virtual corrections, e.g.: (Note that in the Born and virtual contributions no recombination takes place.) There are two different solutions to this problem. In refs. [20,21] it has been advocated to use a modified recombination procedure in which the clustered jets fulfill the Born kinematics: The inequalities in eq. (2.5) become equalities. At the same time the recombination procedure can be chosen such that the phase space factorises: (2.6) Using this factorisation it is straightforward to combine virtual and real corrections to define an event weight at NLO accuracy. The price to pay is however a modification of the recombination prescription which requires a significant re-tuning of existing Monte Carlo programs used in the experimental analysis. In this article we propose a second solution. We stick to the commonly used recombination: the 4-momentum of a jet obtained from the clustering of two partons i, j is given by the sum of the momenta p i and p j : The remaining n − 1 jet momenta are identified with the partonic momenta. (This recombination prescription defines the so-called 'E-scheme', see e.g. ref. [31]). Because the jet masses (J (i j) ) 2 are in general different from the ones occurring in the virtual corrections, the variables {x 1 , . . . , x r } must not depend on the jet masses, since this would prevent the point-wise combination of real and virtual contributions. In fact, this requirement can also be understood as a consequence of infrared safety.
Due to the different parton multiplicities of the phase space integrations in eq. (2.1) it is convenient to split the cross section differential in r variables of n resolved jets into the Born and virtual contribution (BV) and the real contribution (R). The real phase space can be partitioned further into regions R i j where the parton pair i j to be clustered is picked by the respective jet algorithm (see eq. (2.7)), regions R i where one parton i escapes detection and regionsR i where all partons are resolved as jets but (the softest) jet i is not considered in the event definition: The sums run over all i, j which can be clustered/omitted to still end up with the signal signature of the Born process. The last term is absent if additional jet activity is vetoed. The differential cross section can be interpreted as an event weight at NLO accuracy for the event defined by values {y 1 , . . . , y r } for the variables {x 1 , . . . , x r }. The technical implementation of the phase space partitioning is explained in detail in refs. [20,21,25].
In the following, we give explicit parametrisations of the n-and (n + 1)-parton phase space allowing the combination of the virtual and real corrections. Since we have to show that it is possible to factorise dΦ n+1 into a Born-like phase space times some additional contribution due to extra radiation, we start with dΦ n . In the Born and virtual part the n parton momenta are identified with the momenta of the n jets. A useful, although not unique, starting point is given by the following parametrisation (cf. eq. (2.3)): where s 2 = (P a + P b ) 2 and J( x) is the Jacobian of the variable transformation The jet masses, the transverse momentum and the variables {x 1 , . . . , x r } are used as independent integration variables. But on-shell conditions and momentum conservation fix the former, leaving the remaining variables {x 1 , . . . , x r } to be used to define a jet event. This factorisation thereby allows a straightforward integration of the delta-function . . , J (n) n ≡ J n ) in the Born and virtual contributions in eq. (2.1). With this parametrisation the Born and virtual part of the differential cross section reads The parametrisation of the n partonic momenta {p 1 , . . . , p n } in terms of {y 1 , . . . , y r } is fixed by the delta-functions in eq. (2.10) and the identification of n partonic with n jet momenta.
Let us now study the real contribution. We focus on those regions of the (n + 1) parton phase space in which additional radiation is recombined into a jet or associated with the beam resulting in an n-jet final state. The two cases need to be treated separately. (The remaining regions where the additional radiation is resolved as an additional jet do not impose any conceptual problems and can be obtained from a similar parametrisation of the real phase space as presented below.) First we consider the situation in which the additional radiation is clustered with a final-state parton. We denote with R i j a phase space region in which the two final-state partons i and j are combined according to eq. (2.7) to form a jet (i j) with invariant mass squared Note that the momenta of the remaining jets are identified with the underlying partonic momenta as in the Born and virtual contributions. By applying a factorisation of the real phase space Using the parametrisation from the Born and virtual part for the n-particle phase space in eq. (2.1) while the unresolved parton pair is to be integrated over the region R i j . The contribution of the real corrections to the differential cross section from the region R i j reads: 14) The parametrisation of the n + 1 partonic momenta {p 1 , . . . , p n+1 } in terms of {y 1 , . . . , y r }, M i j and the azimuthal (φ i ) and polar angle (θ i ) of parton i is fixed by the delta-functions in eq. (2.10) and the clustering of n + 1 partonic momenta to n jet momenta in region R i j . The above procedure can be applied to all regions R kl in which two partons k, l are recombined into one jet.
Let us now investigate the contribution from regions R i in which the parton i is considered as part of the beam and needs to be integrated out. Since there is no recombination in this case, the momenta of the resolved jets are identified with the underlying partonic momenta. A parametrisation of dΦ n+1 in terms of the same (3n − 2) independent variables as in eq. (2.10) can be used: The jet masses, the transverse momentum and the 3-momentum of the unresolved radiation appear as additional independent integration variables. But on-shell conditions and momentum conservation again fix the former while the 3-momentum p i and the variables {x 1 , . . . , x r } are left as independent integration variables. It is thus straightforward to integrate out the second delta-function in eq. (2.1) leaving the the unresolved radiation to be integrated over the region R i .
The contribution of the real corrections from the region R i to the differential cross section is then given by where the parametrisation of the n + 1 partonic momenta {p 1 , . . . , p n+1 } in terms of {y 1 , . . . , y r } and p i is fixed by the delta-functions in eq. (2.15) and the identification of n resolved partons' momenta with n jet momenta.
It is worth mentioning that with a slight modification of eq. (2.16) the contribution from the real corrections where the extra radiation i can be resolved as an additional jet with momentum J i but does not enter the event definition can be obtained: The integration over J i has to be carried out over the regionR i of the real phase space where parton i is resolved as an additional jet . . , p n+1 ( y, J i ) .

(2.17)
Collecting the Born and virtual contributions (see eq. (2.12)) and the real corrections from the different regions R i j (see eq. (2.14)), R i (see eq. (2.16)) andR i (see eq. (2.17)) we finally obtain the fully differential cross section including the NLO corrections (cf. eq. (2.9)).

Application I: Single top-quark events at the LHC
In this section the calculation of the jet event weight, including NLO corrections as defined in eq. (2.9), is illustrated. As a concrete example we study hadronic single t-channel top-quark production: The event signature is given by a top-tagged jet t (a jet containing a top quark) in association with at least one light jet j (a jet not containing a top-quark). The NLO calculation for t-channel single top-quark production is performed as described in ref. [21] 2 . The jets are defined by the k t -algorithm (as given in ref. [31] with R = 1) employing a 2 → 1 recombination according to the E-scheme (see eq. (2.7)). Resolved jets have to pass the following cuts on the transverse momentum and the pseudo rapidity Assuming that the detector is blind outside these cuts any radiation not satisfying the cuts has to be integrated out according to eq. (2.16). If there is more than one resolved light jet, the hardest is used in the event definition. We assume that it is always possible to identify the jet containing the top quark.
Due to momentum conservation and on-shell conditions each event is fully described by four variables x 1 , . . . , x 4 . The pseudo rapidity of the top-tagged jet and the energy, the pseudo rapidity and the azimuthal angle of the (hardest) light jet are chosen. Each set of measured values x = (η t , E j , η j , φ j ) specifies the 4-momenta of the resolved jets as functions of the squared jet masses J j 2 , J t 2 : From eq. (2.12) the Born and virtual contribution of the NLO event weight follows as with the jet momenta from eq. (3.2). According to eq. (2.10) the squared jet masses are given by J t 2 = m t 2 and J j 2 = 0 and the parton momentum fractions follow as: In the real corrections the top-quark (momentum p t ) and a light parton (momentum p l ) are produced together with additional radiation (momentum p r ) in the final state: The phase space of the real corrections contributing to single top-quark production in association with a light jet can thus be split into the regions R tr , R lr , R r andR r . For the clustering of the extra radiation with the top quark eq. (2.14) yields with the jet momenta from eq. (3.2). Because of the on-shell conditions the squared jet masses follow as J t 2 = M tr 2 = (p t + p r ) 2 and J j 2 = 0. The parton momentum fractions are again given by eq. (3.5). Note that eq. (3.7) already takes into account that either of the massless quarks can be clustered together with the top quark while the other one constitutes the light jet.
In the case where the light jet is obtained from the clustering of the light quark and the additional radiation eq. (2.14) yields with the jet momenta from eq. (3.2). The on-shell conditions result in J t 2 = m t 2 and J j 2 = M lr 2 = (p l + p r ) 2 while the parton momentum fractions are again given by eq. (3.5).
In the regions of the phase space where the extra radiation is associated with the beam, the jet momenta have to be parametrised by the variables (η t , E j , η j , φ j ) and the 3-momentum of the extra radiation in order to ensure momentum conservation: Note that the jet momenta defined in eq. (3.9) fulfill the on-shell conditions in eq. (2.15) by construction: J t 2 = m t 2 and J j 2 = 0. Using this parametrisation and eq. (2.16) together with eq. (3.3) yields with the parton momentum fractions given by Again, either of the massless quarks can escape detection while the other one constitutes the light jet.
The regionR r , corresponding to the extra radiation being resolved as an additional but softer light jet, contributes to the event weight as (cf. eq. (2.17)) with the parametrisation of the jet momenta and parton momentum fractions given in eq. (3.9) and eq. (3.11) with p r → J r .
According to eq. (2.9) the weight including NLO corrections for a jet event t, j defined by η t , E j , η j and φ j is given by the sum of the Born and virtual contribution (see eq. (3.4)), the real corrections from the regions R tr (see eq. (3.7)), R lr (see eq. (3.8)), R r (see eq. (3.10)) andR r (see eq. (3.12)). We stress that a 2 → 1 recombination is used where the momentum of the resulting jet is defined as the 4-momentum sum of the recombined partons. With this event weight unweighted jet events can be generated using a simple acceptance and rejection algorithm. As a consistency check we have generated ≈ 40000 unweighted events. The red dashed histogram in fig. 3.1 shows the E j , n j and η t distributions of these events. These distributions are compared with the result of an independent calculation using a parton level Monte-Carlo integration  (blue solid histogram in fig. 3.1). In the lower part of the plots the respective pull distributions, p-values and reduced χ 2 of the comparison of the two histograms as described in ref. [35] and implemented in refs. [36][37][38] are given. The results in fig. 3.1 show that the unweighted events are indeed distributed according to the cross section calculated at NLO accuracy. Within the statistical uncertainties the results from the parton-level Monte Carlo are in perfect agreement with the results obtained from the sample of unweighted events.

Application II: The MEM at NLO for POWHEG events
So far the matrix element method (MEM) at NLO has been restricted by the requirement of a 3 → 2 clustering in the jet algorithm [20,21]. Therefore, only closure tests of the method have been presented where the pseudo-data has been generated with the modified clustering prescription in order to generate unweighted events. The method presented in this paper and ref. [25] extends the MEM at NLO to incorporate 2 → 1 clusterings and therefore allows to analyse realistic events.
As a proof of concept the extraction of the top-quark mass from t-channel single top events at the LHC is presented. Contrary to the analysis presented in ref. [21,25], where the pseudo-data has been generated using the same NLO calculation as used in the MEM, the pseudo-data is generated using the POWHEG-BOX [39,40]. The sample is generated for the LHC at √ s = 13 TeV and a top-quark mass of m t = 173.2 GeV and a renormalization and factorization scale of µ R = µ F = µ 0 = m t . The events are subsequently showered using the Pythia8 parton shower (PS) program [41]. As in the previous case all light as well as top-tagged jets have to fulfill p ⊥ > 30 GeV and |η| < 3.5. Jets are defined using the k t -algorithm with a separation parameter R = 0.4. The resulting event sample for the fiducial phase space volume consists of N = 28031 events. During the parton shower evolution the top-quark is kept stable and corrections from hadronisation and underlying events have been neglected. Thus the main difference between the pseudo-data and the NLO calculation used for the analysis is the inclusion of additional radiation in the events.
The top-quark mass extraction proceeds as follows. For each showered POWHEG event the set of variables x i = (η t , E j , η j , φ j ) i is calculated from the top-tagged jet and the hardest light jet. Then for each event x i the NLO weight is calculated for a specific value of the top-quark mass and the corresponding Likelihood function 3 L(m t ) is computed according to where only trivial transfer functions (i.e. delta-functions) have been used. Minimizing the negative log-Likelihood function (− log L) with respect to the top-quark mass parameter yields the estimator for the top-quark mass m t . Fig. 4.1 shows the negative log-Likelihood function as a function of the top-quark mass and for various scale choices at LO and NLO in perturbation theory. The leading-order results are shown in blue for three different scales, the central scale µ 0 (blue) and scale variations µ 0 /2 (light blue) as well as 2µ 0 (dark blue). The red curves are obtained by including the full NLO QCD corrections in the definition of the Likelihood functions and are also shown for three different scales: µ 0 (red), µ 0 /2 (light red), 2µ 0 (dark red). The vertical green line denotes the value for the top-quark mass used for the generation of the pseudodata. The minimum of the curves yields the estimator for the top-quark mass m t and the width of the parabola gives the statistical uncertainty ∆ m t on the estimator. The estimated top-quark mass using the Born approximation for the central scale is thus given by m t = (150.19 ± 1.38 stat ) GeV. Based on the results for the two other scale settings the uncertainty due to neglected higherorder corrections is estimated to be of the order of ±5 GeV. We observe a significant difference between the results obtained from the MEM using LO matrix elements and the input value.
In particular, the observed shift of about 22 GeV is not covered by the uncertainties. Using the MEM based on leading-order matrix elements thus requires a significant calibration to reproduce the input value.
Using the full NLO QCD calculation for the determination of the event weights entering the Likelihood function the estimator for the top-quark mass gives m NLO t = (163.75 ± 1.83 stat ) GeV. A significant reduction of the scale uncertainty is observed. The uncertainty goes down from the aforementioned ±5 GeV to +0.85 GeV and −2.65 GeV in NLO. Furthermore the shift compared to the input value is reduced to 10 GeV compared to 22 GeV in leading-order.
In fact, a shift in the extracted top-quark mass to lower mass values is expected due to parton shower effects in the pseudo-data which are not taken into account in the MEM: Multiple parton emissions lead to a modification of the phase space density, which results in shape differences in differential distributions compared to fixed-order NLO computations. It is well known that the MEM is very sensitive to small distortions of the differential distributions. As the parton-shower tends to soften the p t distribution the MEM produces a smaller mass value. (In ref. [25] it has been shown that the input value is reproduced within the uncertainties if the parton shower is omitted.) Note that the analysis based on the Likelihood function given in eq. (4.1) is only sensitive to the normalised multi-differential cross section. As argued in ref. [21,25] the information about the fiducial cross section can improve the parameter determination. To incorporate the information on the total number of events in the sample the Extended Likelihood function, defined by is used. L denotes the integrated luminosity of the experiment. For the pseudo-data used in the analysis it is given by N = L·σ NLO+PS , where σ NLO+PS is the fiducial cross section corresponding to the simulated data. Accordingly, the expected number of events is given by the prediction for the fiducial cross section σ(m t ) as ν(m t ) = L · σ(m t ). Including the information on the number of recorded events should give a significant improvement: While the parton shower tends to soften the p t -distribution leading to the aforementioned small mass values, the total number of events is only mildly affected through acceptance effects. Since a smaller mass value leads While the statistical uncertainty is reduced by more than a factor of 2 the systematic uncertainty is enlarged by more than a factor 2, which amounts to roughly ±10 GeV. Repeating the extraction including the NLO QCD corrections in the Extended Likelihood function improves the analysis significantly with respect to the Born approximation. The top-quark mass estimator and its statistical uncertainty is given by m NLO t = (173.81 ± 0.72 stat ) GeV, which is compatible with the input mass within one standard deviation. Also, in the case of the Extended Likelihood function at NLO accuracy the statistical uncertainty is reduced by more than a factor of 2. In contrast to the Born approximation, the systematic uncertainty, as obtained by scale variations, is also reduced at NLO and amounts to an uncertainty of −0.30 GeV and +1.77 GeV. Furthermore, in the Extended Likelihood approach the NLO estimator is covered by the uncertainty of the leading-order estimator, which justifies to use scale variations as an estimate for the missing higher-order corrections. As mentioned above, the parton shower does not have a significant  impact on the number of accepted events in comparison to the fixed-order NLO calculation. However the predicted number of events depends approximately linearly on the value of topquark mass parameter [21,42]. Therefore, the Extended Likelihood analysis is driven towards top-quark mass values which correspond to a compatible prediction of the number of events with the pseudo-data. It is worth noting that with roughly 30000 events the uncertainty on the extracted top-quark mass is already dominated by systematic uncertainties. In the above analysis, no uncertainty for the luminosity is included. To estimate the impact of this additional uncertainty we assume a relative uncertainty ∆L/L ≈ 2% (see also [43]). Repeating the Extended Likelihood analysis while varying the value of L by ±2% results in additional shifts in both the LO and the NLO estimators for the top-quark mass of about ±2 GeV.

Conclusion
In this article an algorithm to calculate fully differential event weights including NLO QCD corrections is presented. We emphasise that the weight is calculated for events containing jets and not partons. This is a significant extension to existing approaches used in MC tools which allow the generation of jet events together with the corresponding weight but do not allow to evaluate the weight for a given jet event. In difference to previous work [20,21] the conventional 2 → 1 recombination is used-no change of the recombination procedure is required. The method developed in this article relies on properly chosen variables to describe the jet event. As explained in detail in the article, the variables must not depend on the jet masses since this would spoil the point-wise one-to-one correspondence of real and virtual corrections. We give a concrete example for a consistent choice of variables and show how the phase space is rewritten to allow the evaluation of the fully differential event weight.
Finally, the method is illustrated with two example applications. In the first example unweighted events for t-channel single top-quark production are generated using an acceptance and rejection algorithm. Using these unweighted events to calculate distributions, the results are in perfect agreement with results obtained using a traditional parton-level Monte Carlo. In a second example the presented method is used within the context of parameter determination. As a toy application we study the top-quark mass extraction with the Matrix Element Method. Here pseudo-data is generated using the POWHEG-BOX and subsequent showering using the Pythia8 parton shower. These events are analysed with Likelihood functions that include the fully NLO QCD corrected event weight. Within the Extended Likelihood approach perfect agreement with the input value for the top-quark mass of the pseudo-data is found. In particular, we show that in this approach no calibration is required.
We stress, once again, that the method presented here allows to calculate the corresponding weight including NLO corrections for a given event. The approach thus allows for the generation of unweighted events at NLO accuracy. Furthermore, the approach is not restricted to a specific process. For observables where corrections due to the parton shower are negligible, the unweighted events provide a useful approximation. It is also possible to use the unweighted events in a detector simulation. Comparison with generated events including the effects of the parton shower may allow for a more detailed study of parton shower effects. Ultimately, the ability to predict event weights at NLO accuracy for jet events enables the application of the Matrix Element Method at NLO for events recorded by the experiments without a modification of the recombination procedure in conventional jet algorithms.