Azimuthal asymmetries of charged hadrons produced in high-energy muon scattering off longitudinally polarised deuterons

Single hadron azimuthal asymmetries of positive and negative hadrons produced in muon semi-inclusive deep inelastic scattering off longitudinally polarised deuterons are determined using the 2006 COMPASS data and also combined all deuteron COMPASS data. For each hadron charge, the dependence of the azimuthal asymmetry on the hadron azimuthal angle ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document} is obtained by means of a five-parameter fitting function that besides a ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document}-independent term includes four modulations predicted by theory: sinϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sin \phi $$\end{document}, sin2ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sin 2 \phi $$\end{document}, sin3ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\sin 3\phi $$\end{document} and cosϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\cos \phi $$\end{document}. The amplitudes of the five terms have been extracted, first, for the hadrons in the whole available kinematic region. In further fits, performed for hadrons from a restricted kinematic region, the ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document}-dependence is determined as a function of one of three variables (Bjorken-x, fractional energy of virtual photon taken by the outgoing hadron and hadron transverse momentum), while disregarding the others. Except the ϕ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\phi $$\end{document}-independent term, all the modulation amplitudes are very small, and no clear kinematic dependence could be observed within experimental uncertainties.


Introduction
Measurements of Semi-Inclusive Deep-Inelastic Scattering (SIDIS) μ + N → μ + nh + X, n = 1, 2, . . . (1) of high-energy polarised muons μ off nucleons N in the initial state and scattered muons μ , n measured hadrons h and unobserved particles X in the final state are sensitive to the spin-dependent Parton Distribution Functions (PDFs) of nucleons. The SIDIS cross section depends, in particular, on the azimuthal angle of each produced and measured hadron (see e.g. Refs. [1][2][3][4][5]), which leads to azimuthal asymmetries related to convolutions of the nucleon Transversee-mail: Gerhard.Mallot@cern.ch Momentum-Dependent (TMD) PDFs and parton-to-hadron Fragmentation Functions (FFs). These asymmetries can appear in SIDIS off unpolarised, longitudinally or transversely polarised nucleons. The TMD PDFs were studied in a number of experiments. The short overview of earlier results obtained by the HERMES, CLAS and COMPASS collaborations on azimuthal asymmetries in SIDIS production of charged hadrons was given in Ref. [6,7]. The COMPASS collaboration has published results on asymmetries off unpolarised 6 LiD (referred to as "deuteron") target [8], transversely polarised deuterons [9,10] and transversely polarised NH 3 (referred to as "proton") target [11][12][13][14][15][16]. The common analysis of transversely polarised deuteron and proton data is included in Ref. [11][12][13][14][15][16] also. The updated overview of the TMD PDFs including the COMPASS results can be found in Ref. [17]. The COMPASS results on azimuthal asymmetries off longitudinally polarised deuterons based on the data collected in 2002, 2003 and 2004 were published in Ref. [6] for so called "integrated" asymmetries and asymmetries as functions of kinematic variables extracted in the restricted kinematic region. Similar data have been collected in 2006 but not published yet. The results on the integrated asymmetries for 2006 and for the combined 2002-2006 data which are presented in this Paper are extracted using hadrons from the whole available kinematic region, at variance with Ref. [6], while asymmetries as functions of the kinematic variables are extracted using the hadrons from a restricted kinematic region, similar to Ref. [6].
The Paper is organised as follows. The SIDIS kinematics, basic formulae and a brief theoretical overview are given in Sect. 2. The analysis of the 2006 data is described in Sect. 3. The results on the asymmetries of the combined 2002-2006 data are presented in Sect. 4. Systematic uncertainties are discussed in Sect. 5 and conclusions are given in Sect. 6.

Theoretical framework
The SIDIS kinematics is illustrated in Fig. 1.
The 4-momenta of the incident and scattered muons are denoted by l and l , respectively. The 4-momentum of the virtual photon is given by q = l − l with Q 2 = −q 2 . The angle of the momentum vector q of the virtual photon with respect to the incident muon is denoted by θ γ . The vectors p h and P denote the hadron momentum and the longitudinal target deuteron polarisation, respectively. Their transverse components p h T and P T are defined with respect to the virtual-photon momentum. The longitudinal component | P L | = |P | cos θ γ is approximately equal to | P | due to the smallness of the angle θ γ . The small transverse component is equal to | P T | = |P | sin θ γ where sin θ γ ≈ 2(M x/Q) √ 1 − y. Here, M is the nucleon mass and y = (qp)/( pl) is the fractional energy of the virtual photon, where p is the 4-momentum of the target nucleon. The angle φ denotes the azimuthal angle between the lepton scattering plane and the hadron production plane, while φ S denotes the angle of the deuteron polarisation vector with respect to the scattering plane: φ S = 0 • or 180 • for deuteron polarisation parallel or antiparallel to the beam direction, respectively. Furthermore, the Bjorken variable, x B j ≡ x = Q 2 /(2 pq), the fraction of the virtual-photon energy taken by a hadron, z = ( pp h )/( pq), the transverse momentum of a hadron, p h T , and the invariant mass of the photon-nucleon system, W 2 = ( p + q) 2 , that, together with Q 2 > 1 (GeV/c) 2 and 0 < y < 1, characterise SIDIS under study.
The general expression for the differential SIDIS cross section (see Refs. [1][2][3][4][5] and references therein) is a linear function of the incident muon polarisation P μ and of the longitudinal and transverse components P L and P T of the target deuteron polarisation P :  Fig. 1 The SIDIS kinematics shown for target deuteron polarisation P antiparallel to the beam direction Here, the first (second) subscript of the partial cross sections refers to the beam (target) polarisation: 0, L or T denote unpolarised, longitudinally or transversely polarised.
The azimuthal asymmetries of charged hadron production a h ± (φ) are defined as follows: where all cross sections are functions of the angle φ. The Eq.
(3) represents a definition of the experimentally measured asymmetries common for this Paper and for Refs. [6,7,18,19]. The first (second) superscript denotes the beam (target) spin orientation. The symbol ← denotes the incident muon spin orientation that, in the case of a positive charge of the incident muons, is mainly opposite to the beam direction. For the CERN muon beam, the average value of | P μ | is equal to 0.8. The beam polarisation does not enter in the definition of measured asymmetries. The symbols ⇒ and ⇐ denote the target deuteron spin orientations (polarisations) the first of which is parallel, considered further as positive (+), and the second one is antiparallel (−) to the beam direction (see Sect. 3.1). Substituting the general expression for dσ (Eq. 2) in the cross sections of the Eq. (3), one can obtain the expected contributions of the partial cross sections to the azimuthal asymmetries. As the result, when taking into account the signs of the target polarisations, one can see that only four partial cross sections contribute to the numerator of Eq. (3) and two to its denominator. In the numerator, we expect to have contributions from dσ 0L , P μ dσ L L and tan θ γ (dσ 0T + P μ dσ LT ), while in the denominator from dσ 00 and P μ dσ L0 . The explicit expression for these partial cross sections in terms of the PDFs and their dependences on the hadron azimuthal angle have been given in Refs. [6,7] and briefly commented below.
Following phenomenological considerations based on the QCD parton model of the nucleon and SIDIS in one-photon exchange approximation, the squared modulus of the matrix element, defining the cross sections, is represented by a number of diagrams. As an example 1 , the diagram accounting for the contribution to the SIDIS cross section of the chiralodd transversity PDF h 1 (x) convoluted with the chiral-odd Collins FF H ⊥ 1 (z) is shown in Fig. 2. The diagram contributes to the asymmetry Eq. (3) via the term dσ 0T . Other PDFs convoluted with corresponding FFs can also contribute to the cross section of spinless or unpolarised hadron production off longitudinally polarised deuterons and their expected azimuthal modulations. As motivated in Refs. [6,7], out of predicted terms the φ-independent term a 0 h ± and four mod- ulation terms up to the order of M/Q are retained for the analysis: The sign and the amplitude of each modulation is a subject of the a h ± (φ) data analysis (see Sect. 3.6). In Eq. (4), the term a 0 h ± is related to the well known helicity PDF g 1L contributing to asymmetries via dσ L L . The two terms with amplitudes a sin 2φ h ± and a sin φ h ± reported in Ref. [20] are related to the wormgear-L PDF h ⊥ 1L , and to the PDFs h L and f L , respectively, contributing to the asymmetries via dσ 0L . The transversity PDF h 1 and Sivers PDF f ⊥ 1T can also contribute to the term with amplitude a sin φ h ± via dσ 0T with a small factor tan θ γ . Other terms in Eq. (4), not seen yet in experiments with longitudinally polarised deuterons, are the terms with amplitudes a sin 3φ h ± and a cos φ h ± . They are related to the pretzelosity PDF h ⊥ 1T and the worm-gear-T PDF g 1T contributing to asymmetries via dσ 0T and dσ LT , respectively, suppressed by the small factor tan θ γ . The COMPASS results [6] obtained from the 2002-2004 data suggested some indications for a possible xdependence of terms with amplitudes a sin 2φ h ± and a cos φ h ± . The contribution of the term with amplitude a cos 2φ h ± which could have appeared from dσ 00 in the denominator of Eq. (3) is disregarded. This amplitude is expected [21,22] to be of the order of 0.1 and would enter Eq. (4) with the factor a 0 h ± , that is of the order of 10 −3 for integrated asymmetries (see Table  2), or with a 0 h ± (x) ≤ 0.05 for asymmetries as functions of kinematic variables (see Fig. 6). This is beyond our experimental accuracy. The same comments apply to possible contributions of terms with amplitudes a cos φ h ± and a sin φ h ± which also could originate from dσ 00 and dσ L0 of the denominator of Eq. (3), respectively. The negligible impact of the disregarded modulations on the amplitudes in Eq. (4) is confirmed by the 2006 data (see Sect. 3.6). All modulation amplitudes obtained in this Paper refer to the average value of the beam polarisation equal to −0.8.
The aim of this study is to continue searches for possible modulations in a h ± (φ) as manifestation of TMD PDFs describing the nucleons in the deuteron and to investigate the x, z and p h T dependences of the corresponding modulation amplitudes. For these purposes, we used first the 2006 deuteron data and then the combination of all 2002-2006 COMPASS deuteron data with longitudinal target polarisation.

Experimental set-up
The COMPASS set-up is a two-stage forward spectrometer with the world's largest polarised target and various types of tracking and particle identification detectors (PID) in front and behind of two large-aperture magnets SM1 and SM2. These detectors provide data for reconstruction of corresponding tracks. The spectrometer was operated in the high energy (160 GeV) muon beam at CERN. Its initial configuration (see Ref. [23]) was used for data taking in 2002-2004. During the long accelerator shutdown in 2005, the set-up was modified (see Ref. [24]). The major modifications influencing the present analysis were as follows: (i) the replacement of the two 60 cm long target cells (denoted as U and D) by three cells U, M and D of lengths 30 cm, 60 cm and 30 cm, (ii) the replacement of the target solenoid magnet by the new one with a wider aperture and (iii) the installation of the electromagnetic calorimeter ECAL1 in front of the hadron calorimeter HCAL1. The ECAL1 is not used in the analysis because it was not fully operational yet in 2006 and partially acted as a hadron absorber. These modifications of the apparatus where aimed at further reduction of systematic uncertainties, enlargement of the spectrometer acceptance and improvement of e/γ PID capabilities. These modifications have required reconsidering the Refs. [6,7] methods of data stability tests and asymmetry calculations (see Sects. 3.3, 3.6.) The data in 2006 were taken in two groups of periods. Each group is characterised by its initial set of polarisations in the target cells which are obtained by using different frequencies of the microwave field to polarise the target material (deuterons) in different cells at the certain direction of the target magnet solenoid field. The solenoid fields holds the polarisation. The field direction is denoted as f = +, if it coincides with the beam direction, or f = −, if opposite. The first group of the periods is denoted by G1 and other one by G2. Each period includes a certain number of intervals of continuous data taking (referred to as runs). The G1 data taking periods started with the initial setting of positive deuteron polarisation in target cells U and D and the negative one in cell M, both corresponding to f = +. After taking some number of runs, the field was reversed to f = − causing the reversal of the target cell polarisations, so that the were measured for each run in order to make sure that they are stable at the level of about 55%. If polarisations dropped below this limit, they were restored by the microwave field before the beginning of the next period. For G2 periods, the procedure was analogous but the initial setting of polarisation in the cells was opposite to the one in G1 at the same field f = +. The periodic reversal of the cell polarisations within each group of periods was used to estimate a possible timedependent systematics of the data. The change of the initial setting of the cell polarisations was used to estimate a possible systematic change of the spectrometer acceptance due to superposition of the solenoid field and the field of SM1. If there is no such systematic change, the acceptance in G1 and G2 periods must be the same for stable performance of the spectrometer.

Selection of SIDIS events and hadrons
Let us call as "SIDIS event" an event determined by Eq. (1) and reconstructed with tracks using the data recorded by the tracking and PID detectors. The overall statistics of 2006 is about 44.6 × 10 6 of preselected candidates for inclusive DIS and SIDIS event with The sample was obtained after rejection of runs that did not pass the data stability tests (see Sect. 3.3) and events that did not pass the reconstruction tests. The latter ones were rejected if the Z -coordinate (along the beam) of the interaction point (vertex) was determined with an uncertainty larger than 3σ of average which varied within 1.5-2 cm for different target cells.
The selection of SIDIS events from the preselected sample was done as described in Refs. [6,7]. For each SIDIS event, a reconstructed vertex with incident (μ) and scattered (μ ) muons and one or more additional tracks were required. Trajectories of the incident muons were required to traverse all target cells in order to have the same beam intensity for each of them. The track crossing more than 30 radiative lengths along the reconstructed trajectory was associated with μ . The cuts were applied on the quality of the reconstructed tracks forming vertices, the effective lengths of the target cells (28 cm, 56 cm, 28 cm), the momentum of incident muons (140 GeV/c − 180 GeV/c), the fractional energy carried by all tracks from the event (z < 1) and the fractional virtual-photon energy (0.1 < y < 0.9). About 36.6 × 10 6 SIDIS event candidates remained after cuts.
The distribution of track multiplicities per SIDIS candidates peaks at four. These tracks include scattered μ and hadron candidates. For a track to be identified as hadron, it was required that: its transverse momentum was larger than 0.05 GeV/c, it was produced in the current fragmentation region, as defined by the c.m. Feynman variable and it was associated with a cluster in one of the hadron calorimeters HCAL1 or HCAL2 with an energy deposit greater than 5 GeV in HCAL1 or 7 GeV in HCAL2. The efficiencies of the calorimeters above these energies are close to 100%. The energy of hadrons extended up to 120 GeV in the former and up to 140 GeV in the latter. All hadrons of the SIDIS candidates were included in the analysis of asymmetries. For the final selection of the SIDIS events and hadrons, the SIDIS candidates have to pass stability tests, as described in Sects. 3.3 and 3.4. The total number of hadrons in 2006 after afore-mentioned selections is 15.6 × 10 6 including 8.6 × 10 6 h + and 7.0 × 10 6 h − .
To summarise, the SIDIS events and hadrons have been selected from preselected candidates requiring: 140 GeV/c < p μ < 180 GeV/c, In order to test the influence of the stronger p h T and z kinematic cuts on the "integrated" azimuthal asymmetries, they were calculated summing up all selected hadrons (see Sect. 3.6) at variance with Ref. [6]. Azimuthal asymmetries as functions of the kinematic variables x, z or p h T were calculated in a restricted region following Refs. [6,7], i.e. summing up hadrons within the intervals given in Table 1. The number of hadrons within these intervals is reduced by a factor of about two compared to the total number. The distributions of selected SIDIS events as a function of Q 2 and y and of charged hadrons as a function of z and p h T for different data samples are presented in Fig. 3.

Tests of data stability
Taking advantage of the three-cell polarised target, stability tests for the 2006 data were performed by investigating variations of events from run to run for certain observables via ratios R i , where i is the run number, using the combined information from cells U and D denoted by (U + D), and that of cell M. One expects the ratio R i = (U i + D i )/M i to be independent of luminosity, close to unity and stable from run to run. In order to confirm this expectation, the ratios R i per run were obtained for the following observables that are relevant to the selection of SIDIS events and hadrons: number of SIDIS events, number of tracks per SIDIS event, number of clusters in HCAL1 (HCAL2) with E > 5 (7) GeV, average energy of clusters in HCAL1 (HCAL2), average energy of the associated clusters per event in HCAL1 (HCAL2) and average angle φ . The R i values as a function of the run number were fitted by constants R for all runs. It was found that most of these R i were stable within the ±3σ limits around the average values R ≈ 1.05, except for some runs and one of the periods. The stability of the measurement of the hadron azimuthal angle φ in the range ±180 • is essential for determination of asymmetries. Distributions of φ-values in this range were obtained for each run of data taking and average values φ i per run determined. The distribution of φ i had a Gaussian shape around the mean value equal to zero for almost all runs.

Acceptance-cancelling method for calculations of cross section ratios
The modified "acceptance-cancelling" double ratio method was used to calculate the ratios of the SIDIS cross sections for positive and negative target polarisations denoted as σ + /σ − and the 2006 asymmetries. In this Paper, the modified double ratio method is applied in three forms: first, in a form of "acceptance-cancelling", second, in a form of "cross sectioncancelling" and, third (see Sect. 5), again in the "acceptancecancelling" form to test the hadron yield stability.
In order to cancel acceptances, the method utilises double ratios, i.e. the product of two ratios of events. For the threecell target the method was modified as follows. The target cell M was artificially divided in two sub-cells M1 and M2, each 28 cm long, and two pairs of cells (U and M1) and (M2 and D) are considered below. The cells in each pair have equal lengths, i.e. equal densities of deuterons, but opposite polarisations p (+ or −) at a given solenoid field direction f (+ or −). For each pair of cells at a given f , one can construct the double ratio using the number of selected SIDIS events or hadrons. These numbers obtained from the cell i and denoted as N i p f are usually expressed via a product of a cell luminosity (L i f ) given by the beam intensity times the target cell material density a target cell acceptance (A i f ) and the corresponding cross section (σ p ): N i p f = L i f × A i f × σ p , i.e. luminosity, acceptance and cross section are folded in the number of events. Taking this relation into account as well as the COMPASS procedure of measurements divided in two groups of runs G1 and G2, one can construct for each pair of the target cells the "acceptance-cancelling" double ratio of events (hadrons) that provides a way to unfold the (σ + /σ − ) 2 . Particularly, for the polarisation settings at f = +, the two double ratios of event numbers constructed for the (U, M1) and for the (M2, D) pairs have the forms given by Eq. (5), where events for the first (second) ratio of each pair are taken from the G1 (G2) runs.
Substituting in the left parts of Eq. (5) the above expressions for N i p f one can see that, after "cancellations" of L i f and A i f , the double ratios of events are directly related to the cross section ratios squared. Because the luminosities of cells are equal, they contribute equally to the numerators and denominators and their cancellations are expected in each of the Eq. (5) ratios. If the acceptances A i f of the cells are similar at the same f in the G1 and G2 groups of runs (it is subject for tests below), they are also folded equally in the corresponding number of events and "cancel" in the double ratios, i.e. it is not necessary to calculate them (see Sect. 5). In the numerator of each ratio, the number of events (hadrons) are taken from the runs with the positive target cell polarisation, while in the denominator they are taken from the runs with negative target polarisation. Thus under above conditions each ratio of events (hadrons) in the left parts the Eq. (5) is equal to the ratio σ + /σ − , which is known to be close to unity (it is subject for test below). Hence, each double ratio in Eq. (5), which is equal to (σ + /σ − ) 2 , also have to be equal within statistical uncertainties and is expected to be close to unity. The stability of acceptances during the G1 and G2 runs have been checked using the cross sections cancelling double ratios of events (hadrons) in the forms similar to ones given in Eqs. (11,12) of Ref. [6] which are related to the ratios of acceptances.
Similarly, at f = − the two double ratios constructed for the same target pairs are: Because no requirements except time stabilities are imposed on the data, the Eqs. (5, 6) are valid and can be used for calculations of either the cross section ratios in the restricted kinematic regions or in the whole available region (see Sect. 3.6). Thus each of four double ratios of the events (hadrons) in Eqs. (5, 6) calculated with SIDIS events are related to the squared ratio of the SIDIS cross sections for positive and negative target polarisations determined with a part of the data. When statistically averaged, they can be used to calculate asymmetries with the whole data provided that (i) acceptances in the G1 and G2 periods are indeed stable and equal, (ii) the values of the double ratios calculated for SIDIS events and for hadrons with polarisation settings at f = + and f = − are stable and equal within statistical uncertainties. These requirements were checked with SIDIS event candidates and hadrons and final selections of them were determined.
Altogether, the stability tests have shown that (i) acceptances are stable and equal during the G1 and G2 groups of runs, (ii) the double ratios in Eqs. (5, 6) calculated with SIDIS events or hadrons are stable over the data taking peri-ods and contained inside the ±3σ corridors around the average values which are close to unity. In order to be accepted for analysis the average value per a given run of the acceptances, the angles φ i and the double ratio values defined by Eqs. (5,6) have to be within ±3σ limits of the corresponding mean value for all runs. Otherwise the run was rejected. The rejected runs contained about 10% of 2006 hadrons.

Extraction of azimuthal asymmetries in hadron production
For the extraction of the azimuthal asymmetries a h ± (φ) off the cross section ratios, the distributions of the charged hadrons h + and h − were separately analysed as a function of the azimuthal angle φ in the region from −180 • to +180 • divided into 10 φ-bins. For both h + and h − , the double ratios of the hadrons defined by Eqs. (5,6), , were calculated and combined as follows: where the symbol ⊕ means statistically weighted averaging.
As it was shown in Refs. [6,7], in first approximation, the squared ratios of cross sections (σ + /σ − ) 2 h ± (φ) are related to the asymmetries a h ± (φ) multiplied by polarisation terms. For each hadron charge, the polarisation term is given by the sum of the P i p f (x) values, each of them being the product of target cell polarisations |P i p f | and dilution factor f i (x), as defined in Refs. [6,8], where i, p and f are those used to calculate the ratio (σ + /σ − ) 2 k (φ), i.e. four polarisation values at each k. The weight W k is equal to the ratio of the number of hadrons, N k , to the total number of hadrons, N tot . Therefore, the a h ± (φ), referred to as single-hadron asymmetries, are:

The 2006 asymmetries
Following Sect. 3.2, the asymmetries a h ± (φ) were calculated in this Paper (i) as the "integrated" asymmetries using the total number of h + or h − , and (ii) as the asymmetries vs. one of kinematic variables x, p h T or z disregarding the others and using the numbers of h + or h − within the intervals defined in Table 1. In each case, the asymmetries were fitted   by the function from Eq. (4) using the standard least-squaremethod and extracting all asymmetry modulation amplitudes simultaneously.
Integrated asymmetries. These asymmetries for the 2006 data as a function of the azimuthal angle φ are shown in Fig. 4 together with results of the fits given in Table 2. In order to check the impact of the disregarded modulations, which could have appeared from SIDIS partial cross sections dσ 00 and dσ L0 , on the modulation amplitudes in Eq. (4), we have performed fits of the 2006 integrated asymmetries by a new fitting function which contains a numerator and denominator. In the numerator we have used the same modulations as in Eq. (4), but in the denominator we included the disregarded modulation with average amplitudes determined in Ref. [8]. For the asymmetry a h − (φ), the fitting function is as follows: h − cos φ 1 + 0.0412 cos 2φ + 0.0552 cos φ + 0.0008 sin φ .
By comparing results of this fit with results of the standard fit of the 2006 data shown in Fig. 5, it was found that the differences between values of modulation amplitudes in the numerator are smaller than 1% of the fit uncertainties. Similar results are obtained for a h + (φ) replacing amplitudes in the denominator of Eq. (9) by corresponding values from Ref. [8].

Integrated asymmetries
For the integrated asymmetries, the values of the combined 2002-2006 modulation amplitudes, which are shown in Fig. 5 and denoted by AV, were obtained using a statistical combination of four corresponding amplitudes. They are given in Table 2. The combined 2002-2004 modulation amplitudes, calculated as for 2006, are also shown in Table 2 in order to allow comparison to those of Refs. [6,7], which were calculated in the restricted kinematic region. As expected for the iso-scalar deuteron target, consistent results are obtained for the φ-independent terms a 0 h + and a 0 h − . All φmodulation amplitudes are consistent with zero within uncertainties of fits. Comparing results for the 2002-2006 combined data presented in Table 2 to the results for [2002][2003][2004] and to the results of Ref. [6], one can see that they are in agreement between themselves within the quoted uncertainties. This indicates that the integrated asymmetries with and without kinematic cuts of Table 1 Fig. 6.
Except for the a 0 h ± (x), all amplitudes when fitted by constants are found to be consistent with zero within statistical uncertainties (χ 2 /N DF 1). As expected, the a 0 h ± (x) for deuteron target have the same x-dependence for positive and negative hadrons. Additionally, the x-dependence of the a 0 h ± (x)/D 0 (x, y) values are presented in Fig. 7, where D 0 (x, y) is the virtual-photon depolarisation factor for each x interval multiplied by the average beam polarisation |P μ |, as defined in Refs. [6,7]. If the a 0 h ± (x) represent the main contributions to the asymmetries of Eq. (3), the values of a 0 h ± (x)/D 0 (x, y) by definition (see e.g. Ref. [25]) should be equal to the asymmetries A h ± 1d (x). Within experimental uncertainties, there is a good agreement between our data on a 0 h ± (x)/D 0 (x, y) and the data of Ref. [26] on A h ± 1d (x), which confirms the correctness of the results on the asymmetries calculated by the modified acceptance-cancelling method. The values of A h ± 1d (x) were obtained with the 2002-2004 data. A similar x-dependence was also observed with 2002-2006 data for the asymmetries A π ± 1d (x) and A K + 1d (x) obtained with the identified hadrons (see Ref. [24]).

Systematic uncertainties
The global compatibility test of the results on the asymmetries a h± (φ), that were obtained separately for 2002, 2003, 2004 and 2006 years, was performed by building the pull distribution: pull i = (a i − a ) × |σ 2 a i − σ 2 a | −1/2 , where a i is the asymmetry for a given year, hadron charge and kinematic interval, a is the corresponding weighted mean value over four years and σ denotes the corresponding standard deviation. The pull distribution had in total 750 entries compared to 540 for 2002-2004 years. The asymmetries reported in this Paper, in principle, could have the non-cancelled-acceptance or luminosity time-dependent effects folded in the event numbers of one of the Eqs. (5, 6)-ratio and, consequently, in one or several values of a i distorting the pull distribution. In the absence of such effects, as expected, the pull distribution follows the Gaussian distribution with the mean value consistent with zero (− 0.038 ± 0.033 in our case) and σ with unity (0.978 ± 0.024). This indicates that no significant time-dependent systematic effects are present in the data on asymmetries, i.e. systematic uncertainties in values of a i are smaller than statistical ones. Quantitative measures of possible systematic effects have been obtained by estimating additive and multiplicative uncertainties. Main contributions to possible additive systematic uncertainties could come from the instabilities of the hadron yields. The φ-stability of hadron yields in the 2006 data was checked following the procedure described in Refs. [6,7]. For this purpose, the double ratios of hadron numbers as a function of the azimuthal angle φ for different polarisation settings at the field f during the G1 and G2 runs were calculated as follows: Here, N i p f is the number of hadrons per φ-bin from target cell i with polarisation p and field f , as explained in Sect. 3.4. The ratios given by Eqs. (10) are modifications of ratios used in Refs. [6,7] for the case of two target cells. In the above double ratios, we expect the acceptance and luminosity cancellations and, as a result, the φ-stability of the hadron yields. If unstable, they could indicate possible systematics in the acceptance as well. The fits by constants (see Fig. 8) of the weighted sums F(φ) = F + (φ)⊕ F − (φ) in the φ-region from −180 • to +180 • for h + , h − and h + + h − of the 2006 data gave results consistent with unity within statistical uncertain-   [6,7] and hence adopted also for the combined 2002-2006 deuteron data.
Possible sources of multiplicative systematic uncertainties of the asymmetry evaluation are uncertainties in the determination of the beam and target polarisations and estimations of the dilution factor. The multiplicative systematic uncertainties of the extracted asymmetries due to uncertainties in the determination of the beam and target polarisations were estimated to be less than 5% each and those due to uncertainties of the dilution factor to be less than about 2%. When combined in quadrature, overall possible multiplicative systematic uncertainty of less than 6% was obtained.

Conclusions
The searches for possible azimuthal modulations in the single-hadron azimuthal asymmetries a h ± (φ), as a manifestation of TMD PDFs describing the nucleons in the longitudinally polarised deuteron, have been performed using all COMPASS deuteron data and the acceptance-cancelling method of analysis. For each hadron charge, beside the φindependent term, four possible modulations predicted by theory (sin φ, sin 2φ, sin 3φ and cos φ) and their dependence on kinematical variables are considered. The asymmetries have been calculated both for all selected hadrons ("integrated" asymmetries) and for hadrons as functions of kinematic variables within the restricted region.
For the "integrated" asymmetries, it was found that results in the restricted range of kinematic variables are consistent with those of the wider range. In other words, the restricting of the kinematic region reduces the statistics but does not change the values of the asymmetries beyond the sensitivity of this experiment. The same result was obtained in Ref. [6] for the asymmetries as a function of z.
The φ-independent terms a 0 h ± (x) of the asymmetries a h ± (φ), which are expected to originate mostly from the known helicity PDFs g 1L (x) ≡ g 1 (x), are connected to the virtual photon asymmetry A h ± 1d (x) = a 0 h ± (x)/D 0 (x, y). There is good agreement between the COMPASS data on a 0 h ± (x)/D 0 (x, y) and A h ± 1d (x) from Refs. [24,26] which confirms this expectation.
No statistically significant dependences of φ-modulation amplitudes were observed as functions of x, z or p h T when fitted by constants. Still, there are some hints (statistically not confirmed) for a possible x-dependence of the sin 2φ, sin 3φ and cos φ modulation amplitudes. The sin 2φ amplitude for h − is mostly positive and rises with increasing x, while for h + it is mostly negative and decreases with x. This behaviour agrees with that discussed in Refs. [17,20,27], if one takes into account different definitions of asymmetries by the HERMES and COMPASS Collaborations. The increase with x of the modulus of the cos φ amplitudes, related to the Cahn-effect [28,29] and predicted in Ref. [30], was already  [6] and persists for the combined 2002-2006 data. Hints for a possible x-dependence of sin 3φ modulation amplitudes are discussed in Ref. [17]. Quantitative estimates of a possible contribution of the cos φ modulation to the deuteron asymmetries, related to TMD PDFs g ⊥ L and e L , have been obtained in Ref. [31]. They are in agreement with our data.
Altogether, one can conclude that contributions of TMD PDFs convoluted with FFs to the azimuthal asymmetries in the cross sections of hadron production in muon SIDIS off longitudinally polarised deuterons are small. This is either due to possible cancellations of the contributions to the asymmetries by the deuteron's up and down quarks, or/and due to the smallness of the transverse component of the target polarisation and of the suppression factor that behaves as M/Q. Some of these conclusions can be checked studying these asymmetries in muon SIDIS off longitudinally polarised protons.