A search for pair production of new light bosons decaying into muons

: A search for the pair production of new light bosons, each decaying into a pair of muons, is performed with the CMS experiment at the LHC, using a dataset corresponding to an integrated luminosity of 20.7 fb−1 collected in proton–proton collisions at center-of-mass energy of s=8 TeV . No excess is observed in the data relative to standard model background expectation and a model independent upper limit on the product of the cross section, branching fraction, and acceptance is derived. The results are compared with two benchmark models, the first one in the context of the next-to-minimal supersymmetric standard model, and the second one in scenarios containing a hidden sector, including those predicting a nonnegligible light boson lifetime. A search for the pair production of new light bosons, each decaying into a pair of muons, is performed with the CMS experiment at the LHC, using a dataset corresponding to an integrated luminosity of 20 . 7 fb − 1 collected in proton–proton collisions at center-of-mass energy of √ s = 8 TeV. No excess is observed in the data relative to standard model background expectation and a model independent up- per limit on the product of the cross section, branching fraction, and acceptance is derived. The results are compared with two benchmark models, the ﬁrst one in the context of the next-to-minimal super- symmetric standard model, and the second one in scenarios containing a hidden sector, including those predicting a nonnegligible light boson lifetime. © 2015 CERN for the beneﬁt of the CMS Collaboration. Published by Elsevier B.V. This is an open access article under the CC BY license (http://creativecommons.org/licenses/by/4.0/). Funded by SCOAP 3 .


Introduction
In July 2012 the ATLAS and CMS Collaborations at the CERN LHC announced the discovery of a particle [1-3] with properties consistent with the standard model (SM) Higgs boson [4][5][6][7]. Direct measurements of the production and decay rates of the new particle, using SM decay channels, have so far played a key role in determining whether or not it is indeed consistent with the SM predictions. However, substantially increasing the precision of these measurements will require further data. Searches for Higgs bosons through production mechanisms not predicted by the SM, or decay modes involving particles not included in the SM, provide a complementary approach and have the advantage of probing specific types of new physics models with the existing data.
This letter presents a search for the pair production of new light bosons (denoted as 'a') decaying to pairs of isolated, oppositely charged muons (dimuons). One production mechanism for these new bosons is in the decay chain of a Higgs boson h, which can be SM-like or not: h → 2a + X → 4μ + X, where X denotes possible additional particles from cascade decays of the Higgs boson. A range of new physics scenarios predict this decay topology, including the next-to-minimal supersymmetric standard model (NMSSM) [8] and models with hidden (or dark) sectors [9][10][11].
The NMSSM is an extension of the minimal supersymmetric standard model (MSSM) [12,13] that includes an additional gauge singlet field. It resolves the so-called μ problem [14] and significantly reduces the amount of fine tuning required in the MSSM [15]. The NMSSM Higgs sector consists of three CP-even neutral Higgs bosons h 1,2,3 , two CP-odd neutral Higgs bosons a 1,2 and a pair of charged Higgs bosons h ± . The h 1 and h 2 can decay via h 1,2 → 2a 1 , where either the h 1 or h 2 can be the boson observed at 125 GeV. The a 1 boson can be light and couple weakly to SM particles with a coupling to fermions proportional to the fermion mass. Therefore it can have a substantial branching fraction B(a 1 → μ + μ − ) if its mass is within the range 2m μ < m a 1 < 2m τ [16,17] (benchmark model 1 in this letter). A search for final states containing muon pairs provides sensitivity to models of this form.
Supersymmetry (SUSY) models with dark sectors (dark SUSY) offer an explanation for the excess in the ratio of the positron flux to the combined flux of positrons and electrons observed by the satellite experiments [18][19][20] in primary cosmic rays as well as predict cold dark matter with a scale of O(1 TeV). A simple realization of these models includes a new U(1) D symmetry (the subscript "D" stands for "Dark") which is broken and gives rise to massive dark photons (denoted as γ D ). Kinetic mixing of the new U(1) D with the SM hypercharge U(1) Y provides a small mixing between γ D and the SM photon which allows γ D to decay to SM particles [21]. Depending on the value ε of the kinetic mixing, the γ D may also be long-lived. The lack of an antiproton to proton ratio surements of the cosmic ray spectrum constrains the mass of γ D to be less than twice the mass of the proton [22]. If the hidden sector directly or indirectly interacts with the Higgs field, a number of possible scenarios may be realized. One such scenario, denoted in this letter as benchmark model 2, is a model of SUSY where the SM-like Higgs boson can decay via h → 2n 1 , where n 1 is the lightest neutralino in the visible (as opposed to hidden) part of the SUSY spectrum. The n 1 can decay via n 1 → n D + γ D , where n D is a dark neutralino that escapes detection. Assuming that γ D can only decay to SM particles, the branching fraction B(γ D → μ + μ − ) can be as large as 45%, depending on the mass of γ D [11].
Previous searches for pair production of new light bosons decaying into dimuons were performed at the Tevatron [23] and the LHC [24,25]. Searches for associated production of the light CPodd scalar bosons have been performed at e + e − colliders [26,27] and the Tevatron [28]. Direct a 1 production has been studied at the LHC [29], but this is heavily suppressed by the typically very weak couplings of the new bosons to SM particles. The constraints on the allowed NMSSM parameter space are driven by the measurements of relic density by WMAP [30] and more recently by PLANCK [31], while specifically for the Higgs sector the most relevant measurements come from LEP [32][33][34][35][36][37] Results are presented in this Letter in the context of the two benchmark scenarios discussed earlier, one in the context of NMSSM and another one in the framework of dark SUSY scenarios. However, the search has been designed to be independent of the details of these two specific models, and the results can be interpreted in the context of other models predicting the production of the same final states. Compared to the previous version [25], the present analysis has been redesigned to be sensitive to signatures with the intermediate bosons traversing a nonnegligible distance before decaying into a pair of muons. Such signatures can be realized in dark SUSY models if the mixing of the dark photon with its SM counterpart is sufficiently weak. In addition, the present analysis uses a dataset four times larger than the previous analysis, and at a higher center-of-mass energy, further extending the reach for signatures with prompt muons.

The CMS detector
This search is based on a data sample corresponding to an integrated luminosity of 20.7 fb −1 of proton-proton collisions at a center-of-mass energy √ s = 8 TeV, recorded by the CMS detector in 2012. The central feature of the CMS apparatus is a superconducting solenoid of 6 m internal diameter, providing a magnetic field of 3.8 T. Within the solenoid volume are a silicon pixel and strip tracker, a lead tungstate crystal electromagnetic calorimeter, and a brass and scintillator hadron calorimeter, each composed of a barrel and two endcap sections. Muons are measured in gas-ionization detectors embedded in the steel flux-return yoke outside the solenoid. Extensive forward calorimetry complements the coverage provided by the barrel and endcap detectors.
Muons are measured in the pseudorapidity range |η| < 2.4, with detection planes made using three technologies: drift tubes, cathode strip chambers, and resistive plate chambers. Matching muon candidates to tracks measured in the silicon tracker results in an accurate measurement of the transverse momentum (p T ). As an example, for muons with p T < 10 GeV the relative p T resolution is found to be 0.8%-3.0% (depending on |η|) and for muons with 20 < p T < 100 GeV it is 1.3-2.0% in the barrel and better than 6% in the endcaps [57]. A more detailed description of the CMS detector, together with definitions of the coordinate system used and the relevant kinematic variables, can be found in [58].

Data selection
The data were collected with an online trigger selecting events containing at least two muon candidates, one with p T > 17 GeV and another with p T > 8 GeV. In this analysis offline muon candidates are defined as particle-flow (PF) muons [57]. The PF reconstruction algorithm combines information from all CMS subdetectors to identify and reconstruct individual particles, such as electrons, photons, hadrons or muons. Events are further selected by requiring at least four offline muon candidates with p T > 8 GeV and |η| < 2.4 that form two oppositely charged pairs. At least one of these muons must additionally satisfy the requirement of p T > 17 GeV and |η| < 0.9, which ensures that the trigger efficiency is high and independent of the event topology, including effects related to overlaps of nearby muon trajectories. Tracks associated with a pair of opposite-charge muon candidates are fit for a common vertex using a Kalman filter algorithm [59]. If the vertex is reconstructed, a muon pair is combined into a dimuon system if its invariant mass measured at the common vertex m μ + μ − < 5 GeV and the vertex fit prob- Muon pairs failing these requirements are still retained for the analysis if at the point of closest approach of the two trajectories they are within imuthal angles in radians. This recovery step is designed to compensate for the reduced efficiency of the vertex selection for dimuons in which the two muon tracks are nearly parallel to each other, therefore a good efficiency is maintained for dimuon masses down to the 2m μ threshold (0.2114 GeV). For dimuons in which this is the case the point of closest approach is selected as the vertex position, with the additional selection requirement that the distance between the tracks be ≤ 0.5 mm. The dimuon kinematic variables are measured at the dimuon vertex position. There is no restriction on the number of ungrouped additional muons. Both dimuons are required to have at least one hit in the first layer of the barrel or endcaps of the pixel detector, and this defines an effective "fiducial" region. This requirement, along with the muon p T and |η| criteria, ensures high trigger (>96%), reconstruction, and selection efficiencies, with a greatly reduced dependence on the p T , η, or opening angle between the muons.
The projected z coordinate of the dimuon system at the point of the closest approach to the beam line (z μμ ) is reconstructed using the dimuon momentum. The requirement |z 1μμ − z 2μμ | < 1 mm is imposed to ensure that both dimuons are consistent with the same pp interaction; no explicit requirements are made on the impact parameter or the z coordinate at the point of closest approach to the beam line of the individual reconstructed muons to preserve sensitivity to signatures with displaced muons.
To suppress background events in which the muons are produced in the decay of heavy quarks (and thus appear in jets), the dimuons are required to be isolated from other event activity using the criterion I sum < 2 GeV. The isolation parameter I sum is defined as the scalar sum of the p T of charged tracks with p T > 0.5 GeV within a cone of size R = 0.4 centered on the momentum vector of the dimuon system, excluding the tracks corresponding to the two muon candidates. The tracks used in the calculation of I sum must also have a z coordinate at the point of closest approach to the beam line that lies within 1 mm of z μμ . The I sum selection Table 1 Event selection efficiencies sim (m h1 , m a1 ) and sim (m γD , cτ γD ), as obtained from simulation, the geometric and kinematic acceptances α gen (m h1 , m a1 ) and α gen (m γD , cτ γD ), calculated using only generator-level information, and their ratios (with statistical uncertainties), for a few representative NMSSM and dark SUSY benchmark samples. The experimental data-to-simulation scale factor ( data / sim , described later) is not applied. suppresses the contamination from bb production by about a factor of 40, as estimated using a bb enriched control sample with one dimuon recoiling off a jet containing an unpaired muon, while rejecting less than 20% of events with the signal topology. The invariant mass m 1μμ always refers to the dimuon containing a muon with p T > 17 GeV and |η| < 0.9. For events with both dimuon systems containing such a muon, the assignment of m 1μμ and m 2μμ is random for compatibility with the background modeling schema described in Section 5. The invariant masses of both reconstructed dimuons are required to be compatible within the detector resolution, specifically |m 1μμ − m 2μμ | < 0.13 GeV + 0.065 (m 1μμ + m 2μμ )/2, which defines a diagonal signal region in the plane of the invariant masses of the two dimuons. The numerical parameters in the requirement correspond to at least five times the size of the core resolution in the dimuon mass.

Signal modeling
The results from this analysis are designed to be model independent, but are also presented in the context of the two benchmark models introduced earlier. NMSSM simulation samples for benchmark model 1 are generated with pythia 6.4.26 [60] [62] is used to force the Higgs bosons to undergo a non-SM decay to a pair of neutralinos, each of which decays via n 1 → n D + γ D , where m n 1 = 10 GeV, m n D = 1 GeV, which is representative of the type of models considered [42]. Dark photons are generated with m γ D in the range 0.25-2.0 GeV and a decay length cτ γ D in the range of 0-20 mm. Each of the two dark photons is forced to decay to two muons, while both dark neutralinos escape detection. The narrow-width approximation is imposed by setting the widths of the dark photons to a small value (10 −3 GeV).
All benchmark samples are generated using the leading-order CTEQ6.6 [63] set of parton distribution functions (PDF), and are interfaced with pythia using the Z2* tune [64] for the underlying event activity at the LHC and to simulate jet fragmentation.
The signal samples are processed through a detailed simulation of the CMS detector based on Geant4 [65] and are reconstructed with the same algorithms used for data. Table 1 shows the event selection efficiencies sim obtained using the simulated benchmark samples for a few representative choices of (m h 1 , m a 1 ), and (m γ D , cτ γ D ). To provide a simple recipe for future reinterpretations of the results in the context of other models, the variable α gen is separately defined as the geometric and kinematic acceptance of this analysis calculated using only generator-level information. It is defined by selecting events containing at least four muons with p T > 8 GeV and |η| < 2.4, with at least one of these muons having p T > 17 GeV and |η| < 0.9. The new light boson must also decay with transverse decay length L xy < 4.4 cm and longitudinal decay length L z < 34.5 cm (both defined in the detector reference frame), to satisfy the "fiducial" region of the analysis. Table 1 shows α gen along with the ratio sim /α gen .

Background estimation
The SM background for this search is dominated by bb production and has small contributions from the electroweak production of four muons and direct J/ψ pair production. The leading part of the bb contribution is due to b quark decays that result in a pair of muons, via either the semileptonic decays of both the b quark and the resulting c quark, or via resonances, i.e. ω, ρ, φ, J/ψ . A smaller contribution comes from events with one genuine dimuon candidate and a second dimuon candidate containing one muon from a semileptonic b quark decay and a charged hadron misidentified as another muon.
Using data control samples, the bb background is modeled as a in the momentum resolution of the barrel and endcap regions of the tracker. The B 17 and B 8 distributions are fitted with a parametric analytical function using a sum of Bernstein polynomials and Crystal Ball functions [66] describing resonances. These event samples do not overlap with the sample containing two dimuons that is used for the main analysis and they have negligible contributions from non-bb backgrounds. Once the B bb (m 1μμ , m 2μμ ) template is constructed, it is used to provide a description of the bb background shape in the signal region. This technique assumes that each b quark fragments independently and that if the shapes of these distributions are measured using data samples with kinematics very similar to that of the background events then the effects of residual kinematically-induced correlations are small (albeit weakly induced, the shape depends on the b jet p T and the p T of the two b jets in background bb events tend to be similar). The background template is validated in a region where both dimuons fail the I sum < 2 GeV requirement and good agreement with data is observed.
The data events that satisfy all analysis selections but fail the m 1μμ m 2μμ requirement are used to normalize the B bb (m 1μμ , m 2μμ ) template. This selection yields nine events in the offdiagonal sideband region of the (m 1μμ , m 2μμ ) plane, leading in the diagonal signal region to an expected rate of bb background events of 2.0 ± 0.7. This is essentially (9 ± √ 9) × 0.18/0.82, where 0.18 and 0.82 correspond to the integral of the areas under the background template inside and outside the signal diagonal region, respectively. These nine events in the off-diagonal sidebands of the (m 1μμ , m 2μμ ) plane are shown as white circles in Fig. 1.
The contribution from direct J/ψ pair production is estimated using another data control sample. Events are selected with a trigger that requires at least three muon candidates, two of which have a common vertex and an invariant mass consistent with that of the J/ψ particle. Events are further required to contain at least four reconstructed muons with p T > 3.5 GeV, which form dimuon pairs. This control sample does not specifically require that the dimuons satisfy the requirement I sum < 2 GeV since I sum is used to separate the contribution of "prompt" and "nonprompt" (from b quark decays) J/ψ in data. Finally, both dimuons are required to have an invariant mass between 2.8 and 3.3 GeV. Following these requirements the data sample consists of events containing prompt and nonprompt J/ψ . To subtract the nonprompt com-ponent, two independent methods have been studied: the first one divides the control sample based on the values of the isolation variable I sum for each of the two dimuons in each event. The number of events in which both dimuons satisfy the requirement I sum < 2 GeV is extrapolated from the regions in which at least one of the dimuons fails this requirement. The second approach uses the lifetime of the J/ψ candidate, calculated under the hypothesis of it being produced at the beam line, as a discriminating variable. The data distribution is fitted in the isolated region using prompt and nonprompt templates from simulation and nonisolated sideband in data, respectively. Both approaches give consistent results within the associated uncertainties and the results of the isolation-based method are used in the final analysis. There are two mechanisms for the production of prompt double J/ψ events: single-and double-parton scattering (SPS and DPS, respectively), corresponding to whether the two J/ψ mesons are produced from one or two independent parton interactions. The number of prompt events in the control region is further separated into SPS and DPS components using the J/ψ rapidity difference as the discriminating variable. Finally, the data-to-simulation normalization factor and the fraction of SPS and DPS events are extrapolated from the control to the signal region, resulting in a final estimation for the contribution from prompt double J/ψ events of 0.06 ± 0.03 events.
The contribution from other SM processes (low mass Drell-Yan production and pp → Z/γ * → 4μ) is estimated with the CalcHEP 3.6.18 generator [67] using the HEPMDB infrastructure [68], and is found to be 0.15 ± 0.03 events in the entire signal region. The combined expected background contribution to the diagonal signal region is 2.2 ± 0.7 events. This background is represented by the color scale in Fig. 1.

Systematic uncertainties
The selection efficiencies for the offline muon reconstruction, trigger, and dimuon isolation requirements are obtained from simulation, and are corrected with scale factors derived from comparison between data and simulation using Z → μμ and J/ψ → μμ samples. The scale factor per event is found to be data / sim = 0.93 ± 0.07 and it accounts for the differences in the efficiency of the trigger, the efficiency of the muon reconstruction and identification for each of the four muon candidates, and the combined efficiency of the isolation requirement for the two dimuon candidates. The estimate accounts for correlations associated with the presence of multiple muons per event. The main systematic uncertainty is the offline muon reconstruction (4.1%), which includes an uncertainty (1% per muon) to cover variations of the scale factor as a function of the muon p T and η. Other systematic uncertainties include: the uncertainty in dimuon reconstruction effects related to overlaps of muon trajectories in the tracker and in the muon system (3.5%), the trigger efficiency (1.5%), the uncertainty in the efficiency caused by the modeling of the tails in the dimuon invariant mass distribution that arises from the requirement that the two dimuon masses are compatible (1.5%), and the dimuon isolation (negligible). The uncertainty in the integrated luminosity of the data sample (2.6%) [69] is also included. All uncertainties quoted above are related to variations in the signal efficiency due to experimental selection and sum up to 6.3%. The uncertainties related to variations in the signal acceptance due to the model include: the uncertainties related to the PDFs and the knowledge of the strong coupling constant α s , which are estimated by comparing the PDFs in CTEQ6.6 [63] with those in NNPDF2.0 [70] and MSTW2008 [71], following the PDF4LHC recommendations [72,73]. Using the analysis benchmark samples, they are found to be 3% for the signal acceptance. The variation of the renormalization and factorization scales has a negligible effect. In addition a re-weighting procedure is applied to the Higgs boson p T spectra in the benchmark signal samples to reproduce the NNLO+NNLL prediction [74] and account for possible changes to the analysis acceptance. Only a weak sensitivity to this is expected and the result of the re-weighting procedure is limited by the statistical uncertainty of the simulated samples (2%), which is therefore included as a conservative systematic uncertainty on this effect. Thus, the total systematic uncertainty in the signal acceptance and selection efficiency is 7.3%.

Results
After the full analysis selection is applied to the data sample, one event is observed in the diagonal signal region, as shown in Fig. 1. This is consistent with the expected background contribution of 2.2 ± 0.7 events.
For future reinterpretations of this analysis, the results can be presented as a 95% confidence level (CL) upper limit on σ (pp → 2a + X) B 2 (a → 2μ) α gen = N(m μμ )/ (Lr), where α gen is the generator-level kinematic and geometric acceptance defined earlier. The calculation uses the integrated luminosity L = 20.7 fb −1 and central value r of the ratio r = data /α gen = 0.63 ± 0.07. The ratio includes a scale factor correcting for experimental effects not included in the simulation and its uncertainty covers the variation in the ratio over all the benchmark model points used. The limit is calculated as a function of the dimuon mass using the CL S approach [75,76]. The chosen test statistic is based on the profile likelihood ratio and is used to determine how signal-or background-like the data are. Systematic uncertainties are incorporated in the analysis via nuisance parameters with a log-normal probability density function and are treated according to the frequentist paradigm. The overall statistical methodology used in this analysis was developed by the ATLAS and CMS Collaborations in the context of the LHC Higgs Combination Group and is described in [3,77]. The obtained limit as a function of dimuon mass m μμ can be conveniently approximated as a constant everywhere except the vicinity of the observed event, where it follows a Gaussian distribution: where m μμ is measured in GeV and the cross-section limit is expressed in femtobarns. This limit is applicable to models with two pairs of muons coming from light bosons of the same type with a mass in the range 2m μ < m a < 2m τ , where the new light bosons are typically isolated and spatially separated (so as to satisfy the isolation requirements). The weak model dependence of the ratio r allows for a simple reinterpretation of the results in other models. This requires calculating α gen , as defined earlier, and then the full efficiency data can be calculated by multiplying α gen by the ratio r.
There are certain subtleties that must be taken into account when reinterpreting the model-independent results of this analysis in the context of other models, particularly with the isolation requirement. An event should be considered to satisfy the selection requirements if there are at least two well isolated γ D decaying to muon pairs. Experimentally, isolation is based on charged tracks but it may be insufficient to just require the absence of generatorlevel charged particles in the isolation cone. For example a neutral pion decaying to a pair of photons, that convert into electrons, may result in the reconstruction of one or more tracks. This would be particularly relevant for models with more than two dark photons produced in the same event, some of which may decay to hadrons or electrons. In this case the safest approach is to require that there are no particles with p T > 0.5 GeV within the γ D isolation cone. This restriction would result in a more conservative limit but it would be robust against these effects.
The results from this analysis are also interpreted in the context of the NMSSM and the dark SUSY benchmark models, and 95% CL upper limits on the product of the cross section and branching fraction are derived. In these models both the Higgs boson production cross section and the branching fractions can vary significantly, depending on the choice of parameters. In the absence of broadly accepted benchmark scenarios, the production cross sections in these examples are normalized to that of the SM Higgs boson with a mass of 125 GeV [78].
In the case of the NMSSM benchmark scenario, the production cross sections and branching fractions for h 1 and h 2 can vary substantially depending on the chosen parameters. An exact interpretation of these results requires evaluating the experimental acceptance using the generator-level acceptance for each of the two h 1,2 bosons, and then using the measured upper limit on the sum of two contributions to derive limits for any choice of NMSSM model parameters. To present results in a fashion allowing for straightforward interpretation, we note that if one of the two CP-even Higgs bosons is the 125 GeV state observed at the LHC, then the other one is either lighter or heavier. In the NMSSM it is typical that one of the two has approximately the SM production cross section and a small B(h i → 2a), whereas the other one has a suppressed production rate and large B(h i → 2a) due to its large singlet fraction. In Fig. 2 (left) the limit at each mass point is calculated taking the CP-even Higgs boson with the corresponding mass as the only source of signal events; the curve below 125 GeV applies to NMSSM models in which m h 1 < m h 2 = 125 GeV, with h 1 decays dominating the rate of 4μ events. The limit at m h = 125 GeV corresponds to the case where 125 GeV = m h 1 < m h 2 , with h 1 decays still responsible for the vast majority of signal-like events. The points above 125 GeV correspond to model points for which only h 2 (m h 2 > m h 1 = 125 GeV) is allowed to have a sizeable rate of observable 4μ events. Finally, for models with m h 2 > 150 GeV, the limit at 150 GeV can be used as a conservative estimate of the production rate limit. In each of these scenarios it is possible that the other Higgs boson also contributes some fraction of the 4μ signal events, in which case the limit shown is more conservative than would be given by an exact evaluation.
In the case of the dark SUSY scenario, a 95% CL limit on the product of the Higgs boson production cross section and the branching fractions of the Higgs boson (cascade) decay to a pair of dark photons is determined. The limit set in the (m γ D , ε) plane from this analysis is shown in Fig. 2 (right), along with limits from other experimental searches, where the lifetime is directly related to the kinetic mixing parameter ε and the mass of the dark pho- that depends only on the mass of the dark photon [79]. The significant vertical structures in the limits visible in Fig. 2 (right) arise because the total width of the dark photon varies rapidly in those mass regions due to resonant decays to hadrons. This search constrains a large, previously unconstrained area of the parameter space. Unlike the other results in the figure, the CMS and ATLAS limits are model-dependent and only valid under the assumption that B(h → 2n 1 → 4μ + X) = 0. The recent ATLAS analysis [42]  and m a1 = 3.55 GeV (dotted curve) on σ (pp → h 1/2 → 2a 1 ) B 2 (a 1 → 2μ) as a function of m h1 in the range 86 < m h1 < 125 GeV and of m h2 for m h2 > 125 GeV. As an illustration, the limits are compared to the predicted rate (solid curve) obtained using a simplified scenario with σ (pp → h i → 2a 1 ) = 0.008 σ SM , which yields predictions for the rates of dimuon pair events comparable to the obtained experimental limits, and B(a 1 → 2μ) = 7.7%. The chosen B(a 1 → 2μ) is taken from [17] for m a1 = 2 GeV and tan β = 20. Right for benchmark model 2: 95% CL upper limits (black solid curves) from this search on σ (pp → h → 2γ D + X) B(h → 2γ D + X) (with m n1 = 10 GeV, m nD = 1 GeV) in the plane of two of the parameters (ε and m γD ) for the dark SUSY scenarios, along with constraints from other experiments [42][43][44][45][46][47][48][49][50][51][52][53][54][55][56] showing the 90% CL exclusion contours. The colored contours represent different values of B(h → 2γ D + X) in the range 0.1-40%.
focused on highly displaced objects and these searches therefore probe different regions of the available parameter space.

Summary
A search for pairs of new light bosons produced in the decay of a Higgs boson, that subsequently decay to pairs of oppositely charged muons, is presented. One event is observed in the signal region, with 2.2 ± 0.7 events expected from the SM backgrounds.
A model independent upper limit at 95% CL on the product of the cross section, branching fraction, and acceptance is obtained. This limit is valid for light boson masses in the range 2m μ < m a < 2m τ .
The obtained results allow a straightforward interpretation within a broad range of physics models that predict the same type of signature. The results are compared with two benchmark models in the context of the NMSSM and dark SUSY, including scenarios predicting a nonnegligible light boson lifetime.

Acknowledgments
We would like to thank J. Ruderman for his guidance with the theoretically motivated benchmark samples of dark SUSY. We congratulate our colleagues in the CERN accelerator departments for the excellent performance of the LHC and thank the technical and administrative staffs at CERN and at other CMS institutes for their contributions to the success of the CMS effort. In addition, we gratefully acknowledge the computing centres and personnel of the Worldwide LHC Computing Grid for delivering so effectively the computing infrastructure essential to our analyses.