Search for invisible decays of $\omega$ and $\phi$ with $J/\psi$ data at BESIII

Using a data sample of $(1310.6 \pm 7.0)\times 10^6$ $J/\psi$ events collected with the BESIII detector operating at the BEPCII collider, we perform the first experimental search for invisible decays of a light vector meson ($V=\omega,\phi$) via $J/\psi \to V \eta$ decays. The decay of $\eta \to \pi^+\pi^-\pi^0$ is utilized to tag the $V$ meson decaying into the invisible final state. No evidence for a significant invisible signal is observed, and the upper limits on the ratio of branching fractions at the 90\% confidence level are determined to be $\frac{\mathcal{B}(\omega \rightarrow \rm{invisible})}{\mathcal{B}(\omega \to \pi^+\pi^-\pi^0)}<8.1 \times 10^{-5}$ and $\frac{\mathcal{B}(\phi \to \rm{invisible})}{\mathcal{B}(\phi \to K^+K^-)}<3.4\times 10^{-4}$. By using the world average values of $\mathcal{B}(\omega \to \pi^+\pi^-\pi^0)$ and $\mathcal{B}(\phi\to K^+K^-)$, the upper limits on the decay branching fractions at the $90\%$ confidence level are set as $\mathcal{B}(\omega \to \rm{invisible})<7.3 \times 10^{-5}$ and $\mathcal{B}(\phi \to \rm{invisible})<1.7\times 10^{-4}$, respectively.


I. INTRODUCTION
Although there is strong evidence from many astrophysical observations for the existence of dark matter, its nature is still mysterious. Dark matter is invisible in the entire electromagnetic spectrum, and its existence is inferred via gravitational effects only. Any information about its interactions with a Standard Model (SM) particle would shed light on the nature of dark matter. Quarkonium states, whose constituents are a quark and its own anti-quark, are expected to annihilate into a neutrino-pair (νν) via a virtual Z 0 boson. However, the process is very rare in the SM [1]. The branching fraction of the invisible decays might be enhanced by several orders of magnitude in the presence of light dark matter (LDM) particles χ [2][3][4] as described in Refs. [5,6].
The LDM particles, which are in the kinematic reach of BESIII, may provide one possible explanation of the feature of the 511 keV gamma ray excess from the galactic center observed by the INTEGRAL satellite [7]. The smooth symmetric morphology of 511 keV gamma emission is believed to originate from the annihilation of LDM particles into e + e − pairs [2,8]. The LDM particles can have adequate relic abundance to account for the nonbaryonic dark matter [9] in the universe, if they couple with the SM particles via a new light gauge boson U [10] or the exchange of heavy fermions in the case of scalar dark matter [2,3]. One of the most popular LDM candidates is the neutralino predicted by the Next-to-Minimal Supersymmetric Standard Model [11], which is stable due to the conserved R-parity [12].
The BESII [13] and BaBar [14] experiments have set the most stringent upper limits on the invisible decays of J/ψ and Υ(1S), respectively, which are still above the SM predictions [1]. The experimental exploration of invisible decays for other quarkonium states (qq, q = u, d or s) may help to constrain the masses of the LDM particles and the coupling of the U boson to light quarks [15,16]. The branching fraction B(V → χχ) (V = ω, φ) is predicted to be up to the level of 10 −8 by assuming the same cross section for the time reversed processes, σ(qq → χχ) ≃ σ(χχ → qq) [6]. The search for these decays can be performed via a two-body decay process of J/ψ → V η. In this paper, we report the first experimental search for the invisible decays of ω and φ mesons via J/ψ → V η using (1310.6±7.0)×10 6 J/ψ events collected with the BESIII detector in 2009 and 2012 [17].

II. THE BESIII EXPERIMENT AND MONTE
CARLO SIMULATION BESIII is a cylindrical particle physics detector located at the BEPCII facility, a double-ring e + e − collider with a peak luminosity of 10 33 cm −2 s −1 at the center-of-mass (CM) energy of 3.773 GeV. It has four detector subcomponents with a coverage of 93% of the total solid angle as described in Ref. [18]. Charged particle momenta are measured in a 43-layer helium based main drift chamber (MDC) operating with a 1.0 T (0.9 T) solenoidal magnetic field during 2009 (2012) J/ψ runs. Charged particle identification (PID) is performed using the energy loss (dE/dx) measured in the MDC with a resolution better than 6%, and a time-of-flight (TOF) system consisting of 5 cm thick plastic scintillators with a time resolution of 80 ps in the barrel region and 110 ps in the end-cap region, respectively. Photon and electron energies are measured in a CsI(Tl) electromagnetic calorimeter (EMC). The energy (position) resolution of the EMC for 1 GeV electrons and photons is 2.5% (6 mm) in the barrel and 5.0% (9 mm) in the end-cap regions. The muons are identified in a muon counter (MUC) containing nine (eight) layers of resistive plate chamber counters interleaved with steel in the barrel (end-caps) region. The MUC provides a spatial resolution better than 2 cm.
A large number of Monte Carlo (MC) events are produced to optimize the event selection criteria, to study the potential backgrounds and to determine the recon-struction efficiencies. The MC simulation includes the detector response and signal digitization models simulated by Geant4 [19] and takes into account time-dependent detector effects, such as beam related backgrounds and detector running conditions during the data-taking period. An MC sample of 1225 × 10 6 inclusive J/ψ events is generated for background studies. The known J/ψ decay modes are generated by the EvtGen generator package [20] with the branching fractions taken from the Particle Data Group (PDG) [21], while the remaining unknown J/ψ decay modes are generated by the LUNDCHARM [22] generator. The production of the J/ψ resonance via e + e − annihilation is simulated by the KKMC [23] including the effects of the beam energy spread and initial state radiation (ISR). We use a helicity amplitude model for the J/ψ → V η decay, an ω Dalitz plot distribution model for the ω → π + π − π 0 decay [24], an η Dalitz plot distribution model for the η → π + π − π 0 decay [25], a vector meson decaying to a pair of scalar particles model for the φ → K + K − decay, and a phase space model for V → νν decays [20].

III. ANALYSIS STRATEGY
The search for invisible decays of ω and φ mesons is performed by using the two-body J/ψ → V η decay process. The candidate events are tagged with the η reconstructed from its π + π − π 0 decay mode, and the mass distribution of the system recoiling against the η candidate is used to investigate invisible decays of ω and φ mesons. The more prominent decay mode of η → γγ is not used for the tagging due to the huge background contamination.
In order to minimize the systematic uncertainty, the decays ω → π + π − π 0 and φ → K + K − from J/ψ → V η decays are reconstructed as reference channels. The ratio of the branching fraction of the invisible decay to that of the visible decay of V mesons is measured by where N invisible sig and N visible sig are the numbers of signal events for the invisible and visible decays, respectively, the ǫ invisible and ǫ visible are the corresponding detection efficiencies. By applying this method, the systematic uncertainties associated with the total number of J/ψ events, the branching fractions B(J/ψ → V η) and B(η → π + π − π 0 ), and the reconstruction of η candidates (such as tracking, PID and photon detection efficiency etc.) are canceled.

IV. DATA ANALYSIS
The charged tracks are measured in the MDC with the polar angle θ satisfying | cos θ| < 0.93. They must have the points of closest approach to the beam line within ±10.0 cm from the interaction point along the beam direction and ±1.0 cm in the plane perpendicular to the beam. PID for the charged tracks is accomplished by combining the measured energy loss (dE/dx) in the MDC, the flight time obtained from the TOF and the electromagnetic cluster shower information from the EMC to form the likelihoods for electron, kaon and pion hypotheses. A charged pion is identified by requiring the PID probability of its pion hypothesis to be larger than the kaon and electron hypotheses.
The photon candidates, reconstructed using the clusters of energy deposited in the EMC, are selected with a minimum energy of 25 MeV in the barrel region (| cos θ| < 0.8) or 50 MeV in the end-cap region (0.86 < | cos θ| < 0.92). To improve the reconstruction efficiency and energy resolution, the energy deposited in the nearby TOF counters is included. The angle between a photon and the nearest extrapolated track in the EMC is required to be greater than 10 degrees to avoid any overlap between charged and neutral tracks. In order to suppress electronic noise and energy deposits unrelated to the signal events, the EMC timing of the photon candidate is required to be within 700 ns relative to the event start time. A π 0 candidate is reconstructed from a photon pair candidate, and the two-photon invariant mass is constrained to the nominal value of the π 0 meson [21] by performing a kinematic fit.

A. The invisible decays of ω and φ mesons
For studies of the invisible decays of a V meson using the decay chain J/ψ → V η, η → π + π − π 0 , the event candidate is required to have two oppositely charged tracks identified as pions. A vertex fit is performed to these two charged tracks to insure that they originate from a common vertex. The π 0 candidate for which the π + π − π 0 invariant mass (M π + π − π 0 ) is closest to the nominal mass of the η meson [21] is considered as originating from the η decay. An η candidate is required to have M π + π − π 0 within [0.52,0.57] GeV/c 2 . With the above requirements, MC studies indicate that the dominant backgrounds are from J/ψ → V η with the V meson decaying into purely neutral final states, such as ω → γπ 0 and φ → K S K L , K S → π 0 π 0 . Thus, E Extra γ is required to be less than 0.2 GeV, where E Extra γ is the sum of energies of the extra photons, which are not used in the η reconstruction. Furthermore, the polar angle of the system recoiling against the selected η candidate, θ recoil , is required to satisfy | cos θ recoil | < 0.7 to further eliminate the background contributions from J/ψ → Xη, where X can be any final state emitted in the region which is not covered by the acceptance of the detector.
The signals of the invisible decays of ω and φ mesons are inferred from the invariant mass of the system recoiling against the selected η candidate, defined as where E cm is the CM energy, and E π + π − π 0 and P π + π − π 0 are the energy and momentum of the π + π − π 0 system in the CM frame, respectively. The M V recoil distribution of the event candidates for the data range [0.40, 1.35] GeV/c 2 is shown in Fig. 1. The expected distributions for ω and φ invisible decay signals by MC simulation are also depicted in the plot. Detailed studies of the inclusive J/ψ decay sample indicate that the non-peaking backgrounds are dominated by processes with non-η mesons in the final state, which can be evaluated with the normalized events in the η mass sideband regions, as shown by a cyan histogram in Fig. 1. The non-peaking background from J/ψ → γη, which has a large branching fraction, is evaluated to be 1.8 events with negligible uncertainties by using an exclusive MC sample normalized according to the branching fractions quoted in the PDG [21], and is ignored in the following analysis. The possible peaking background is from the decay J/ψ → V η with the V meson decaying visibly. The numbers of peaking backgrounds are evaluated to be 0.1 for J/ψ → ωη and 2.0 for J/ψ → φη with negligible uncertainty using the simulated MC samples normalized according to the measured branching fractions of J/ψ → V η described in Sec. IV B and IV C, respectively, and the corresponding distributions are presented in Fig. 1. The backgrounds from other sources are negligible. The M V recoil distributions of simulated signal MC events for invisible decays of ω and φ mesons are observed to be well consistent with the data and MC simulations of their visible decays described in Sec. IV B and IV C, respectively.  An extended maximum likelihood (ML) fit to the M V recoil distribution is performed to obtain the signal yield (N sig ). The probability density function (PDF) of the V meson invisible decay signal and peaking background is described by their MC simulated shapes, while that of the non-peaking background is represented by an increasing exponential function. In the fit, the number of peaking background events is fixed, while the parameters of the non-peaking background PDF and the yields for signal and non-peaking background events are free parameters in the fit. The ML fit yields N sig = 1.4 ± 3.6 events for the ω → invisible decay and N sig = −0.6 ± 4.5 for the φ → invisible decay, respectively. The obtained N sig events for both decay modes are consistent with zero, and no evidence of invisible decays of ω and φ mesons is observed. The fitted M V recoil are shown in Fig. 2. The corresponding signal detection efficiencies, estimated with the MC simulation, are 20.5% and 21.3% for ω and φ invisible decays, respectively. B. The visible decay mode ω → π + π − π 0 The candidate events of J/ψ → ωη with subsequent decays ω → π + π − π 0 and η → π + π − π 0 are required to have four charged tracks with net charge zero and at least two independent π 0 candidates without sharing the same photon. The four charged tracks are assumed to be pions and required to originate from a common vertex by performing a vertex fit. For an event with multiple π 0 π 0 pair candidates, the one with the least value of p tot is selected, where p tot is the total momentum of the 2(π + π − π 0 ) candidates. The total energy (E tot ) of the selected candidate is also required to satisfy E tot > 2.95 GeV. For a selected 2(π + π − π 0 ) final state, the combinations of π + π − π 0 for ω and η signals are determined by where M X π + π − π 0 (X = ω, η) is the invariant mass of the π + π − π 0 combination for the X candidate, M X is the nominal X meson mass quoted in the PDG [21], and σ X is the corresponding mass resolution determined from the signal MC simulation. All eight combinations of (π + π − π 0 ) ω versus (π + π − π 0 ) η are explored, and the one with the least χ 2 ωη is selected. In order to improve the purity of pions in the η → π + π − π 0 decay and to minimize the systematic uncertainty in the analysis, PID for charged pions from the η decay is performed, but no PID requirement for those from ω decay due to tiny expected background contribution from ω → l + l − π 0 (l = e, µ) in the full J/ψ data sample. Similarly to the invisible decay, the polar angle of the system recoiling against the η candidate θ recoil is required to satisfy | cos θ recoil | < 0.7 to minimize the systematic uncertainty. The selected candidate events are further required to have M ω π + π − π 0 and M η π + π − π 0 in the ranges [0.65, 0.98] and [0.41, 0.65] GeV/c 2 , respectively. Figure 3 shows the two-dimensional (2D) histogram of M ω π + π − π 0 versus M η π + π − π 0 for data.
The remaining backgrounds are dominated by those with the same final state as the signal, but neither ω nor η intermediate states included (named BKGI thereafter) or without either ω or η intermediate state (named BKGII thereafter). In addition, there is a small peaking background for both M ω π + π − π 0 and M η π + π − π 0 simultaneously (named BKGIII thereafter), dominated by J/ψ → ωη with the subsequent decays ω → π + π − π 0 and η → γπ + π − . Consequently, the contributions of BKGI and BKGII are determined by performing a 2D ML fit to M ω π + π − π 0 and M η π + π − π 0 , while BKGIII is determined by using a corresponding exclusive MC sample normalized according to the branching fractions quoted in PDG [21]. The BKGIII yield, estimated to be 1085.8 ± 126.6 events, which uncertainty includes the uncertainties of both the total number of J/ψ events and the branching fractions of the corresponding decay process, is subtracted from the signal yield obtained from the 2D ML fit, eventually. Detailed MC studies indicate that the distributions of M ω π + π − π 0 and M η π + π − π 0 are uncorrelated as no kinematic fit is performed. Therefore, in the 2D ML fit, a 2D PDF can be the direct product of the two one-dimensional (1D) PDFs for the two variables. Furthermore, MC studies validate that the 1D line shapes for the ω (η) meson are identical between the signal and the peaking background (BKGII). Consequently, the full 2D PDF used in the ML fit is constructed by where N data sig is the number of signal events including the contribution from BKGIII, N bkg , N ω bkg and N η bkg are the numbers of events for BKGI, and BKGII with ω and η intermediate states, respectively. F ω sig and F η sig are the PDFs of ω and η signals in the π + π − π 0 invariant mass spectrum, respectively, described by the sum of two crystal ball (CB) functions [26] with common mean and sigma values, but opposite side and different parameters for tails. F ω bkg represents the non-ω component in the M ω π + π − π 0 distribution, and is described by a second order Chebyshev polynomial function. F η bkg is the non-η component in the M η π + π − π 0 distribution, and is represented by a reversed ARGUS function [27], defined as The projections of the ML fit to the M ω π + π − π 0 and M η π + π − π 0 distributions are shown in Fig. 4. The fit yields N data sig = 32528 ± 283. After subtracting the contribution of BKGIII, the net number of signal events is N sig = 31442 ± 314. By taking into account the signal yield N sig , the detection efficiency 6.2% obtained from the corresponding MC sample, and the decay branching fractions of ω → π + π − π 0 and η → π + π − π 0 quoted in PDG [21], the branching fraction of J/ψ → ωη is measured to be larger by 12.0% with respect to its world average value quoted in the PDG [21], but consistent within the uncertainty.
C. The visible decay mode of φ → K + K − For the candidate events of J/ψ → φη with subsequent decays φ → K + K − and η → π + π − π 0 , the η candidate is reconstructed with exactly same process as described in Sec. IV A, and the φ candidate is reconstructed by two additional oppositely charged tracks, which are assumed to be kaons without any PID requirement. The total energy of the selected K + K − π + π − π 0 final state must satisfy E tot > 2.95 GeV. Similarly, the polar angle of the system recoiling against the η candidate θ recoil is required to satisfy | cos θ recoil | < 0.7 to minimize the systematic uncertainty. The candidate events with invariant mass of K + K − in the range [0.987, 1.10] GeV/c 2 are kept for further studies ( Figure 5). The remaining backgrounds are analogous to BKGI, BKGII and BKGIII in the J/ψ → ωη visible decay obtained by replacing ω with φ signal, and the corresponding π + π − π 0 with K + K − . Similarly, the contributions of BKGI and BKGII are determined by a 2D ML fit, and the BKGIII of 238.6 ± 26.0 events, estimated with an exclusive MC sample of J/ψ → φη with subsequent decays φ → K + K − and η → γπ + π − and normalized according to the branching fractions quoted in the PDG [21], is subtracted from the signal yield obtained from the 2D ML fit.
is carried out to obtain the signal yield. The parameterizations of the 1D PDF for η and non-η components in M η π + π − π 0 are the same as those used in the case of the J/ψ → ωη visible decay. The 1D PDF for the φ signal in the M φ K + K − distribution is described by a relativistic Breit-Wigner (BW) [21] function convolved with a Gaussian function representing the mass resolution, where the mass and width of the φ signal are taken from the PDG [21], and the parameters of the Gaussian function are left free in the fit. For the non-φ components in M φ K + K − , its 1D PDF is represented by a reversed ARGUS function as described in Sec. IV B by fixing the threshold parameter t to the upper limit of the fit range. The ML fit yields N data sig = 19534 ± 186, and the M φ K + K − and M η π + π − π 0 projections are shown in Fig. 6. The net number of signal events after subtracting the contribution of BKGIII from N data sig is N sig = 19295 ± 188. We compute the branching fraction of J/ψ → φη by taking into account N sig , the detection efficiency 15.8% from the corresponding MC sample, and the decay branching fractions of φ → K + K − and η → π + π − π 0 from the PDG [21]. The measured branching fraction of J/ψ → φη is also larger by 12.0% over its world average value [21], but consistent within the uncertainty.

V. SYSTEMATIC UNCERTAINTY
Since we measure the relative ratios of the branching fractions of invisible decay to that of corresponding visible decay, the systematic uncertainties associated with the number of J/ψ events, the reconstruction efficiency of η → π + π − π 0 , the requirement on cos θ recoil , the branching fractions of η → π + π − π 0 and J/ψ → V η cancel. The remaining sources of systematic uncertainties are associated with the fit procedure of the invisible and visible decays, the E Extra γ requirement in the invisible decay, charged track reconstruction, trigger efficiency, photon detection and the E tot requirement for the visible decay and the branching fractions of ω → π + π − π 0 and φ → K + K − decays. The details of the evaluation of individual uncertainties are described below and summarized in Table I.
The systematic uncertainty associated with the fit procedure in the invisible decays, which can reduce the significance of any observation, but does not scale with the reconstructed signal yields, is considered to be an additive systematic uncertainty. The other remaining sources of systematic uncertainties, which don't not affect the significance of any observation, but scale with the number of reconstructed signal yield, are considered multiplicative systematic uncertainties.
The systematic uncertainties associated with the tracking efficiency of kaon and pion are 1.0% for each track, obtained by investigating the control samples of J/ψ → K 0 S K ± π ∓ and J/ψ → π + π − pp, respectively. The systematic uncertainty of the photon reconstruction efficiency is less than 1.0% per photon, investigated with a control sample sample of e + e − → γµ + µ − in which the four-momenta of two muons are used to obtain the ISR photon momentum [28]. The uncertainty associated with the E Extra γ requirement in the invisible decay is determined to be 1.1% by comparing the corresponding detection efficiencies between data and MC simulation with a control sample of J/ψ → π + π − π 0 . The systematic uncertainty due to π 0 reconstruction efficiency is determined to be 1.0% using a control sample of J/ψ → ppπ 0 . The uncertainty associated with the E tot requirement in the visible decay processes is explored with the relative efficiency with respect to an alternative requirement E tot > 2.6 GeV, where the signal loss is expected to be negligible. The relative differences in efficiency between data and MC simulation, 2.1% and 1.0% for ω and φ visible decays, respectively, are considered as the uncertainties.
The BESIII trigger system combines the information from the sub-detectors of EMC, MDC and TOF to select the events of interest for readout. We study the trigger efficiency with a control sample of J/ψ → π + π − π 0 , and found the efficiency is almost 100% for an event with two charged and two photons by considering the different kinematics of the final state. We assign 0.1% as a systematic uncertainty related with the trigger efficiency.
The uncertainty associated with the fit procedure in invisible decays originates from the signal PDF, nonpeaking background modelling and the fit bias. The uncertainty due to the signal PDF is estimated by an alternative fit with the sum of the two CB functions for the signal PDF, where the corresponding parameters of the CB functions are obtained by fitting the simulated MC samples and fixed in the fit. The uncertainty due to the non-peaking background shape is estimated by using an alternative PDF of a second order Chebyshev polynomial function in the fit. The relative changes in the results are taken as the uncertainties. A large number of pseudo-experiments with fixed amount of signal, peaking and non-peaking background events expected from the data are generated to examine the bias of the fit procedure. The same fit procedure is repeated for each MC set, and the average shift of resultant signal yields can be taken as systematic uncertainty, and is found to be negligible.
The uncertainty related to the fit procedure for the reference decay is obtained in an analogous way. The uncertainty due to the fixed PDF parameters of the BW in the φ → K + K − decay is evaluated by varying each parameter within its statistical uncertainty while taking the correlations between the parameters into account. The uncertainty associated with the PDFs for the non-ω component on the M ω π + π − π 0 distribution is estimated by changing the order of the Chebyshev polynomial function. The uncertainties associated with the PDFs for the non-η component in M η π + π − π 0 and the non-φ component in M φ K + K − are estimated by modifying the formula of Eq. (4) as F η bkg (m) = m 2 /t · (1 − (X − m) 2 /t 2 ) a · exp(−b · (1 − (X − m) 2 /t 2 )). The largest relative change of the signal yields of individual alternative fits is considered as systematic uncertainty, and the total uncertainty associ- ated with the fit procedure is the quadrature sum of the individual values. The uncertainties associated with the branching fractions of ω or φ visible decays are taken from the PDG [21], and included only in the results for the branching fractions.

VI. RESULTS
No obvious signal for ω and φ invisible decays is observed. We compute the upper limits on the ratio of branching fractions of the invisible decay to that of the corresponding visible decay, B(ω→invisible) B(ω→π + π − π 0 ) and B(φ→invisible) B(φ→K + K − ) , at the 90% confidence level (C.L.) using the Bayesian approach [21], individually. The branching fraction ratios of B(ω→invisible) B(ω→π + π − π 0 ) and B(φ→invisible) B(φ→K + K − ) are calculated using the formula of Eq. (1) after incorporating obtained signal yields and the corresponding detection efficiencies for the visible and invisible decays as presented above. The systematic uncertainty is included by convolving the likelihood versus the branching fraction ratio curve with a Gaussian function with a width equal to the systematic uncertainty. The upper limits on the branching fraction ratios are measured to be B(ω→invisible) B(ω→π + π − π 0 ) < 8.1×10 −5 and B(φ→invisible) B(φ→K + K − ) < 3.4×10 −4 for ω and φ mesons, respectively, at the 90% C.L. after integrating their likelihood versus branching fraction ratio curves from zero to 90% of the total curve. By using the branching fractions of ω → π + π − π 0 and φ → K + K − quoted in the PDG [21], the upper limits on the invisible decay branching fractions at the 90% C.L. are calculated to be B(ω → invisible) < 7.3 × 10 −5 and B(φ → invisible) < 1.7 × 10 −4 , individually.

VII. SUMMARY
Using a data sample of (1310.6 ± 7.0) × 10 6 J/ψ events collected by the BESIII experiment at the BEPCII collider, a search for the invisible decays of ω and φ mesons in J/ψ → V η decays is performed for the first time. We find no significant signal for these invisible decays and set 90% C.L. upper limits on the ratio of branching fractions of invisible decays to that of the corresponding visible decays to be B(ω→invisible) B(ω→π + π − π 0 ) < 8.1 × 10 −5 and B(φ→invisible) B(φ→K + K − ) < 3.4 × 10 −4 , respectively. The upper limits on the branching fractions B(ω → invisible) and B(φ → invisible) are also determined to be less than 7.3 × 10 −5 and 1.7 × 10 −4 , respectively, at the 90% C.L. by using B(ω → π + π − π 0 ) and B(φ → K + K − ) from the PDG [21]. These results can provide a complementary information to study the nature of dark matter and constrain the parameters of phenomenological models [15,16].