Resonant heavy Higgs searches at the HL-LHC

TL;DR

This paper evaluates the potential of the High Luminosity LHC to discover or constrain heavy Higgs bosons in extended scalar sectors, focusing on key decay channels and their implications for supersymmetric models.

Contribution

It provides detailed projections for heavy Higgs searches at the HL-LHC, analyzing multiple channels and their sensitivities, and explores constraints on the MSSM parameter space.

Findings

$bar{b}\gamma\gamma$ channel constrains heavy Higgs production cross-section to 14.45 fb at 600 GeV

$bar{b}bar{b}$ channel limits cross-section to 2.5 fb at 1000 GeV

Strong constraints on MSSM parameters, excluding certain $ an{eta}$ and $m_A$ ranges at 95\% CL

Abstract

In this work, we show the importance of searches for heavy resonant scalars ($H$) and pseudoscalars ($A$). Taking cue from the present searches, we make projections for searches in an extended scalar sector at the high luminosity run of the Large Hadron Collider. We study the three most relevant search channels, \textit{i.e.}, $H \to hh, \; H/A \to t\bar{t}$ and $b\bar{b}H/A$. Upon studying multifarious final states for the resonant double Higgs production, we find that the $b\bar{b}\gamma\gamma$ ($\sigma(pp \to H \to hh) \in [81.27,14.45]$ fb for $m_H \in [300,600]$ GeV at 95\% C.L.) and $b\bar{b}b\bar{b}$ ($[5.4,2.5]$ fb for $m_H \in [800,1000]$ GeV at 95\% C.L.) channels are the most constraining. For the $b\bar{b}H$ channel, we can exclude $\sigma(pp \to b\bar{b}H) \in [22.2,3.7]$ fb for $m_H \in [300,500]$ GeV. Finally, we consider the phenomenological Minimal Supersymmetric Standard Model as an example and impose various present constraints and our future direct search-limits and obtain strong constraints on the $m_A-\tan{\beta}$ parameter space, where $m_A$ and $\tan{\beta}$ are respectively the mass of the pseudoscalar and the ratio of the vacuum expectation values of the two Higgs doublets. Assuming that the heavy Higgs boson decays only to Standard Model (SM) states, we find that the $H \to hh \to b\bar{b}\gamma\gamma$ ($H \to t\bar{t}$) channel excludes $\tan{\beta}$ as low as 4 ($m_A \in [400,800]$ GeV) at 95\% CL. This weakens up to $\sim 5.5$ when the $b\bar{b}H$ channel dominates. Upon allowing for non-SM decay modes, the limits weaken.

Figures

Tables

Table 1. Table 1 : Bounds obtained on the di-Higgs cross-sections (in fb) from CMS and ATLAS studies dedicated to the search for non-resonant (NR) and resonant (R) double Higgs production in various channels. The numbers in brackets show the range of the heavy scalar mass considered in that particular study.

Channel	CMS (NR)	CMS (R)	ATLAS (NR)	ATLAS (R)
	($\times$SM)	[fb, (GeV)]	($\times$SM)	[fb, (GeV)]
$b\overline{b}b\overline{b}$	$75$	$1500-45$	$13$	$2000-2$
Sirunyan:2018tki ; CMS-PAS-HIG-17-009 ; Aaboud:2018knk		$\left(260-1200\right)$		$\left(260-3000\right)$
$b\overline{b}\gamma \gamma$	$24$	$240-290$	$22$	$1100-120$
CMS-PAS-HIG-17-008 ; Sirunyan:2018iwt ; Aaboud:2018ftw		$\left(250-900\right)$		$\left(260-1000\right)$
$b\overline{b}{\tau }^{+}{\tau }^{-}$	$30$	$3110-70$	$12.7$	$1780-100$
Sirunyan:2017djm ; Aaboud:2018sfw		$\left(250-900\right)$		$\left(260-1000\right)$
$\gamma \gamma W{W}^{\ast }$ Aaboud:2018ewm			$200$	$40000-6100$
($\gamma \gamma \ell \nu jj$)				$\left(260-500\right)$
$b\overline{b}\ell \nu \ell \nu$	$79$	$20500-800$	$300$	$6000-170$
Sirunyan:2017guj ; Aaboud:2018zhh		$\left(300-900\right)$		$\left(500-3000\right)$
$W{W}^{\ast }W{W}^{\ast }$			$160$	$9300-2800$
Aaboud:2018ksn				$\left(260-500\right)$

Table 2. Table 2 : Applied fixed cuts for the cut-based analysis in the b ​ b ¯ ​ γ ​ γ 𝑏 ¯ 𝑏 𝛾 𝛾 b\bar{b}\gamma\gamma channel.

Fixed cuts

$N_{\ell }=0$

$p_{T,b}>40\left(30\right)\mathrm{GeV},p_{T,\gamma }>30\left(30\right)\mathrm{GeV}$

, , $\Delta R_{\gamma b}>0.4$

Table 3. Table 3 : Details of the final optimised cuts with signal efficiency and background yields after applying all cuts.

Heavy Higgs mass,	Optimised cuts		After all cuts
$m_H$ (GeV)	$m_{b\overline{b}\gamma \gamma }$ (GeV)	$p_{T,\gamma \gamma }>$ (GeV)	Signal Efficiency ($ϵ$)	Background yield at $3000{\text{fb}}^{-1}$
$275$	[$235$ , $275$]	$50$	$0.012$	$30.01$
$300$	[$255$ , $305$]	$50$	$0.024$	$55.62$
$350$	[$300$ , $355$]	$100$	$0.024$	$23.33$
$400$	[$345$ , $405$]	$100$	$0.032$	$15.80$
$450$	[$395$ , $455$]	$100$	$0.042$	$13.75$
$500$	[$445$ , $510$]	$100$	$0.051$	$10.87$
$550$	[$460$ , $570$]	$100$	$0.068$	$14.39$
$600$	[$460$ , $615$]	$100$	$0.076$	$18.11$
$800$	[$560$ , $830$]	$100$	$0.091$	$9.54$
$1000$	[$780$ , $1030$]	$100$	$0.090$	$2.31$

Table 4. Table 4 : The cut-flow table for heavy Higgs of mass 400 400 400 GeV. The table also shows the various perturbative orders at which the cross-sections have been used.

	Signal Efficiency	Event rates with $ℒ=3$ ab-1
Cut flow	$pp\to H$	SM Backgrounds
	$\to hh\to 2b2\gamma$)	$hh\to 2b2\gamma$	$hb\overline{b}$	$t\overline{t}h$	$Zh$	$Z\gamma \gamma$	$b\overline{b}\gamma \gamma \ast$ 111$b\overline{b}\gamma \gamma +c\overline{c}\gamma \gamma +jj\gamma \gamma$.	Fake 1 222$b\overline{b}j\gamma +c\overline{c}j\gamma$.	Fake 2 333$b\overline{b}jj$.	$hj{j}^{\ast }$ 444$\left(gg\to hjj\right)+\left(gg\to hc\overline{c}\right)$
Order		NNLO hhtwiki	NNLO (5FS) +	NLO bkg_twiki_cs	NNLO (QCD) +	LO	LO	LO	LO	LO
			NLO (4FS) bkg_twiki_cs		NLO EW bkg_twiki_cs
$m_{\gamma \gamma }$	$0.123$	$39.71$	$36.68$	$397.97$	$62.21$	$32.86$	$1071.38$	$837.45$	$403.98$	$9.60$
$N_{\ell }$	$0.122$	$39.70$	$36.68$	$290.10$	$62.21$	$32.86$	$1071.34$	$837.40$	$403.98$	$9.60$
$p_{T,b/\gamma }$	$0.081$	$27.65$	$16.34$	$197.83$	$35.87$	$14.00$	$510.73$	$361.01$	$183.70$	$4.91$
$\Delta R$ cuts	$0.052$	$20.56$	$5.09$	$36.73$	$22.32$	$4.86$	$56.24$	$35.60$	$27.05$	$1.53$
$m_{bb}$	$0.036$	$14.19$	$1.41$	$12.74$	$4.43$	$1.02$	$16.44$	$11.47$	$7.47$	$0.41$
$p_{T,\gamma \gamma }$	$0.035$	$14.01$	$1.36$	$12.29$	$4.28$	$0.98$	$15.53$	$10.90$	$6.70$	$0.40$
$m_{b\overline{b}\gamma \gamma }$	$0.032$	$2.96$	$0.29$	$3.31$	$0.87$	$0.21$	$3.84$	$3.18$	$1.03$	$0.08$

Table 5. Table 5 : Applied fixed cuts for the cut-based analysis.

Fixed cuts

$p_{T,b}>35\left(25\right)\mathrm{GeV},p_{T,\gamma }>30\left(30\right)\mathrm{GeV}$

, , $\Delta R_{\gamma b}>0.4$

Table 6. Table 6 : Details of the final optimised cuts with signal efficiency and background yields after applying all cuts.

Heavy Pseudoscalar mass,	Optimised cuts		After all cuts
$m_A$ (GeV)	$m_{b\overline{b}\gamma \gamma }$ (GeV)	$p_{T,\gamma \gamma }>$ (GeV)	Signal Efficiency ($\times {10}^{-2}$)	Background yield at $3000{\text{fb}}^{-1}$
$220$	[$170$ , $235$]	$50$	$0.48$	$30.40$
$300$	[$255$ , $305$]	$50$	$3.24$	$91.87$
$350$	[$290$ , $360$]	$100$	$4.01$	$51.08$
$400$	[$345$ , $420$]	$100$	$5.00$	$34.10$
$600$	[$470$ , $625$]	$100$	$10.19$	$24.28$
$800$	[$590$ , $830$]	$100$	$10.68$	$10.69$
$1000$	[$780$ , $1040$]	$100$	$7.50$	$4.15$

Table 7. Table 7 : The cut-flow table for a pseudoscalar mass of 400 GeV. The various perturbative orders used in the calculations are the same as in Table  4 .

	Signal Efficiency	Event rates with $3000{\text{fb}}^{-1}$ of integrated luminosity
Cut flow	for $pp\to A\to Zh$	SM Backgrounds
	$\to 2b2\gamma$	$hh\to 2b2\gamma$	$hb\overline{b}$	$t\overline{t}h$	$Zh$	$Z\gamma \gamma$	$b\overline{b}\gamma \gamma \ast$ 777$b\overline{b}\gamma \gamma +c\overline{c}\gamma \gamma +jj\gamma \gamma$.	Fake 1 888$b\overline{b}j\gamma +c\overline{c}j\gamma$.	Fake 2 999$b\overline{b}jj$.	$hj{j}^{\ast }$ 101010$\left(gg\to hjj\right)+\left(gg\to hc\overline{c}\right)$
$m_{\gamma \gamma }$	$0.115$	$39.71$	$36.68$	$397.97$	$62.21$	$32.86$	$1071.38$	$837.45$	$403.98$	$9.60$
$p_{T,b/\gamma }$	$0.091$	$32.19$	$26.29$	$314.39$	$47.00$	$19.82$	$670.56$	$483.95$	$241.41$	$7.12$
$\Delta R$ cuts	$0.071$	$23.08$	$8.28$	$62.03$	$25.88$	$5.86$	$72.77$	$48.40$	$34.27$	$2.32$
$m_{bb}$	$0.059$	$9.02$	$3.03$	$23.85$	$22.56$	$5.03$	$33.45$	$21.94$	$15.72$	$0.88$
$p_{T,\gamma \gamma }$	$0.058$	$8.89$	$2.86$	$22.77$	$21.50$	$4.80$	$31.30$	$20.48$	$13.91$	$0.84$
$m_{b\overline{b}\gamma \gamma }$	$0.050$	$2.94$	$0.64$	$6.39$	$5.97$	$1.28$	$8.16$	$5.68$	$2.83$	$0.21$

Table 8. Table 8 : Cut-flow table before performing the multivariate analysis, in the b ​ b ¯ ​ b ​ b ¯ 𝑏 ¯ 𝑏 𝑏 ¯ 𝑏 b\bar{b}b\bar{b} channel.

	Signal Efficiency ($\times {10}^{-3}$)				Cross section [fb]
Cut flow	($pp\to H\to hh\to b\overline{b}b\overline{b}$)				SM Backgrounds
	For $m_H$ of				$b\overline{b}b\overline{b}$	$b\overline{b}c\overline{c}$	$t\overline{t}$	$b\overline{b}jj$	$hh\to 4b$	$t\overline{t}b\overline{b}$
	$400$ GeV	$600$ GeV	$800$ GeV	$1000$ GeV
Order	-	-	-	-	LO	LO	NNLO + NNLL ttbarNNLO	LO	NNLO hhtwiki	LO
$4$ $b$-jets	$5.40$	$30.79$	$56.36$	$67.29$	$278.15$	$26.85$	$2.83$	$66.74$	$0.16$	$8.29$
$H_T$	$5.05$	$30.69$	$56.32$	$67.26$	$263.80$	$25.48$	$2.78$	$64.15$	$0.16$	$8.27$
$\Delta R_{bb}$	$1.67$	$24.80$	$48.00$	$57.34$	$93.67$	$7.94$	$1.72$	$11.90$	$0.12$	$1.87$
$p_{T,\text{di-jet}}$	$0.41$	$21.19$	$46.44$	$56.52$	$54.25$	$4.69$	$1.28$	$7.54$	$0.09$	$1.61$
$m_t$	$0.33$	$18.25$	$38.34$	$45.83$	$46.43$	$3.80$	$0.20$	$6.36$	$0.08$	$0.66$

Table 9. Table 9 : Respective background yields for the b ​ b ¯ ​ b ​ b ¯ 𝑏 ¯ 𝑏 𝑏 ¯ 𝑏 b\bar{b}b\bar{b} channel after the BDT analyses optimised for a heavy Higgs mass of ( a ) ​ 400 𝑎 400 (a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The tables also list the perturbative order at which the cross-sections are considered.

(a)	Process	Order	Events
Background	$b\overline{b}b\overline{b}$	LO	$203.60$
	$b\overline{b}c\overline{c}$	LO	$121.51$
	$b\overline{b}jj$	LO	$46.11$
	$t\overline{t}$	NNLO ttbarNNLO	$10.65$
	$pp\to hh$	NNLO hhtwiki	$6.77$
	$t\overline{t}b\overline{b}$	LO	$77.28$
	Total		$465.92$

Table 10. Table 10 : Background yields after the BDT analysis for heavy Higgs mass of ( a ) ​ 400 𝑎 400 (a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1000 𝑑 1000 (d)~{}1000 GeV for the b ​ b ¯ ​ τ h ​ τ h 𝑏 ¯ 𝑏 subscript 𝜏 ℎ subscript 𝜏 ℎ b\bar{b}\tau_{h}\tau_{h} channel.

(a)	Process	Order	Events
Background	$t\overline{t}$ had	NNLO ttbarNNLO	$662.20$
	$t\overline{t}$ semi-lep	NNLO ttbarNNLO	$5366.58$
	$\tau \tau b\overline{b}$	LO	$3143.59$
	$t\overline{t}h$	NLO bkg_twiki_cs	$296.20$
	$t\overline{t}Z$	NLO Lazopoulos:2008de	$141.56$
	$t\overline{t}W$	NLO Campbell:2012dh	$33.50$
	$pp\to hh$	NNLO hhtwiki	$50.46$
	$b\overline{b}h$	NNLO	$2.36$
	$Zh$	NNLO	$132.88$
	$b\overline{b}jj$	LO	$9558.83$
	Total		$19388.16$

Table 11. Table 11 : Respective background yields for the 2 ​ ℓ + 2 ​ b + E̸ T 2 ℓ 2 𝑏 subscript italic-E̸ 𝑇 2\ell+2b+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The tables also list the perturbative order at which the cross-sections are considered.

(a)	Process	Order	Events
Background	$t\overline{t}$ lep	NNLO ttbarNNLO	$356309.30$
	$t\overline{t}h$	NLO bkg_twiki_cs	$1310.44$
	$t\overline{t}Z$	NLO Lazopoulos:2008de	$1264.20$
	$t\overline{t}W$	NLO Campbell:2012dh	$627.97$
	$pp\to hh$	NNLO hhtwiki	$90.72$
	$\ell \ell b\overline{b}$	LO	$5013.31$
	Total		$364615.94$

Table 12. Table 12 : Respective background yields for the 1 ​ ℓ + 2 ​ j + 2 ​ b + E̸ T 1 ℓ 2 𝑗 2 𝑏 subscript italic-E̸ 𝑇 1\ell+2j+2b+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The various orders of the signal and backgrounds are same as in Table  11 .

(a)	Process	Events
Background	$t\overline{t}$ semi-lep	$9740640.28$
	$t\overline{t}$ lep	$1614225.57$
	$Wb\overline{b}+\text{jets}$ [LO]	$569181.53$
	$t\overline{t}h$	$28364.75$
	$t\overline{t}Z$	$24846.06$
	$t\overline{t}W$	$16935.36$
	$pp\to hh$	$318.61$
	$\ell \ell b\overline{b}$	$20252.87$
	Total	$12014765.03$

Table 13. Table 13 : Respective background yields for the γ ​ γ ​ 1 ​ ℓ ​ 2 ​ j + E̸ T 𝛾 𝛾 1 ℓ 2 𝑗 subscript italic-E̸ 𝑇 \gamma\gamma 1\ell 2j+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The tables also list the perturbative order at which the cross-sections are considered.

(a)	Process	Order	Events
Background	$t\overline{t}h$	NLO bkg_twiki_cs	$6.16$
	$Zh+$ jets	NNLO (QCD) + NLO (EW) bkg_twiki_cs	$1.28$
	$Wh+$ jets	NNLO (QCD) + NLO (EW) bkg_twiki_cs	$11.27$
	$pp\to hh$	NNLO hhtwiki	$1.35$
	$\ell \nu \gamma \gamma +$ jets	LO	$3.33$
	$\ell \ell \gamma \gamma +$ jets	LO	$\sim 1.00$
	Total		$24.39$

Table 14. Table 14 : Respective background yields for the γ ​ γ ​ 2 ​ ℓ + E̸ T 𝛾 𝛾 2 ℓ subscript italic-E̸ 𝑇 \gamma\gamma 2\ell+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The various perturbative orders for the backgrounds are the same as in Table  13 .

(a)	Process	Events
Background	$t\overline{t}h$	$4.78$
	$Zh+$ jets	$1.03$
	$pp\to hh$	$0.74$
	$\ell \ell \gamma \gamma +$ jets	$2.44$
	Total	$8.99$

Table 15. Table 15 : Respective background yields for the 2 ​ ℓ + 2 ​ b + E̸ T 2 ℓ 2 𝑏 subscript italic-E̸ 𝑇 2\ell+2b+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The various orders of the signal and backgrounds are same as in Table  11 .

(a)	Process	Events
Background	$t\overline{t}$ lep	$4979032.27$
	$t\overline{t}h$	$6211.98$
	$t\overline{t}Z$	$6769.81$
	$t\overline{t}W$	$4018.71$
	$pp\to hh$	$111.55$
	$\ell \ell b\overline{b}$	$38875.30$
	Total	$5035019.62$

Table 16. Table 16 : Respective background yields for the 1 ​ ℓ ​ 2 ​ j ​ 2 ​ b + E̸ T 1 ℓ 2 𝑗 2 𝑏 subscript italic-E̸ 𝑇 1\ell 2j2b+\not{E}_{T} channel after the BDT analyses optimised for m H = ( a ) ​ 400 subscript 𝑚 𝐻 𝑎 400 m_{H}=(a)~{}400 GeV, ( b ) ​ 600 𝑏 600 (b)~{}600 GeV, ( c ) ​ 800 𝑐 800 (c)~{}800 GeV and ( d ) ​ 1 𝑑 1 (d)~{}1 TeV. The various orders of the signal and backgrounds are same as in Table  11 .

(a)	Process	Events
Background	$t\overline{t}$ semi-lep	$15257053.17$
	$t\overline{t}$ lep	$2037363.62$
	$Wb\overline{b}+\text{jets}$	$513737.45$
	$t\overline{t}h$	$23963.91$
	$t\overline{t}Z$	$20628.78$
	$t\overline{t}W$	$14852.43$
	$pp\to hh$	$232.54$
	$\ell \ell b\overline{b}$	$25865.11$
	Total	$17893697.01$

Table 17. Table 17 : The details of final optimised cuts with signal efficiency and background yields after all the applied cuts.

Heavy Higgs mass,	Optimised cuts (GeV)				After all cuts
$m_H$ (GeV)	$p_{T,{\tau }_1}>$	$\sum _{\tau }\cos \Delta \varphi >$	$M_T$	${\mathrm{E̸}}_T>$	Signal Efficiency ($\times {10}^{-4}$)	Background yield at $3000{\text{fb}}^{-1}$
$200$	$70$	$-0.10$	[$80$ , $200$]	$0$	$7.94$	$3725.90$
$300$	$75$	$-0.06$	[$160$ , $320$]	$0$	$49.40$	$17172.83$
$400$	$180$	$-0.04$	[$380$ , $400$]	$80$	$1.39$	$2.22$
$500$	$180$	$-0.02$	[$380$ , $420$]	$80$	$5.14$	$5.16$
$600$	$240$	$-0.10$	[$500$ , $580$]	$140$	$5.38$	$0.62$
$800$	$260$	$0.00$	[$400$ , $840$]	$220$	$7.91$	$0.62$
$1000$	$260$	$-0.02$	[$780$ , $1020$]	$200$	$23.03$	$0.36$

Table 18. Table 18 : Best-fit value of signal strength variables, along with the associated errors, derived by ATLAS and CMS using LHC Run-II data, which have been imposed on the parameter space region.

Production mode	CMS Sirunyan:2018koj		ATLAS
	Decay channel	Best fit value	Decay channel	Best fit value
$ggh$	$b\overline{b}$	${2.51}_{-2.01}^{+2.43}$	$WW$	${1.21}_{-0.21}^{+0.22}$ Aaboud:2018jqu
	${\tau }^{+}{\tau }^{-}$	${1.05}_{-0.47}^{+0.53}$	$ZZ$	${1.17}_{-0.50}^{+0.41}$ ATLAS-CONF-2017-047
	$W{W}^{\ast }$	${1.35}_{-0.19}^{+0.21}$	$\gamma \gamma$	${0.81}_{-0.18}^{+0.19}$ Aaboud:2018xdt
	$Z{Z}^{\ast }$	${1.22}_{-0.21}^{+0.23}$
	$\gamma \gamma$	${1.16}_{-0.18}^{+0.21}$
$VBF$	$\tau \tau$	${1.12}_{-0.43}^{+0.45}$	$b\overline{b}$	${3.00}_{-1.60}^{+1.70}$ Aaboud:2018gay
	$W{W}^{\ast }$	${0.28}_{-0.60}^{+0.64}$	$WW$	${0.62}_{-0.36}^{+0.37}$ Aaboud:2018jqu
	$Z{Z}^{\ast }$	$-{0.09}_{-0.76}^{+1.02}$	$\gamma \gamma$	${2.00}_{-0.50}^{+0.60}$ Aaboud:2018xdt
	$\gamma \gamma$	${0.67}_{-0.46}^{+0.59}$
$Wh$	$b\overline{b}$	${1.73}_{-0.68}^{+0.70}$	$b\overline{b}$	${1.08}_{-0.43}^{+0.47}$ Aaboud:2018zhk
	$W{W}^{\ast }$	${3.91}_{-2.01}^{+2.26}$	$b\overline{b}$	${1.21}_{-0.42}^{+0.45}$ Aaboud:2017xsd
	$Z{Z}^{\ast }$	${0.00}_{-0.00}^{+2.33}$
	$\gamma \gamma$	${3.76}_{-1.35}^{+1.48}$
$Zh$	$b\overline{b}$	${0.99}_{-0.45}^{+0.47}$	$b\overline{b}$	${1.20}_{-0.31}^{+0.33}$ Aaboud:2018zhk
	$W{W}^{\ast }$	${0.96}_{-1.46}^{+1.81}$	$b\overline{b}$	${0.69}_{-0.33}^{+0.35}$ Aaboud:2017xsd
	$Z{Z}^{\ast }$	${0.00}_{-0.00}^{+4.26}$
	$\gamma \gamma$	${0.00}_{-0.00}^{+1.14}$
$t\overline{t}h$	$b\overline{b}$	${0.91}_{-0.43}^{+0.45}$	$\frac{{\sigma }_{t\overline{t}h}}{{\sigma }_{t\overline{t}{h}_{SM}}}$	${0.84}_{-0.61}^{+0.64}$ Aaboud:2017rss
	${\tau }^{+}{\tau }^{-}$	${0.23}_{-0.88}^{+1.03}$	$b\overline{b}$	${0.80}_{-0.60}^{+0.60}$ Aaboud:2017jvq
	$W{W}^{\ast }$	${1.60}_{-0.59}^{+0.65}$	$\gamma \gamma$	${0.60}_{-0.60}^{+0.70}$ Aaboud:2017jvq
	$Z{Z}^{\ast }$	${0.00}_{-0.00}^{+1.50}$
	$\gamma \gamma$	${2.18}_{-0.75}^{+0.88}$

Table 19. Table 19 : Generation level cuts and cross-sections for the various backgrounds used in the analyses. The backgrounds labelled with ∗ * are multiplied by the fake rates before doing the analysis. The fake rates used are 0.05 % percent 0.05 0.05\%   ATL-PHYS-PUB-2017-001 for j → γ → 𝑗 𝛾 j\to\gamma , ∼ 1.75 % similar-to absent percent 1.75 \sim 1.75\% (average from the fake rate function) for j → b → 𝑗 𝑏 j\to b and 0.35 % percent 0.35 0.35\%   CMS-PAS-TAU-16-002 for j → τ → 𝑗 𝜏 j\to\tau .

Process	Backgrounds	Generation-level cuts ($\ell =e^{\pm },{\mu }^{\pm }$) (NA : Not Applied)	Cross section (fb)
$pp\to H\to hh$, $pp\to A\to Zh$ and $pp\to H\to t\overline{t}$ final states
$b\overline{b}\gamma \gamma$	$hh\to b\overline{b}\gamma \gamma$	NA	$0.10$
	$b\overline{b}\gamma \gamma +$ jets	$p_{T,j/b/\gamma }>20\text{GeV}$, , , $\Delta R_{b,j,\gamma }$212121$\Delta R_{b,j,\gamma }$ means $\Delta R$ between all possible combination of $b,j$ and $\gamma$.$>0.2$, $m_{bb}>50$ GeV, GeV	$18.78$
	$c\overline{c}\gamma \gamma$	$p_{T,j/b/\gamma }>20\text{GeV}$, , $\Delta R_{b,j,\gamma }$$>0.2$, GeV	$162.22$
	$jj\gamma \gamma$	same as $c\overline{c}\gamma \gamma$	${2770.67}^{\ast }$
	$t\overline{t}h$, $h\to \gamma \gamma$	NA	$1.39$
	$hb\overline{b}$, $h\to \gamma \gamma$	NA	$1.32$
	$Zh$, $h\to \gamma \gamma$, $Z\to b\overline{b}$	NA	$0.33$
	$b\overline{b}jj$	$p_{T,j}>10\text{GeV}$, $p_{T,b}>20\text{GeV}$, , $m_{jj}>50$ GeV, $m_{bb}>50$ GeV	${549583730.00}^{\ast }$
	$b\overline{b}j\gamma$	$p_{T,j/b/\gamma }>20\text{GeV}$, , , $\Delta R_{b/b/\gamma /\gamma ,b/j/j/b}$222222$\Delta R_{a/b,c/d}$ signifies $\Delta R_{ac}$ and $\Delta R_{bd}.$$>0.2$, $m_{bb}>50$ GeV	${201800}^{\ast }$
	$c\overline{c}j\gamma$	$p_{T,j/\gamma }>20\text{GeV}$, , , $\Delta R_{j,\gamma }$$>0.2$	${1132709.63}^{\ast }$
	$Z\gamma \gamma +$ jet, $Z\to b\overline{b}$	$p_{T,j/b/\gamma }>20\text{GeV}$, , , $\Delta R_{b,j,\gamma }$$>0.2$, $m_{bb}>50$ GeV, GeV	$0.87$
	$gg\to h+c\overline{c}$, $h\to \gamma \gamma$	$p_{T,j}>20\text{GeV}$, , $\Delta R_{j,j}$$>0.2$	$0.31$
	$gg\to h+jj$, $h\to \gamma \gamma$	same as $gg\to h+c\overline{c}$	${27.89}^{\ast }$
$b\overline{b}b\overline{b}$	$hh\to b\overline{b}b\overline{b}$	NA	$13.42$
	multijet $b\overline{b}b\overline{b}$	$p_{T,j/b}>50\text{GeV}$, , $\Delta R_{b,j}>0.3$, $H_T>250$ GeV	$14541.30$
	multijet $b\overline{b}c\overline{c}$	$same$ as multijet $b\overline{b}b\overline{b}$	$28633.60$
	multijet $b\overline{b}jj$	$same$ as multijet $b\overline{b}b\overline{b}$	${3602560.00}^{\ast }$
	$t\overline{t}$, $W^{\pm }\to c\left(\overline{c}\right)\overline{s}\left(s\right)$	$same$ as multijet $b\overline{b}b\overline{b}$	$860.17$
	$t\overline{t}b\overline{b}$	$p_{T,b}>50\text{GeV}$, , $\Delta R_{b,b}>0.3$, $H_T>250$ GeV	$170.58$
$b\overline{b}{\tau }^{+}{\tau }^{-}$	$hh\to b\overline{b}{\tau }^{+}{\tau }^{-}$	NA	$2.89$
	$t\overline{t}$ hadronic	$p_{T,j/b}>20\text{GeV}$, $p_{T,l}>8\text{GeV}$, , , $\Delta R_{b,j,\ell }>0.2$, $m_{bb}>50$ GeV	$135623.50$
	$t\overline{t}$ semi-leptonic	$same$ as $t\overline{t}$ hadronic	$173409.88$
	$\tau \tau b\overline{b}$	$p_{T,b}>20\text{GeV}$, $p_{T,\tau }>8\text{GeV}$, , $\Delta R_{b,\tau }>0.2$, $m_{bb}>50$ GeV, $m_{\tau \tau }>30$ GeV	$2128.56$
	$b\overline{b}h$, $h\to \tau \tau$	$p_{T,b}>20\text{GeV}$, $p_{T,\tau }>10\text{GeV}$, , , $\Delta R_{b,\tau }>0.2$, $m_{bb}>50$ GeV	$1.23$
	Zh, $h\to \left(b\overline{b}+\tau \tau \right)$, $Z\to \left(\tau \tau +b\overline{b}\right)$	NA	$28.21$
	$t\overline{t}h$	NA	$611.30$
	$t\overline{t}Z$	NA	$731.54$
	$t\overline{t}W$	NA	$437.87$
	$b\overline{b}jj$	$p_{T,j}>10\text{GeV}$, $p_{T,b}>20\text{GeV}$, , $m_{jj}>50$ GeV, $m_{bb}>50$ GeV	${549583730.00}^{\ast }$

Table 20. Table 20 : Generation level cuts and cross-sections for the signals and various backgrounds used in the analyses. The backgrounds labelled with ∗ * are multiplied by the fake rates before doing the analysis. The fake rates used are 0.35 % percent 0.35 0.35\%   CMS-PAS-TAU-16-002 for j → τ → 𝑗 𝜏 j\to\tau .

Process	Backgrounds	Generation-level cuts ($\ell =e^{\pm },{\mu }^{\pm }$) (NA : Not Applied)	Cross section (fb)
$pp\to H\to hh$, $pp\to A\to Zh$ and $pp\to H\to t\overline{t}$ final states
$b\overline{b}W{W}^{\ast }$	$hh\to b\overline{b}{W}^{+}{W}^{-}$	NA	$9.85$
	$t\overline{t}$ semi-leptonic	$p_{T,j/b}>20\text{GeV}$, $p_{T,l}>8\text{GeV}$, , , $\Delta R_{b,j,\ell }>0.2$, $m_{bb}>50$ GeV	$173409.88$
	$t\overline{t}$ leptonic	$same$ as $t\overline{t}$semileptonic	$55319.44$
	$\ell \ell b\overline{b}$	$p_{T,b}>20\text{GeV}$, $p_{T,l}>8\text{GeV}$, , $\Delta R_{b,\ell }>0.2$, $m_{bb}>50$ GeV	$7393.72$
	$Wbb+jets$, $W\to \ell \nu$, $\ell$ also includes $\tau$	$p_{T,j/b}>20\text{GeV}$, $p_{T,l}>8\text{GeV}$, , , $\Delta R_{j,b,l}>0.2$	$32576.60$
	$t\overline{t}h$	NA	$611.30$
	$t\overline{t}Z$	NA	$731.54$
	$t\overline{t}W$	NA	$437.87$
$\gamma \gamma W{W}^{\ast }$	$hh\to \gamma \gamma W^{+}{W}^{-}$	NA	$0.04$
	$t\overline{t}h$, $h\to \gamma \gamma$	NA	$1.39$
	$Zh$ + jets, $h\to \gamma \gamma$, $Z\to \ell \ell$($\ell$ includes $\tau$ also)	$p_{T,\gamma /\ell }>10\text{GeV}$, , , $\Delta R_{\gamma ,\ell ,j}>0.2$, 120 GeV 130 GeV	$0.12$
	$Wh$ + jets, $h\to \gamma \gamma$, $W\to \ell \nu$($\ell$ includes $\tau$ also)	same as $Zh$ + jets	$0.70$
	$\ell \nu \gamma \gamma$ + jets, $\ell$ also includes $\tau$	$p_{T,\gamma /\ell }>10\text{GeV}$, , , $\Delta R_{\gamma \gamma }>0.2$, $\Delta R_{\gamma \ell }>0.2$, $\Delta R_{\gamma j}>0.4$, 120 GeV 130 GeV	$3.17$
	$\ell \ell \gamma \gamma$ + jets, $\ell$ also includes $\tau$	same as $\ell \nu \gamma \gamma$ + jets, with $m_{\ell \ell }>$ 20 GeV	$1.00$
$pp\to b\overline{b}H$ final state
$b\overline{b}{\tau }^{+}{\tau }^{-}$	$t\overline{t}$	$p_{T,b}>20\text{GeV}$, , $\Delta R_{bb}>0.2$	$633946.81$
	single top s-channel	NA	$11390.00$
	single top t-channel	NA	$248090.00$
	single top Wt-channel	NA	$84400.00$
	$\tau \tau$ + jets, via $Z/{\gamma }^{\ast }$	$p_{T,j/b}>20\text{GeV}$, $p_{T,\ell }>60\text{GeV}$, , $\Delta R_{j,b,\ell }>0.2$, $m_{\ell \ell }>50$ GeV	$884370.24$
	$W$ + jets, $W\to \ell \nu$($\ell$ includes $\tau$ also)	$p_{T,j/b/\ell }>20\text{GeV}$, , $\Delta R_{j,b,\ell }>0.2$	$112358.64$
	$VV$ ($V$ includes $W^{\pm }$ and $Z$)	NA	$106510.72$
	$b\overline{b}jj$	$p_{T,j}>65\text{GeV}$, $p_{T,b}>20\text{GeV}$, , $\Delta R_{j,b}>0.2$, $m_{jj}>50$ GeV	${12091572.60}^{\ast }$

Table 21. Table 21 : Comparison table for the τ h ​ τ h subscript 𝜏 ℎ subscript 𝜏 ℎ \tau_{h}\tau_{h} channel in the b 𝑏 b -tag category.

Process	Event rates at $13$ TeV with $36.1{\text{fb}}^{-1}$ of integrated luminosity	Total
	b-tag category
multijet	$97.06$	$215.26$
$Z/{\gamma }^{\ast }+$ jets	$11.03$
$W+$ jets	$2.82$
$t\overline{t}$	$83.66$
$VV+$ jets	$1.87$
$500$ GeV Signal	$18.82$
ATLAS numbers Aaboud:2017sjh
multijet	$106\pm 32$	$180\pm 60$
$Z/{\gamma }^{\ast }\to \tau \tau$	$7.5\pm 2.9$
$W(\to \tau \nu )+$ jets	$4.0\pm 1.0$
$t\overline{t}+$ single top	$60\pm 50$
Others	$1.0\pm 0.5$
$500$ GeV Signal	$28\pm 12$

Equations

$$\begin{array}{r}{m}_{bb},\Delta R_{bb},p_{T,\gamma \gamma },\Delta R_{\gamma \gamma },m_{b\bar{b}\gamma \gamma },\Delta R_{b_i{\gamma }_j},\\ \Delta R_{bb,\gamma \gamma },p_{T,b_1},p_{T,b_2},p_{T,{\gamma }_1},p_{T,{\gamma }_2}\end{array}$$

$$p_{T,\text{di-jet},k},m_{\text{di-jet},k},\Delta {\varphi }_{bb,\text{di-jet},k},\Delta {\eta }_{\text{di-jets}},\Delta R_{\text{di-jets}},m_{b\bar{b}b\bar{b}}.$$

$${\chi }_{hh}^2$$

$$\begin{array}{r}{p}_{T,bb},m_{bb},p_{T,{\tau }_h{\tau }_h},\Delta R_{{\tau }_h{\tau }_h},\Delta {\varphi }_{{\tau }_{h1}/E_T},\Delta {\varphi }_{{\tau }_{h2}/E_T},\\ M_T,m_{T2},p_{T,\text{tot}},m_{\text{tot}},m_{\text{eff}},\Delta R_{b_1{\tau }_{h1}},\Delta R_{bb,{\tau }_h{\tau }_h}\end{array}$$

$$\begin{array}{r}{p}_{T,bb},{\eta }_{bb},{\varphi }_{bb},m_{bb},\Delta R_{bb},\Delta {\varphi }_{bb},p_{T,\ell \ell },{\eta }_{\ell \ell },{\varphi }_{\ell \ell },m_{\ell \ell },\Delta R_{\ell \ell },\\ M_T,m_{T2},m_{\text{tot}},p_{T,\text{tot}},{\varphi }_{\text{tot}},m_{\text{eff}},\Delta R_{b_1{\ell }_2},\Delta R_{bb,\ell \ell },p_{T,{\ell }_2},\end{array}$$

$$\begin{array}{r}{\eta }_{bb},m_{bb},m_t,m_{jj},\Delta R_{jj},\Delta R_{\ell ,jj},M_T,m_{T2},\\ m_{bb{j}_1},m_{t11},m_{t12},p_{T,\ell \nu },p_{T,b_1},p_{T,{\ell }_1},p_{T,j_1},\end{array}$$

$$p_{T,\gamma \gamma },\Delta R_{\gamma \gamma },\Delta {\varphi }_{\ell /E_T},mt,\Delta R_{\ell jj},m_{tot},m_{eff},M_T,\Delta R_{{\gamma }_1\ell },\Delta R_{{\gamma }_2\ell },p_{T,{\gamma }_2},p_{T,\ell },/E_T,$$

$$p_{T,\gamma \gamma },\Delta R_{\gamma \gamma },\Delta {\varphi }_{\gamma \gamma },m_{\ell \ell },\Delta R_{\ell \ell },M_T,m_{\text{tot}},m_{\text{eff}},\Delta R_{\gamma \gamma ,\ell \ell },p_{T,\ell 1},$$

$$p_{T,bb},{\eta }_{bb},{\varphi }_{bb},m_{bb},M_T,m_{\text{tot}},m_{\text{eff}},\Delta R_{b_1{\ell }_1},p_{T,b_1},p_{T,b_2},p_{T,{\ell }_1},p_{T,{\ell }_2},$$

$$\begin{array}{r}{p}_{T,bb},m_{bb},\Delta R_{bb},m_{jj},m_{\text{eff}},M_T,m_{t11},\\ m_{t12},m_{t\bar{t}1},m_{t\bar{t}2},p_{T,\ell \nu },\Delta R_{bb,jj},p_{T,{\ell }_1},p_{T,j_1},\end{array}$$

$${\sigma }^{\text{matched}}=\frac{{\sigma }^{\text{4FS}}+w{\sigma }^{\text{5FS}}}{1+w},$$

$$w=\text{ln}\frac{m_H}{m_b}-2$$

$$M_T=\sqrt{(p_{T,{\tau }_1}+p_{T,{\tau }_2}+/E_T{)}^2-({\vec{p}}_{T,{\tau }_1}+{\vec{p}}_{T,{\tau }_2}+\overrightarrow{/E_T}{)}^2},$$

$${\mu }_{if}=\frac{{\sigma }_i\times {\text{BR}}_f}{{\sigma }_i^{SM}\times {\text{BR}}_f^{SM}},$$

$$1<\tan \beta <60,200\mathrm{GeV}<m_A<1\mathrm{TeV},1\mathrm{TeV}<M_3<10\mathrm{TeV}$$

$$1\mathrm{TeV}<M_{{\tilde{Q}}_3,{\tilde{u}}_3,{\tilde{d}}_3}<20\mathrm{TeV},-10\mathrm{TeV}<A_{t,b}<10\mathrm{TeV}$$

$$1\mathrm{TeV}<M_{{\tilde{Q}}_1,{\tilde{u}}_1,{\tilde{d}}_1}<20\mathrm{TeV},M_{{\tilde{Q}}_2}=M_{{\tilde{Q}}_1},M_{{\tilde{u}}_2}=M_{{\tilde{u}}_1},M_{{\tilde{d}}_2}=M_{{\tilde{d}}_1}$$

$$A_{e,\mu ,\tau ,u,d,c,s}=0,M_{{\tilde{e}}_{1_L},{\tilde{e}}_{1_R},{\tilde{e}}_{2_L},{\tilde{e}}_{2_R},{\tilde{e}}_{3_L},{\tilde{e}}_{3_R}}=3\mathrm{TeV},600\mathrm{GeV}<M_{1,2},\mu <5\mathrm{TeV}$$

Full text

††preprint: IPPP/18/107

1 Introduction

The Higgs boson discovered in 2012, was the last missing piece in the Standard Model of particle physics (SM). The SM, however, is inadequate to explain the nature and existence of dark matter, the small but non-negligible masses of the neutrinos, the asymmetry between baryons and anti-baryons, to name a few. Besides, SM can not explain the hierarchy problem which is inherent in the theory. Well motivated theories including supersymmetry have the potential to solve some of these limitations. There are additional fundamental theoretical requirements that the SM can not satisfy. The aforementioned experimental observations and theoretical requirements compel us to look for physics beyond the Standard Model (BSM). However, the possibilities being innumerable, it is extremely difficult to ascertain the nature of such new physics. Since the discovery of the Higgs boson, and a growing convergence of its properties with the SM expectations Khachatryan:2016vau ; CMS-PAS-HIG-16-020 ; ATLAS-CONF-2016-067 ; ATLAS-CONF-2016-079 ; ATLAS-CONF-2016-081 ; CMS:2017jkd ; ATLAS-CONF-2016-112 ; CMS-PAS-HIG-15-003 ; Sirunyan:2017khh ; ATLAS-CONF-2018-021 ; Sirunyan:2018hoz ; Aaboud:2018urx ; Aaboud:2018zhk ; Sirunyan:2018kst ; CMS-PAS-HIG-16-006 ; ATLAS-CONF-2017-014 , the new physics possibilities are gradually getting strongly constrained. Searches for BSM are being performed at the Large Hadron Collider (LHC) by the CMS, ATLAS, ALICE and LHCb collaborations. Except for some excitement in the flavour physics sector, there have not been any strong hints for new physics in the form of new particles or significant deviations in couplings with respect to the SM. Even though supersymmetry is perhaps one of the most elegant theories of our time, it comes with additional new particles, which need to be discovered at some stage. Even though searches are being conducted for a considerable region of parameter space for the Minimal Supersymmetric Standard Model (MSSM), there are more non-traditional searches which need to be performed. The MSSM parameter space has been extensively studied in light of the constraints from cosmology, flavour physics, and Run-I plus Run-II data from LHC Carena:2011aa ; Arbey:2011ab ; Baer:2012mv ; Arbey:2012dq ; Altmannshofer:2012ks ; Cheung:2013kla ; Chowdhury:2015yja ; Bhattacherjee:2015sga ; Bechtle:2016kui ; Bechtle:2012jw ; Djouadi:2013lra ; Bechtle:2015pma ; Buchmueller:2013rsa ; Scopel:2013bba ; deVries:2015hva ; Barr:2016sho ; Kowalska:2016ent ; Han:2016xet ; Buckley:2016kvr ; Zhao:2017qpe ; Bagnaschi:2016afc ; Bagnaschi:2016xfg ; Barman:2016jov ; Bagnaschi:2017tru ; Costa2018 . However, there are simple extensions of the MSSM that can weaken the present bounds considerably. On the positive side, the LHC can potentially pin down the Higgs couplings to weak bosons and most of the fermions at the level of $\mathcal{O}(5-10\%)$ CMS-PAS-FTR-18-011 ; ATL-PHYS-PUB-2018-054 ; Peskin:2012we ; Peskin:2013xra . However, as has been shown in numerous experimental Aad:2015xja ; Aad:2015uka ; ATLAS-CONF-2016-049 ; ATLAS-CONF-2016-071 ; ATLAS-CONF-2016-004 ; CMS-PAS-HIG-16-028 ; CMS-PAS-HIG-16-002 ; CMS-PAS-B2G-16-008 ; CMS-PAS-HIG-16-032 ; CMS-PAS-HIG-16-011 ; CMS-PAS-HIG-17-002 ; CMS-PAS-HIG-17-006 ; Aaboud:2016xco ; Aaboud:2018sfw (including future extrapolations ATL-PHYS-PUB-2014-019 ; ATL-PHYS-PUB-2017-001 ; CMS-DP-2016-064 ) and phenomenological studies Osland:1998hv ; Baur:2003gp ; Liu:2004pv ; Dib:2005re ; Wang:2007zx ; Pierce:2006dh ; Kanemura:2008ub ; Contino:2010mh ; Grober:2010yv ; Dolan:2012ac ; Dolan:2012rv ; Kribs:2012kz ; Contino:2012xk ; Dawson:2012mk ; Dolan:2012ac ; Ellwanger:2013ova ; Nishiwaki:2013cma ; Barr:2013tda ; Chen:2013emb ; deLima:2014dta ; Chen:2014xra ; Baglio:2014nea ; Hespel:2014sla ; Slawinska:2014vpa ; Englert:2014uqa ; Goertz:2014qta ; Maltoni:2014eza ; Baglio:2014nea ; Hespel:2014sla ; Azatov:2015oxa ; Lu:2015jza ; Lu:2015qqa ; Carvalho:2015ttv ; Bian:2016awe ; Banerjee:2016nzb ; Gorbahn:2016uoy ; Carvalho:2016rys ; Cao:2016zob ; Crivellin:2016ihg ; Grober:2016wmf ; Huang:2017jws ; Kribs:2017znd ; DiLuzio:2017tfn ; Nakamura:2017irk ; Alves:2017ued ; Adhikary:2017jtu ; Chala:2018ari ; Basler:2018dac , the measurement of the elusive triple Higgs coupling ($\lambda_{hhh}$) is a difficult feat at the LHC. Studies show that future colliders are expected to be more adept in constraining or even measuring this coupling to a great precision Plehn:2005nk ; Binoth:2006ym ; Yao:2013ika ; Liu:2014rba ; Barr:2014sga ; He:2015spf ; Papaefstathiou:2015iba ; Azatov:2015oxa ; Kotwal:2016tex ; Fuks:2017zkg ; Banerjee:2018yxy . In order to be completely sure whether or not there is any extended Higgs sector, it is of utmost importance to measure the Higgs quartic coupling, $\lambda_{hhhh}$ and the Higgs trilinear coupling, $\lambda_{hhh}=\lambda_{hhhh}v$, where $v$ is the vacuum expectation value of the SM Higgs boson. Now, independent measurements of the Higgs couplings to the gauge bosons and fermions will constrain $v$ and we already have a precise Higgs mass measurement. To confirm this sector of the SM, one needs to measure $\lambda_{hhh}$ or $\lambda_{hhhh}$ Glover:1987nx ; Boudjema:1995cb ; Plehn:1996wb ; Djouadi:1999rca ; Baur:2002qd ; Baur:2003gp .

In the following sections, we focus on the various production and decay processes of a resonant scalar, viz., a resonant decay to a pair of SM-like Higgs bosons, to a pair of top quarks, a heavy pseudoscalar, $A$, decaying to an SM-like Higgs boson and a $Z$-boson and the production of a heavy scalar in association with a pair of bottom quarks. The final theme of this work is in the context of the phenomenological MSSM (pMSSM). However, our results are presented in such a way that they can be mapped onto most models with an extended scalar sector. Table 1 summarises the various bounds set on the double-Higgs production cross-section by CMS and ATLAS in the non-resonant and resonant categories. The resonant production results are mostly interpreted in terms of spin-0 and spin-2 hypotheses. Besides, there are many supersymmetric interpretations for the resonant scalar searches. As an example, for the $b\bar{b}\tau^{+}\tau^{-}$ resonant search performed by CMS Sirunyan:2017djm , the MSSM parameters $m_{A}$ (mass of the $CP$-odd scalar, $A$) and $\tan{\beta}$ (ratio of the vacuum expectation values of the two Higgs doublets in the model, viz., $H_{u}$ and $H_{d}$) are excluded in the range 230 GeV $<m_{A}<$ 360 GeV and $\tan{\beta}\lesssim 2$, at 95% CL. Thus, besides measuring deviations to the Higgs self-coupling, there are other possible channels to look for in order to establish an extended scalar sector. Some of these new channels include the production of the SM-like Higgs in association with a $Z$-boson reconstructing a resonance. Another possible channel is the production of a pair of top-quarks. Now, the first of these channels can be via a heavy pseudoscalar resonance Khachatryan:2015tha ; Aaboud:2017cxo , whereas the $t\bar{t}$ production can be either through a heavy scalar or pseudoscalar Carena:2016npr ; Aaboud:2017hnm . However, both these channels can also come about from a heavy $Z^{\prime}$ CMS:2016zte ; Aaboud:2017ahz ; Aaboud:2017cxo ; Aaboud:2018mjh ; Aaboud:2018bun . Now, in order to be sure whether the $Zh$ or $t\bar{t}$ production is via a spin-0 or spin-1 resonance, one needs to delve deeper into the angular observables. Lastly, we also study the effects of the high $\tan{\beta}$ regime for a heavy scalar produced in association with a pair of $b$-quarks and decaying to a pair of $\tau$-leptons Aaboud:2017sjh .

Our paper is organised as follows. We study the reach of the HL-LHC for the $H\to hh$ channel in various final states, in section 2, by showing the 95% CL bounds on $\sigma(pp\to H\to hh)$ as functions of the heavy Higgs mass, $m_{H}$. Following the Higgs pair production, we address the couplings of the $CP$-even heavy Higgs to a pair of top quarks and to a pair of bottom quarks in sections 3 and 4 respectively. In section 5, we use the previous results to recast our limits in the purview of the phenomenological MSSM (pMSSM) and show the future reach of these searches in the $m_{A}-\tan{\beta}$ parameter space. We finally summarise our results and present our future outlook in section 6.

2 The $pp\to H\to hh$ Channel

As discussed in the introduction, the objective of this work is to scrutinise the viable scalar extensions of the Standard Model (SM). In this section, we focus on a heavy $CP$-even scalar produced via gluon fusion and subsequently decaying to a pair of SM-like Higgs bosons. The decay width of the heavy Higgs boson is chosen to be $\sim 1~{}{\rm GeV}$ which is within the resolution of the detector. We would like to mention that the heavy Higgs search limits derived in the course of this analysis would stand valid only if the detector resolution is greater than the chosen heavy Higgs decay width. In the next five subsections, we study the reach of the HL-LHC in constraining the resonant Higgs pair production cross-section, $\sigma(pp\to H\to hh)$, upon studying multifarious channels, viz., $b\bar{b}\gamma\gamma,b\bar{b}b\bar{b},b\bar{b}\tau^{+}\tau^{-},b\bar{b}W^{+}W^{-}$ and $\gamma\gamma W^{+}W^{-}$. Many of these channels with $\tau$-leptons and $W$-bosons give different signatures upon considering leptonic or hadronic modes. We study all possible final states giving importance to the total rate as well the cleanliness. Unless otherwise stated, we generate the signal samples with Pythia 6 and for the background samples, we use MG5_aMC@NLO Alwall:2014hca . The showering and hadronisation is performed within the Pythia 6 Sjostrand:2001yu framework. The $b$-tagging efficiency and mistag efficiencies of a $c$-jet or a light jet posing as a $b$-tagged jet are employed as functions of the transverse momentum of the jet Sirunyan:2017ezt . Unless explicitly mentioned, the CTEQ6l PDF set has been used throughout this work. Also, to take into account the detector effects, we use the fast-detector simulation package, Delphes-3.4.1 deFavereau:2013fsa .

2.1 The $b\bar{b}\gamma\gamma$ Channel

The $b\bar{b}\gamma\gamma$ final state is the golden channel when it comes to studying the non-resonant double Higgs production. The cleanliness of this channel, owing to smaller backgrounds, triumphs over the reduced rate ($Br(h\to\gamma\gamma)\sim 0.2\%$). Here however, we turn to a resonant scalar production which decays to a pair of SM-like Higgs bosons. Our goal is to ascertain the reach of the HL-LHC in measuring $\sigma(pp\to H\to hh)$ for a range of heavy Higgs masses ($m_{H}$). One of the reasons for this final state being a favourite is that the reconstruction and identification precision of photons at the LHC is very high.

Even though the signal seems to have a clean final state, there are several backgrounds at play which need to be dealt with carefully. The major backgrounds (Fig. 2) typically have the form of $hh+X$ which includes the SM double Higgs production, $h+X$ which includes $Zh$, $hb\bar{b}$ and $t\bar{t}h$, and null-Higgs processes like $t\bar{t}+t\bar{t}\gamma$ where leptons may fake as photons, $b\bar{b}\gamma\gamma+c\bar{c}\gamma\gamma+jj\gamma\gamma$ (henceforth termed as $b\bar{b}\gamma\gamma^{*}$) where for the latter two, the light-jets may fake $b$-jets. Other fake backgrounds include $b\bar{b}j\gamma+c\bar{c}j\gamma$ (we will refer to it as Fake 1 category), $b\bar{b}jj$ (referred to as the Fake 2 category), where the $c$-jets may pose as $b$-jets and the light jet may mimic a photon, and the single Higgs processes, viz., $hjj+hc\bar{c}$ (classified henceforth as the $hjj^{*}$ category), where the light-jets and $c$-jets may mimic $b$-jets. One of the major differences between most of the backgrounds and the signal lies in the invariant mass distribution of the $b$-jets, $m_{b\bar{b}}$. However, even when the $m_{b\bar{b}}$ distribution of the signal (as well as the non-resonant SM di-Higgs background) peaks around the SM-like Higgs mass, $m_{h}$, it is broad and can have considerable overlap with the $m_{b\bar{b}}$ distribution either ensuing from a $Z$-boson or from a continuum. It should be noted that the most dominant backgrounds come from the QCD-QED $b\bar{b}\gamma\gamma^{*}$, $t\bar{t}h$ and SM-like di-Higgs processes. The former being a continuum, covers a large part of the kinematic variable space with the signal. The SM-like di-Higgs background also has very large overlap with the signal. One of the easiest way to break this degeneracy is to utilise the $m_{b\bar{b}\gamma\gamma}$ or reconstructed $m_{hh}$ distribution which has a clear peak around the heavy scalar mass, $m_{H}$, for the signal.

We generate the QCD-QED $b\bar{b}\gamma\gamma$ and $Z\gamma\gamma\to b\bar{b}\gamma\gamma$ backgrounds upon merging with one additional jet. We employ the MLM merging scheme Mangano:2006rw where the extra jet contains gluon, light quarks, $c$- as well as $b$-quarks. Among the $h+X$ category, the $Zh$ is generated with the Higgs boson decaying to a pair of photons and the $Z$-boson decaying to a pair of bottom quarks. Furthermore, the $t\bar{t}h$ and $b\bar{b}h$ backgrounds are generated with $h\to\gamma\gamma$. The major fake backgrounds with jets in the final state are generated with the aforementioned jet definition with one exception. We define both of the jets in the $jj\gamma\gamma$ channel in a way as to have no overlap with the $b\bar{b}\gamma\gamma$ background. In case of the $t\bar{t}+X$ backgrounds, we generate the $t\bar{t}$ events with both of the top quarks decaying leptonically which ultimately fakes as photons. However, for the $t\bar{t}\gamma$ background, we require one of the tops to decay leptonically and the other hadronically. Finally, we generate separate single Higgs backgrounds via gluon fusion in association with $c$-quarks and also with light jets. The separation between the $hc\bar{c}$ and $hjj$ backgrounds are necessary in order to appropriately take into account the different fake rates for $c\to b$ and $j\to b$. All of these backgrounds are generated with specific cuts at the generation level which we summarise in Appendix A.

The idea of this section is to understand the reach of the HL-LHC in constraining models with extended scalar sectors. We thus employ optimised search strategies for a varied range of scalar masses. We vary $m_{H}$ in the mass range 275 GeV and 1 TeV. Specifically, we consider the following benchmark points, viz., $m_{H}=275,300,350,400,450,500,550,600,800$ GeV and 1 TeV. In line with our previous work Adhikary:2017jtu , we first perform a classical cut and count analysis to gauge the sensitivity of various benchmark points. We closely follow various cuts from the ATLAS projection study ATL-PHYS-PUB-2017-001 . Namely, we require exactly two $b$-tagged jets and two photons in the final state. The photons are required to have transverse momenta, $p_{T}>10$ GeV and a pseudorapidity coverage of $|\eta|<2.5$. Moreover, the two photons are also required to lie within the pseudorapidity range, $|\eta_{\gamma}|<1.37$ (barrel region) or $1.52<|\eta_{\gamma}|<2.37$ (endcap region). After imposing these basic requirements, we apply some stronger selection cuts in order to enhance the signal to background ratio, $S/B$. We require the invariant mass of the pair of photons, $m_{\gamma\gamma}$, to reconstruct sharply about the SM-like Higgs mass in the range (122,128) GeV. Furthermore, we veto events containing lepton(s) in the final state in order to reduce the impact of the $t\bar{t}h$ background when it decays semi-leptonically or leptonically. We also impose lower bounds on the transverse momenta of the leading and sub-leading $b$-jets and photons. Moreover, upon inspecting the distribution of $\Delta R_{\gamma\gamma}$ and $\Delta R_{b\bar{b}}$ (Fig. 3), we find that with larger values of $m_{H}$, the SM-like Higgs bosons are more boosted yielding more collimated final states.

We thus require $\Delta R_{\gamma\gamma}$ and $\Delta R_{b\bar{b}}$ to lie in the range (0.4,3.0), (0.4,2.0) and (0.4,1.5) for $m_{H}=275,300$ and 350 GeV, $m_{H}=400,450,500,550$ and 600 GeV, and $m_{H}=800$ GeV and 1 TeV respectively. Besides, we require the invariant mass of the $b$-jets, $m_{b\bar{b}}$ to lie in the range (90,130) GeV. This choice is related to account for the jet-energy correction and has been described in Ref. Adhikary:2017jtu . We summarise these cuts in Table 2. As the next logical step, we delve deeper into the kinematics. We reconstruct the invariant mass of the $b\bar{b}\gamma\gamma$ system, $m_{b\bar{b}\gamma\gamma}$ and its total visible energy, $H_{T}$. These two variables are intrinsically correlated. Also, from Fig. 4, it is evident that the $H_{T}$ distribution is broader leading to more background contamination as compared to the $m_{b\bar{b}\gamma\gamma}$ distribution which we show in Fig. 6 (left). Thus, we proceed with $m_{b\bar{b}\gamma\gamma}$ in order to further optimise our analysis.

We also reconstruct the transverse momentum of the SM-like Higgs decaying to a pair of photons, $p_{T,\gamma\gamma}$. As can be seen from the $p_{T,\gamma\gamma}$ distribution in Fig. 6, the spectrum is harder for heavier values of $m_{H}$. We choose $p_{T,\gamma\gamma}>50$ GeV for $m_{H}=275$ and 300 GeV and for all other masses, we choose the transverse momentum of this reconstructed Higgs to be larger than 100 GeV. Thus, after these fixed cuts, we perform a simplified optimisation with the $m_{b\bar{b}\gamma\gamma}$ variable in order to enhance $S/B$. These cuts are finally tabulated in Table 3 where we also present the signal efficiency and the background yield at an integrated luminosity, $\mathcal{L}=3$ ab*-1* with $m_{H}$ being varied. The signal efficiency, $\epsilon$, here points to the ratio of the total number of signal events remaining after all the cuts applied in sequence to the generated number of events. The second column refers to the range of $m_{b\bar{b}\gamma\gamma}$ that optimises the signal and the third column denotes the minimum $p_{T}$ for the diphoton system. This optimisation is different as compared to the SM di-Higgs production scenario as shown in Ref. Adhikary:2017jtu . This difference is related to the kinematics of the event topology. As an example, the $p_{T,\gamma\gamma}$ distribution changes with the heavy Higgs mass as shown in Fig. 5. From this figure, it is evident that different optimisation is required for each mass point and also for the SM scenario. We must note that in Table 3, the choices of the upper and lower cuts on $m_{b\bar{b}\gamma\gamma}$ for the different heavy Higgs masses can be understood from Fig. 7. These ranges are obtained after optimising for each value of $m_{H}$. Finally, we provide a detailed cut-flow table for $m_{H}=400$ GeV in Table 4 with a complete information for the backgrounds including the perturbative order at which the rates are considered as well as the number of events remaining at the HL-LHC.

Utilising these results, we derive the projected upper limits on the production cross section of the heavy Higgs in a model independent manner 555We consider the cut-based optimisations as final as we did not obtain any observable improvement with a multivariate analysis.. We calculate the cross-section reach by using the significance formula: $S/\sqrt{B}=N$, where $N$ denotes the number of confidence intervals. Here, the signal yield, $S$, is defined as $\sigma(pp\to H\to hh\to b\bar{b}\gamma\gamma)\times\mathcal{L}\times\epsilon$ and $B$ represents the total background yield after the cut-based analysis. With this, we derive $\sigma(pp\to H\to hh)$ at the $N\sigma$ level, with $N=2~{}\textrm{and}~{}5$, respectively, corresponding to a $95\%$ and $99.7\%$ confidence level (CL) upper limit, also referred to as the exclusion and discovery limits, respectively. We show the final results in Fig. 8 with the upper limit on $\sigma(pp\to H\to hh)$ as a function of $m_{H}$ and we display the $2\sigma$ and $5\sigma$ lines. The $2\sigma$ upper limit is strong between $400$ GeV and 1 TeV, varying between $31.74$ fb and $4.24$ fb. Upon adding $5\%$ systematic uncertainty 666The significance formula for a systematic uncertainty of $N\%$ has the form: $S/\sqrt{B+\left(N*B/100\right)^{2}}$, the upper limit becomes $32.35$ fb and $4.25$ fb respectively within the previously mentioned mass range. It must be noted that the upper limit is mildly affected by incorporating a systematic uncertainty of $\sim$5%. The reason can be attributed to the fact that the signal over background ratio ($S/B$) is high.

Next, we perform a multivariate analysis in order to improve upon the cut-based analysis. We use the following variables:

[TABLE]

Here, the variable names have their usual meaning. The $\Delta R_{b_{i}\gamma_{j}}$ is the distance in the $\eta-\phi$ plane between the $b$-jets and photons with $i~{}{\rm and}~{}j=1,2$. Also, $\Delta R_{bb,\gamma\gamma}$ is the $\Delta R$ separation between the system of $b$-jets and the two photon system. However, after performing this analysis we obtain comparable results. Thus, we do not show the results of the multivariate analysis in this section.

2.1.1 The $pp\to A\to Zh$ Channel

With the accumulation of more data, we are on the brink of accepting the fact that in MSSM or in generic two Higgs doublet models, the SM-like Higgs is in the decoupling regime with its coupling to the SM gauge bosons being proportional to $\sin(\beta-\alpha)$, with $\alpha$ and $\tan{\beta}$ being the mixing angle in the neutral $CP$-even sector and the ratio of the two vacuum expectation values of the two doublets, respectively. For the up and down type fermions, the Yukawa couplings for the SM-like Higgs boson are proportional to $\cos{\alpha}/\sin{\beta}$ and $\sin{\alpha}/\cos{\beta}$ respectively. In the decoupling regime, $\sin(\beta-\alpha)\sim 1$ and hence the decay width of $A\to Zh$ which is proportional to the coupling $\cos(\beta-\alpha)$, is small. In a non-decoupling regime, $pp\to A\to Zh$ can give us deep insight into two scalars simultaneously.

The CMS CMS-PAS-FTR-13-024 collaboration has derived projected upper limits on $\sigma(pp\to A\to Zh)$ from searches in the $\ell\ell b\bar{b}$ ($Z\to\ell\ell,~{}h\to b\bar{b}$) final state for HL-LHC. In the present study, keeping in continuation to the analysis prescribed in the previous section (Section 2.1), we explore the prospects of directly probing $A$ in the $A\to Zh\to b\bar{b}\gamma\gamma$ channel in the context of HL-LHC.

We will remain agnostic to the fact that the prospects of observing $A\to Zh$ in the decoupling regime may be extremely small. The main difference here is the fact that in the previous analysis both the diphoton and the $b\bar{b}$ pairs are required to peak around the SM-like Higgs boson mass, whereas in the present case, the $b$-jets are required to peak around the $Z$-boson mass. We follow a cut-based analysis as before and optimise the $m_{b\bar{b}\gamma\gamma}$ and $p_{T,\gamma\gamma}$ cuts for different values of $m_{A}$. These variables are shown to have substantial discriminatory power and are shown in Fig. 9 for $m_{A}=400$ GeV and $600$ GeV. Details of the fixed and optimised cuts are presented in Table 5 and Table 6 respectively. In Table 7, we show the cut-flow table for $m_{A}=400$ GeV. After a full optimisation, we show the 95% and 99.7% CL exclusion for $\sigma(pp\to A\to Zh)$ in Fig. 10. The bounds are weaker than their $H\to hh$ counterpart mainly because of a larger overlap with the $Zh$ background. It is to be noted that the projected upper limits on $\sigma(pp\to A\to Zh)$ derived from searches in the $b\bar{b}\gamma\gamma$ channel (see Figure 10) are weaker than the projected limit obtained by CMS from searches in the $\ell\ell b\bar{b}$ CMS-PAS-FTR-13-024 final state.

2.2 The $b\bar{b}b\bar{b}$ Channel

After having studied the cleanest possible di-Higgs final state, we turn our attention to the one with the largest rate, viz., $pp\to H\to hh\to 4b$. Several searches have already been conducted in this channel Aaboud:2016xco ; Aaboud:2018knk ; CMS-PAS-HIG-17-009 ; CMS-PAS-HIG-16-026 and provide some of the strongest bounds in both the non-resonant and resonant sectors. ATLAS Aaboud:2018knk has computed the observed (expected) upper bound on $\sigma(pp\to hh\to 4b)$ to be 13 (20.7) times that of the SM expectation with an integrated luminosity of 27.5 fb*-1*. This channel has been further combined in ATLAS’ HL-LHC projection ATL-PHYS-PUB-2018-053 alongside $pp\to hh\to b\bar{b}\gamma\gamma/b\bar{b}\tau\tau$. The above channel will gain an improvement in sensitivity according to the Ref. ATL-PHYS-PUB-2018-053 because of a projected 8% improvement in $b$-tagging efficiency, besides having larger data sets. In this work, we closely follow the search strategy used by the ATLAS collaboration in Ref. Aaboud:2016xco . Even though this channel has the highest signal rate, the enormous multijet and $t\bar{t}$ backgrounds may considerably overwhelm the signal yield.

The dominant backgrounds (Fig. 11) include the multijet production from QCD processes and the top pair production. For the multijet backgrounds, we dissect our sample generation into three different categories each having at least two $b$-quarks, viz., $b\bar{b}b\bar{b}$, $b\bar{b}c\bar{c}$ and $b\bar{b}jj$, in order to have sufficient statistics to take into account the different tagging efficiencies and fake rates. We do not generate the $h$ + jets and $Z$ + jets backgrounds separately but we include their tree-level diagrams while generating the multijet backgrounds as they have negligible contributions Aaboud:2016xco . We must mention here that we do not consider other possible sources of multijet production viz., $c\bar{c}c\bar{c}$, $c\bar{c}jj$ etc. as these processes will be highly suppressed (with respect to $b\bar{b}b\bar{b}$) upon multiplying by the fake efficiency factors, in succession. We generate the $t\bar{t}$ background with the top quark decaying to a $b$-quark and a $W$-boson. The $W$-bosons are then further decayed to $c\bar{s}$ or $\bar{c}s$. We avoid the $W\to u\bar{d}$ mode as the probability of a light jet faking a $b$-jet is $\sim 10$ times smaller than that of a $c$-jet. Lastly, we also consider the subdominant background coming from the non-resonant di-Higgs production ($gg\to hh$) and also from $t\bar{t}b\bar{b}$ (including $t\bar{t}Z/t\bar{t}h$).

We select events containing exactly 4 $b$-tagged jets with the requirement of $p_{T,b}>60~{}\text{GeV}$ and $|\eta_{b}|<2.5$. The scalar sum, $H_{T}$, of the transverse momenta of all the visible particles in an event must fulfil, $H_{T}>300$ GeV. Finally, we form two di-jet systems from these four $b$-jets. The two jets within a dijet system must satisfy $0.4<\Delta R_{bb}<1.5$. We choose the leading (sub-leading) di-jet system to have $p_{T}>200\;(150)$ GeV 111111These are preliminary cuts before performing the multivariate analysis. The rationale behind these cuts are (a) some of these cuts have been applied at the generation level on some of the backgrounds in order to have better control over event statistics owing to large production cross-sections (see Appendix A) and (b) some other cuts are applied by observing the kinematic distributions of these observables. However, we apply stronger cuts on these variables in the following where we optimise them alongside other correlated variables through a more sophisticated multivariate analysis.. Furthermore, to reduce the contamination from the $t\bar{t}$ background, we reconstruct the top by combining extra jets in an event with the di-jet systems. These jets must be within $\Delta R<1.5$ in the $\eta-\phi$ plane with the di-jet system. If an event contains exactly one extra jet, then we choose the di-jet system which is closest to it and combine to form a top quark system, $m_{t_{1}}$. However, when there are two such jets, we compute the minimum of all possible $\Delta R$ combinations between these two jets and the two di-jets before reconstructing two other top masses, $m_{t_{2}}$ and $m_{t_{3}}$. Because for our signal, we do not expect any proper top quark reconstruction, we thus veto events if the reconstructed mass of any of these possible choices for the top quark exceeds 120 GeV. After imposing this cut the $t\bar{t}$ background reduces to half with more than $80\%$ of the signal events still to spare. We detail these cuts one by one alongside the signal efficiency and cross-sections for the background processes in Table 8.

Finally, after all the aforementioned cuts are applied in succession, we check for any possible improvement upon performing a multivariate analysis. We utilise the BDT algorithm for our purposes and choose the following nine kinematic variables with maximal potency,

[TABLE]

Here we use the kinematic variables reconstructed from the two di-jet systems viz., invariant mass ($m_{\textrm{di-jet}}$), transverse momentum ($p_{T,\textrm{di-jet}}$) and azimuthal angle separation between the $b$-jets forming the dijet systems ($\Delta\phi_{bb,\textrm{di-jet}}$). The subscript $k=1,2$ refers to the $p_{T}$ ordering of the di-jets. We also take the separation in the $\eta$ and $\eta$-$\phi$ plane between the two di-jets, viz., $\Delta\eta_{\textrm{di-jets}}$ and $\Delta R_{\textrm{di-jets}}$ respectively. $m_{b\bar{b}b\bar{b}}$ is the invariant mass of the four $b$-jet system. The top four variables with the highest discriminatory power are shown in Fig. 12. We can see that the lower masses have significantly longer tails while performing the mass reconstructions.

Finally, in Table 9, we present the background yields after the BDT optimisation has been completed. Like in the previous section, we translate these results into an exclusion diagram showing the upper limits on $\sigma(pp\to H\to hh)$ as a function of the heavy Higgs mass. We show these in Fig. 13. The limit is very strong between $600$ GeV and $1$ TeV with the 95% CL upper limit varying between $15.26$ fb and $2.51$ fb. The upper limit becomes between $[82.70,~{}5.77]$ fb within the aforementioned range by adding $5\%$ systematic uncertainty,.

Discussion about $m_{H}=400$ GeV

The $95\%$ and $99.7\%$ CL upper limits on $\sigma(pp\to H\to hh)$ for the heavy Higgs with a mass around $400$ GeV is very large (256.66 fb and 641.65 fb respectively) as compared to the other mass points, even after the BDT optimisation. The reason for this is the following. The signal efficiency for $m_{H}=400$ GeV reduces by $\sim$ 67% after imposing the $\Delta R_{bb}$ selection as can be seen from Table 8. Since the heavy Higgs mass ($400$ GeV) is near the threshold of the non-resonant di-Higgs production, the SM-like Higgs bosons for the resonant case are produced with low $p_{T}$. This further leads to the Higgs decay products being widely separated in the $\eta-\phi$ plane and thus obviously does not satisfy our di-jet selection criteria of $\Delta R_{bb}<1.5$ within each di-jet system. With the sole intention of improving the sensitivity, we adopt a $\chi^{2}$ minimisation technique as described below. We define a new kinematic variable $\chi^{2}_{hh}$ for the events which do not satisfy the $\Delta R_{bb}<1.5$ selection criteria as follows

[TABLE]

where $m_{h}=125$ GeV and $\sigma_{hj}=0.1\times m_{\textrm{di-jet},j}$ with $j=1,2$ marks the $p_{T}$ ordering. Thus, in addition to the events satisfying $\Delta R_{bb}<1.5$, we also consider those events which contain di-jet pairs with $\Delta R$ separation between the $b$-jets to be more than $1.5$. Following this, we construct the aforementioned $\chi^{2}_{hh}$ variable for each possible pair of reconstructed di-jet. The event is finally selected if the non-zero minimum value of the $\chi^{2}_{hh}$ variable is less than 50 121212We checked our results upon choosing $\chi^{2}_{hh}$ both higher and lower than 50. For higher values, the signal yield increases but the background increases at a higher rate, generating a weaker limit. On the other hand, upon lowering the $\chi^{2}_{hh}$ value below 50, the signal yield decreases substantially. This makes the value 50 an ideal and optimal choice.. Upon using this modification, the signal efficiency increases by $\sim 26\%$ at the di-jet selection level while simultaneously increasing the dominant backgrounds like $b\bar{b}b\bar{b}$ by $\sim 5\%$ and $t\bar{t}$ by $\sim 7\%$. However, the limit on the upper limit of the cross-section improves to $245.75$ fb and $614.37$ fb at 95% and 99.7% CL respectively.

2.3 The $b\bar{b}\tau\tau$ Channel

Next, we turn our attention to one of the best probes for the di-Higgs searches, viz., the $b\bar{b}\tau\tau$ channel. The intricacy and potential of this channel lies in our ability to reconstruct the $\tau$-leptons as these come with neutrinos which show up as missing transverse energy in the detector. This channel gives rise to three phenomenologically different final states, viz., $b\bar{b}\ell^{+}\ell^{-}+\not{E}_{T}$, $b\bar{b}\ell\tau_{h}+\not{E}_{T}$ and $b\bar{b}\tau_{h}\tau_{h}+\not{E}_{T}$. In this work, we will only consider the last category, i.e., the one with the fully hadronic $\tau$ decays. The hadronically decayed $\tau$-leptons are termed as $\tau$-hadrons or $\tau$-jets which may either contain one (one-pronged) or three (three-pronged) charged particle(s) inside the jet cone. Thus, it is essential to tag these $\tau$-jets in order to segregate them from regular QCD jets ensuing from the various backgrounds that we will discuss below. We will not discuss the fully leptonic case here as from our previous analysis Adhikary:2017jtu we know that the sensitivity is extremely low even at the HL-LHC.

We generate two different samples for the dominant $t\bar{t}$ background (Fig. 14), where either both the $W$-bosons decay to jets or where one decays into a lepton ($viz.~{}e^{\pm},\mu^{\pm}~{}or~{}\tau^{\pm}$) and the other to a pair of jets. The QCD-QED background $\tau\tau b\bar{b}$ also contributes significantly. Besides, we also generate the subdominant backgrounds which include $t\bar{t}h$, $t\bar{t}W$, $t\bar{t}Z$, $b\bar{b}h$, $Zh$ and the non-resonant Higgs pair production i.e., $gg\to hh$. We simulate the $Zh$ background upon considering two processes where in one case the $Z$-boson decays to a pair of bottom quarks and the Higgs boson decays to a pair of $\tau$-leptons and in the other the decays are reversed. Moreover, we also generate the dominant fake background for the hadronic channel in the form of $b\bar{b}jj$ (Fig. 14), where the light jets can be fake $\tau$-tagged jets131313The $\tau$-leptons decay hadronically (each with a branching fraction of $\sim 65\%$) leading to jets in the final state. In our analyses, we use the $\tau$-tagging method as discussed in Ref. Bagliesi:2007qx .. We detail the generation level cuts for these various backgrounds in Appendix A. Following the generation level cuts, we further apply some basic cuts on the signal and background samples in order to ensure a common kinematic phase space. The $b/\tau$-jets and the leptons ($e,\mu$) are required to lie within $|\eta|<2.5$ and have $p_{T,b/\tau(\ell)}>20\;(10)$ GeV. The light jets must satisfy $p_{T,j}>20$ GeV and $|\eta_{j}|<4.5$. The minimum distance in the $\eta-\phi$ plane between the $b$-jets and the leptons, and also among themselves is required to be $\Delta R>0.2$. The reconstructed invariant mass of the bottom pair and the visible $\tau$ pair must obey $m_{bb/\tau\tau}>50$ GeV.

2.3.1 The $b\bar{b}\tau_{h}\tau_{h}$ Channel

In this sub-section we briefly outline the prospects of searching for the heavy Higgs in the $b\bar{b}\tau_{h}\tau_{h}$ final state. In doing so, we select events containing exactly two $b$-tagged jets and two $\tau$-tagged jets alongside the cuts described above. Having seen the strength of the multivariate analyses for this channel, in Ref. Adhikary:2017jtu , rather than opting for the classical cut-based analysis, we perform a BDT analysis with the following 13 variables with the maximal discerning capability:

[TABLE]

where, $M_{T}$ is the transverse mass of the $h\to\tau\tau$ system 141414For this whole section, we use the conventional definition Han:2005mu of $M_{T}^{2}=(\sum\limits_{i}E_{T,i})^{2}-(\sum\limits_{i}\vec{p}_{T,i})^{2}$, where $i$ runs over the relevant objects. In Section 4, we use a modified definition of $M_{T}$., $p_{T,\textrm{tot}}$ and $m_{\textrm{tot}}$ are respectively the transverse momenta and mass of the full visible system and $m_{\textrm{eff}}$ is the scalar sum of the transverse mass of all the visible products plus $\not{E}_{T}$. The rest of the variables have usual definitions. The top five variables are shown in fig. 15. As can be seen, the $m_{T2}$ variable is particularly useful for heavier Higgs masses as it can be used to completely eradicate the $t\bar{t}$ background. We train the signal and background samples and they are optimised for each benchmark signal point. We list the background events after optimising the BDT and imposing the cut for four values of $m_{H}$, in Table 10. Finally, we show the upper limit on the heavy Higgs production cross-section (assuming BR($H\to hh$) $=100\%$) in Fig. 16. The $95\%$ CL upper limit on the cross-section between $m_{H}=$ 600 GeV and 1 TeV varies between $67.34$ fb and $39.56$ fb. With $5\%$ systematic uncertainty, the limits become $251.76$ fb and $116.30$ fb respectively.

2.4 The $b\bar{b}WW^{*}$ Channel

In this section, we consider the situation where a heavy scalar decays to a pair of SM-like Higgs bosons with one of them decaying to a pair of $b$-quarks and the other to $WW^{*}$, leading to three possible final states depending on the decays of the $W$-bosons. We perform our analyses for the fully leptonic (leptons at this stage include $e,\mu,\tau$) and the semi-leptonic channels. We avoid studying the fully hadronic mode as the signal will be overwhelmed by the huge QCD background.

The dominant contribution to the background (Fig. 17) for both the channels mentioned above comes from top pair production because of its large production cross-section. We generate this background where either or both the $W$-bosons decay leptonically. The fully hadronic $t\bar{t}$ mode is not considered as a potential background as the fake rate for $j\to\ell$ is negligible for all practical purposes. The fully leptonic $t\bar{t}$ background contributes to the fully leptonic channel final state whereas for the semi-leptonic scenario, the contribution comes from both the fully-leptonic as well as the semi-leptonic $t\bar{t}$. The second most dominant background for the semi-leptonic channel is $Wb\bar{b}~{}+$ jets, where the $W$-boson decays leptonically ($e,\mu,\tau$). We generate this background upon merging with two additional jets by exploiting the MLM merging scheme Mangano:2006rw . While generating the $Wb\bar{b}~{}+$ jets background we ensure that there is no double counting ensuing from the semi-leptonic $t\bar{t}$ background. Besides the aforementioned backgrounds, a significant contribution also comes from the $\ell^{+}\ell^{-}b\bar{b}$ production where $\ell$ refers to $e,\mu$ and $\tau$. Finally, we also consider the subdominant backgrounds viz., $t\bar{t}h$, $t\bar{t}Z$, $t\bar{t}W$ and the non-resonant $gg\to hh$.

In this subsection and the following section (section 3), the top-pair production is the dominant background. Thus, the reconstruction of the top quarks is a very powerful tool in order to reduce the contribution from this background. For the semi-leptonic case, the only source of missing transverse energy, $\not{E}_{T}$ 151515To incorporate the $\not{E}_{T}$ smearing, we use the standard module of the Delphes ATLAS card., arises from the neutrino of the leptonically decaying $W$-boson from the top decay. We reconstruct the top from its decay products 161616First the $W$-boson mass is reconstructed in order to attain the $p_{Z}$ component of the neutrino.. The quadratic equation gives two possible solutions for the neutrino $p_{z}$. Besides, because there are two $b$-jets in the final state, we get four possible choices for the mass of the leptonically decaying top. We use these variables during our analysis. After reconstructing both the tops, we reconstruct the total system from all the final state particles. We also use this later in section 3 which exhibits the same final state. These variables help us greatly in reducing the semi-leptonic $t\bar{t}$ background for high values of $m_{H}$.

Before embarking on the final analysis, we impose a common set of trigger cuts for both the leptonic and the semi-leptonic channels. The $p_{T}$, $|\eta|$ and $\Delta R$ cuts for the various objects are discussed in subsection 2.3 and also in Appendix A. Furthermore, we require generation-level cuts on the invariant mass of b-jets, viz., $m_{bb}>50$ GeV. The selected events are also require to have $\not{E}_{T}>40$ GeV upon scrutinising the distribution. The $\not{E}_{T}$ distribution for the 1$\ell$ and 2$\ell$ cases are shown in Fig. 18. Finally, we perform separate multivariate analyses for the two final states upon using the BDTD algorithm. While training samples for both the leptonic and semi-leptonic analyses, we only consider the $t\bar{t}$ background since it constitutes the bulk of the total background. This training is used for testing all other backgrounds which are subdominant in front of $t\bar{t}$.

2.4.1 The $2\ell 2b+\not{E}_{T}$ Channel

For the fully leptonic final state, we select events with exactly two $b$-tagged jets, and two isolated leptons having opposite charge meeting the trigger criteria as mentioned above. We choose the following set of kinematic variables in order to perform the multivariate analysis:

[TABLE]

where, $M_{T}$ is the transverse mass of the SM-like Higgs decaying to $W$-bosons. The rest of the variables have either been defined before or have usual meaning. The top four variables are shown in Fig. 19. The signal distributions are significantly different from the various backgrounds.

Finally, in Table 11, we summarise the number of background events after imposing the optimised cut on the BDT variable. Like in the other channels, we impose an upper limit on $\sigma(pp\to H\to hh)$ as a function of the heavy Higgs mass. This is shown in Fig. 20. The $95\%$ CL upper limit varies between $67.41$ fb and $26.18$ fb ($357.51$ fb and $82.21$ fb with $5\%$ systematic uncertainty) within 600 GeV $<m_{H}<$ 1 TeV and is somewhat weaker than the channels discussed earlier owing to smaller $S/B$.

2.4.2 The $1\ell 2b2j+\not{E}_{T}$ Channel

Finally, we discuss the potential of the semi-leptonic final state as well. We require events with exactly two $b$-tagged jets, one isolated lepton and at least two light jets satisfying the trigger criteria discussed earlier. Besides, we consider the same set of cuts as for the dileptonic channel before performing the multivariate analysis. We find the following kinematic variables to have the best discriminatory power and use them for our multivariate analysis:

[TABLE]

where, $m_{t}$ is the transverse mass of the leptonically decaying $W$-boson. $\Delta R_{\ell,jj}$ is the distance in the $\eta-\phi$ plane between the system comprising of the two hardest jets and the lepton. $m_{bbj_{1}}$ refers to the invariant mass of the two $b$-tagged jets and the hardest $p_{T}$ jet. The reconstructed transverse momentum of the leptonically decaying $W$-boson is denoted as $p_{T,\ell\nu}$. $m_{tij}$ is the mass of the leptonically decaying top quark with the reconstruction procedure outlined before. The first index $i=1,2$ indicates the $p_{T}$ ordering of the $b$-jet. The second index $j=1,2$ refers to the choice of the $z$-component of the neutrino momentum. The other variables have usual definitions. The best discriminatory variables are listed in Fig. 21. However, we can see that the separation power for the 1$\ell$ category is significantly less compared to its 2$\ell$ counterpart.

Coming to the results, table 12 summarises the background yields after the BDT cut. The upper limit on $\sigma(pp\to H\to hh)$ as a function of $m_{H}$ are shown in Fig. 22. The limits are considerably weak in this channel.

2.5 The $\gamma\gamma WW^{*}$ Channel

After the $b\bar{b}\gamma\gamma$ channel this is the second most cleanest channel in terms of the final state particles but with the pitfall of having very low event rate. In this channel, one of the SM-like Higgs decays to a pair of photons and the other to lepton(s) through $h\to WW^{*}$. Similar to the $b\bar{b}WW^{*}$ analysis in subsection 2.4, here also we divide the channel into the leptonic and semi-leptonic category. Because of the relatively clean final states, these channels have low contaminations due to backgrounds. We simulate the $Zh$ and $Wh$ backgrounds upon merging with two additional jets (the definition of jet is given in subsection 2.1). Here we decay the $Z$- and the $W$-bosons leptonically ($e,\mu,\tau$). The $Wh$ channel contributes only to the semi-leptonic category. Besides, there are $\ell\nu\gamma\gamma$ and $\ell\ell\gamma\gamma$ productions coming from $\gamma^{*}$, diagrams containing triple and quartic gauge boson interactions and various other diagrams not involving the Higgs. We generate these two backgrounds upon merging with an additional jet and using the same scheme as before. Next, we also consider the $t\bar{t}h$ background with Higgs-boson decayed to a pair of photons. Finally, we also consider the SM Higgs pair production which is subdominant. We show the Feynman diagram of the dominant backgrounds in Fig 23.

Before performing the multivariate analyses, we impose the generic trigger cuts. The $p_{T}$, $|\eta|$ and $\Delta R$ 171717$\Delta R_{\gamma\gamma/\gamma\ell}>0.4$ and $\Delta R_{\ell\ell>0.2}$. cuts are the same as has been defined in subsection 2.3. The above cuts for the photons are the same as those on the leptons. Owing to an excellent resolution for the diphoton invariant mass, we require $122~{}\textrm{GeV}<m_{\gamma\gamma}<128~{}\textrm{GeV}$. Finally, we also require $m_{\ell\ell}>20$ GeV because we generate the $\ell\ell\gamma\gamma$ background with this invariant mass cut at the generation level (the details of these cuts are mentioned in Appendix A). We now describe the results of the multivariate analyses for these two channels in the following two subsections.

2.5.1 The $\gamma\gamma 1\ell 2j+\not{E}_{T}$ Channel

Before performing the BDT analysis, we select events with exactly two isolated photons, one isolated lepton and at least two jets in the final state, which fulfils all the aforementioned trigger requirements. Like all the other channels, we consider the following variables to train our signal and background samples for the multivariate analysis:

[TABLE]

where the variables carry their usual meaning. The five best variables are shown in Fig. 24. The background yields after the BDT optimisation are shown in Table 13. In Fig. 25, we show the upper limit on $\sigma(pp\to H\to hh)$ as a function of $m_{H}$. The $95\%$ CL upper limit changes from $220.11$ fb ($226.74$ fb with $5\%$ systematic) for $m_{H}=400$ GeV to $112.20$ fb ($113.24$ fb with $5\%$ systematic) for $m_{H}=1$ TeV.

2.5.2 The $\gamma\gamma 2\ell+\not{E}_{T}$ Channel

This is the final channel that we study for the $pp\to H\to hh$ case. We choose events with exactly two isolated photons, and two isolated leptons with opposite charge, following the trigger cuts mentioned above. Finally, we choose the following kinematic variables for the multivariate analysis:

[TABLE]

with the usual definitions for the variables. Some of the variables of interest are plotted in Fig. 26. The background yields after the BDT cut are tabulated in Table 14 whereas the upper limit on $\sigma(pp\to H\to hh)$ as a function of heavy Higgs mass is shown in Fig. 27. The 95% CL upper limit for the leptonic scenario is stronger than its semi-leptonic counterpart in the heavy Higgs mass range of 600 GeV and 1 TeV. The upper limit varies in between $109.80$ fb and $56.30$ fb ($110.24$ fb and $56.30$ fb with $5\%$ systematic uncertainty) at 95% CL, in the aforementioned range.

2.6 Summarising the $H\to hh$ channel

Having studied five different channels with more than one sub-processes in three instances, we summarise the results in this subsection. The 95% CL upper limits on $\sigma(pp\to H\to hh)$ for all these channels is shown in Fig. 28. We find that the strongest limits come from the $b\bar{b}\gamma\gamma$ and $4b$ channels. The $b\bar{b}\gamma\gamma$ is strongest up to $m_{H}\sim$ 600 GeV. From 600 GeV onward, the $4b$ channel is more constraining owing to its larger cross-section. The present limits on $\sigma(pp\to H\to hh)$ from the 13 TeV analysis vary between $\sim$ 970 fb (300 GeV) and $\sim$ 225 fb (600 GeV) from the $b\bar{b}\gamma\gamma$ channel Aaboud:2018ftw and between $\sim$ 70 fb (800 GeV) and $\sim$ 25 fb (1 TeV) from the $4b$ analysis Aaboud:2018knk . Our projected limits vary between $79.03$ fb (300 GeV) and $14.10$ fb (600 GeV) from the $b\bar{b}\gamma\gamma$ channel and between $5.36$ fb (800 GeV) and $2.51$ fb (1 TeV) from the $4b$ channel. We find an order of magnitude improvement in the sensitivity. We must note in passing that from the results obtained in Ref. Aaboud:2018ksn , the $H\to hh\to 4W$ category does not show much promise even at the HL-LHC, owing to very small signal yields in all categories. The maximum sensitivity is expected in the $e\mu$ 2-lepton category with 5 expected events.

3 The $pp\to H\to t\bar{t}$ Channel

After having studied the $H\to hh$ in multifarious channels in detail, we now turn our attention to a heavy scalar (or pseudoscalar) resonance being produced predominantly by gluon fusion and decaying to a pair of top quarks (Fig 29). This channel has already gained some attention in the experimental community Aaboud:2017hnm ; Aaboud:2018mjh . Searches for resonant scalars, pseudoscalars, $Z^{\prime}$-bosons, Kaluza-Klein gluons and Kaluza-Klein gravitons have been performed. The aim of this section is to try and improve upon these existing searches and provide potential reach of the $\sigma(H\to t\bar{t})$ by studying the fully leptonic and the semi-leptonic final states. The branching ratio of $t\to bW$ being close to 100% makes the channel essentially become a search for $H\to b\bar{b}W^{+}W^{-}$. However, unlike the $H\to hh\to b\bar{b}WW^{*}$ channel studied in subsection 2.4, where one of the $W$-bosons is off-shell, here both of them are on-shell. This is the first essential difference between the two channels and the reason why one requires a completely different search strategy for the two cases. In the previous section 2, we required BR$(H\to hh)=100\%$. However, in realistic scenarios, if the heavy scalar is produced predominantly via gluon fusion (top/bottom loops), it should also decay to a pair of top quarks (and also bottom quarks) if it is above the $t\bar{t}$ threshold. Similar to the $H\to hh\to b\bar{b}WW^{*}$ channel, here also we divide the analysis into two parts, viz., the leptonic and the semi-leptonic channels. We apply the same trigger-level cuts to the various objects as sketched in subsection 2.4. The backgrounds are the same as before. As before, we implement the production and decay of the heavy scalar in the Pythia 6 framework.

3.1 The leptonic Channel

Like in section 2.4.1, here also we select events with two oppositely charged isolated leptons and two $b$-tagged jets. Without performing a classical cut-based analysis, we optimise our results to obtain the best-possible sensitivity by employing a boosted decision tree analysis. The set of variables which discriminate the signal from the backgrounds are as follows:

[TABLE]

where all the variables have their usual meaning as mentioned earlier. We would like to mention here that we also consider the $m_{T2}$ variable during our analysis. However, this variable is $\sim 80\%$ correlated with $p_{T,bb}$. Moreover, $m_{T2}$ has a lower BDT ranking as compared to $p_{T,bb}$. We explicitly checked that adding this correlated variable does not improve our BDT sensitivity. Thus for this analysis, we do not use $m_{T2}$ ($m_{T2}$ was used in the $b\bar{b}\tau\tau$ (section 2.3) and $b\bar{b}WW^{*}$ (section 2.4) analyses.). The top four discriminatory variables are shown in Fig. 30. In Table 15, the number of background events at an integrated luminosity of 3000 fb*-1*, optimised to maximise the sensitivity for various values of $m_{H}$ and after imposing cuts on the BDT observable, are presented. Like in all the other channels, we present the 95% and 99.7% upper limit on $\sigma(pp\to H\to t\bar{t})$ as a function of $m_{H}$, in Fig. 31. We find that the 95% upper limit on the cross-section lies between $380.43$ fb and $135.25$ fb ($42683.56$ fb and $3940.56$ fb with $5\%$ systematic uncertainty) for $m_{H}$ varying between 400 GeV and 1 TeV.

3.2 The semi-leptonic Channel

We end this section by analysing the semi-leptonic final state ensuing from the semi-leptonic decays of $t\bar{t}$. We select events which contain a single isolated lepton, two $b$-tagged jets and at least two light jets after applying the same set of trigger cuts as discussed in section 2.4. Finally, we perform a multivariate analysis with the following set of kinematic variables:

[TABLE]

where, $m_{t\bar{t}i}$ are the possible combinations for the invariant mass of the heavy Higgs reconstructed from the top pair. The four most sensitive variables are listed in Fig 33. The reconstruction procedure is discussed at the beginning of section 2.4. We show the reconstructed $t\bar{t}$ invariant masses in Fig. 32. Finally, we summarise the boosted decision tree results in Table 16. For heavy Higgs mass ranging between 400 GeV and 1 TeV, we show the upper limit on $\sigma(pp\to H\to t\bar{t})$ in Fig. 34. The $95\%$ CL upper limit varies between $186.57$ fb ($39460.45$ fb) and $32.81$ fb ($2021.51$ fb) for $m_{H}$ varying between 400 GeV and 1 TeV with zero ($5\%$) systematic uncertainty. The $H\to t\bar{t}$ channel has a small $S/B$ ratio. Hence, adding systematic uncertainty will drastically change the upper limit on the cross-section.

4 The $(H/A)b\bar{b}$ channel

Finally, we study the process where the resonant (pseudo)scalar is produced in association with a pair of bottom quarks, viz., $pp\rightarrow(b\bar{b})H/A$. The need to study this process lies in the fact that one can probe and impose strong limits on the lower part in the $m_{A}-\tan{\beta}$ plane, as will be discussed in section 5. The cross-section of the inclusive $(b\bar{b})H$ process receives contribution from both the 4-flavour (4F) and the 5-flavour (5F) processes. There are two QCD processes (Fig. 35) contributing to the 4F scheme at LO where the heavy Higgs is produced in association with two $b$-quarks, one is via the gluon fusion process ($gg\rightarrow b\bar{b}H$) and the other is via quarks ($q\bar{q}\rightarrow b\bar{b}H$). The 4F inclusive cross-section suffers from large logarithms due to an almost collinear splitting of a gluon into a pair of bottom quarks. This is of the form $\textrm{ln}(\frac{\mu_{F}}{m_{b}})$ ($\mu_{F}\equiv$ factorisation scale) and may lead to a breakdown of the perturbative theory. However, these logarithms can be absorbed inside the bottom quark parton distribution function (PDF) by re-summing at all orders in the perturbation theory. This forms the basis of the 5F scheme. At leading order (LO), the 5F scheme is dominated by the QCD process $b\bar{b}\rightarrow H$ (Fig. 35). However, for scenarios involving $b$-jet in the final state, the processes where the resonant scalar is produced in association with a $b$-quark or a gluon, becomes important, viz., $gb~{}(b\bar{b})\rightarrow bH~{}(gH)$ (Fig. 35). Owing to different perturbative expansions, these two schemes give different results when truncated at any finite order. Thus, higher order calculations become important to match these two results. The 4F scheme calculation is available up to NLO in QCD Dittmaier:2003ej ; Dawson:2005vi ; Dawson:2004wq ; Dawson:2003kb , while the 5F scheme is known up to NNLO accuracy in QCD Harlander:2003ai . Here the LO process in the 4F scheme i.e. $gg\rightarrow b\bar{b}H$ appears at the NNLO order in 5F scheme. The resonant scalar production in association with a $b$-quark or gluon i.e. $gb~{}(b\bar{b})\rightarrow bH~{}(gH)$ has been calculated up to NLO in QCD Dawson:2004sh and electroweak (EW) Dawson:2010yz .

It has been argued that with a proper choice of factorisation scale $\sim\frac{m_{H}}{4}$, the inclusive cross-section in the 4F and 5F schemes agree very well Maltoni:2003pn ; Boos:2003yi ; Plehn:2002vy . There is a proposed way to combine these two approaches. This is known as the Santander matching Harlander:2011 scheme. The total inclusive cross-section is obtained by matching the 4F and 5F scheme numbers in which both these cross-sections are multiplied by their proper weight factors. These weight factors change logarithmically with the heavy scalar mass ($m_{H}$) because of logarithmic difference between these two scheme approaches. The matched cross-section is computed as follows:

[TABLE]

where

[TABLE]

is the weight factor181818Here, $m_{b}$ is the bottom quark pole mass which enters in the re-summed logarithms.. The analysis can be subdivided according to the number of $b$-tagged jets. However, we specifically focus on the category with $\geq 1$ b-jets upon following a recent study performed by the ATLAS collaboration Aaboud:2017sjh . Furthermore, we consider the heavy Higgs decaying to a pair of $\tau$-leptons and we specifically focus on the scenario where both the $\tau$s decay hadronically. The $H$ and $A$ masses are varied between $200$ GeV and $1$ TeV.

The various backgrounds at play are $Z/\gamma^{*}+$ jets, multijets, $W+$ jets, $VV$ ($V=W^{\pm},Z$), $t\bar{t}$ and single top. The $Z/\gamma^{*}+$ jets with the $Z$-boson decaying to a pair of leptons ($e$, $\mu$ and $\tau$) is the dominant background for the $\tau_{h}\tau_{\ell}$ category (a category that we will not address in the present work) but also gives significant contribution to the $\tau_{h}\tau_{h}$ category. We simulate this background merged with three additional partons and some specific generation level cuts which are described in Appendix-A. Similarly, the $W+$ jets is also generated with up to three additional partons and the $W$-boson is then decayed leptonically. In order to include the dominant multijets background in the $\tau_{h}\tau_{h}$ category, we generate an exclusive $b\bar{b}jj$ sample where $j$ includes light quarks and gluon. These light jets can fake hadronically decaying $\tau$s. Finally, we include the top-quark related backgrounds viz., $t\bar{t}$ and single top. Next, we describe our analysis for the $\geq 1$ $b$-tagged jets category upon closely following Ref. Aaboud:2017sjh .

4.1 The $\tau_{h}\tau_{h}$ Channel : $b$-tag category

We select events containing at least one $b$-tagged jet with $p_{T}>20$ GeV and two $\tau$-tagged jets with $p_{T}>65$ GeV 191919Before performing this analysis, we validated our setup with the ATLAS analysis at 13 TeV. The validation is shown in Appendix B.. These two $\tau$-tagged jets must have opposite electric charge (from their track reconstruction). We also veto events having leptons ($e,\mu$) or $\tau$-tagged jets with $1.37<|\eta_{\tau}|<1.52$, in the final state. The azimuthal angle separation between the two $\tau$-tagged jets has to fulfil the condition, $|\Delta\phi(\tau,\tau)|>2.7$. The $b$- and the $\tau$-jets must have an angular separation in the $\eta-\phi$ plane, viz., $\Delta R(b,\tau)>0.2$. Besides, we also impose a minimum bound on the visible invariant mass of the two hadronically decaying $\tau$ leptons to be $m_{\tau\tau}^{vis.}>50$ GeV. For the fake $b\bar{b}jj$ background, we demand the two light jets to satisfy the $\tau$ jet configuration during the analysis and we later multiply the event yield with the $j\to\tau$ fake rate. Similarly, for the $W(\to\tau\nu\;\textrm{or}\;\ell\nu)+$ jets background, we demand at least one extra light jet satisfying the $\tau$ jet requirement on top of the $\tau$ jet ensuing from $W$-boson decay. After imposing the aforementioned cuts, we improve our analysis by optimising over some other kinematic variables viz., the transverse momentum of the hardest $\tau$-tagged jet ($p_{T,\tau_{1}}$), sum of the cosine of the azimuthal angle separation between the $\tau$-jets and $\not{E}_{T}$ ($\sum\limits_{\tau_{1,2}}\cos\Delta\phi$) and the transverse mass of the total system which is defined below,

[TABLE]

where the symbols have their usual meanings. These four kinematic variables are shown in Fig. 36. The optimised cuts along with the signal efficiencies and the background yields for each benchmark point are shown in Table 17. We perform our analysis upon considering both the 4F and 5F signal samples separately. Finally, we add them by multiplying these cross-sections with the aforementioned weight factor in order to obtain the upper limit on the matched $b\bar{b}H$ production cross section. We show the 95% and 99.7% exclusion for $\sigma(pp\to b\bar{b}H)\times\textrm{BR}(H\to\tau_{h}\tau_{h})$ in Fig. 37. The 95% upper limit varies between $22.16$ fb and $3.68$ fb (within $[146.91,~{}3.70]$ fb with $5\%$ systematic uncertainty) for $m_{H}$ varying between 300 GeV and 500 GeV. This is close to an order of magnitude improvement over the existing bounds at 13 TeV Aaboud:2017sjh . The effects of systematic uncertainties for this $b\bar{b}H$ channel become negligible for $m_{H}>400$ GeV. However, for lower masses, the backgrounds are larger and thus the inclusion of uncertainties weaken the limits.

5 The future of the pMSSM parameter space

The Higgs sector in the MSSM comprises two Higgs doublets which give rise to five massive Higgs states after the electroweak symmetry breaking. The Higgs spectrum is thus composed of two $CP$-even scalars, $h$ and $H$, one $CP$-odd scalar, $A$, and two charged scalars, $H^{\pm}$ (detailed studies on the Higgs sector of MSSM can be found in Ref. Djouadi:2015jea ; Zhao:2017qpe ). In addition to the extended Higgs sector, the SUSY particle spectrum boasts a multitude of particles, viz., the sleptons, squarks, gluinos and electroweakinos. A majority of direct searches at the LHC have excluded stops and gluinos below the TeV scale (Refs. Moortgat:2001pp ; Belanger:2015vwa ; Aaboud:2018kya ; Aaboud:2018mna ; Aaboud:2018ujj ; CMS-PAS-SUS-18-002 ; CMS-PAS-SUS-16-017 ; CMS-PAS-SUS-17-012 ; Sirunyan:2017bsh show some such limits in various supersymmetric interpretations). This more or less nullifies the prospect of observing these particles unless the luminosity is enhanced significantly. The electroweakino sector has also been probed in numerous studies Zhang:2002fu ; Ibrahim:2008rq ; ATLAS:2009zmv ; Heinemeyer:2015pfa ; Medina:2017bke ; Gori:2018pmk and bounds have been obtained on their masses within simplified scenarios Kulkarni:2017xtf ; Aaboud:2018htj ; Aaboud:2018sua ; Aaboud:2018zeb ; Aaboud:2018jiw ; Aaboud:2017leg ; Aaboud:2017mpt ; Aaboud:2017nhr ; Sirunyan:2018lul ; Sirunyan:2018ubx . In many of these studies, the electroweakino masses are excluded from between a few hundred GeVs to about half a TeV and are comparatively weakly coupled compared to the gluinos and stops. Within a generic SUSY parameter space without any correlation between the choice of the electroweakino mass parameters, these bounds can become considerably weaker. The ATLAS and CMS collaborations have also performed several studies to search for resonant Higgs through their decay into SM final states Aaboud:2017rel ; Aaboud:2017gsl ; Aaboud:2017sjh ; CMS:2017epy ; Aaboud:2017yyg ; CMS:2017ihs ; CMS-PAS-HIG-17-009 ; Sirunyan:2017djm ; Aaboud:2018cwk ; Aaboud:2018sfw . However, none of these searches could find any significant excess over the SM expectations and thus only imposed upper limits on the production cross-section of the heavy Higgs bosons times their branching ratio into various SM final states. In this section, we present a brief discussion on the relevant constraints and discuss the parameter scan for the 14 TeV HL-LHC. We follow this up with an analysis to capture the present status of the MSSM parameter space in light of the latest results from the Run-II data ($13~{}{\rm TeV},36~{}{\rm fb^{-1}}$) of LHC. Finally, we study the implications of the projected heavy Higgs direct search limits derived in sections 2, 3 and 4, on the MSSM parameter space.

The initial constraint on the parameter space ensues from the allowed mass of the 125 GeV SM-like Higgs boson. A combined measurement by the ATLAS and CMS collaborations constrains $m_{h}$ within the range [124.4,125.8] GeV at $3\sigma$. It is to be duly noted that in the context of MSSM, the available calculation of the Higgs mass is not exact. Thus, in order to correctly account for the existing uncertainties, we allow a window of $\pm$3 GeV about 125 GeV and restrict the light Higgs mass in our parameter space to lie within [122, 128] GeV. Furthermore, both collaborations have performed numerous measurements on the coupling strengths of the SM-like Higgs bosons. These results are presented through bounds on the signal strength variable ($\mu^{if}$) which is defined as follows:

[TABLE]

where, $\sigma_{i}$ represents the MSSM (or any specific model in question) Higgs production cross-section in the $i^{th}$ production mode ($i~{}=~{}ggF$, $VBF$, $t\bar{t}h$ or $Vh$) at the LHC and $\sigma_{i}^{SM}$ denotes the corresponding SM cross-section. BRf corresponds to the branching fraction of the SM-like Higgs into a particular SM final state ($f$ = $WW,ZZ,b\bar{b},\gamma\gamma,\tau\tau$) and BR${}_{f}^{SM}$ is the corresponding SM value. We apply all these constraints over our parameter space by demanding that all our signal strengths simultaneously lie within 2$\sigma$ of their experimental counterparts. The latest Higgs signal strengths (13 TeV, $36~{}{\rm fb^{-1}}$) measured by both the CMS and ATLAS collaborations are listed in Table 18.

Additionally, the flavour physics bounds also potentially constrain the MSSM parameter space, as shown in Barman:2016jov . In this regard, the bounds on the branching fraction of rare $B$-decay processes: $B\to X_{s}\gamma$, $B_{s}\to\mu^{+}\mu^{-}$, $B^{+}\to\tau^{+}\nu_{\tau}$, are among the most sensitive probes of new physics searches. As shown in Barman:2016jov , constraints from $B\to X_{s}\gamma$ disfavours the low $M_{A}$ regime while bounds on $Br(B_{s}\to\mu^{+}\mu^{-})$ constrains the low $M_{A}$ and high $\tan\beta$ regions. The low $M_{A}$ and high $\tan\beta$ region gets further constrained by $B^{+}\to\tau^{+}\nu_{\tau}$. On the other hand, the current limits from direct heavy Higgs searches in $b\bar{b}H/A,~{}H/A\to\tau^{+}\tau^{-}$, imposes much stringent constraint in the $M_{A}\gtrsim 300~{}{\rm GeV}$ region and excludes $\tan\beta$ up to $\sim 18$ for $M_{A}\sim 1~{}{\rm TeV}$ (a detailed discussion concerning this can be found in Barman:2016jov ). Consequently, within the scope of this section where our major emphasis is on exploring the future reach of direct heavy Higgs searches at the HL-LHC, we do not consider the implications from flavour physics bounds and impose only the light Higgs mass constraint and Higgs signal strength constraints in order to obtain the allowed parameter space region relevant for studying the current and future reach of direct heavy Higgs searches on the MSSM parameter space.

In order to evaluate the current allowed region in the parameter space of the phenomenological MSSM (pMSSM), we perform a random scan over a wide range of pMSSM input parameters, as described below. The parameters relevant to our study are the pseudo-scalar mass variable ($m_{A}$), ratio of the vacuum expectation values of the two Higgs doublets ($\tan\beta$), the third generation soft squark mass parameters ($M_{\tilde{Q}_{3}},~{}M_{\tilde{u}_{3}},~{}M_{\tilde{d}_{3}}$), the trilinear coupling of the stop ($A_{t}$) and sbottom ($A_{b}$) and the gluino mass parameter ($M_{3}$). These parameters are varied in the following range:

[TABLE]

The bino, wino and higgsino mass parameters, viz., $M_{1},~{}M_{2}$ and $\mu$ respectively, are varied from 600 GeV in order to prevent the heavier Higgs bosons from having any decays to the electroweakinos. This choice ensures only SM final states for the heavy Higgs boson decays. The second generation soft squark mass parameters ($M_{\tilde{Q}_{2},\tilde{u}_{2},\tilde{d}_{2}}$) are taken to be equal to their corresponding first generation counterparts ($M_{\tilde{Q}_{1},\tilde{u}_{1},\tilde{d}_{1}}$). The slepton mass parameters ($\tilde{e}_{1_{L},1_{R},2_{L},2_{R},3_{L},3_{R}}$) are fixed at $3~{}{\rm TeV}$ while the trilinear couplings of the first and the second generation squarks ($A_{u,d,c,s}$) and all three generations of sleptons ($A_{e,\mu,\tau}$) are taken to be zero.

The particle spectra and the branching fractions of the SM and SUSY particles are obtained using FeynHiggs Heinemeyer:1998yj . We consider only those parameter points which satisfy the light Higgs mass constraint defined above. Furthermore, we allow only those points which lie within 2$\sigma$ uncertainty of each of the signal strength variables listed in Table 18. The parameter space points which are allowed by the aforementioned light Higgs mass constraint and the Higgs signal strength constraints are referred to as the allowed parameter space points in the remainder of this section and are shown in grey in Fig. 39.

The ATLAS and CMS collaborations have also performed numerous searches for the heavy Higgs bosons through their decay into the SM final states, however, none of these searches have been able to observe any significant excess over the SM expectation. Consequently, upper limits have been set on the production cross-section of the heavy Higgs boson ($\sigma_{H/A}$) times its branching ratio into SM states. In this analysis, we consider the latest search limits on $\sigma_{ggH}\times Br(H\to ZZ,WW,\tau\tau)$ Aaboud:2017rel ; Aaboud:2017gsl ; Aaboud:2017sjh ; CMS:2017epy , $\sigma_{bbH/A}\times Br(H/A\to\tau\tau)$ Aaboud:2017sjh ; CMS:2017epy , $pp\to H\to\gamma\gamma$ Aaboud:2017yyg and $pp\to H\to hh\to 4b,2b2\gamma,2b2\tau$ CMS:2017ihs ; CMS-PAS-HIG-17-009 ; Sirunyan:2017djm derived by the CMS and ATLAS collaborations upon using the Run-II dataset with an integrated luminosity of $\sim 36~{}\rm{fb}^{-1}$. The gluon fusion channel is undoubtedly the dominant Higgs production mode at the LHC for low values of $\tan\beta$. However, it gets overrun by the $b\bar{b}H/A$ production channel at high $\tan\beta$ values. In the current analysis, while evaluating the impact of the existing upper limits on $pp\to H\to hh\to 4b,2b2\gamma,2b2\tau$, only the contributions from the gluon fusion production are taken into account. This choice is motivated by the fact that the $H\to hh$ decay modes gain dominance only in the low and intermediate $\tan\beta$ values where the gluon fusion mode overshadows the $b\bar{b}H/A$ channel. Although the current search limits on $H\to hh$ do not impose any constraints on our parameter space, the future runs have the potential to probe the low $m_{A}$ and low $\tan\beta$ regime. The impact of these future limits are discussed in the later part of this section. The $H\to ZZ/WW$ limits also turn out to be ineffective in constraining our parameter space and will require improvements of about three orders of magnitude for making any impact. We would like to mention that the upper limits derived by ATLAS in the $H\to\gamma\gamma$ search channel is on the fiducial cross-section times BR$(H\to\gamma\gamma)$. We compare these upper limits against a combination of the $ggF+b\bar{b}H/A$ production cross-sections and observe that an improvement of around two orders of magnitude will be required in order to affect our parameter space. Limits from searches in the $H/A\to\tau\tau$ channel impose the strongest constraints on the parameter space. Constraints from $\sigma_{b\bar{b}H/A}\times\textrm{BR}(H/A\to\tau\tau)$ yield stronger limits compared to their gluon fusion counterparts and exclude the low $m_{A}$ and high $\tan\beta$ region. The current search limits from ATLAS and CMS furnish roughly equivalent impact and rule out $\tan\beta\gtrsim 16$ for $m_{A}\sim 1~{}{\rm TeV}$. Before presenting the results in the $m_{A}-\tan{\beta}$ plane, we show the current allowed branching fractions, viz., $H\to hh,H\to t\bar{t}$ and $H\to\tau^{+}\tau^{-}$ in Fig. 38. The $H\to hh$ branching ratio dominates for $\tan{\beta}\lesssim 8$ and for $m_{H}\leq 2m_{t}$. All points are allowed by the Higgs mass and Higgs signal strength constraints. However, the grey regions are excluded by the present direct searches of the heavy Higgs. In Fig. 39, we show the impact of the latest direct search limits from $\sigma_{b\bar{b}H/A}\times Br(H/A\to\tau\tau)$ in the $m_{A}-\tan\beta$ plane. The parameter space points shown in Fig. 39 (grey and orange) are obtained after implementing the light Higgs mass constraints and the Higgs signal strength measurements. The grey points are excluded upon imposing the aforementioned direct search limits.

Our main concern in this section is to quantify the impact of the projected direct search limits for the HL-LHC which were derived in the previous sections. In this regard, we consider the projected direct search limits for the HL-LHC in the $H\to hh$ (Sec. 2), $H\to t\bar{t}$ (Sec. 3) and $b\bar{b}H/A\to b\bar{b}\tau_{h}\tau_{h}$ (Sec. 4) channels. Among the various final states of the $H\to hh$ channel, the $b\bar{b}\gamma\gamma$ final state furnishes the strongest limit in the $m_{A}\lesssim 600~{}{\rm GeV}$ regime, while the $4b$ final state imposes the strongest upper limits in the $m_{A}\gtrsim 600~{}{\rm GeV}$ region. The $H\to hh$ decay mode gains dominance in the low $\tan\beta$ region and especially before the $t\bar{t}$ mass threshold is attained. The same is reflected in the left panel of Fig. 40 where the brown points represent the region excluded at 95% CL by the projected $2\sigma$ reach from the $H\to hh\to b\bar{b}\gamma\gamma$ channel. The $4b$ final state, on the other hand, is rendered ineffective on account of reduced production cross-section at high values of $m_{A}$. The upper limits derived from searches in the remaining $H\to hh$ channels furnish much weaker bounds and will not be able to probe the pMSSM parameter even at the HL-LHC. The couplings of the heavy Higgs bosons with the up-type quarks have an inverse dependence on $\tan\beta$ and thus consequently the $H\to t\bar{t}$ channel has the potential to play an important role in the low $\tan\beta$ regime. The parameter space points excluded at 95% CL by the $H\to t\bar{t}$ HL-LHC search limits derived in Sec. 3 are shown in green in Fig. 40. The strongest future limits are obtained by the $b\bar{b}H\to b\bar{b}\tau\tau$ channel (derived in Sec. 4). This will be able to exclude (at 2$\sigma$) until $\sim\tan\beta\sim 5.5$ at $m_{A}\sim 1~{}{\rm TeV}$ as shown in Fig. 40, where the orange points are excluded by the same. The blue points in Fig. 40 denotes the parameter space which will evade the direct searches at the HL-LHC as well. The right panel in Fig. 40, however, shows the discovery potential at 5$\sigma$.

At this point, we would like to briefly discuss the implications from direct charged Higgs search limits on the parameter space of our interest. A detailed analysis of the exclusion reach of current limits from direct charged Higgs searches in the $pp\to H^{\pm}\to\tau^{\pm}\nu_{\tau}$ and $pp\to H^{+}\to t\bar{b}$ can be found in Barman:2016jov ; Djouadi2015 . Fig. 11 and Fig. 12 of Barman:2016jov shows that the allowed MSSM parameter space points (obtained by imposing the light Higgs mass constraint, Higgs signal strength limits and flavour physics constraints) are outside the current reach of charged Higgs searches in the $\tau\nu_{\tau}$ and $t\bar{b}$ channels, respectively, and, the direct charged Higgs search limits would require an improvement of roughly an order of magnitude in order to be capable of probing some of the MSSM parameter space points considered in Barman:2016jov . The implications for future direct charged Higgs search limits for HL-LHC, obtained by scaling the current limits, has been analysed in Djouadi2015 , where the projected reach of direct searches in the $H^{\pm}\to\tau^{\pm}\nu_{\tau}$ and $H^{\pm}\to t\bar{b}$ channels has been translated to the MSSM parameter space and presented in the $M_{A}-\tan\beta$ plane (see Fig. 20 of Djouadi2015 ). A comparison with the analysis in Djouadi2015 indicates that the future reach of direct charged Higgs searches at the HL-LHC is weaker than the future reach of direct heavy Higgs searches in the $b\bar{b}H/A,~{}H/A\to\tau^{+}\tau^{-}$ channel derived in this work.

The results discussed till now assume that the heavy Higgs bosons underwent decays only into SM final states. The branching fractions of the heavy Higgs bosons into SM final states can, in principle, undergo significant modifications in the presence of light SUSY particles202020SUSY particles with their masses less than $M_{H/A}/2$, such that it is kinematically possible for the heavy Higgs bosons to decay into them, are referred to as the light SUSY particles. For example, for intermediate values of $\tan\beta\sim 7-10$, the branching fraction of the heavy Higgs bosons into charginos and neutralinos may attain significantly large values ($\gtrsim 50\%$) Ananthanarayan:2015fwa ; Barman:2016kgt . In the remainder of this section we study the impact on the pMSSM parameter space in the presence of non-SM decay modes of the heavy Higgs bosons at the HL-LHC. Here we will restrict ourselves to the case of light electroweakinos. These electroweakinos are required to be an admixture of gauginos (bino and wino) and higgsinos in order to have couplings with the Higgs bosons. The pure gauginos and higgsinos do not couple with the Higgs states. LHC searches in the chargino-neutralino pair production mode furnishes the most stringent constraints on the electroweakino sector and excludes degenerate wino-like $\chi_{2}^{0}$ and $\chi_{1}^{\pm}$ of mass $\lesssim 450~{}{\rm GeV}$ Sirunyan:2018ubx , for an LSP neutralino of mass $\sim 100~{}{\rm GeV}$. However, such constraints do not apply to scenarios where the LSP and NLSP are almost degenerate in mass. We explore this fact and vary $M_{2}$ and $\mu$ in such a way that $|M_{2}-\mu|<10~{}{\rm GeV}$, $M_{H/A}>(\chi_{1}^{0}+\chi_{2}^{0})$ and $M_{2},\mu>200~{}{\rm GeV}$ Sirunyan:2018ubx . Closeness between $M_{2}$ and $\mu$ ensures that the $\chi_{1}^{0},~{}\chi_{2}^{0},~{}\chi_{3}^{0}$ and $\chi_{1}^{\pm}$ have significant admixtures from both winos and higgsinos. The remaining input parameters are randomly varied within the range specified in Eqn. 3 except for $M_{1}$ which we fix at $1~{}{\rm TeV}$. In presence of these $H/A\to ino$ decay modes, the branching fraction of the heavy Higgs bosons to SM final states undergoes modifications and manifests in weaker limits on the parameter space, as shown in Fig. 41. Correspondingly, the orange and blue regions shift upward. The brown and green regions shrink further down. For $m_{A}$ varying between 400 GeV and 700 GeV, $\tan{\beta}$ as low as 3 is excluded at 95% CL. In presence of these non-SM decay modes, the current ($13~{}{\rm TeV}$, $36~{}\rm{fb}^{-1}$) limits on $\sigma_{b\bar{b}H/A}\times Br(H/A\to\tau\tau)$ exclude $\tan\beta\gtrsim 22$ for $m_{A}\sim 1~{}{\rm TeV}$. The HL-LHC reach weakens out till $\tan\beta\sim 10$ at $m_{A}\sim 1~{}{\rm TeV}$. In the current scenario, the projected limits from $H\to t\bar{t}$ lose sensitivity on the pMSSM parameter space under study. The HL-LHC projections from $H\to hh\to b\bar{b}\gamma\gamma$ also imposes weaker constraints and excludes $\tan\beta<8$ at $m_{A}\sim 400~{}{\rm GeV}$.

6 Summary

In this work, we have studied the prospects of observing or excluding a resonant heavy Higgs or pseudoscalar in the purview of the HL-LHC at 14 TeV. Various experimental observations and theoretical motivations necessitate physics beyond the Standard Model (BSM). Several searches are performed, either in the context of specific models or in a generic model-independent fashion, to gauge the type of new physics. In this work, we specifically focus on neutral heavy Higgs bosons (both $CP$ odd and even). Run-II data at the LHC has already constrained strongly interacting BSM particles like gluinos and stops to $\mathcal{O}(\geq 1)$ TeV. However, the LHC still hasn’t imposed such strong constraints on extended Higgs sectors. The Standard Model Higgs self-coupling being still unknown, we are yet to fully understand the scalar sector of new physics. Here, we studied three major search channels for such a heavy Higgs (or pseudoscalar). Specific to the $CP$-even heavy Higgs, we studied the prospects of constraining $\sigma(pp\to H\to hh)$ in multifarious channels, viz., $b\bar{b}\gamma\gamma$, $b\bar{b}b\bar{b}$, $b\bar{b}\tau^{+}\tau^{-}$, $b\bar{b}WW^{*}$ and $\gamma\gamma WW^{*}$. We took guidance from the present searches and optimised each channel carefully to obtain upper limits on $\sigma(pp\to H\to hh)$ from each channel. Corroborating the present searches, we find that the $b\bar{b}\gamma\gamma$ and $b\bar{b}b\bar{b}$ final states serve as the golden channels for $m_{H}$ in the range $[300,600]$ GeV and $[600,1000]$ GeV, respectively. The $b\bar{b}\gamma\gamma$ sets a 95% CL upper limit on $\sigma(pp\to H\to hh)$ between $[79.03,14.10]$ fb in the aforementioned mass range. The $4b$ channel on the other hand sets a corresponding cross-section limit between $[5.36,2.51]$ fb for $m_{H}\in[800,1000]$ GeV. The limits from the remaining three channels are not so promising. On the other hand, if the mass of the scalar or the pseudoscalar Higgs is above the $t\bar{t}$ threshold and the Higgs is dominantly produced via gluon fusion, then it can also have a dominant decay into $t\bar{t}$. Upon studying the fully leptonic as well as the semi-leptonic final states, we find the strongest limits on $\sigma(pp\to H\to t\bar{t})$ from the semi-leptonic category. The 95% CL upper limits lie between $186.57$ fb and $32.81$ fb for $m_{H}\in[400,1000]$ GeV. Finally, we studied the $b\bar{b}H/A$ production in the $b\bar{b}\tau\tau$ final state upon demanding at least one $b$-tagged jet in the final state and demanding two hadronic $\tau$s. The 95% CL upper limit on $\sigma(pp\to b\bar{b}H)$ varies between $22.16$ fb and $3.68$ fb for $m_{H}$ lying between 300 GeV and 500 GeV. All our searches for the scalars are mostly model independent and can be translated to models with multiple Higgs bosons with narrow widths.

In this work, we considered the specific example of supersymmetry, more specifically the pMSSM. We apply present constraints from the SM-like Higgs boson mass measurement and all its signal strengths into multiple final states. The future limits obtained in this work constrain different regimes of the parameter space. The $H\to hh$ search, mostly in the $b\bar{b}\gamma\gamma$ channel excludes $\tan{\beta}$ to as low as $4$, at $95\%$ CL, for $m_{A}\sim 2m_{t}$ GeV. The $H\to t\bar{t}$ has a similar exclusion on $\tan{\beta}$ for $m_{A}$ varying between $[400,800]$ GeV. The $b\bar{b}H$ channel in the di-$\tau+\geq 1b$-tagged jet final state excludes $\tan{\beta}$ as low as $5.5$ for $m_{A}=1$ TeV. The blue region in Fig. 40 will not be probed even by direct searches if the heavy Higgs bosons decay only to SM particles. We will require higher energy colliders in order to be able to probe this region. This scenario might change if there are light electroweakinos, below $m_{H/A}/2$. In that situation, the $m_{A}-\tan{\beta}$ parameter region changes. Upon considering a scenario where $|M_{2}-\mu|<10$ GeV, $M_{H/A}>(\chi_{1}^{0}+\chi_{2}^{0})$ and $M_{2}/\mu>200$ GeV, one finds that the $H\to hh$ and $H\to t\bar{t}$ excludes $\tan{\beta}$ down to 3 for $m_{A}\in[400,700]$ GeV. The exclusion bound on $\tan{\beta}$ from the $b\bar{b}H$ search decreases to $10$ for $m_{A}=1$ TeV at $95\%$ CL.

Appendix A Detailing the cross section with generation cuts for the signal and backgrounds

Appendix B Validation of the $b\bar{b}H$ analysis

Before performing our analysis for the $b\bar{b}H$ category, we validate our setup with the existing analysis in this channel Aaboud:2017sjh ; Pickering:2258362 ; Wahrmund:2279938 . We generate the signal events (at LO in SM) with MG5_aMC@NLO and shower them via Pythia-8 Sjostrand:2014zea . We use different parton distribution functions (PDFs) for the sample generations. Specifically, we use the CT10nlo_nf4 Dulat:2015mca for the 4F $b\bar{b}H$ process, MSTW2008nnlo68cl Martin:2009iq for the 5F $b\bar{b}H$ process and CT10 Lai:2010vv for the ggF process. Next, we impose the following cuts in sequence. For the $b$-tag category, we demand at least one $b$-tagged jet in the final state. The events are required to have at least two $\tau$ jets with opposite charge (from their reconstructed charged tracks). The leading and sub-leading ($p_{T}$ ordered) $\tau$-tagged jets are required to have $p_{T}>65$ GeV. The $\tau$ jets lying inside the transition region viz., $1.37<|\eta|<1.52$, are removed. Furthermore, the azimuthal angle separation between the leading and the sub-leading $\tau$ jets is required to be $\Delta\phi(\tau,\tau)>2.7$. Furthermore, we require $m_{\tau\tau}^{vis.}>50$ GeV, $\not{E}_{T}>20$ GeV and $\sum\limits_{\tau}\cos(\Delta\phi)>0$. We show the validation in Table 21.

Acknowledgements.

We thank Mikael Chala, Arghya Choudhury and Marius Wiesemann for helpful discussions at various stages of the work. The work of BB is supported by the Department of Science and Technology, Government of India, under the Grant Agreement number IFA13-PH-75 (INSPIRE Faculty Award). SB is supported by a Durham Junior Research Fellowship COFUNDed by Durham University and the European Union, under grant agreement number 609412.