Amplitude Analysis and Branching Fraction Measurement of $D^0\rightarrow K^-\pi^+\pi^0\pi^0$

TL;DR

This paper presents the first amplitude analysis and precise branching fraction measurement of the decay D^0→K^-π^+π^0π^0 using BESIII data, revealing sub-structures and intermediate process contributions.

Contribution

It introduces the first amplitude analysis of this decay mode and provides an accurate measurement of its branching fraction, enhancing understanding of D^0 decay dynamics.

Findings

Branching fraction measured as (8.86 ± 0.13(stat) ± 0.19(syst))%

Identified sub-structures and intermediate processes in the decay

Provided relative fractions and phases among decay components

Abstract

Utilizing the dataset corresponding to an integrated luminosity of $2.93$ fb$^{-1}$ at $\sqrt{s}=3.773$ GeV collected by the BESIII detector, we report the first amplitude analysis and branching fraction measurement of the $D^0\rightarrow K^-\pi^+\pi^0\pi^0$ decay. We investigate the sub-structures and determine the relative fractions and the phases among the different intermediate processes. Our results are used to provide an accurate detection efficiency and allow measurement of ${\cal B}(D^0\rightarrow K^-\pi^+\pi^0\pi^0) \,=\, (8.86 \pm 0.13(\text{stat}) \pm 0.19(\text{syst}))\%$.

Figures

Tables

Table 1. Table 1: Spin factor for each decay chain. All operators, i.e.  t ~ ~ 𝑡 \tilde{t} , have the same definitions as Ref.  Zou2003 . Scalar, pseudo-scalar, vector, axial-vector, and tensor states are denoted by S 𝑆 S , P 𝑃 P , V 𝑉 V , A 𝐴 A , and T 𝑇 T , respectively.

Decay chain	$S\left(p\right)$
$D\left[S\right]\to V_1{V}_2$	${\tilde{t}}^{\left(1\right)\mu }\left(V_1\right){\tilde{t}}_{\mu }^{\left(1\right)}\left(V_2\right)$
$D\left[P\right]\to V_1{V}_2$	${ϵ}_{\mu \nu \lambda \sigma }{p}^{\mu }\left(D\right){\tilde{T}}^{\left(1\right)\nu }\left(D\right){\tilde{t}}^{\left(1\right)\lambda }\left(V_1\right){\tilde{t}}^{\left(1\right)\sigma }\left(V_2\right)$
$D\left[D\right]\to V_1{V}_2$	${\tilde{T}}^{\left(2\right)\mu \nu }\left(D\right){\tilde{t}}_{\mu }^{\left(1\right)}\left(V_1\right){\tilde{t}}_{\nu }^{\left(1\right)}\left(V_2\right)$
$D\to A{P}_1,A\left[S\right]\to V{P}_2$	${\tilde{T}}^{\left(1\right)\mu }\left(D\right)P_{\mu \nu }^{\left(1\right)}\left(A\right){\tilde{t}}^{\left(1\right)\nu }\left(V\right)$
$D\to A{P}_1,A\left[D\right]\to V{P}_2$	${\tilde{T}}^{\left(1\right)\mu }\left(D\right){\tilde{t}}_{\mu \nu }^{\left(2\right)}\left(A\right){\tilde{t}}^{\left(1\right)\nu }\left(V\right)$
$D\to A{P}_1,A\to S{P}_2$	${\tilde{T}}^{\left(1\right)\mu }\left(D\right){\tilde{t}}_{\mu }^{\left(1\right)}\left(A\right)$
$D\to VS$	${\tilde{T}}^{\left(1\right)\mu }\left(D\right){\tilde{t}}_{\mu }^{\left(1\right)}\left(V\right)$
$D\to V_1{P}_1,V_1\to V_2{P}_2$	${ϵ}_{\mu \nu \lambda \sigma }{p}_{V_1}^{\mu }{r}_{V_1}^{\nu }{p}_{P_1}^{\lambda }{r}_{V_2}^{\sigma }$
$D\to P{P}_1,P\to V{P}_2$	$p^{\mu }\left(P_2\right){\tilde{t}}_{\mu }^{\left(1\right)}\left(V\right)$
$D\to TS$	${\tilde{T}}^{\left(2\right)\mu \nu }\left(D\right){\tilde{t}}_{\mu \nu }^{\left(2\right)}\left(T\right)$

Table 2. Table 2: Parameters of K ​ π 𝐾 𝜋 K\pi S 𝑆 S -wave, by BABAR  PhysRevD.78.034023 , where the uncertainties include the statistical and systematic uncertainties.

$M$(GeV/$c^2$)	$1.463\pm 0.002$
$\Gamma$(GeV)	$0.233\pm 0.005$
$F$	$0.80\pm 0.09$
${\varphi }_F$	$2.33\pm 0.13$
$R$	1(fixed)
${\varphi }_R$	$-5.31\pm 0.04$
$a$	$1.07\pm 0.11$
$r$	$-1.8\pm 0.3$

Table 3. Table 3: FFs, phases, and significances of the optimal set of amplitude modes. The first and second uncertainties are statistical and systematic, respectively. The details of systematic uncertainties are discussed in Section  VI.1 .

Amplitude mode	FF $(\%)$	Phase $\left(\varphi \right)$	Significance $\left(\sigma \right)$
$D\to SS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_{S\text{-wave}}{\left({\pi }^0{\pi }^0\right)}_S$	$6.92\pm 1.44\pm 2.86$	$-0.75\pm 0.15\pm 0.47$	$>10$
$D\to {\left(K^{-}{\pi }^0\right)}_{S\text{-wave}}{\left({\pi }^{+}{\pi }^0\right)}_S$	$4.18\pm 1.02\pm 1.77$	$-2.90\pm 0.19\pm 0.47$	$6.0$
$D\to AP,A\to VP$
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[S\right]$	$28.36\pm 2.50\pm 3.53$	$0$ (fixed)	$>10$
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[D\right]$	$0.68\pm 0.29\pm 0.30$	$-2.05\pm 0.17\pm 0.25$	$6.1$
$D\to K_1{\left(1270\right)}^{-}{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	$0.15\pm 0.09\pm 0.15$	$1.84\pm 0.34\pm 0.43$	$4.9$
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	$0.39\pm 0.18\pm 0.30$	$-1.55\pm 0.20\pm 0.26$	$4.8$
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	$0.11\pm 0.11\pm 0.11$	$-1.35\pm 0.43\pm 0.48$	$4.0$
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{-}{\rho }^{+}\left[S\right]$	$2.71\pm 0.38\pm 0.29$	$-2.07\pm 0.09\pm 0.20$	$>10$
$D\to {\left(K^{\ast -}{\pi }^0\right)}_A{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	$1.85\pm 0.62\pm 1.11$	$1.93\pm 0.10\pm 0.15$	$7.8$
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	$3.13\pm 0.45\pm 0.58$	$0.44\pm 0.12\pm 0.21$	$>10$
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	$0.46\pm 0.17\pm 0.29$	$-1.84\pm 0.26\pm 0.42$	$5.9$
$D\to {\left({\rho }^{+}{K}^{-}\right)}_A{\pi }^0,K^{-}{\rho }^{+}\left[D\right]$	$0.75\pm 0.40\pm 0.60$	$0.64\pm 0.36\pm 0.53$	$5.1$
$D\to AP,A\to SP$
$D\to {\left({\left(K^{-}{\pi }^{+}\right)}_{S\text{-wave}}{\pi }^0\right)}_A{\pi }^0$	$1.99\pm 1.08\pm 1.55$	$-0.02\pm 0.25\pm 0.53$	$7.0$
$D\to VS$
$D\to {\left(K^{-}{\pi }^0\right)}_{S\text{-wave}}{\rho }^{+}$	$14.63\pm 1.70\pm 2.41$	$-2.39\pm 0.11\pm 0.35$	$>10$
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_S$	$0.80\pm 0.38\pm 0.26$	$1.59\pm 0.19\pm 0.24$	$4.1$
$D\to K^{\ast 0}{\left({\pi }^0{\pi }^0\right)}_S$	$0.12\pm 0.12\pm 0.12$	$1.45\pm 0.48\pm 0.51$	$4.1$
$D\to VP,V\to VP$
$D\to {\left(K^{\ast -}{\pi }^{+}\right)}_V{\pi }^0$	$2.25\pm 0.43\pm 0.45$	$0.52\pm 0.12\pm 0.17$	$>10$
$D\to VV$
$D\to K^{\ast -}{\rho }^{+}\left[S\right]$	$5.15\pm 0.75\pm 1.28$	$1.24\pm 0.11\pm 0.23$	$>10$
$D\to K^{\ast -}{\rho }^{+}\left[P\right]$	$3.25\pm 0.55\pm 0.41$	$-2.89\pm 0.10\pm 0.18$	$>10$
$D\to K^{\ast -}{\rho }^{+}\left[D\right]$	$10.90\pm 1.53\pm 2.36$	$2.41\pm 0.08\pm 0.16$	$>10$
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[P\right]$	$0.36\pm 0.19\pm 0.27$	$-0.94\pm 0.19\pm 0.28$	$5.7$
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[D\right]$	$2.13\pm 0.56\pm 0.92$	$-1.93\pm 0.22\pm 0.25$	$>10$
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_V\left[D\right]$	$1.66\pm 0.52\pm 0.61$	$-1.17\pm 0.20\pm 0.39$	$7.6$
$D\to {\left(K^{-}{\pi }^0\right)}_V{\left({\pi }^{+}{\pi }^0\right)}_V\left[S\right]$	$5.17\pm 1.91\pm 1.82$	$-1.74\pm 0.20\pm 0.31$	$7.6$
$D\to TS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_{S\text{-wave}}{\left({\pi }^0{\pi }^0\right)}_T$	$0.30\pm 0.21\pm 0.30$	$-2.93\pm 0.31\pm 0.82$	$5.8$
$D\to {\left(K^{-}{\pi }^0\right)}_{S\text{-wave}}{\left({\pi }^{+}{\pi }^0\right)}_T$	$0.14\pm 0.12\pm 0.10$	$2.23\pm 0.38\pm 0.65$	$4.0$
TOTAL	$98.54$

Table 4. Table 4: FF systematic uncertainties (in units of statistical standard deviations) for: (I) the amplitude model, (II) background, (III) experimental effects, and (IV) fit bias. The total uncertainty is obtained by adding all contributions in quadrature.

Amplitude mode	I	II	III	IV	Total
$D\to SS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_S{\left({\pi }^0{\pi }^0\right)}_S$	1.518	1.258	0.072	0.235	1.987
$D\to {\left(K^{-}{\pi }^0\right)}_S{\left({\pi }^{+}{\pi }^0\right)}_S$	1.524	0.835	0.078	0.004	1.740
$D\to AP,A\to VP$
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[S\right]$	1.293	0.436	0.030	0.363	1.412
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[D\right]$	0.938	0.368	0.024	0.284	1.046
$D\to K_1{\left(1270\right)}^{-}{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	1.643	1.175	0.160	0.182	2.035
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	1.562	0.567	0.034	0.036	1.662
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	0.989	0.541	0.035	0.068	1.201
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{-}{\rho }^{+}\left[S\right]$	0.713	0.221	0.098	0.172	0.772
$D\to {\left(K^{\ast -}{\pi }^0\right)}_A{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	1.253	1.254	0.076	0.237	1.790
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	1.145	0.524	0.022	0.162	1.278
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	0.865	1.468	0.052	0.106	1.708
$D\to {\left({\rho }^{+}{K}^{-}\right)}_A{\pi }^0,K^{-}{\rho }^{+}\left[D\right]$	1.249	0.812	0.084	0.186	1.504
$D\to AP,A\to SP$
$D\to {\left({\left(K^{-}{\pi }^{+}\right)}_S{\pi }^0\right)}_A{\pi }^0$	1.377	0.372	0.102	0.164	1.439
$D\to VS$
$D\to {\left(K^{-}{\pi }^0\right)}_S{\rho }^{+}$	1.308	0.252	0.070	0.476	1.416
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_S$	0.381	0.549	0.023	0.166	0.689
$D\to K^{\ast 0}{\left({\pi }^0{\pi }^0\right)}_S$	0.880	0.417	0.078	0.232	1.005
$D\to VP,V\to VP$
$D\to {\left(K^{\ast -}{\pi }^{+}\right)}_V{\pi }^0$	0.688	0.752	0.033	0.273	1.056
$D\to VV$
$D\to K^{\ast -}{\rho }^{+}\left[S\right]$	0.980	1.354	0.059	0.371	1.713
$D\to K^{\ast -}{\rho }^{+}\left[P\right]$	0.425	0.506	0.031	0.348	0.747
$D\to K^{\ast -}{\rho }^{+}\left[D\right]$	1.365	0.598	0.049	0.398	1.543
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[P\right]$	0.695	1.223	0.027	0.140	1.414
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[D\right]$	1.335	0.848	0.237	0.401	1.649
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_V\left[D\right]$	0.751	0.894	0.049	0.074	1.171
$D\to {\left(K^{-}{\pi }^0\right)}_V{\left({\pi }^{+}{\pi }^0\right)}_V\left[S\right]$	0.818	0.443	0.046	0.211	0.955
$D\to TS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_S{\left({\pi }^0{\pi }^0\right)}_T$	1.171	0.936	0.084	0.273	1.528
$D\to {\left(K^{-}{\pi }^0\right)}_S{\left({\pi }^{+}{\pi }^0\right)}_T$	0.803	0.188	0.068	0.018	0.828

Table 5. Table 5: Phase, ϕ italic-ϕ \phi , systematic uncertainties (in units of statistical standard deviations) for: (I) the amplitude model, (II) background, (III) experimental effects, and (IV) fit bias. The total uncertainty is obtained by adding all contributions in quadrature.

Amplitude mode	I	II	III	IV	Total
$D\to SS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_S{\left({\pi }^0{\pi }^0\right)}_S$	3.137	0.093	0.043	0.030	3.139
$D\to {\left(K^{-}{\pi }^0\right)}_S{\left({\pi }^{+}{\pi }^0\right)}_S$	2.330	0.850	0.044	0.109	2.483
$D\to AP,A\to VP$
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[S\right]$	0.000	0.000	0.000	0.000	0.000
$D\to K^{-}{a}_1{\left(1260\right)}^{+},{\rho }^{+}{\pi }^0\left[D\right]$	1.194	0.761	0.081	0.479	1.497
$D\to K_1{\left(1270\right)}^{-}{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	0.953	0.820	0.054	0.124	1.264
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	1.051	0.556	0.029	0.565	1.316
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	1.002	0.483	0.045	0.121	1.120
$D\to K_1{\left(1270\right)}^0{\pi }^0,K^{-}{\rho }^{+}\left[S\right]$	2.007	0.188	0.079	0.847	2.188
$D\to {\left(K^{\ast -}{\pi }^0\right)}_A{\pi }^{+},K^{\ast -}{\pi }^0\left[S\right]$	1.208	0.706	0.048	0.455	1.472
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[S\right]$	1.711	0.365	0.053	0.214	1.750
$D\to {\left(K^{\ast 0}{\pi }^0\right)}_A{\pi }^0,K^{\ast 0}{\pi }^0\left[D\right]$	1.501	0.605	0.051	0.187	1.630
$D\to {\left({\rho }^{+}{K}^{-}\right)}_A{\pi }^0,K^{-}{\rho }^{+}\left[D\right]$	1.195	0.613	0.133	0.611	1.482
$D\to AP,A\to SP$
$D\to {\left({\left(K^{-}{\pi }^{+}\right)}_S{\pi }^0\right)}_A{\pi }^0$	2.039	0.410	0.045	0.446	2.127
$D\to VS$
$D\to {\left(K^{-}{\pi }^0\right)}_S{\rho }^{+}$	3.159	0.471	0.053	0.216	3.201
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_S$	1.207	0.258	0.045	0.156	1.245
$D\to K^{\ast 0}{\left({\pi }^0{\pi }^0\right)}_S$	0.938	0.476	0.062	0.116	1.060
$D\to VP,V\to VP$
$D\to {\left(K^{\ast -}{\pi }^{+}\right)}_V{\pi }^0$	1.260	0.471	0.032	0.490	1.432
$D\to VV$
$D\to K^{\ast -}{\rho }^{+}\left[S\right]$	1.995	0.154	0.070	0.712	2.125
$D\to K^{\ast -}{\rho }^{+}\left[P\right]$	1.612	0.214	0.035	0.864	1.842
$D\to K^{\ast -}{\rho }^{+}\left[D\right]$	1.586	1.108	0.051	0.588	2.022
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[P\right]$	1.429	0.324	0.023	0.128	1.471
$D\to {\left(K^{-}{\pi }^0\right)}_V{\rho }^{+}\left[D\right]$	0.401	0.832	0.133	0.666	1.146
$D\to K^{\ast -}{\left({\pi }^{+}{\pi }^0\right)}_V\left[D\right]$	1.445	1.313	0.040	0.190	1.962
$D\to {\left(K^{-}{\pi }^0\right)}_V{\left({\pi }^{+}{\pi }^0\right)}_V\left[S\right]$	1.354	0.213	0.041	0.726	1.551
$D\to TS$
$D\to {\left(K^{-}{\pi }^{+}\right)}_S{\left({\pi }^0{\pi }^0\right)}_T$	2.544	0.724	0.057	0.189	2.653
$D\to {\left(K^{-}{\pi }^0\right)}_S{\left({\pi }^{+}{\pi }^0\right)}_T$	1.533	0.718	0.050	0.135	1.699

Table 6. Table 6: D 0 → K − ​ π + ​ π 0 ​ π 0 → superscript 𝐷 0 superscript 𝐾 superscript 𝜋 superscript 𝜋 0 superscript 𝜋 0 D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0} BF systematic uncertainties. The total uncertainty is obtained by adding all contributions in quadrature.

Source	Systematic (%)
Tracking efficiency	0.8
PID efficiency	0.4
${\pi }^0$ efficiency	1.2
Decay model	0.5
Yield fits (ST)	0.6
Yield fits (DT)	1.2
Peaking background	0.2
DCS decay correction	0.6
Total	2.3

Equations

$$\begin{array}{rl}{M}_{\text{BC}}& =\sqrt{E_{\text{beam}}^2-|{\vec{p}}_D{|}^2},\\ \Delta E& =E_D-E_{\text{beam}},\end{array}$$

$$\begin{array}{r}S(a,p)=\frac{ϵ(p)|A(a,p){|}^2{R}_4(p)}{\int ϵ(p)|A(a,p){|}^2{R}_4(p)dp},\end{array}$$

$$\begin{array}{r}{R}_4(p)dp=(2\pi {)}^4{\delta }^4(p_D-\sum _{\alpha }^4{p}_{\alpha })\prod _{\alpha }^4\frac{d^3{p}_{\alpha }}{(2\pi {)}^{3}2E_{\alpha }},\end{array}$$

$$\begin{array}{r}A(a,p)=\sum _i{a}_i{A}_i(p),\end{array}$$

$$\begin{array}{r}{A}_i(p)=P_i^1(p)P_i^2(p)S_i(p)F_i^1(p)F_i^2(p)F_i^D(p),\end{array}$$

$$\begin{array}{r}L=\prod _{k=1}^{N_s}S(a,p^k),\end{array}$$

$$\begin{array}{rl}\ln L& =\sum _{k=1}^{N_s}\ln S(a,p^k)\\ & =\sum _{k=1}^{N_s}\ln \left(\frac{|A(a,p^k){|}^2}{\int ϵ(p)|A(a,p){|}^2{R}_4(p)dp}\right)\\ & +\sum _{k=1}^{N_s}\ln R_4(p^k)+\sum _{k=1}^{N_s}\ln ϵ(p^k).\end{array}$$

$$\begin{array}{rl}& \int ϵ(p)|A(a,p){|}^2{R}_4(p)dp\\ & \approx \frac{1}{N_{\text{g,ph}}}\sum _{l=1}^{N_{\text{s,ph}}}|A(a,p^l){|}^2,\end{array}$$

$$\begin{array}{rl}& \int ϵ(p)|A(a,p){|}^2{R}_4(p)dp\\ & \approx \frac{1}{N_{\text{MC}}}\sum _{l=1}^{N_{\text{MC}}}\frac{|A(a,p^l){|}^2}{|A(a^{\text{gen}},p^l){|}^2},\end{array}$$

$$\begin{array}{r}{\gamma }_{ϵ}(p)=\prod _j\frac{{ϵ}_{j,\text{data}}(p)}{{ϵ}_{j,\text{MC}}(p)},\end{array}$$

$$\begin{array}{rl}& \int ϵ(p)|A(a,p){|}^2{R}_4(p)dp\\ & \approx \frac{1}{N_{\text{MC}}}\sum _l^{N_{\text{MC}}}\frac{|A(a,p^l){|}^2{\gamma }_{ϵ}(p^l)}{|A(a^{\text{gen}},p^l){|}^2}.\end{array}$$

$$\begin{array}{rl}{P}_{\mu {\mu }^{\prime }}^{(1)}(a)& =-g_{\mu {\mu }^{\prime }}+\frac{p_{a,\mu }{p}_{a,{\mu }^{\prime }}}{p_a^2},\\ P_{\mu \nu {\mu }^{\prime }{\nu }^{\prime }}^{(2)}(a)& =\frac{1}{2}(P_{\mu {\mu }^{\prime }}^{(1)}(a)P_{\nu {\nu }^{\prime }}^{(1)}(a)+P_{\mu {\nu }^{\prime }}^{(1)}(a)P_{\nu {\mu }^{\prime }}^{(1)}(a))\\ & +\frac{1}{3}{P}_{\mu \nu }^{(1)}(a)P_{{\mu }^{\prime }{\nu }^{\prime }}^{(1)}(a).\end{array}$$

$$\begin{array}{rl}{\tilde{t}}_{\mu }^{(1)}(a)& =-P_{\mu {\mu }^{\prime }}^{(1)}(a)r_a^{{\mu }^{\prime }},\\ {\tilde{t}}_{\mu \nu }^{(2)}(a)& =P_{\mu \nu {\mu }^{\prime }{\nu }^{\prime }}^{(2)}(a)r_a^{{\mu }^{\prime }}{r}_a^{{\nu }^{\prime }}.\end{array}$$

$$\begin{array}{r}{F}_L(q)=z^L{X}_L(q),\end{array}$$

$$\begin{array}{r}{q}^2=\frac{(s_a+s_b-s_c{)}^2}{4s_a}-s_b,\end{array}$$

$$\begin{array}{rl}{X}_{L=0}(q)& =1,\\ X_{L=1}(q)& =\sqrt{\frac{2}{z^2+1}},\\ X_{L=2}(q)& =\sqrt{\frac{13}{9z^4+3z^2+1}}.\end{array}$$

$$\begin{array}{r}P(m)=\frac{1}{(m_0^2-m^2)-i{m}_0\Gamma (m)},\end{array}$$

$$\begin{array}{r}\Gamma (m)={\Gamma }_0{\left(\frac{q}{q_0}\right)}^{2L+1}\left(\frac{m_0}{m}\right){\left(\frac{X_L(q)}{X_L(q_0)}\right)}^2,\end{array}$$

$$\begin{array}{r}{P}_{\text{GS}}(m)=\frac{1+d\frac{{\Gamma }_0}{m_0}}{(m_0^2-m^2)+f(m)-i{m}_0\Gamma (m)}.\end{array}$$

$$\begin{array}{rl}& f(m)={\Gamma }_0\frac{m_0^2}{q_0^3}\times \\ & \left[q^2(h(m)-h(m_0))+(m_0^2-m^2)q_0^2{\frac{dh}{d(m^2)}|}_{m_0^2}\right],\end{array}$$

$$\begin{array}{r}h(m)=\frac{2q}{\pi m}\ln \left(\frac{m+2q}{2m_{\pi }}\right),\end{array}$$

$${\frac{dh}{d(m^2)}|}_{m_0^2}$$

$$=h(m_0)\left[(8q_0^2{)}^{-1}-(2m_0^2{)}^{-1}\right]+(2\pi m_0^2{)}^{-1}.$$

$$\begin{array}{r}d=\frac{3m_{\pi }^2}{\pi q_0^2}\ln \left(\frac{m_0+2q_0}{2m_{\pi }}\right)+\frac{m_0}{2\pi q_0}-\frac{m_{\pi }^2{m}_0}{\pi q_0^3}.\end{array}$$

$$\begin{array}{r}A(m)=F\sin {\delta }_F{e}^{i{\delta }_F}+R\sin {\delta }_R{e}^{i{\delta }_R}{e}^{i2{\delta }_F},\end{array}$$

$$\begin{array}{rl}{\delta }_F& ={\varphi }_F+{\cot }^{-1}\left[\frac{1}{aq}+\frac{rq}{2}\right]\\ {\delta }_R& ={\varphi }_R+{\tan }^{-1}\left[\frac{M\Gamma (m_{K\pi })}{M^2-m_{K\pi }^2}\right],\end{array}$$

$$\begin{array}{rl}{\text{FF}}_i& =\frac{\int |a_i{A}_i(p){|}^2{R}_4(p)dp}{\int |{\sum }_k{a}_k{A}_k(p){|}^2{R}_4(p)dp}\\ & \approx \frac{{\sum }_{l=1}^{N_{g,\text{ph}}}|a_i{A}_i(p^l){|}^2}{{\sum }_{l=1}^{N_{g,\text{ph}}}|{\sum }_k{a}_k{A}_k(p^l){|}^2},\end{array}$$

$$\begin{array}{r}{N}_{\text{tag}}^{\text{ST}}=2N_{D^0{\bar{D}}^0}{\mathcal{B}}_{\text{tag}}{\epsilon }_{\text{tag}},\end{array}$$

$$\begin{array}{r}{N}_{\text{tag,sig}}^{\text{DT}}=2N_{D^0{\bar{D}}^0}{\mathcal{B}}_{\text{tag}}{\mathcal{B}}_{\text{sig}}{\epsilon }_{\text{tag,sig}},\end{array}$$

$$\begin{array}{r}{\mathcal{B}}_{\text{sig}}=\frac{N_{\text{tag,sig}}^{\text{DT}}}{N_{\text{tag}}^{\text{ST}}}\frac{{\epsilon }_{\text{tag}}}{{\epsilon }_{\text{tag,sig}}}.\end{array}$$

$$\begin{array}{r}f\times \text{CB}(x;\mu ,\sigma ,\alpha ,n)+(1-f)\text{Gaussian}({\mu }_G,{\sigma }_G),\end{array}$$

Full text

I Introduction

Many measurements of $D$ meson decays have been performed since the $D$ mesons were discovered in 1976 by Mark I PhysRevLett.37.255 ; PhysRevLett.37.569 . Today, most of the low-multiplicity $D$ decay mode branching fractions (BFs) are well-measured. The largest decay modes are Cabibbo-favored (CF) hadronic (semileptonic) decay modes resulting from $c\rightarrow sW^{+},W^{+}\rightarrow u\bar{d}\,(l^{+}\nu_{l})$ transitions, but some of these decays are still unmeasured, in which the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay should be the largest unmeasured mode. Charge-conjugate states are implied throughout this paper.

The $D^{0}/D^{+}$ meson is the lightest meson containing a single charm quark. No strong decays are allowed, which makes the $D^{0}/D^{+}$ meson a perfect place to study the weak decay of the charm quark. The CF $\bar{K}\pi$, $\bar{K}2\pi$, and $\bar{K}3\pi$ modes are the most common hadronic decay modes of $D^{0}/D^{+}$ mesons. All $\bar{K}\pi$ and $\bar{K}2\pi$ branching fractions have been measured, but only four of the seven $\bar{K}3\pi$ KLModes have been determined. The Mark III and The E691 collaboration performed amplitude analyses of all four $D\rightarrow\bar{K}\pi\pi\pi$ decay modes, $K^{-}\pi^{+}\pi^{+}\pi^{-}$, $K^{0}_{S}\pi^{+}\pi^{+}\pi^{-}$, $K^{-}\pi^{+}\pi^{+}\pi^{0}$, and $K^{0}_{S}\pi^{+}\pi^{-}\pi^{0}$ PhysRevD.45.2196 ; PhysRevD.46.1941 . Recently, BESIII has remeasured the structure of the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{+}\pi^{-}$ decay with better precision PhysRevD.95.072010 . However, $\bar{K}3\pi$ modes with two or more $\pi^{0}$’s remain unmeasured.

Furthermore, the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay has a large BF and is often used as a $D^{0}$ meson “tag mode” in experiment, such as in the CLEO and BESIII studies of $D^{0}$ semileptonic decays PhysRevD.79.052010 ; PhysRevD.92.072012 . This mode contributes up to $10\%$ of the total reconstructed tags. Therefore, the accurate measurement of its sub-structures and branching fraction is essential to reduce systematic uncertainties of such analyses. While it is true that tag-mode BFs and sub-structure effects cancel to first order, higher-order systematic effects are increasingly important as statistics and precision increase.

The BESIII detector collected a $2.93$ fb*-1* dataset in 2010 and 2011 at $\sqrt{s}=3.773$ GeV 1674-1137-37-12-123001 ; Ablikim2016629 , which corresponds to the mass of the $\psi(3770)$ resonance. The $\psi(3770)$ decays predominantly to $D^{0}\bar{D}^{0}$ or $D^{+}D^{-}$ without any additional hadrons. The excellent tracking, precision calorimetry, and a large $D\bar{D}$ threshold data sample at BESIII provide an excellent opportunity for study of the unmeasured $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay mode. The knowledge of intermediate structure will be crucial for determining the detection efficiency and useful for future usage as a tagging mode. We report here the first partial wave analysis (PWA) and BF measurement of the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay.

II Detection and Data Sets

The BESIII detector is described in detail in Ref. ABLIKIM2010345 . The geometrical acceptance of the BESIII detector is 93% of the full solid angle. Starting from the interaction point (IP), it consists of a main drift chamber (MDC), a time-of-flight (TOF) system, a CsI(Tl) electromagnetic calorimeter (EMC) and a muon system (MUC) with layers of resistive plate chambers (RPC) in the iron return yoke of a 1.0 T superconducting solenoid. The momentum resolution for charged tracks in the MDC is 0.5% at a transverse momentum of 1 GeV/$c$.

Monte Carlo (MC) simulations of the BESIII detector are based on geant4 sim . The production of $\psi(3770)$ is simulated with the kkmc KKMC package, taking into account the beam energy spread and the initial-state radiation (ISR). The photos FSR package is used to simulate the final-state radiation of charged particles. The evtgen EvtGen package is used to simulate the known decay modes with BFs taken from the Particle Data Group (PDG)  PDG , and the remaining unknown decays are generated with the LundCharm model LundCharm . The MC sample referred to as “generic MC”, including the processes of $\psi(3770)$ decays to $D\bar{D}$, non-$D\bar{D}$, ISR production of low mass charmonium states and continuum ($e^{+}e^{-}\rightarrow e^{+}e^{-}$, $\mu^{+}\mu^{-}$, $\gamma\gamma$ and $q\bar{q}$) processes, is used to study the background contribution. The effective luminosities of the generic MC sample correspond to at least 5 times the data luminosity. The signal MC sample includes $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ versus $\bar{D}^{0}\rightarrow K^{+}\pi^{-}$ events generated according to the results of the fit to data.

III Event Selection

Photons are reconstructed as energy clusters in the EMC. The shower time is required be less than 700 ns from the event start time in order to suppress fake photons due to electronic noise or $e^{+}e^{-}$ beam background. Photon candidates within $|\cos\theta|<0.80$ (barrel) are required to have larger than $25$ MeV energy deposition and those with $0.86<|\cos\theta|<0.92$ (endcap) must have larger than $50$ MeV energy deposition. To suppress noise from hadronic shower splitoffs, the calorimeter positions of photon candidates must be at least $10\degree$ away from all charged tracks.

Charged track candidates from the MDC must satisfy $|\cos\theta|<0.93$, where $\theta$ is the polar angle with respect to the direction of the positron beam. The closest approach to the interaction point is required to be less than 10 cm in the beam direction and less than 1 cm in the plane perpendicular to the beam.

Charged tracks are identified as pions or kaons with particle identification (PID), which is implemented by combining the information of $dE/dx$ in the MDC and the time-of-flight from the TOF system. For charged kaon candidates, the probability of the kaon hypothesis is required to be larger than that for a pion. For charged pion candidates, the probability for the pion hypothesis is required to be larger than that for a kaon.

The $\pi^{0}$ candidates are reconstructed through $\pi^{0}\rightarrow\gamma\gamma$ decays, with at least one barrel photon. The diphoton invariant mass is required to be in the range of $0.115<M_{\gamma\gamma}<0.150$ GeV/$c^{2}$.

Two variables, beam constrained mass $M_{\text{BC}}$ and energy difference $\Delta E$, are used to identify the $D$ meson,

[TABLE]

where $|\vec{p}_{D}|^{2}$ and $E_{D}$ are the total reconstructed momentum and energy of the $D$ candidate in the center-of-mass frame of the $\psi(3770)$, respectively, and $E_{\text{beam}}$ is the calibrated beam energy. The $D$ signals will be consistent with the nominal $D$ mass in $M_{\text{BC}}$ and with zero in $\Delta E$.

After charged kaons and charged pions are identified, and neutral pions are reconstructed, hadronic $D$ decays can be reconstructed with a DTag technique. There are two types of samples used in the DTag technique: single tag (ST) and double tag (DT) samples. In the ST sample, only one $D$ or $\bar{D}$ meson is reconstructed through a chosen hadronic decay without any requirement on the remaining measured tracks and showers. For multiple ST candidates, only the candidate with the smallest $|\Delta E|$ is kept. In the DT sample, both $D$ and $\bar{D}$ are reconstructed, where the meson reconstructed through the hadronic decay of interest is called the “signal side”, and the other meson is called the “tag side”. For multiple DT candidates, only the candidate with the smallest summation of $|\Delta E|$s in the signal side and the tag side is kept.

In this amplitude analysis, the DT candidates used are required to have the $D^{0}$ meson decaying to $K^{-}\pi^{+}\pi^{0}\pi^{0}$ as the signal, and the $\bar{D}^{0}$ meson decaying to $K^{+}\pi^{-}$ as the tag. For charged tracks of the signal side, a vertex fit is performed and the $\chi^{2}$ must be less than 100. To improve the resolution and ensure that all events fall within the phase-space boundary, we perform a three-constraint kinematic fit in which the invariant masses of the signal $D$ candidate and the two $\pi^{0}$’s are constrained to their PDG values PDG . The events with kinematic fit $\chi^{2}$ $>$ 80 are discarded.

The tag side is required to satisfy $1.8575<M_{\text{BC}}<1.8775$ GeV/$c^{2}$ and $-0.03<\Delta E<0.02$ GeV. The signal side is required to satisfy $1.8600<M_{\text{BC}}<1.8730$ GeV/$c^{2}$ and $-0.04<\Delta E<0.02$ GeV. A $K^{0}_{S}\rightarrow\pi^{0}\pi^{0}$ mass veto, $M_{\pi^{0}\pi^{0}}\notin(0.458,0.520)$ GeV/$c^{2}$, is also applied on the signal side to remove the dominant peaking background, $D^{0}\rightarrow K^{-}K^{0}_{S}\pi^{+}$. The $M_{\text{BC}}$ and $\Delta E$ distributions of the data and generic MC samples are given in Fig. 1, where the generic MC sample is normalized to the size of data. Note that we always apply the $\Delta E$ requirements before plotting $M_{\text{BC}}$, and vice-versa.

The generic MC sample is used to estimate the background of the DT candidates in the amplitude analysis. The dominant peaking background arises from $D^{0}\rightarrow K^{-}K^{0}_{S}\pi^{0}$, which is suppressed by the $K^{0}_{S}$ mass veto from $2.2\%$ to $0.07\%$. The remaining non-peaking background is about 1.0$\%$. With all selection criteria applied, 5,950 candidate events are obtained with a purity of 98.9%.

IV Amplitude Analysis

This analysis aims to determine the intermediate-state composition of the $D^{0}\to K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay. This four-body decay spans a five-dimensional space. The daughter particle momenta are used as inputs to the probability density function (PDF) which describes the distribution of signal events. This is then used in an unbinned maximum likelihood fit to determinate the intermediate-state composition.

IV.1 Likelihood Function Construction

The PDF is used to construct the likelihood of the amplitude mode and it is given by

[TABLE]

where $p$ is the set of the four daughter particles’ four momenta and $a$ is the set of the complex coefficients for amplitude modes. The $\epsilon(p)$ is the efficiency parameterized in terms of the daughter particles’ four momenta. The four-body phase-space, $R_{4}$, is defined as

[TABLE]

where $\alpha$ indicates the four daughter particles. This analysis uses an isobar model formulation, where the signal decay amplitude, $A(a,p)$, is represented as a coherent sum of a number of two-body amplitude modes:

[TABLE]

where $a_{i}$ is written in the polar form as $\rho_{i}e^{i\phi_{i}}$ ($\rho_{i}$ is the magnitude and $\phi_{i}$ is the phase), and $A_{i}(p)$ is the amplitude for the $i^{\text{th}}$ amplitude mode modeled as

[TABLE]

where the indexes 1 and 2 correspond to the two intermediate resonances. Here, $F^{D}_{i}(p)$ is the Blatt-Weisskopf barrier factor for the $D$ meson, while $P^{1,2}_{i}(p)$ and $F^{1,2}_{i}(p)$ are propagators and Blatt-Weisskopf barrier factors, respectively. The spin factor $S_{i}(p)$ is constructed with the covariant tensor formalism Zou2003 . Finally, the likelihood is defined as

[TABLE]

where $k$ sums over the selected events and $N_{s}$ is the number of candidate events. Consequently, the log likelihood is given by

[TABLE]

Since the second term of Eq. (7) is independent of $a$ and the normalization integration in the denominator of the first term can be approximated by a phase-space MC integration, one can execute an amplitude analysis without knowing efficiency in advance. The phase-space MC integration is obtained by summing over a phase-space MC sample,

[TABLE]

where $N_{\text{g,ph}}$ is the number of generated phase-space events and $N_{\text{s,ph}}$ is the number of selected phase-space events. This holds since the generated sample is uniform in phase space, while the nonuniform distribution after selection reflects the efficiency.

For signal MC samples, the amplitude squared for each event should be normalized by the PDF which generates the sample. The normalization integration using signal MC samples is given by

[TABLE]

where $N_{\text{MC}}$ is the number of the signal MC sample and $a^{\text{gen}}$ is the set of the parameters used to generate the signal MC sample, which is obtained from the preliminary results using the phase-space MC integration. We allow for possible biases caused by tracking, PID, and $\pi^{0}$ data versus MC sample efficiency differences by introducing the correction factors $\gamma_{\epsilon}$,

[TABLE]

where $\epsilon_{j,\text{data}}$ and $\epsilon_{j,\text{MC}}$ are the $\pi^{0}$ reconstruction, the PID, or the tracking efficiencies as a function of $p$ for the data and the MC sample, respectively. By weighting each signal MC event with $\gamma_{\epsilon}$, the MC integration is given by

[TABLE]

IV.1.1 Spin Factor

For a decay process of the form $a\rightarrow bc$, we use $p_{a}$, $p_{b}$, $p_{c}$ to denote the momenta of the particles $a$, $b$, $c$, respectively, and $r_{a}\,=\,p_{b}-p_{c}$. The spin projection operators Zou2003 are defined as

[TABLE]

The covariant tensors are given by

[TABLE]

We list the 10 kinds of spin factors used in this analysis in Table 1, where scalar, pseudo-scalar, vector, axial-vector, and tensor states are denoted by $S$, $P$, $V$, $A$, and $T$, respectively.

IV.1.2 Blatt-Weisskopf Barrier Factors

The Blatt-Weisskopf barrier $F_{i}(p_{j})$ is a barrier function for a two-body decay process, $a\rightarrow bc$. The Blatt-Weisskopf barrier depends on angular momenta and the magnitudes of the momenta of daughter particles in the rest system of the mother particle. The definition is given by

[TABLE]

where $L$ denotes the angular momenta, and $z\,=\,qR$ with $q$ the magnitudes of the momenta of daughter particles in the rest system of the mother particle and $R$ the effective radius of the barrier. For a process $a\rightarrow bc$, we define $s_{i}=E^{2}_{i}-p_{i}^{2}$, $i=a,b,c$, such that

[TABLE]

while the values of $R$ used in this analysis, $3.0$ GeV*-1* and $5.0$ GeV*-1* for intermediate resonances and the $D$ meson, respectively, are used in the BESIII MC generator (based on evtgen). However, these values will also be varied as a source of systematic uncertainties. The $X_{L}(q)$ are given by

[TABLE]

IV.1.3 Propagators

We use the relativistic Breit-Wigner function as the propagator for the resonances $\bar{K}^{*0}$, $K^{*-}$, and $a_{1}(1260)^{+}$, and fix their widths and masses to their PDG values PDG . The relativistic Breit-Wigner function is given by

[TABLE]

where $m=\sqrt{E^{2}-p^{2}}$ and $m_{0}$ is the rest mass of the resonance. $\Gamma(m)$ is given by

[TABLE]

where $q_{0}$ indicates the value of $q$ when $s_{a}=m^{2}_{0}$. Resonances $\bar{K}_{1}(1270)^{0}$ and $K_{1}(1270)^{-}$ are also parameterized by the relativistic Breit-Wigner function but with constant width $\Gamma(m)=\Gamma_{0}$ since these two resonances are very close to the threshold of $\rho K$ and $\Gamma(m)$ vary very rapidly as $m$ changes. We parameterize the $\rho$ with the Gounaris-Sakurai lineshape PhysRevLett.21.244 , which is given by

[TABLE]

The function $f(m)$ is given by

[TABLE]

where

[TABLE]

and

[TABLE]

The normalization condition at $P_{\text{GS}}(0)$ fixes the parameter $d=f(0)/(\Gamma_{0}m_{0})$. It is found to be

[TABLE]

IV.1.4 $K\pi$ $S$-Wave

The kinematic modifications associated with the $K\pi$ $S$-wave are modeled by a parameterization from scattering data ASTON1988493 ; PhysRevD.78.034023 , which are described by a $K^{*0}$ Breit-Wigner along with an effective range non-resonant component with a phase shift,

[TABLE]

with

[TABLE]

where $a$ and $r$ are the scattering length and effective interaction length, respectively. The parameters $F(\phi_{F})$ and $R(\phi_{R})$ are the magnitude (phase) for non-resonant state and resonance terms, respectively. The parameters $M$, $F$, $\phi_{F}$, $R$, $\phi_{R}$, $a$, $r$ are fixed to the results of the $D^{0}\rightarrow K^{0}_{S}\pi^{+}\pi^{-}$ analysis by BABAR PhysRevD.78.034023 , given in Table 2.

Note that we have also tested different parametrizations of the $\pi\pi$ $S$-wave, but no significant improvement is observed. We decide to use phase-space for the $\pi\pi$ $S$-wave.

IV.2 Fit Fraction

The fit fraction (FF) is independent of the normalization and phase conventions in the amplitude formalism, and hence provides a more meaningful summary of amplitude strengths than the raw amplitudes, $\rho_{i}$ in Eq. (4), alone. The definition of the FF for the $i^{\text{th}}$ amplitude is

[TABLE]

where the integration is approximated by a MC integration with a phase-space MC sample. Since the FF does not involve efficiency, the MC sample used here is at the generator level instead of at the reconstruction level, as shown previously in Eq. (LABEL:eq:MCintegration).

As for the statistical uncertainty of the FF, it is not practical to analytically propagate the uncertainties of the FFs from that of the amplitudes and phases. Instead, we randomly perturb the variables determined in our fit (by a Gaussian-distributed amount controlled by the fit uncertainty and the covariance matrix) and calculate the FFs to determine the statistical uncertainties. We fit the distribution of each FF with a Gaussian function and the width is reported as the uncertainty of the FF.

IV.3 Results of Amplitude Analysis

We perform an unbinned likelihood fit using the likelihood described in Section IV.1, where only the complex $a_{i}$ are floating. Starting with amplitude modes with significant contribution, we add (remove) amplitude modes into (from) the fit one by one based on their statistical significances, which are obtained by the change of the log-likelihood value $\Delta\ln L$ with or without the amplitude mode under study. There are 26 amplitudes each with a significance larger than $4\sigma$ chosen as the optimal set, listed in Table 3 and the uncertainties are discussed in Section VI.1. There are more than 40 amplitudes tested but not used in the optimal set ($<4\sigma$ significance), listed in Appendix A.

The amplitude $D\rightarrow K^{-}a_{1}(1260)^{+}$, $a_{1}(1260)^{+}\rightarrow\rho^{+}\pi^{0}[S]$ is expected to have the largest FF. Thus, we choose this amplitude as the reference (phase is fixed to 0) in the PWA. Other important amplitudes are $D\rightarrow(K^{-}\pi^{0})_{S}\rho^{+}$, $D\rightarrow K^{-}a_{1}(1260)^{+}$ with $a_{1}(1260)^{+}[S]\rightarrow\rho^{+}\pi^{0}$, and $D\rightarrow K^{-}a_{1}(1260)^{+}$ with $a_{1}(1260)^{+}[S]\rightarrow\rho^{+}\pi^{0}$. The notation $[S]$ denotes a relative $S$-wave between daughters in a decay, and similarly for $[P],[D]$. A MC sample is generated based on the PWA results, called the PWA signal MC sample. The projections of the data sample and the PWA signal MC sample on the invariant masses squared and the cosines of helicity angles for the $K^{-}\pi^{+}$, $K^{-}\pi^{0}$, $\pi^{+}\pi^{0}$ and $\pi^{0}\pi^{0}$ systems are shown in Fig. 2. The helicity angle $\theta_{ij}$ ($i$ or $j$ is $K^{-}$, $\pi^{+}$ and $\pi^{0}$) is defined as the angle between the momentum vector of the particle $i$ in the $ij$ rest frame and the direction of the $ij$ system in the $D$ rest frame. There are clear $K^{*}(892)^{0}$ and $K^{*}(892)^{-}$ resonances around 0.796 GeV2$/c^{4}$ in the $M^{2}_{K^{-}\pi^{+}}$ and $M^{2}_{K^{-}\pi^{0}}$ projections, respectively, and a $\rho^{+}(770)$ resonance around 0.593 GeV2$/c^{4}$ in the $M^{2}_{\pi^{+}\pi^{0}}$ projection. The gap in the $M^{2}_{\pi^{0}\pi^{0}}$ projection is due to the $K^{0}_{S}$ mass veto. A more detailed goodness-of-fit study is presented in the next section. The PWA signal MC sample improves the accuracy of the DT efficiency (needed to determine the BF), which is discussed in more detail in Section V.3.

IV.4 Goodness-of-Fit

While the one-dimensional projections of the data sample and the PWA signal MC sample shown in Fig. 2 look quite good, much information is lost in projecting down from the full five-dimensional phase space. It is thus desirable to have a more rigorous test of the fit quality. We have programmed a “mixed-sample method” for determining the goodness of our unbinned likelihood fit 1748-0221-5-09-P09004 . According to the method, we can calculate the “T” value of the mixing of two samples, the expectation mean, $\mu_{\text{T}}$, and the variance, $\sigma^{2}_{\text{T}}$. From these values, we can calculate a “pull”, $(\text{T}-\mu_{\text{T}})/\sigma_{\text{T}}$, which should distribute as a normal Gaussian function due to statistical fluctuations. The pull is expected to center at zero if the two samples come from the same parent PDF, and be biased toward larger values otherwise. In the case of our PWA fit, the pull is expected to be a little larger than zero because some amplitudes with small significance are dropped. In other words, adding more amplitudes into the model is expected to decrease the pull.

To check the goodness-of-fit of our PWA results, we calculate the pull of the T value of the mixing of the data sample and the PWA signal MC sample, and it is determined to be 0.97, which indicates good fit quality.

V Branching Fraction

We determinate the BF of $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ using the efficiency based on the results of our amplitude analysis.

V.1 Tagging Technique and Branching Fraction

To extract the absolute BF of the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay, we obtain the ST sample by reconstructing the $\bar{D}^{0}$ meson through the $\bar{D}^{0}\rightarrow K^{+}\pi^{-}$ decay, and the DT sample by fully reconstructing both $D^{0}$ and $\bar{D}^{0}$ through the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay and the $\bar{D}^{0}\rightarrow K^{+}\pi^{-}$ decay as the signal side and the tag side, respectively. The ST yield is given by

[TABLE]

and the DT yield is given by

[TABLE]

where $N_{D^{0}\bar{D}^{0}}$ is the total number of produced $D^{0}\bar{D}^{0}$ pairs, ${\cal B}_{\text{tag(sig)}}$ is the BF of the tag (signal) side, and $\varepsilon$ are the corresponding efficiencies.

The BF of the signal side is determined by isolating ${\cal B}_{\text{sig}}$ such that

[TABLE]

V.2 Fitting Model

The ST yield, $N^{\text{ST}}_{\text{tag}}$, is obtained by a maximum-likelihood fit to the $M_{\text{BC}}$ ($K^{+}\pi^{-}$) distribution. A Crystal Ball (CB) function CBcitation , along with a Gaussian, is used to model the signal while an ARGUS function ARGUScitation is used to model the background. The signal shape is

[TABLE]

where $f$ is a fraction ranging from 0 to 1, $\mu_{G}$ and $\sigma_{G}$ are the mean and the width of the Gaussian function, respectively. The CB function has a Gaussian core transitioning to a power-law tail at a certain point, and is given by

[TABLE]

where $N$ is the normalization and $\alpha$ controls the start of the tail. The beam energy (end point of the ARGUS function) is fixed to be 1.8865 GeV, while all other parameters are floating.

The DT yield, $N^{\text{DT}}_{\text{tag,sig}}$, is obtained by a maximum-likelihood fit to the two dimensional (2-D) $M_{\text{BC}}$ ($K^{-}\pi^{+}\pi^{0}\pi^{0}$) versus $M_{\text{BC}}$ ($K^{+}\pi^{-}$) distribution for the signal and the tag side with a 2-D fitting technique introduced by CLEO PhysRevD.76.112001 . This technique analytically models the signal peak, and considers ISR and mispartition (i.e., where one or more daughter particles are associated with the incorrect $D^{0}$ or $\bar{D}^{0}$ parent) effects, which are non-factorizable in the 2-D plane. In this fitting, the mass of $\psi(3770)$ is fixed to be 3.773 GeV and the beam energy is fixed to be 1.8865 GeV.

V.3 Efficiency and Data Yields

An updated MC sample based on our PWA results, called the PWA MC sample, is used to determine the efficiency. The PWA MC sample is the generic MC sample with the $K^{-}\pi^{+}\pi^{0}\pi^{0}$ versus $K^{+}\pi^{-}$ events replaced by the PWA signal MC sample. All event selection criteria mentioned in Section III are applied except the $M_{\text{BC}}$ requirements. The projections to the signal and the tag side of the fit to the $M_{\text{BC}}$ distributions of the DT of data are shown in Figs. 3(a) and (b), respectively. The background peak in the projection to the signal (tag) side axis is caused by events with a correct signal (tag) and a fake tag (signal). The fit to the $M_{\text{BC}}$ distribution of the ST of data is shown in Fig. 3(c), where both the mean values of the Gaussian function and the CB function agree well with our expectation for the $D^{0}$ mass. The ST and DT data yields are determined to be $534,581\pm 769$ and $6,101\pm 83$, respectively. The ST and DT efficiencies based on the PWA MC sample are $(66.01\pm 0.03)\%$ and $(8.39\pm 0.04)\%$, respectively.

We further take the differences in efficiencies for $\pi^{0}$ reconstruction, tracking, and PID between the data and the PWA MC sample into account. For these differences, we obtain weighted-average efficiency differences $(\varepsilon_{\text{data}}/\varepsilon_{\text{MC}}-1)$ of $-0.69\%$, $1.83\%$, and $0.22\%$ for $\pi^{0}$ reconstruction, kaon tracking, and pion tracking, respectively, while that for PID is negligible. More details are discussed in Section VI.2. This correction is applied to obtain the corrected DT efficiency to be $(8.50\pm 0.04)\%$.

V.4 Result of Branching Fraction

Inserting the values of the DT and ST data yields, the ST efficiency, and the corrected DT efficiency into Eq. (28), we determine the BF of the $K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay, ${\cal B}(D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0})\,=\,(8.86\pm 0.13(\text{stat})\pm 0.19(\text{syst}))\%$. The systematic uncertainties are discussed in Section VI.2.

VI Systematic Uncertainties

The systematic uncertainties of the PWA and BF measurement are discussed in Sections VI.1 and VI.2, respectively.

VI.1 Uncertainties for Amplitude Analysis

The systematic uncertainties for our amplitude analysis are studied in four categories: amplitude model, background, experimental effects, and fit bias. The contributions from the different categories to the systematic uncertainties for the FFs and phases are given in Tables 4 and 5, respectively. The uncertainties of these categories are added in quadrature to obtain the total systematic uncertainties.

The effects of the amplitude model arise from three possible sources: the $K\pi$ $S$-wave model, the effective barrier radii, and the masses and widths of intermediate particles. To determine the systematic uncertainties due to the $K\pi$ $S$-wave model, the fixed parameters of the model are varied according to the BABAR measurement uncertainties ASTON1988493 ; PhysRevD.78.034023 , listed in Table 2. The effective barrier radius $R$ is varied from $1.5$ to $4.5$ GeV*-1* for intermediate resonances, and from $3.0$ to $7.0$ GeV*-1* for the $D^{0}$. The masses and widths of intermediate particles are perturbed according to their published uncertainties in the PDG. The consequent changes of fitting results are considered as the systematic uncertainties inherent in the amplitude model.

The effects of background estimation are separated into non-peaking background, and peaking background. The uncertainties associated with non-peaking background are studied by widening the $M_{\text{BC}}$ and $\Delta E$ windows on the signal side to increase the non-peaking background. The peaking background can be mostly removed by the $K^{0}_{S}$ mass veto. However, this veto is also a source of uncertainties. Its uncertainty is studied by widening this veto from the nominal $M_{\pi^{0}\pi^{0}}\notin(0.458,0.520)$ GeV/$c^{2}$ to $M_{\pi^{0}\pi^{0}}\notin(0.418,0.542)$ GeV/$c^{2}$.

The experimental effects are related to the acceptance difference between data and MC sample caused by $\pi^{0}$ reconstruction, tracking, and PID efficiencies, which weight the normalization of the signal PDF, Eq. (10). To estimate the uncertainties associated with the experimental effects, the amplitude fit is performed varying $\pi^{0}$ reconstruction, tracking and PID efficiencies according to their uncertainties, and the changes of the nominal results are taken as the systematic uncertainties.

The fit bias is tested with 200 pseudoexperiment samples generated based on the PWA model. The distribution of each FF or phase is fitted with a Gaussian function. The difference of the mean and the nominal value is considered as the uncertainty associated with fit bias.

VI.2 Uncertainties for Branching Fraction

We examine the systematic uncertainties for the BF from the following sources: tag side efficiency, tracking, PID, and $\pi^{0}$ efficiencies for signal, $K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay (PWA) model, yield fits, $K^{0}_{S}$ peaking background and the $K^{0}_{S}$ mass veto, and doubly Cabibbo-suppressed (DCS) decay.

The efficiency for reconstructing the tag side, $\bar{D}^{0}\rightarrow K^{+}\pi^{-}$, should almost cancel, and any residual effects caused by the tag side are expected to be negligible. Unlike the case of the tag side, the reconstruction efficiency of the signal side does not cancel in the DT to ST ratio. This efficiency of the signal side is determined with the PWA signal MC sample. The mismatches of tracking, PID, and $\pi^{0}$ reconstruction between the data and MC samples, therefore, bring in systematic uncertainties.

One possible source of those uncertainties is that the momentum spectra simulated in the MC sample do not match those in data, if there are any variations in efficiency versus momentum. This effect, however, is expected to be small due to the PWA MC sample’s successful modeling of the momentum spectra in data, as shown in Fig. 2. The major possible source of the $\pi^{0}$ reconstruction, tracking, and PID systematic uncertainties is that, although the momentum spectra in the MC sample and data follow each other well, the efficiency of the MC sample disagrees with that of data as a function of momentum. This disagreement results in taking a correctly weighted average of incorrect efficiencies. We have performed an efficiency correction and choose $0.6\%$, $0.5\%$, $0.3\%$, and $0.2\%$ as the systematic uncertainties for the $\pi^{0}$ reconstruction, kaon tracking, pion tracking, and kaon/pion PID, respectively. The uncertainty of the $\pi^{0}$ reconstruction efficiency is investigated with the control sample of $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}$ decays and the uncertainties for charged tracks and PID are determined using the control sample of $D^{+}\rightarrow K^{-}\pi^{+}\pi^{+}$ decays, $D^{0}\rightarrow K^{-}\pi^{+}$ decays, and $D^{0}\rightarrow K^{-}\pi^{+}\pi^{+}\pi^{-}$ decays.

To estimate the systematic uncertainty caused by the imperfections of the decay model, we compare our PWA model to another PWA model which only includes amplitudes with significance larger than $5\sigma$. The relative shift on efficiency is less than $0.5\%$. We therefore assign $0.5\%$ as the systematic for the effect of any remaining decay modeling imperfections on efficiency.

To get the uncertainty of the yield fit, we change the nominal $\Delta E$ window to a wider one, $-0.05<\Delta E<0.03$ GeV, and the change of the BF is considered as the associated uncertainty.

The $K^{0}_{S}$ mass veto can veto most $K^{0}_{S}$ peaking background and reduce it to be only $0.07\%$ of the total events. However, the peaking background simulation is not perfect and the $K^{0}_{S}$ mass veto also removes about 13$\%$ of the signal events. Thus, we estimate the uncertainty by narrowing the veto from $M_{\pi^{0}\pi^{0}}\notin(0.458,0.520)$ GeV/$c^{2}$ to $M_{\pi^{0}\pi^{0}}\notin(0.470,0.510)$ GeV/$c^{2}$, while the $K^{0}_{S}$ peaking background increases from $0.07\%$ to $0.15\%$ and the BF change is 0.18$\%$ of itself. We take this full shift as the corresponding uncertainty.

The smooth ARGUS background level is about $1.0\%$ in the signal region. In addition, the 2-D $M_{\text{BC}}$ ($K^{-}\pi^{+}\pi^{0}\pi^{0}$) versus $M_{\text{BC}}$ ($K^{+}\pi^{-}$) fit works well for the background determination. Thus, we believe the uncertainty of the background with such small size will be very small and neglect it.

Our tag and signal sides are required to have opposite-sign kaons. This means that our DT decays are either both CF or both DCS. These contributions can interfere with each other, with amplitude ratios that are approximately known, but with a priori unknown phase. The fractional size of the interference term varies between approximately $\pm 2|A_{\text{DCS}}/A_{\text{CF}}|^{2}\simeq\pm 2\tan^{4}\theta_{C}$, where $\theta_{C}$ is the Cabibbo angle (the square in the first term arises as one power from each decay mode in the cross-term). The two amplitude ratios are not exactly equal to $\tan^{2}\theta_{C}$, due to differing structure in the CF and DCS decay modes, but nonetheless we believe $2\tan^{4}\theta_{C}$ is a conservative uncertainty to set as an approximate “$1\sigma$” scale to combine with other uncertainties.

Systematic uncertainties on the BF are summarized in Table 6, where the total uncertainty is calculated by a quadrature sum of individual contributions.

VII Conclusion

Based on the $2.93$ fb*-1* sample of $e^{+}e^{-}$ annihilation data near $D\bar{D}$ threshold collected by the BESIII detectors, we report the first amplitude analysis of the $D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0}$ decay and the first measurement of its decay branching fraction. We find that the $D^{0}\rightarrow K^{-}a_{1}(1260)^{+}$ decay is the dominant amplitude occupying $28\%$ of total FF (98.54%) and other important amplitudes are $D\rightarrow K_{1}(1270)^{-}\pi^{+}$, $D\rightarrow(K^{-}\pi^{0})_{S\text{-wave}}\rho^{+}$, and $D\rightarrow K^{*-}\rho^{+}$, which is similar, in general, with the results of the BESIII $D^{0}\rightarrow K^{-}\pi^{-}\pi^{+}\pi^{+}$ amplitude analysis PhysRevD.95.072010 . Our PWA results are given in Table 3. With these results in hand, which provide access to an accurate efficiency for the $K^{-}\pi^{+}\pi^{0}\pi^{0}$ data sample, we obtain ${\cal B}(D^{0}\rightarrow K^{-}\pi^{+}\pi^{0}\pi^{0})\,=\,(8.86\pm 0.13(\text{stat})\pm 0.19(\text{syst}))\%$.

Acknowledgements.

The BESIII collaboration thanks the staff of BEPCII and the IHEP computing center for their strong support. This work is supported in part by National Key Basic Research Program of China under Contract No. 2015CB856700; National Natural Science Foundation of China (NSFC) under Contract No. 11835012; National Natural Science Foundation of China (NSFC) under Contracts Nos. 11475185, 11625523, 11635010, 11735014; the Chinese Academy of Sciences (CAS) Large-Scale Scientific Facility Program; Joint Large-Scale Scientific Facility Funds of the NSFC and CAS under Contracts Nos. U1532257, U1532258, U1732263, U1832207; CAS Key Research Program of Frontier Sciences under Contracts Nos. QYZDJ-SSW-SLH003, QYZDJ-SSW-SLH040; 100 Talents Program of CAS; INPAC and Shanghai Key Laboratory for Particle Physics and Cosmology; German Research Foundation DFG under Contract No. Collaborative Research Center CRC 1044; Istituto Nazionale di Fisica Nucleare, Italy; Koninklijke Nederlandse Akademie van Wetenschappen (KNAW) under Contract No. 530-4CDP03; Ministry of Development of Turkey under Contract No. DPT2006K-120470; National Science and Technology fund; The Knut and Alice Wallenberg Foundation (Sweden) under Contract No. 2016.0157; The Swedish Research Council; U. S. Department of Energy under Contracts Nos. DE-FG02-05ER41374, DE-SC-0010118, DE-SC-0012069; University of Groningen (RuG) and the Helmholtzzentrum fuer Schwerionenforschung GmbH (GSI), Darmstadt.

VIII Appendix A: Amplitudes Tested

The following is a list of all amplitude modes tested and found to have a significance smaller than $4\sigma$. These are not included in the final fit set.

 

$D\rightarrow PP,P\rightarrow VP$

$D\rightarrow(K^{*-}\pi^{0})_{P}\pi^{+}$

$D\rightarrow K^{-}(\rho^{+}\pi^{0})_{P}$

 

$D\rightarrow AP,A\rightarrow VP$

$D\rightarrow K_{1}(1270)^{-}\pi^{+},K^{*-}\pi^{0}[D]$

$D\rightarrow K_{1}(1270)^{0}\pi^{0},K^{*-}\pi^{+}[S]$

$D\rightarrow K_{1}(1270)^{0}\pi^{0},K^{*-}\pi^{+}[D]$

$D\rightarrow K_{1}(1270)^{0}\pi^{0},K^{-}\rho^{+}[D]$

$D\rightarrow K^{-}(\rho^{+}\pi^{0})_{A},\rho^{+}\pi^{0}[D]$

$D\rightarrow K^{-}(\rho^{+}\pi^{0})_{A}[S]$

$D\rightarrow(K^{*-}\pi^{+})_{A}\pi^{0},K^{*-}\pi^{+}[S]$

$D\rightarrow(K^{*-}\pi^{0})_{A}\pi^{+},K^{*-}\pi^{0}[D]$

$D\rightarrow(K^{*-}\pi^{+})_{A}\pi^{0},K^{*-}\pi^{+}[D]$

$D\rightarrow(\rho^{+}K^{-})_{A}\pi^{0},K^{-}\rho^{+}[S]$

 

$D\rightarrow AP,A\rightarrow SP$

$D\rightarrow K^{-}((\pi^{+}\pi^{0})_{S}\pi^{0})_{A}$

$D\rightarrow K^{-}((\pi^{0}\pi^{0})_{S}\pi^{+})_{A}$

$D\rightarrow((K^{-}\pi^{0})_{S}\pi^{+})_{A}\pi^{0}$

$D\rightarrow((K^{-}\pi^{0})_{S}\pi^{0})_{A}\pi^{+}$

$D\rightarrow(K^{-}(\pi^{+}\pi^{0})_{S})_{A}\pi^{0}$

$D\rightarrow(K^{-}(\pi^{0}\pi^{0})_{S})_{A}\pi^{+}$

 

$D\rightarrow VS$

$D\rightarrow(K^{-}\pi^{0})_{S}(\pi^{+}\pi^{0})_{V}$

$D\rightarrow(K^{-}\pi^{+})_{V}(\pi^{0}\pi^{0})_{S}$

$D\rightarrow(K^{-}\pi^{0})_{V}(\pi^{+}\pi^{0})_{S}$

 

$D\rightarrow VP,V\rightarrow VP$

$D\rightarrow(K^{*0}\pi^{0})_{V}\pi^{0}$

$D\rightarrow(K^{-}\rho^{+})_{V}\pi^{0}$

 

$D\rightarrow VV$

$D[S]\rightarrow(K^{-}\pi^{0})_{V}\rho^{+}$

$D[S]\rightarrow K^{*-}(\pi^{+}\pi^{0})_{V}$

$D[P]\rightarrow K^{*-}(\pi^{+}\pi^{0})_{V}$

$D[P]\rightarrow(K^{-}\pi^{0})_{V}(\pi^{+}\pi^{0})_{V}$

$D[D]\rightarrow(K^{-}\pi^{0})_{V}(\pi^{+}\pi^{0})_{V}$

 

$D\rightarrow TS$

$D\rightarrow(K^{-}\pi^{+})_{T}(\pi^{0}\pi^{0})_{S}$

$D\rightarrow(K^{-}\pi^{0})_{T}(\pi^{+}\pi^{0})_{S}$

 

**Other

$K^{*-}(1410)\pi^{+}$, $K^{*0}(1410)\pi^{0}$, $K^{*-}(1680)\pi^{+}$, $K^{*0}(1680)\pi^{0}$

$K^{*-}_{2}(1430)\pi^{+}$, $K^{*0}_{2}(1430)\pi^{+}$, $K^{*-}_{2}(1770)\pi^{+}$, $K^{*0}_{2}(1770)\pi^{+}$

$K^{-}a^{+}_{2}(1320)$

$K^{-}\pi^{+}(1300)$

$K^{-}\omega^{+}(1420)$

$K^{-}a^{+}_{1}(1260)$

$K^{*0}f_{0}(980)$

$K^{*0}_{2}(1430)(\pi^{+}\pi^{-})_{S}$, $K^{*-}_{2}(1430)(\pi^{+}\pi^{0})_{S}$

$K^{*-}_{2}(1430)\rho^{+}$

$K^{*0}_{2}(1430)f_{2}(1270)$

$(K^{-}\pi^{+})_{S\text{-wave}}f_{2}(1270)$

$(K^{*-}\pi^{+})_{T}\pi^{0}$, $(K^{*-}\pi^{0})_{T}\pi^{+}$, $(K^{*0}\pi^{0})_{T}\pi^{0}$

$(K^{-}\rho^{+})_{T}\pi^{0}$

**