Connecting the ISM to TeV PWNe and PWNe candidates

TL;DR

This study explores the interstellar medium around seven TeV gamma-ray sources linked to pulsar wind nebulae, revealing dense molecular clouds and potential cosmic-ray interactions, with implications for understanding their origins and contributions to gamma-ray emissions.

Contribution

It provides new molecular line observations and distance estimates for several PWNe and SNRs, and discusses their potential roles in cosmic-ray acceleration and gamma-ray production.

Findings

Detection of dense molecular clouds near TeV sources.

Identification of possible SNR-PWN interactions.

Assessment of PWNe contributions to gamma-ray emission.

Abstract

We investigate the interstellar medium (ISM) towards seven TeV gamma-ray sources thought to be pulsar wind nebulae (PWNe) using Mopra molecular line observations at 7mm [CS(1-0), SiO(1-0,v=0)], Nanten CO(1-0) data and the SGPS/GASS HI survey. We have discovered several dense molecular clouds co-located to these TeV gamma-ray sources , which allows us to search for cosmic-rays (CRs) coming from progenitor SNRs or, potentially, from PWNe. We notably found SiO(1-0,v=0) emission towards HESS J1809-193, highlighting possible interaction between the adjacent supernova remnant SNR G011.0-0.0 and the molecular cloud at d $\sim$ 3.7 kpc. Using morphological features, and comparative studies of our column densities with those obtained from X-ray measurements, we claim a distance d $\sim$ 8.6 - 9.7 kpc for SNR G292.2-00.5, d $\sim$ 3.5 - 5.6 kpc for PSR J1418-6058 and d $\sim$ 1.5 kpc for the new SNR candidate found towards HESS J1303-631. From our mass and density estimates of selected molecular clouds, we discuss signatures of hadronic/leptonic components from PWNe and their progenitor SNRs. Interestingly, the molecular gas, which overlaps HESS J1026-582 at d $\sim$ 5 kpc, may support a hadronic origin. We find however that this scenario requires an undetected cosmic-ray accelerator to be located at d $\lt$ 10 pc from the molecular cloud. For HESS J1809-193, the cosmic-rays which have escaped SNR G011.0-0.0 could contribute to the TeV gamma-ray emission. Finally, from the hypothesis that at most 20% the pulsar spin down power could be converted into CRs, we find that, among the studied PWNe, only those from PSR J1809-1917 could potentially contribute to the TeV emission.

Tables

Table 1. Table 1: Mopra 7mm coverage of the TeV sources studied in this paper. The central position in (RA, Dec) and size is given for all observations undertaken towards the TeV sources.

TeV sources	central position	Area
	(${\alpha }^{\circ }$, ${\delta }^{\circ }$)	$\Delta l\times \Delta b$
	(J2000.0)
HESS J1018–589B	(154.7∘,-58.9∘)	${20}^{\prime }\times {20}^{\prime }$
HESS J1026–582	(156.8∘,-57.8∘)	${20}^{\prime }\times {20}^{\prime }$
	(156.4∘,-58.1∘)	${20}^{\prime }\times {20}^{\prime }$
	(156.1∘,-58.4∘)	${20}^{\prime }\times {20}^{\prime }$
	(155.9∘,-58.0∘)	${10}^{\prime }\times {10}^{\prime }$
HESS J1119–614	(169.7∘,-61.4∘)	${20}^{\prime }\times {20}^{\prime }$
HESS J1303–631	(195.4∘,-63.1∘)	${20}^{\prime }\times {20}^{\prime }$
HESS J1418–609	(215.1∘,-60.8∘)	${20}^{\prime }\times {20}^{\prime }$
HESS J1420–607	(214.4∘,-60.9∘)	${20}^{\prime }\times {20}^{\prime }$
HESS J1809–193	(272.4∘,-19.2∘)	${20}^{\prime }\times {20}^{\prime }$
	(272.2∘,-19.5∘)	${20}^{\prime }\times {20}^{\prime }$
	(272.5∘,-19.7∘)	${20}^{\prime }\times {20}^{\prime }$
	(272.7∘,-19.4∘)	${20}^{\prime }\times {20}^{\prime }$

Table 2. Table 2: Cosmic-ray enhancement factors k cr subscript 𝑘 cr k_{\text{cr}} derived using eq. 10 from Aharonian ( 1991 ) , required to reproduce the observed TeV emission above 1 TeV F ( > 1 ​  TeV ) = ∫ N 0 ​ E γ − Γ ​ d ​ E γ annotated 𝐹 absent 1  TeV subscript 𝑁 0 superscript subscript 𝐸 𝛾 Γ d subscript 𝐸 𝛾 F\left(>1\text{ TeV}\right)=\int{N_{0}E_{\gamma}^{-\Gamma}\text{d}E_{\gamma}} via CR-ISM interaction.

	Molecular
	Regions	$k_{\text{cr}}$	$F\left(>1\text{TeV}\right)$
			[ph cm-2 s-1]
HESS J1809$-$193	1a	32	${}_{}{}^{\ast }6.2\times {10}^{-13}$
	2a	12	${}_{}{}^{\ast }6.2\times {10}^{-13}$
	3a	22	${}_{}{}^{\ast }4.7\times {10}^{-13}$
HESS J1026$-$582	B	60	${}_{}{}^{†}1.1\times {10}^{-12}$
HESS J1119$-$193	D	143	${}_{}{}^{‡}9.2\times {10}^{-13}$
HESS J1420$-$607	6	53	${}_{}{}^{\coprod }1.3\times {10}^{-13}$
	7	88	${}_{}{}^{\coprod }6.7\times {10}^{-14}$
	8	81	${}_{}{}^{\coprod }4.2\times {10}^{-14}$
	9	81	${}_{}{}^{\coprod }3.2\times {10}^{-14}$
	10	97	${}_{}{}^{\coprod }1.1\times {10}^{-13}$
HESS J1303$-$631	Ab	257	${}_{}{}^{\pm }1.7\times {10}^{-12}$
	Cc	194	${}_{}{}^{\pm }2.1\times {10}^{-13}$
	Dc	237	${}_{}{}^{\pm }2.1\times {10}^{-13}$
HESS J1018$-$589A	A	91	${}_{}{}^{\star }2.4\times {10}^{-13}$
∗Aharonian et al. (2007), †Abramowski et al. (2011), ‡Abdalla et al. (2018a), ±, ∐Aharonian et al. (2006), Abramowski et al. (2012b),     ⋆Abramowski et al. (2015)
${}^{\text{a}}$:We scaled down the HESS J1809$-$193 photon flux by 16%, 16% and 12% for regions ‘1’, ‘2’, ‘3’ respectively, corresponding to the ratio between the molecular regions and the TeV emission sizes.
${}^{\text{b}}$:We scaled down the HESS J1420$-$607 photon flux by 4%,$\sim 1$%, %, % and 3% for the regions 6 to 10 respectively.
${}^{\text{c}}$: We scaled down the HESS J1303$-$631 photon flux by 44%, 7% and 7% for region ‘A’,‘B’ and ‘C’ respectively.

Table 3. Table 3: Table showing the diffusion time-scale t diff = ( ⟨ r ⟩ + d source ) 2 / 6 ​ D ​ ( E , B ) subscript 𝑡 diff superscript delimited-⟨⟩ 𝑟 subscript 𝑑 source 2 6 𝐷 𝐸 𝐵 t_{\text{diff}}=\left(\langle r\rangle+d_{\text{source}}\right)^{2}/6D\left(E,B\right) for a particle with energy E = 10 𝐸 10 E=10  TeV to cross the listed molecular regions with mean radius ⟨ r ⟩ delimited-⟨⟩ 𝑟 \langle r\rangle and positioned at d source subscript 𝑑 source d_{\text{source}} from the centre of the listed sources. The age of the source t age subscript 𝑡 age t_{\text{age}} , the CR-ISM interaction time-scale τ pp subscript 𝜏 pp \tau_{\text{pp}} for protons, and the synchrotron time scale τ sync subscript 𝜏 sync \tau_{\text{sync}} for electrons with energy E = 10 𝐸 10 E=10  TeV are also displayed as means of comparison.

		Molecular region properties					Source particle properties
	Reg.	$⟨r⟩$	$B^{\text{a}}$	${\tau }_{\text{pp}}^{\text{b}}$	${\tau }_{\text{sync}}^{\text{c}}$		Source	$t_{\text{age}}$	$d_{\text{source}}$	$t_{\text{diff}}{\left(10\text{ TeV}\right)}^{\text{d}}$
		[pc]	[$\mu$G]	[kyr]	[kyr]			[kyr]	[pc]	[kyr/${\chi }_{0.1}$]
HESS J1809$-$193	1	9	29	69	1.4	$\{$	SNR G011.0$-$0.0	51	8	2
							PSR J1809$-$5158	51	13	3
	2	10	41	35	0.6	$\{$	SNR G011.0$-$0.0	51	20	6
							PSR J1809$-$5158	51	28	9
HESS J1026$-$582	A	11	10	1200	12		PSR J1028$-$5819	89	5	1
HESS J1303$-$631	A	12	18	26	3		PSR J1301$-$6305	11	45	14
	B	3	32	56	1		SNRce	11	8	1
	C	11	19	167	3		PSR J1301$-$6305	11	48	15
	D	24	10	555	11		PSR J1301$-$6305	11	26	4
HESS J1418$-$609	A	6	33	41	1		PSR J1418$-$6058	1.6	21	4
	1	4	26	85	2		PSR J1418$-$6058	1.6	18	2
	2	4	12	400	8		PSR J1418$-$6058	1.6	19	2
	3	3	22	115	2		PSR J1418$-$6058	1.6	26	4
	4	5	24	107	2		PSR J1418$-$6058	1.6	30	6
	5	4	32	60	1		PSR J1418$-$6058	1.6	27	5
	6	3	25	90	2		PSR J1418$-$6058	1.6	22	3
	7	2	24	107	2		PSR J1418$-$6058	1.6	15	1
	8	2	35	49	1		PSR J1418$-$6058	1.6	12	1
	9	1	22	125	3		PSR J1418$-$6058	1.6	16	1
	10	3	30	41	1		PSR J1418$-$6058	1.6	32	4
HESS J1018$-$589	A	4	26	85	2		SNR G284.3$-$1.8	10	8	1
${}^{\text{a}}$:Obtained from eq. 21 from Crutcher et al. (2010) (see Section 4.2).
${}^{\text{b}}$ ${\tau }_{\text{pp}}=6\times {10}^7/2n_{{\text{H}}_2}$ where $n_{{\text{H}}_2}$ can be found in Tables D.1 to D.6
${}^{\text{c}}$: ${\tau }_{\text{sync}}=2.45\times {10}^7{\left(B_{\text{ mG}}^{2}E/m_{\text{e}}\right)}^{-1}$ yr (Ginzburg & Syrovatskii 1964).
${}^{\text{d}}$: See eqs. 2 and 3 from Gabici et al. (2007) to obtain the diffusion time-scale. The results shown assume a diffusion coefficient suppression factor $\chi =0.1$
${}^{\text{e}}$ SNR candidate, see text and Sushch et al. (2017) for further details and Figure 8 for location. We used standard SNR expansion parameters to provide an age estimate of the SNRc if located at this distance (see text).

Table 4. Table C.4: Position, size and fitting parameters of the emission traced by CS(1–0), its isotopologue C 34 S(1–0), CO(1–0) and HC 3 N(5–4,F=4–3) towards HESS J1809–193 . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using a main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 and text). Finally, the { { \big{\{} indicates that the two components with little velocity separation are physically connected.

Molecular tracer	Region	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	(arcsec2)	(K)	(km/s)	(km/s)	(K km/s)
		(J2000.0)
Nanten CO(1–0)	1	(272.61∘, -19.38∘)	819$\times$373	6.5	30.7	8.7	67.6
	2	(272.36∘, -19.71∘)	663$\times$437	6.9	31.0	13.6	111.8
	3	(272.37∘, -19.34∘)	600$\times$360	7.1	32.0	13.5	114.0
Mopra CS(1–0)	1	(272.43∘, -19.74∘)	819$\times$373	0.2	30.2	4.6	2.1
	1-1	(272.61∘, -19.38∘)	146$\times$45	0.6	29.4	3.1	4.6
	1-2	(272.62∘, -19.46∘)	30$\times$30	0.3	30.0	2.3	1.9
	1-3	(272.57∘, -19.46∘)	30$\times$30	0.3	29.5	2.6	2.2
	1-4	(272.54∘, -19.47∘)	30$\times$30	0.5	29.1	3.6	4.3
	2	(272.36∘, -19.71∘)	663$\times$437	0.2	29.2	3.9	2.0
	2-1	(272.57∘, -19.46∘)	30$\times$30	$\{$ 0.5	31.5	2.7	3.3
				0.3	31.7	1.1	0.8
	2-2	(272.53∘, -19.72∘)	30$\times$30	$\{$ 0.5	29.4	3.9	4.4
				0.1	31.7	1.2	0.4
	2-3	(272.41∘, -19.77∘)	30$\times$30	$\{$ 0.5	30.8	2.6	3.0
				0.3	28.6	2.2	1.3
	3	(272.37∘, -19.34∘)	600$\times$360	$\{$ 0.1	34.2	3.6	1.3
				0.1	31.3	1.5	0.4
Mopra C34S(1–0)a	1-1	(272.61∘, -19.38∘)	146$\times$45	0.1	30.0	3.7	0.9
	1-2	(272.62∘, -19.46∘)	30$\times$30	0.2	30.0	1.2	0.5
	1-3	(272.57∘, -19.46∘)	30$\times$30	0.2	29.4	0.8	0.4
	1-4	(272.54∘, -19.47∘)	30$\times$30	0.2	29.3	0.8	0.4
	2-1	(272.57∘, -19.46∘)	30$\times$30	$\{$ 0.2	31.6	1.0	0.4
				0.1	30.1	1.0	0.3
	2-2	(272.53∘, -19.72∘)	30$\times$30	0.1	30.5	1.6	0.6
	2-3	(272.41∘, -19.77∘)	30$\times$30	0.2	30.8	1.0	0.4
Mopra HC3N(5–4,F=4–3)	HC1	(272.69∘, -19.27∘)	187$\times$60	0.2	31.7	1.8	1.0
	HC2	(272.61∘, -19.38∘)	205$\times$67	0.2	29.3	2.1	1.1
	HC3	(272.53∘, -19.46∘)	88$\times$88	0.2	29.2	2.1	1.2
	HC4	(272.49∘, -19.74∘)	257$\times$214	0.1	29.1	5.1	1.3
Mopra SiO(1–0,v=0)	S1	(272.62∘, -19.37∘)	30$\times$30	0.06	30.2	4.0	0.6
	S2	(272.54∘, -19.46∘)	40$\times$40	0.05	26.9	3.6	0.4
	S3	(272.44∘, -19.70∘)	66$\times$66	0.06	29.4	4.6	0.7
a:The values represent the semi-major and semi-minor axes respectively.

Table 5. Table C.5: Position, size and fitting parameters of the emission traced by CS(1–0) and CO(1–0) towards HESS J1026–582 . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using a main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 for 7mm tracers)). Finally, the { { \big{\{} indicates that the two components with little velocity separation are physically connected.

Molecular tracer	Region	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	(arcsec2)	(K)	(km/s)	(km/s)	K km/s
		(J2000.0)
Nanten CO(1–0)	A	(157.47∘, -58.58∘)	1096$\times$1096	1.1	-15.2	4.3	5.7
	B	(156.35∘, -58.20∘)	773$\times$773	2.3	3.2	4.3	11.9
Mopra CS(1–0)	1	(156.36∘, -58.15∘)	87$\times$87	0.3	-19.2	1.4	1.0
	2	(156.15∘, -58.10∘)	87$\times$87	0.2	-16.57	1.2	0.7
	3	(156.85∘, -57.77∘)	85$\times$260	0.2	11.4	3.1	1.3
	4	(156.66∘, -57.79∘)	173$\times$272	0.3	11.3	3.6	2.3
	5	(156.57∘, -57.89∘)	153$\times$61	0.2	1.2	2.2	1.3
	6	(156.46∘, -58.20∘)	57$\times$234	$\{$ 0.2	3.7	1.8	0.8
				0.2	2.0	1.0	0.5
	7	(156.33∘, -58.24∘)	98$\times$98	0.5	3.6	1.6	1.9
a:The values represent the semi-major and semi-minor axes respectively.

Table 6. Table C.6: Position, size and fitting parameters of the emission traced by CO(1–0) towards HESS J1119–582 . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using the main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 for 7mm tracers). Finally, the { { \big{\{} indicates that two the two components with little velocity separation are physically connected.

Molecular tracer	Region	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	(arcsec2)	(K)	(km/s)	(km/s)	(K km/s)
		(J2000.0)
Nanten CO(1–0)	A	(170.20∘, -61.35∘)	522$\times$522	1.4	-12.3	5.9	9.6
	B	(169.40∘, -61.48∘)	495$\times$159	0.6	-19.1	5.1	3.8
	C	170.08∘, -61.19∘)	219$\times$219	$\{$ 1.1	20.4	8.2	4.1
				1.0	35.1	2.2	1.5
	D	(169.38∘, -61.44∘)	705$\times$480	1.7	30.2	4.7	4.4
a:The values represent the semi-major and semi-minor axes respectively.

Table 7. Table C.7: Position, size and fitting parameters of the emission traced by CS(1–0) and CO(1–0) towards HESS J1418–609 and HESS J1420–607 . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using the main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 for tracers). Finally, the { { \big{\{} indicates that the two components with little velocity separation are physically connected.

Molecular tracer tracer	Regions	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	(arcsec2)	(K)	(km/s)	(km/s)	(K km/s)
		(J2000.0)
Nanten CO(1–0)	A	(214.0∘, -61.1∘)	364$\times$364	5.9	-45.0	10.2	71.3
	B	(214.47∘, -61.11∘)	352$\times$352	1.9	-3.8	6.9	15.3
	6	(215.33∘, -60.84∘)	110$\times$288	$\{$ 4.8	-50.1	8.8	50.4
				3.6	-43.3	6.3	27.3
	7	(215.08∘, -60.85∘)	94$\times$171	5.3	-49.1	7.2	45.9
	8	(215.11∘, -60.93∘)	58$\times$173	4.7	-48.2	9.0	50.5
	9	(215.03∘, -60.79∘)	87$\times$87	4.8	-49.0	8.7	49.7
	10	(215.05∘, -60.69∘)	179$\times$143	5.0	-47.2	7.0	41.9
Mopra CS(1–0)	1	(214.52∘, -60.80∘)	71$\times$217	0.3	-59.1	1.4	1.1
	2	(214.37∘, -60.85∘)	85$\times$207	0.1	-56.7	0.9	0.3
	3	(214.21∘, -60.84∘)	61$\times$140	0.3	-56.7	1.1	0.7
	4	(214.08∘, -60.88∘)	90$\times$250	0.2	-57.5	2.6	1.1
	5	(214.11∘, -60.99∘)	82$\times$201	$\{$ 0.2	-62.2	2.2	1.2
				0.2	-60.3	1.4	0.6
	6	(215.33∘, -60.84∘)	110$\times$288	0.3	-44.8	1.3	1.0
	7	(215.08∘, -60.85∘)	94$\times$171	0.2	-48.8	1.6	0.8
	8	(215.11∘, -60.93∘)	58$\times$173	0.2	-48.6	1.6	0.9
	9	(215.03∘, -60.79∘)	87$\times$87	0.2	-48.6	1.2	0.6
	10	(215.05∘, -60.69∘)	179$\times$143	0.3	-46.7	2.6	1.8
	11	(214.75∘, -61.00∘)	87$\times$87	$\{$ 0.2	-51.3	1.7	0.6
				0.1	-42.4	0.7	0.3
	12	(214.46∘, -61.01∘)	87$\times$87	0.3	-47.2	2.8	1.8
	13	(214.31∘, -61.09∘)	87$\times$87	0.2	-50.5	1.0	0.6
	14	(214.88∘, -60.85∘)	109$\times$178	$\{$ 0.3	-3.4	4.3	3.1
				0.1	-4.9	0.7	0.3
	15	(214.99∘, -60.80∘)	87$\times$87	0.2	0.9	1.6	0.7
a:The values represent the semi-major and semi-minor axes respectively.

Table 8. Table C.8: Position, size and fitting parameters of the emission traced by CS(1–0) and CO(1–0) transitions towards HESS J1303–631 . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using the main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 for 7mm tracers). Finally, the { { \big{\{} indicates that two the two components with little velocity separation are physically connected.

Molecular tracer	Region	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	arcsec2	(K)	(km/s)	(km/s)	K km/s
		(J2000.0)
Nanten CO(1–0)	A	(195.63∘, -63.46∘)	512$\times$300	2.8	-30.0	10.3	34.7
	B	(196.06∘, -63.21∘)	438$\times$438	3.8	-19.7	7.9	35.6
	C	(195.89∘, -63.17∘)	130$\times$131	$\{$ 5.6	33.4	5.4	36.6
				1.0	24.7	3.1	3.8
	D	(195.22∘, -63.15∘)	109$\times$208	$\{$ 2.2	32.7	7.4	19.4
				1.6	26.3	4.2	8.2
Mopra CS(1–0)	1	(195.52∘, -63.23∘)	87$\times$87	0.3	-30.8	2.8	1.8
	2	(195.62∘, -63.18∘)	87$\times$87	0.1	-30.1	3.1	1.0
	3	(195.71∘, -62.99∘)	87$\times$87	0.1	-29.7	3.0	0.7
	4	(195.78∘, -63.22∘)	87$\times$87	0.3	-23.0	2.1	1.5
	5	(195.60∘, -63.00∘)	87$\times$87	0.2	-18.1	1.8	0.9
	6	(195.46∘, -63.08∘)	87$\times$87	0.3	-30.6	2.1	1.3
	7	(195.49∘, -62.96∘)	87$\times$87	0.2	-23.7	2.0	1.1
	8	(195.46∘, -63.07∘)	92$\times$92	0.1	30.7	2.7	0.9
	9	(195.48∘, -62.95∘)	47$\times$47	0.2	24.1	3.1	1.5
	10	(195.20∘, -63.14∘)	36$\times$36	0.2	25.9	1.6	0.9
a:The values represent the semi-major and semi-minor axes respectively.

Table 9. Table C.9: Position, size and fitting parameters of the emission traced by CO(1–0) towards HESS J1018–589B . T ∗ A superscript subscript absent A {}_{\text{A}}^{*} denotes the peak temperature at v lsr = v cent subscript 𝑣 lsr subscript 𝑣 cent v_{\text{lsr}}=v_{\text{cent}} , Δ ​ v Δ 𝑣 \Delta v indicates the full width half maximum (FWHM) of a Gaussian fit. W 𝑊 W represents the main beam integrated intensity using a main beam efficiency factor η mb subscript 𝜂 mb \eta_{\text{mb}} (see Urquhart et al. 2010 for 7mm tracers).

Molecular tracer	Region	Position	Areaa	$T_{\text{A}}^{\ast }$	$v_{\text{cent}}$	$\Delta v$	$W=\int T_{\text{mb}}\text{d}v$
		(RA, Dec)	arcsec2	(K)	(km/s)	(km/s)	K km/s
		(J2000.0)
Nanten CO(1–0)	A	(154.68∘, $-{58.82}^{\circ }$)	400$\times$400	0.9	$-15.9$	12.0	13.8
a:The values represent the semi-major and semi-minor axes respectively.

Table 10. Table D.1: Physical parameters obtained from our CS and CO analyses for the different selected regions located towards HESS J1809–193 . In the case where C 34 S(1–0) is detected, the derived optical depth τ CS(1–0) subscript 𝜏 CS(1–0) \tau_{\text{CS(1--0)}} is shown as superscript next to the CS column density N CS subscript 𝑁 CS N_{\text{CS}} . Otherwise, an optical thin scenario is assumed and the derived column densities N H 2 subscript 𝑁 subscript H 2 N_{\text{H}_{2}} and N CS subscript 𝑁 CS N_{\text{CS}} , mass M H 2 subscript 𝑀 subscript H 2 M_{\text{H}_{2}} (CS) and H 2 averaged density n H 2 subscript 𝑛 subscript H 2 n_{\text{H}_{2}} (CS) act as lower limits.

HESS J1809–193
CS(1–0)
Reg.	Distance	$N_{\text{CS}}{\left[{10}^{12}\right]}^a$	$N_{{\text{H}}_2}{\left[{10}^{20}\right]}^{ab}$	$M_{{\text{H}}_2}$(CS)abc	$n_{{\text{H}}_2}$(CS)abc
	(kpc)	(cm-2)	(cm-2)	(${\text{M}}_{\odot }$)	(cm-3)
1	3.7	20	52	$3.2\times {10}^4$	$1.7\times {10}^2$
1-1	3.7	210(4)	520	$7.5\times {10}^3$	$1.8\times {10}^4$
1-2	3.7	270(7)	670	$1.3\times {10}^3$	$3.2\times {10}^4$
1-3	3.7	180(4)	460	$6.0\times {10}^3$	$2.2\times {10}^4$
1-4	3.7	190(2)	460	$9.0\times {10}^4$	$2.2\times {10}^4$
2	3.7	19	48	$3.2\times {10}^4$	$1.5\times {10}^2$
2-1	3.7	380(5)	950	$1.7\times {10}^3$	$4.4\times {10}^4$
2-2	3.7	290(3)	720	$1.3\times {10}^3$	$3.3\times {10}^4$
2-3	3.7	200(2)	510	$1.0\times {10}^3$	$2.3\times {10}^4$
3	2.7	15	30	$6.2\times {10}^3$	$2.5\times {10}^2$
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}$(CO)cd	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
1	3.7	$8.1\times {10}^4$	$4.4\times {10}^2$	$5.5\times {10}^3$
2	3.7	$2.3\times {10}^5$	$1.7\times {10}^2$	$6.0\times {10}^3$
3	2.7	$5.0\times {10}^4$	$1.0\times {10}^3$	$2.6\times {10}^3$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.
d:The $X_{\text{CO}}=2.0\times {10}^{20}$ cm-2/(K km/s) have been used to convert the integrated intensity $W_{\text{CO}}$ into H2 column density $N_{{\text{H}}_2}$

Table 11. Table D.2: Physical parameters obtained from our CS and CO analyses for the different selected regions located towards HESS J1026–582 . In the case of our CS analysis, we assumed a optically thin scenario and the derived column densities N H 2 subscript 𝑁 subscript H 2 N_{\text{H}_{2}} and N CS subscript 𝑁 CS N_{\text{CS}} , mass M H 2 subscript 𝑀 subscript H 2 M_{\text{H}_{2}} (CS) and H 2 averaged density n H 2 subscript 𝑛 subscript H 2 n_{\text{H}_{2}} (CS) act as lower limits.

HESS J1026–582
CS(1–0)
Reg.	Distance	$N_{\text{CS}}{\left[{10}^{12}\right]}^a$	$N_{{\text{H}}_2}{\left[{10}^{20}\right]}^{ab}$	$M_{{\text{H}}_2}$(CS)abc	$n_{{\text{H}}_2}$(CS)abc
	(kpc)	(cm-2)	(cm-2)	(${\text{M}}_{\odot }$)	(cm-3)
1	2.3	10	26	$1.7\times {10}^2$	$6.3\times {10}^2$
2	2.3	7	18	$1.2\times {10}^2$	$4.6\times {10}^2$
3	6.1	13	33	$4.5\times {10}^3$	$3.2\times {10}^2$
4	6.1	23	58	$1.7\times {10}^4$	$2.8\times {10}^2$
5	4.9	14	35	$1.3\times {10}^3$	$5.9\times {10}^2$
6	4.9	14	35	$1.8\times {10}^3$	$7.2\times {10}^2$
7	4.9	19	49	$1.9\times {10}^3$	$5.1\times {10}^2$
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}$(CO)cd	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
A	2.3	$1.1\times {10}^4$	$2.2\times {10}^1$	$5.0\times {10}^3$
B	4.9	$5.4\times {10}^4$	$3.2\times {10}^1$	$3.2\times {10}^4$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.
d:The $X_{\text{CO}}=2.0\times {10}^{20}$ cm-2/(K km/s) have been used to convert the integrated intensity $W_{\text{CO}}$ into H2 column density $N_{{\text{H}}_2}$

Table 12. Table D.3: Physical parameters obtained from our CO analysis for the different selected regions located towards HESS J1119–582 .

HESS J1119–582
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}$(CO)cd	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
A	5.0	$2.2\times {10}^4$	$5.0\times {10}^1$	$6.3\times {10}^3$
B	4.2	$7.1\times {10}^3$	$5.0\times {10}^1$	$3.6\times {10}^3$
C	8.6	$1.3\times {10}^4$	$6.1\times {10}^1$	$5.8\times {10}^3$
D	9.7	$1.3\times {10}^5$	$1.4\times {10}^1$	$7.8\times {10}^4$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.
d:The $X_{\text{CO}}=2.0\times {10}^{20}$ cm-2/(K km/s) have been used to convert the integrated intensity $W_{\text{CO}}$ into H2 column density $N_{{\text{H}}_2}$

Table 13. Table D.4: Physical parameters obtained from our CS and CO analyses for the different selected regions located towards HESS J1303–631 . In the case of our CS analysis, we assumed a optically thin scenario and the derived column densities N H 2 subscript 𝑁 subscript H 2 N_{\text{H}_{2}} and N CS subscript 𝑁 CS N_{\text{CS}} , mass M H 2 subscript 𝑀 subscript H 2 M_{\text{H}_{2}} (CS) and H 2 averaged density n H 2 subscript 𝑛 subscript H 2 n_{\text{H}_{2}} (CS) act as lower limits.

HESS J1303–631
CS(1–0)
Reg.	Distance	$N_{\text{CS}}{\left[{10}^{12}\right]}^a$	$N_{{\text{H}}_2}{\left[{10}^{20}\right]}^{ab}$	$M_{{\text{H}}_2}$(CS)abc	$n_{{\text{H}}_2}$(CS)abc
	(kpc)	(cm-2)	(cm-2)	(${\text{M}}_{\odot }$)	(cm-3)
1	6.6	18	45	$2.4\times {10}^3$	$3.8\times {10}^2$
2	6.6	10	25	$1.4\times {10}^3$	$2.2\times {10}^2$
3	6.6	7	19	$1.0\times {10}^3$	$1.6\times {10}^2$
4	1.5	15	37	$1.1\times {10}^2$	$1.4\times {10}^3$
5	1.5	9	22	$6.2\times {10}^2$	$8.6\times {10}^2$
6	12.6	13	33	$3.1\times {10}^3$	$1.6\times {10}^2$
7	12.6	11	28	$5.1\times {10}^3$	$1.3\times {10}^2$
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}$(CO)cd	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
A	6.6	$7.0\times {10}^4$	$1.8\times {10}^2$	$8.4\times {10}^3$
B	1.5	$4.8\times {10}^3$	$5.4\times {10}^2$	$7.4\times {10}^2$
C	12.6	$6.8\times {10}^4$	$1.8\times {10}^2$	$9.4\times {10}^3$
D	12.6	$4.1\times {10}^4$	$2.5\times {10}^2$	$8.5\times {10}^3$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.

Table 14. Table D.5: Physical parameters obtained from our CS and CO analyses for the different selected regions located towards HESS J1420–607 and HESS J1418–609 . In the case of our CS analysis, we assumed a optically thin scenario and the derived column densities N H 2 subscript 𝑁 subscript H 2 N_{\text{H}_{2}} and N CS subscript 𝑁 CS N_{\text{CS}} , mass M H 2 subscript 𝑀 subscript H 2 M_{\text{H}_{2}} (CS) and H 2 averaged density n H 2 subscript 𝑛 subscript H 2 n_{\text{H}_{2}} (CS) act as lower limits.

HESS J1420-607+HESS J1418-609
CS(1–0)
Reg.	Distance	$N_{\text{CS}}{\left[{10}^{12}\right]}^a$	$N_{{\text{H}}_2}{\left[{10}^{20}\right]}^{ab}$	$M_{{\text{H}}_2}$(CS)abc	$n_{{\text{H}}_2}$(CS)abc
	(kpc)	(cm-2)	(cm-2)	(${\text{M}}_{\odot }$)	(cm-3)
1	5.6	11	27	$2.2\times {10}^3$	$3.5\times {10}^2$
2	5.6	3	7	$6.4\times {10}^2$	$7.5\times {10}^1$
3	5.6	7	17	$7.8\times {10}^2$	$2.6\times {10}^2$
4	5.6	11	28	$3.3\times {10}^3$	$2.8\times {10}^2$
5	5.6	17	45	$3.8\times {10}^3$	$5.0\times {10}^2$
6	3.5	10	25	$1.6\times {10}^3$	$3.3\times {10}^2$
7	3.5	8	19	$6.2\times {10}^2$	$2.8\times {10}^2$
8	3.5	10	25	$5.1\times {10}^2$	$6.2\times {10}^2$
9	3.5	6	15	$2.1\times {10}^2$	$2.4\times {10}^2$
10	3.5	18	46	$2.3\times {10}^3$	$4.5\times {10}^2$
11	3.5	6	16	$2.5\times {10}^2$	$2.6\times {10}^2$
12	3.5	18	45	$7.0\times {10}^2$	$7.4\times {10}^2$
13	3.5	6	14	$2.1\times {10}^2$	$2.3\times {10}^2$
14	0.1	33	83	$2.7\times {10}^0$	$3.8\times {10}^4$
15	0.1	7	18	$2.3\times {10}^{-1}$	$1.1\times {10}^4$
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}^{cd}$	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
A	3.5	$3.8\times {10}^4$	$5.6\times {10}^2$	$2.2\times {10}^3$
6	3.5	$1.0\times {10}^4$	$2.0\times {10}^3$	$4.6\times {10}^2$
7	3.5	$3.0\times {10}^3$	$1.4\times {10}^3$	$2.4\times {10}^2$
8	3.5	$2.1\times {10}^3$	$2.5\times {10}^3$	$1.5\times {10}^2$
9	3.5	$1.5\times {10}^3$	$1.6\times {10}^3$	$1.1\times {10}^2$
10	3.5	$4.4\times {10}^3$	$8.3\times {10}^2$	$3.6\times {10}^2$
B	0.1	$6.8\times {10}^0$	$1.9\times {10}^3$	$1.0\times {10}^0$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.
d:The $X_{\text{CO}}=2.0\times {10}^{20}$ cm-2/(K km/s) have been used to convert the integrated intensity $W_{\text{CO}}$ into H2 column density $N_{{\text{H}}_2}$

Table 15. Table D.6: Physical parameters obtained from our CO analysis for the different selected regions located towards HESS J1018–589B .

HESS J1018–589B
CO(1–0)				Hi
Reg.	Distance	$M_{{\text{H}}_2}$(CO)cd	$n_{{\text{H}}_2}$(CO)cd	$M_{{\text{H}}_I}$
	(kpc)	(${\text{M}}_{\odot }$)	(cm-3)	(${\text{M}}_{\odot }$)
A	2.0	$2.9\times {10}^3$	$1.8\times {10}^2$	$8.0\times {10}^2$
a:Parameters have been derived using the LTE assumption.
b:The H2 physical parameters derived using a CS abundance ratio ${\chi }_{\text{CS}}=4\times {10}^{-9}$
c:A prolate geometry has been used in order to derive the mass and density.
d:The $X_{\text{CO}}=2.0\times {10}^{20}$ cm-2/(K km/s) have been used to convert the integrated intensity $W_{\text{CO}}$ into H2 column density $N_{{\text{H}}_2}$

Equations

$$n\left(E,R,t\right)$$

$$R_{\text{d}}$$

$$w_{\text{CR}}\left(R,t\right)={\int }_{{ϵ}_0}n\left(E,R,t\right)E\text{d}E$$

$${\dot{E}}_{\text{SD}}\left(t\right)$$

$$\Gamma$$

$$n\left(E,R,t\right)=\frac{E^{-\alpha }}{{\pi }^{3/2}}{\int }_{t_{\text{age}}}^0-\frac{S\left(\xi -t_{\text{age}}\right)}{R_d^3\left(\xi \right)}\exp \left(-\frac{R^2}{R_d^2\left(\xi \right)}\right)\text{d}\xi$$

$$B=10\mu \text{G}$$

$$B=10{\left(\frac{n_{\text{H}}}{300{\text{ cm}}^{-3}}\right)}^{0.6}\mu \text{G}$$

$$\frac{W_{CS(1-0)}}{W_{C^{34}S(1-0)}}=\frac{1-e^{-{\tau }_{CS(1-0)}}}{1-e^{\alpha {\tau }_{CS(1-0)}}}$$

$$N_{{\text{CS}}_1}=\frac{8k\pi {\nu }_{10}^2}{A_{10}h{c}^3}\left(\frac{\Delta {\Omega }_{\text{A}}}{\Delta {\Omega }_{\text{S}}}\right)\left(\frac{{\tau }_{\text{CS(1-0)}}}{1-e^{-{\tau }_{\text{CS(1-0)}}}}\right)\int T_{\text{mb}}\left(v\right)dv$$

$$N_{CS}=N_{C{S}_1}\left(1+\frac{1}{3}{e}^{2.35/T_{\text{kin}}}+\frac{5}{3}{e}^{-4.7/T_{\text{kin}}}+\text{...}\right)$$

$$M_{{\text{H}}_2}=\mu m_{\text{H}}\pi ab{N}_{{\text{H}}_2}$$

$$n_{{\text{H}}_2}=\frac{M}{4/3\left(\mu m_{\text{H}}\right)\pi a{b}^2}$$

Full text

Keywords: gamma-rays: ISM - ISM: individual objects (HESS J1809$-$193, HESS J1026$-$589, HESS J1119$-$614, HESS J1418$-$609, HESS J1420$-$607, HESS J1303$-$631, HESS J1018$-$589) - ISM: clouds - ISM: cosmic-rays - ISM: supernova remnants - ISM: pulsar wind nebulae - molecular data

1 Introduction

Pulsar wind nebulae (PWNe) represent the majority of the identified TeV gamma-ray sources in the Galactic plane (Abdalla et al. 2018a, b). The TeV emission is generally believed to be of leptonic origin, where high energy electrons are accelerated after crossing the pulsar termination shock, and scatter soft photons to produce inverse-Compton (IC) radiation at TeV gamma-ray energies. The interstellar medium (ISM) greatly affects the morphology of the PWN observed from radio up to gamma-ray energies. For instance, the interaction between the progenitor SNR and a nearby molecular cloud (MC) leads to an offset position of the pulsar with respect to the TeV peak along the pulsar-MC axis (Blondin et al. 2001).

In this paper, we investigate the ISM towards seven TeV PWNe and PWNe candidates (Acero et al. 2013; Abdalla et al. 2018a). We make use of the Nanten CO(1–0) survey (Mizuno & Fukui 2004) to illustrate the wide-field morphology of the diffuse molecular gas towards the PWNe and 7mm Mopra spectral line observations to probe the dense molecular, and possibly shocked, gas along with star forming regions. An accurate description of the ISM can help explain the morphology of the PWN (e.g Vela X, Moriguchi et al. 2001; Slane et al. 2018). Linking each TeV source to its local ISM can also provide additional constraints on its distance. We can identify target material for hadronic components (i.e. cosmic-rays) escaping the TeV source, for example from a progenitor supernova remnant (e.g. Voisin et al. 2016). Lastly, by combining ISM mapping with the improved sensitivity and angular resolution of the next generation Cerenkov Telescope Array (CTA, Acharya et al. 2013, 2017), we will be able to study the diffusion process of high energy particles towards and inside the ISM clouds.

In section 2, we briefly outline the technical properties of the Mopra and Nanten telescopes as well as the analysis used for our 7mm Mopra data reduction. We provide the results towards the different individual sources in section 3 and discuss the nature of the TeV source in section 4.

2 ISM Data and Analysis Procedure

2.1 7mm Mopra data and Analysis

We conducted $20^{\prime}\times 20^{\prime}$ and $10^{\prime}\times 10^{\prime}$ observations of the ISM towards several HESS TeV sources as part of our ‘MopraGam’ survey111http://www.physics.adelaide.edu.au/astrophysics/MopraGam/. The 7mm observations towards these studied sources were carried out between 2012 and 2015. Table 1 indicates the position and size of the different observations towards the seven regions which were mapped.

We used the Mopra spectrometer MOPS in ‘zoom’ mode, which records 16 sub-bands, each consisting of 4096 channels over a 137.5 MHz bandwidth. The ‘on the fly mapping’ (OTF) Nyquist sampled these regions with a 1*′* beam size and a velocity resolution of $\sim 0.2$ km/s. We can thus simultaneously observe tracers that can help understand the structure and morphology of molecular clouds. Among the tracers observed in this paper, both the carbon monosulfide and cyanoacetylene transitions CS(1$-$0) and HC3N(5$-$4,F=4$-$3) can be found in dense molecular gas and star forming regions (Irvine et al. 1987). Silicon monoxide SiO emission can also be detected inside shocked dense molecular clouds (Schilke et al. 1997; Gusdorf et al. 2008) and is a good signpost to claim physical association between a SNR and a nearby molecular cloud (e.g see Nicholas et al. 2012). Finally, the 44 GHz methanol class I maser CH3OH(I) generally indicates nearby active star forming regions (Voronkov et al. 2014).

Data reduction

We used Livedata222http://www.atnf.csiro.au/computing/software/Livedata/ to produce the spectra of each observation, calibrated by an OFF position, and subtracted the baseline using a linear fit. Then, we used Gridzilla333http://www.atnf.csiro.au/computing/software/Gridzilla/ to obtain a 3D cube showing the variation of the antenna temperature $T_{\text{A}}^{*}$ as a function of position and line of sight velocity $v_{\text{lsr}}$. Each region was mapped with a pixel spacing of $15^{\prime\prime}$. After a careful look at the data, we recursively performed linear and sinusoidal fits on each pixel in order to remove baseline ripples. Notably, each sinusoidal fit was performed with different initial wavelengths to account for ripples of various sizes. Velocity ranges with significant emission were masked during these fits. The cleaned data cubes were then smoothed via a Gaussian with a Full-Width-Half-Maximum (FWHM) of $1.25^{\prime}$ so as to mitigate spatial fluctuations. Finally, the cubes were Hanning smoothed using 5 channels in order to remove spikes, reduce the $T_{\text{rms}}$, and highlight broader CS(1–0) emission. An overview of the reduction can be found in Appendix A.

Producing CS(1–0) integrated intensity maps

We used the 7mm CS(1–0) transition to probe the denser gas. The standard method to produce CS(1–0) integrated intensity maps is to sum the temperature $T_{\text{A}}^{*}$ of each channel within a given velocity range. However, in the case where detections span a large velocity range ($\Delta v_{\text{lsr}}\sim 10\text{ to }30$ km/s), any narrow CS(1–0) detections are likely to fall below the noise level, as a large portion of the noise may be included. In order to deal with this issue, we adopted the method described in Appendix A. This method, although more complex, provides a lower $T_{\text{rms}}$ and cleaner integrated intensity maps, which could be used for multi-wavelength studies towards TeV sources.

Physical properties of CS(1–0) regions

From the CS(1–0) integrated intensity maps, we have selected clumps whose angular diameter were equal or greater than our Mopra CS beam size $\theta_{\text{FWHM}}\sim 1.6^{\prime}$. At the component velocity range, we fitted the CS(1–0) emission with Gaussian distributions.

We then used the Galactic rotation curve model from Brand & Blitz (1993) to obtain near/far kinematic distance estimates based on the peak velocity of each detection. In the case where the isotopologue C34S(1–0) is detected, we use Eq. D.9 to derive the averaged optical depth $\tau_{\text{CS(1--0)}}$, using the isotopologue ratio $\alpha=^{32}$S/${}^{34}S\sim 24$, based on terrestrial measurements (Fink 1981). Otherwise, an optically thin scenario $\tau_{\text{CS(1--0)}}=0$ was adopted. We obtained the column density of the upper state $N_{\text{CS}_{1}}$ using Eq. D.10. Assuming the gas to be in local thermal equilibrium, we thus obtained the total CS column density $N_{\text{CS}}$ using Eq. D.11, assuming a kinetic temperature $T_{\text{kin}}$=10 K typical of cold dense molecular clouds. The CS to H2 abundance ratio $X_{\text{CS}}$ inside dense molecular clumps varies between $10^{-9}$ and $10^{-8}$ as suggested by Irvine et al. (1987). In this work we chose $X_{\text{CS}}=4\times 10^{-9}$ as per Zinchenko et al. (1994) who studied molecular clouds associated with star-forming regions. As a result, we expect our H2 column density estimates to systematically vary by a factor of 2. Finally, we use Eqs. D.12 and D.13 to determine the total mass $M_{\text{H}_{2}}\left(\text{CS}\right)$ and H2 density $n_{\text{H}_{2}}\left(\text{CS}\right)$. The Gaussian fits for the CS(1–0) components towards the regions of the individual TeV sources can be found in Tables C.4 to C.9 and their derived physical properties in Tables D.1 to D.6.

2.2 CO(1–0) Data and Analysis

While the Mopra CO(1–0) survey is well under way (see Braiding et al. 2015, 2018), we used for this work the 4 metre Nanten CO(1–0) survey (Mizuno & Fukui 2004), as it encompassed all our sources. The Nanten telescope CO(1–0) survey covered the Galactic plane with a $4^{\prime}$ sampling grid, a velocity resolution $\Delta v=0.625$ km/s, and an averaged noise temperature per channel of $\sim 0.4$ K.

The CO(1–0) emission probes the more diffuse molecular gas surrounding PWNe. As PWNe expand inside a low density medium, we seek velocity ranges where the CO(1–0) emission spatially anti-corresponds with the TeV emission as that would support the PWN scenario (see Section 4.2). Molecular gas may also provide sufficient target material to produce hadronic TeV emission. Therefore, probing extended CO(1–0) emission overlapping the TeV emission is also a powerful means to test the hadronic scenario in the vicinity of a cosmic-ray (CR) accelerator.

We arbitrarily selected CO regions based on the prominence of the CO(1–0) integrated emission compared to the rest of the maps and their proximity with the TeV gamma-ray source, as these can be helpful to test the hadronic/leptonic scenario (see section 4). In the case where the pulsar is offset from the observed TeV gamma-ray source (i.e HESS J1026$-$583), we highlight molecular regions which could produce such asymmetry along the pulsar-molecular cloud axis. We also derive physical parameters towards some extended regions defined from our CS detections (e.g see HESS J1809$-$193 and HESS J1418$-$609). Assuming all the gas traced by CO is embedded within this CS region, we are able to derive an upper-limit on the averaged density $n_{\text{H}_{2}}$.

We fitted the CO components at the velocity range of interest, with Gaussian distributions. We then used the X-factor $X_{\text{CO}}=2\times 10^{20}$ cm*-2*/(K km/s) to convert the CO integrated intensity into H2 column density. Bolatto et al. (2013) have argued that this value is correct to $\sim 30$% across the Galactic plane. We then used Eq. D.12 to obtain the total molecular mass $M_{\text{H}_{2}}\left(\text{CO}\right)$ of the cloud, accounting for a 20% contribution from Helium. We finally assumed a prolate geometry to obtain the H2 volume occupied by the molecular gas, and the averaged particle density $n_{\text{H}_{2}}\left(\text{CO}\right)$ using Eq. D.13 (see Appendix D). The Gaussian fits for the CO(1–0) components towards the regions of the individual TeV source can be found in Tables C.4 to C.9 and their derived physical properties in Tables D.1 to D.6.

2.3 HI analysis with SGPS/GASS surveys

From the Southern Galactic Plane Survey (SGPS) and GASS Hi surveys (McClure-Griffiths et al. 2005; McClure-Griffiths et al. 2009), we also obtained the atomic Hi column density $N_{\text{H}\textsc{i}}$ assuming the (optically thin) conversion factor $X_{\text{H}\textsc{i}}=1.8\times 10^{18}$ cm*-2*/(K km/s) (Dickey & Lockman 1990). Hi self absorption may however occur towards colder Hi regions and, as a result, the optically thin assumption may underestimate $N_{\text{H}\textsc{i}}$ by a factor of 2 (Fukui et al. 2015). Comparing the total column density $N_{\text{H}}=N_{\text{H}\textsc{i}}+2N_{\text{H}_{2}}$ to the absorbed X-ray column density from X-ray counterparts can provide further constraints on the distance to the TeV source. The images showcasing the evolution of the total column density as a function of $v_{\text{lsr}}$ and distance towards several TeV sources can be found in Appendix E. Atomic gas should also be accounted for while testing the hadronic scenario towards TeV sources. Consequently, as per the CO(1–0) emission, we derived the atomic masses $M_{\text{H}_{I}}$ towards the selected CO regions. These masses can be found in Tables D.1 to D.6.

3 The ISM found towards the TeV sources

We list here various detections from our 7mm observations. For each source (see Figure 1 to 11) we labelled the regions with CS(1–0) detections in numerical order (e.g. ‘1’) while the regions with SiO(1–0, v=0), HC3N(5–4, F=4–3), CH3OH(I) detections are labelled ‘S’, ‘HC’ and ‘CH’ respectively. In this section however, we focus on combining our CS(1–0) detections with the Nanten CO(1–0) and primarily focus on the morphology of the various molecular regions surrounding our TeV PWNe and PWNe candidates, which is useful to understand their nature and morphology.

3.1 HESS J1809-193

HESS J1809–193 is a bright and extended TeV source whose position is coincident with several potential cosmic-ray (CR) accelerators (Aharonian et al. 2007). HAWC has also recently detected the TeV source 2HWC J1809$-$190, associated with HESS J1809$-$193 (Abeysekara et al. 2017). Araya (2018) recently found extended GeV emission associated with HESS J1809$-$193.

ASCA (Bamba et al. 2003) and Suzaku observations (Anada et al. 2010) revealed non-thermal X-rays likely associated with the pulsar PSR J1809$-$1917 (shown as a cyan diamond in Figure 1 and pink in Figure 2), with spin down power $\dot{E}_{\text{SD}}=1.8\times 10^{36}$ erg s*-1*, a characteristic age $\tau_{c}=51$ kyr and a dispersion measure distance $d\sim 3.7$ kpc (Cordes et al. 2002). The presence of two SNR shells G011.0$-$0.0 and G011.1$+$0.1 (shown as blue circles in Figure 1) both observed at 330 MHz and 1465 MHz (Brogan et al. 2004; Castelletti et al. 2016) adds more complexity to the picture. Additionally, the $\sim 2\,$kyr old millisecond pulsar PSR J1811$-$1925 with spin down energy $\dot{E}_{\text{SD}}=6.4\times 10^{36}$ erg s*-1* (Torii et al. 1999) and its progenitor SNR G011.2$-$0.3, located at $d\sim 4.4\text{ to }7$ kpc, are also positioned adjacent to HESS J1809$-$193 (see Figures 1 and 2) and thus might also contribute to the TeV emission. Finally, the anomalous X-ray magnetar XTE J1810$-$197 with period $P=5.54$ s and period derivative $\dot{P}=1.8\times 10^{-12}$ s s*-1* is located $\sim 0.35^{\circ}$ south of HESS J1809$-$193.

Based on the kinematic distance of PSR J1809$-$1917 ($d\sim 3.7$ kpc), we here focus on components at $v_{\text{lsr}}=10\text{ to }38\,$km/s (see blue and pink regions in Figure 3), although our Nanten CO(1–0) data have revealed several other components in the line of sight (see CO(1–0) emission in Figure 3).

At this velocity range, we remark that the molecular gas generally overlaps the TeV emission shown in black contours. Notably, at $v_{\text{lsr}}=25\text{ to }38$ km/s ($d\sim 3.7$ kpc), prominent CO emission is found south and east of the TeV emission, while the prominent CO emission spatially overlaps the TeV emission at $v_{\text{lsr}}=10\text{ to }25$ km/s ($d\sim 2.7$ kpc). Our Mopra CS(1–0) integrated intensity maps (see Figure 1 panels b., d., f.) however provides a clear view of the dense gas.

At v${}_{\text{lsr}}=25\text{ to }38$ km/s (Figure. 1, panels a., b.), the molecular clouds in the regions labelled ‘1’ and ‘2’ and located east and south of SNR G011.0$+$00.0 respectively appear very extended. From our CO and CS analyses, the masses derived in region ‘1’ attain $M_{\text{H}_{2}}\left(\text{CO}\right)=8.1\times 10^{4}\textrm{M}_{\odot}$ and $M_{\text{H}_{2}}\left(\text{CS}\right)=3.2\times 10^{4}\textrm{M}_{\odot}$ while we obtain $M_{\text{H}_{2}}\left(\text{CO}\right)=2.3\times 10^{5}\textrm{M}_{\odot}$ and $M_{\text{H}_{2}}\left(\text{CS}\right)=3.2\times 10^{4}\textrm{M}_{\odot}$ in region ‘2’. A significant fraction of the molecular gas in region ‘1’ and ‘2’ is therefore concentrated in clumps. The CS(1–0) emission in region ‘1’ also appears to anti-correspond with the ASCA X-ray diffuse emission shown in blue, supposedly produced by high energy electrons from the pulsar PSR J1809$-$1317 (see Anada et al. 2010).

Interestingly, we have found embedded dense filaments in region ‘1’ from C34S(1–0) and HC3N(5–4, F=4–3) detections (see Figure F.1 in Appendix F). In fact, the positions of the HC3N detections labelled ‘HC1 to HC3’ in Figure 2 coincide with the Spitzer infra-red (IR) dark cloud IRDC G011.11$-$00.11a$-$e (Parsons et al. 2009, see dotted white lines in Figure 2), confirming the molecular gas is foreground to the IR emission.

We have also identified two weak but broad (FWHM $\sim 4$ km/s) SiO(1–0) features in region ‘1’ (see Figure 2 and Figure 3 for spectra), labelled ‘S1’ and ‘S2’. The absence of overlapping IR continuum emission indicates the lack of active star-forming regions and could indicate a possible interaction between the molecular cloud and a non star-forming shock (Schilke et al. 1997; Gusdorf et al. 2008), which could here come from the adjacent SNR G011.0$-$0.0.

Based on the potential interaction between the SNR G011.0$-$0.0 and the dense molecular cloud at $v_{\text{lsr}}\sim 30\,$km/s, we thus suggest an alternate SNR distance $d\sim 3.7$ kpc compared to the distance $d\sim 3.0$ kpc claimed by Castelletti et al. (2016). If SNR G011.0$-$0.0 indeed interacted with the molecular cloud at $d\sim$3.7 kpc, it would then be located at the pulsar PSR J1809$-$1917 distance and quite likely be its progenitor SNR. Using eq. 3.33a from Cioffi et al. (1988), we derive the ambient density required to match the small projected radius $r_{\text{SNR}}\sim 5$ pc with the pulsar characteristic age $\tau_{c}=51$ kyr, is $n_{\text{amb}}\sim 370$ cm*-3*. The averaged density $n_{\text{H}_{2}}\left(\text{CO}\right)=440\,$cm*-3* found towards region ‘1’ appears somewhat consistent with this estimated averaged density.

In region ‘2’, we have also detected extended HC3N(5–4,F=4–3) emission and C34S(1–0) emission overlapping the infra-red dark cloud IRDC G010.71$-$00.16a$-$h (Parsons et al. 2009), highlighting another dense region. The morphology of the IR dark cloud and our HC3N detection appear elliptic and embed the anomalous X-ray magnetar XTE J1810$-$197 (see Gotthelf et al. 2004 and references therein), suggesting their potential physical association. Additionally, a prominent SiO(1–0) detection (see region ‘S3’ in Figure 2) was also found inside this dark cloud. This molecular cloud may thus be disrupted by another shock, perhaps caused by the progenitor SNR of XTE J1810$-$197.

At $v_{\text{lsr}}=10-25$ km/s (Figure 1 panel d.), we observe several CS(1–0) components inside the region here labelled ‘3’. We notably find that the CS(1–0) prominent emission corresponds with the JCMT CO(3–2) peaks found by Castelletti et al. (2016) and the Hii regions (purple contours and blue crosses in Figure 1 panel f. respectively).

We observe that the CS(1–0) emission averaged over the grid regions (Figure 1 panels e-f.) exhibits considerable variation of the peak velocity ranging between v${}_{\text{lsr}}=10-22$ km/s. For example, the two peaks at $v_{\text{lsr}}\sim 12$ and $18$ km/s inside ‘boxes 9 and 10’ merge to a single peaked emission with $v_{\text{lsr}}\sim 15$ km/s in ‘box 15’. This dense molecular region appears to host several Hii regions (see ‘H1 to H3’ in Figure 1.f). Consequently, the two spectral components may actually probe the same molecular gas, disrupted by the driving motion forces from various Hii regions. We also note that the molecular cloud anti-corresponds with the two SNRs. The disrupted gas could consequently be caused by one of these SNRs. As a result, we cannot rule out the SNR G011.0$-$0.0 distance $d\sim 3.0$ kpc suggested by Castelletti et al. (2016).

3.2 HESS J1026–583

The TeV source HESS J1026$-$583 was discovered from energy dependent morphology studies towards HESS J1023$-$591 (Abramowski et al. 2011). The latter source is thought to be powered by the colliding winds from Wolf-Rayet stars within the massive stellar cluster Westerlund 2 at $d=5.4^{+1.1}_{-1.4}$ kpc (Furukawa et al. 2009). Acero et al. (2013) catalogued HESS J1026$-$583 as a PWN candidate based on the detection of a nearby radio quiet gamma-ray pulsar PSR J1028$-$5810 (Ray et al. 2011) responsible for the GeV emission towards 3FGL J1028.4$-$5819 (shown as a red circle in Figure 4). PSR J1028$-$5819 has a spin down power $\dot{E}_{\text{SD}}=8.3\times 10^{35}$ erg s*-1*, characteristic age $\tau=89$ kyr, and a dispersion measure distance $d=2.3\pm 0.7$ kpc. However, HESS J1026$-$583 shows a hard VHE spectral index $\Gamma_{\gamma}=1.94$. It also does not exhibit any X-rays that are spatially coincident with the TeV source. Additionally, diffuse GeV gamma-ray emission has been detected with Fermi-LAT towards Westerlund 2 (Yang et al. 2017). The authors have argued for a hadronic origin based on its 200 pc extension of the 1 to 250 GeV emission. Consequently, a clear identification remains to be seen. Because of its proximity to HESS J1023$-$575, several ISM features have already extensively been studied (Dame 2007; Fukui et al. 2009; Furukawa et al. 2009; Furukawa et al. 2014; Hawkes et al. 2014).

Figure 4 shows the Nanten CO(1–0) integrated intensity at $v_{\text{lsr}}=-23\text{ to }-13$ km/s (panel a.) and $0\text{ to }20$ km/s (panel c.), inferring a distance $d\sim 2.3$ kpc and $d\sim 4.9$ kpc respectively, positioning the detections along the Carina arm (see Figure B.1 in Appendix B). At $v_{\text{lsr}}=-23\text{ to }-13$  km/s, we note that the CO emission does not overlap the TeV emission, nor the pulsar’s position. In fact, the molecular region with prominent CO emission, that we labelled ‘A’, is located east of the pulsar position, and could support the crushed PWN scenario (see Blondin et al. 2001). Assuming a kinematic distance $d\sim 2.3$ kpc, we obtain a total mass $M_{\text{H}_{2}}\left(\text{CO}\right)=9.7\times 10^{3}\textrm{M}_{\odot}$. Lastly, from our CS survey, we have also identified two unresolved CS(1–0) features labelled ‘1’ and ‘2’ in Figure 4 panel b.

At $v_{\text{lsr}}=0\text{ to }20$ km/s ($d\sim 4.9$ kpc, see Figure 4 panel c.), we observe prominent CO emission at $v_{\text{lsr}}\sim 4$ km/s, labelled ‘B’, spatially coincident with the HESS J1026$-$582 TeV peak. The molecular cloud in region ‘B’, with total mass $M_{\text{H}_{2}}\left(\text{CO}\right)=4.9\times 10^{4}\textrm{M}_{\odot}$ appears next to a possible shell like structure (see dashed purple and cyan ellipses in Figure 4 panels c. and e. respectively) which overlaps the HESS J1023$-$575 TeV emission. The molecular region is however located at the tangent of the Sagittarius arm (see Figure B.1), thus it is possible that the various molecular clouds found in Figure 4 panel c.) may be unrelated. From Figure 4 panel d., it nonetheless appears that the CO(1–0) emission at $v_{\text{lsr}}=2\text{ to }6$ km/s fills in some missing segment of the shell-like structure and might form a ring (indicated as a dashed-cyan ellipse) whose centre is located south-west of HESS J1026$-$582 observed at $v_{\text{lsr}}=6\text{ to }16$ km/s. The molecular structure at $v_{\text{lsr}}=6\text{ to }16$ km/s could therefore be physically connected to the molecular cloud in region ‘B’.

Among the CS(1–0) features found at $v_{\text{lsr}}=0\text{ to }20$ km/s, the gas clumps labelled ‘5’ to ‘7’, embedded in the molecular gas in region ‘B’ form a partial shell structure overlapping the TeV peak emission with a combined mass $M_{\text{H}_{2}}\left(\text{CO}\right)\sim 5.0\times 10^{3}\textrm{M}_{\odot}$, very similar to the CO mass obtained in region ‘B’. No massive stars have been catalogued towards the centre of this structure. However, we note the presence of the Hii region GAL 284.65$-$00.48 towards region ‘5’ and the reflective nebula GV 10.23.6 at the centre of the shell (see purple box and cross in Figure 4 panel d. respectively).

3.3 HESS J1119–614

Djannati-Ataï (2009) and Abdalla et al. (2018a) reported the detection of the TeV source HESS J1119$-$614 (see solid black contours in Figure 5). It is thought to be associated with either the PWN, powered by the radio pulsar PSR J1119$-$6127 (see cyan diamond in Figure 5), with spin down period $P=400$ ms (Camilo et al. 2000) and characteristic age $\tau_{c}=1.9$ kyr (Weltevrede et al. 2011); or its progenitor SNR G292.2$-$0.5 (see dashed blue circle in Figure 5). The $3^{\prime\prime}\times 3^{\prime\prime}$ PWN has been resolved in X-ray with Chandra (Gonzalez & Safi-Harb 2005; Safi-Harb & Kumar 2008). The shell of progenitor SNR has been observed both in the radio band with ATCA (Crawford et al. 2001), and in X-ray with ASCA between 0.4-10 keV (Pivovaroff et al. 2001). Follow-up studies with Chandra and XMM-Newton (Kumar et al. 2012; Ng et al. 2012) revealed additional information about the nature of the X-ray emission. Indeed, from the four regions studied by (Kumar et al. 2012) (here labelled ‘X1’ to ‘X4’ in Figure 5), the three regions ‘X1’ to ‘X3’ shown as blue ellipses, have a non-thermal component. However, the X-ray emission in region ‘X4’ shown as a blue-pink ellipse, is prominent and is thought to be of thermal origin.

Caswell et al. (2004) estimated a distance $d\sim 8.4\,$kpc for the pulsar and its progenitor SNR, based on the Hi and magnetic field studies. However, reconciling the characteristic age $\tau_{c}=1.9\,$kyr of PSR J1119$-$6127 and the $\sim 25$ pc diameter of the progenitor SNR requires a low density medium (Ng et al. 2012) and a massive progenitor star. Gonzalez & Safi-Harb (2005) argued for a SNR distance at $d=3.6\text{ to }6.3$ kpc by modelling the X-ray spectrum. Our gas study thus aims to provide additional constraints on the distance and how it could affect the nature of the TeV source.

No extended CS emission have been found towards this region. From our Nanten CO data, several components have been detected along the line of sight at $v_{\text{lsr}}\sim-30\,$km/s, $v_{\text{lsr}}\sim-10\,$km/s (near/far distance $d\sim 2.6/4.5$ kpc, near Carina-arm, see Figure B.1), $v_{\text{lsr}}\sim 20\,$km/s (distance $d\sim 8.6$ kpc, far Carina-arm), $v_{\text{lsr}}\sim 30\,$km/s (kinematic distance $d\sim 9.7$ kpc, far Carina-arm).

At $v_{\text{lsr}}=-20\text{ to }-10\,$km/s, we observe two molecular clouds (labelled ‘A’ and ‘B’ in Figure 5 panel a.) positioned north-east and west to the SNR respectively, and with respective masses $M_{\text{H}_{2}}\left(\text{CO}\right)=2.3\times 10^{4}\textrm{M}_{\odot}$ and $M_{\text{H}_{2}}\left(\text{CO}\right)=7.1\times 10^{4}\textrm{M}_{\odot}$. We also note that the molecular ISM anti-corresponds with all X-ray regions.

At $v_{\text{lsr}}\sim 20\text{ to }40\,$km/s we observe that the bulk of the molecular gas is found in two regions labelled ‘C’ and ‘D’, with masses $M_{\text{H}_{2}}\left(\text{CO}\right)=1.3\times 10^{4}\textrm{M}_{\odot}$ and $M_{\text{H}_{2}}\left(\text{CO}\right)=2.3\times 10^{5}\textrm{M}_{\odot}$ respectively. The morphology of the gas overlaps the TeV gamma-ray detection. It also corresponds with the thermal X-ray in region ‘X1’ both at $v_{\text{lsr}}\sim 20\,$km/s and $v_{\text{lsr}}\sim 30\,$km/s (see Figure 5 panel d.), highlighting potential SNR-MC interaction. The molecular gas also anti-corresponds with the other X-ray regions ‘X2’ to ‘X4’, thought to have non-thermal components likely produced by the high energy electrons from both the SNR and the PWN. Therefore, we argue that the morphology of the CO(1–0) both at $v_{\text{lsr}}\sim 20\,$km/s and $v_{\text{lsr}}\sim 30\,$km/s appears consistent with the X-ray results, inferring a source kinematic distance $d=8.6\text{ to }9.7$ kpc. We finally remark that these two components may be physically connected and highlight motion caused by the progenitor star.

We also compare the column densities, shown in Figure E.1, towards the four X-ray regions. Towards the regions ‘X1’ and ‘X2’, we match the X-ray modelled column density $N_{\text{H}}$ (grey regions in Figure E.1 top panels) at $v_{\text{lsr}}=22\text{ to }30$ km/s, also suggesting a distance $d=8.6-9.6$ kpc, positioning the source in the far-Carina arm. Towards the X-ray regions ‘X3’ and ‘X4’ however, our column densities do not reach the X-ray one, or yield to unreasonable large distance $d\sim 15.5$ kpc. Clumps unresolved by the Nanten CO(1–0) may account for the missing components. Our column density studies thus provide a lower-limit kinematic distance $d>8.6$ kpc, consistent with our morphological studies and with Caswell et al. (2004). To explain both the large SNR radius and the enhanced thermal X-ray emission in region ‘X4’ , we argue that the SNR once propagated in a low density medium (Ng et al. 2012), until it reached the denser gas in region ‘D’.

3.4 Kookaburra and Rabbit

The two TeV sources HESS J1418$-$609 and HESS J1420$-$607 have been classified as PWNe based on their spatial coincidence with the X-ray (Roberts & Romani 1998) and GeV gamma-ray counterparts (Acero et al. 2013). Ng et al. (2005) indicated that two diffuse non-thermal X-ray sources were associated with the pulsar PSR J1420$-$6048 (labelled P1 in Figure 6), a 68.2 ms period pulsar with spin down power $\dot{E}_{\text{SD}}=1.0\times 10^{37}$ erg s*-1*, characteristic age $\tau_{c}=13$ kyr and dispersion measure distance $d\sim 5.6$ kpc ; and the 108 ms radio-quiet gamma-ray pulsar PSR J1418$-$6058 (labelled P2), with a characteristic age $\tau_{c}=1.6$ kyr. It has also been noted that P2’s location is offset ($\sim 8.4^{\prime}$) from the HESS J1418$-$609 TeV peak.

Roberts et al. (2005) argues that the expected cometary shape of the PWN produced by a pulsar with a high space velocity, would explain the offset position of the gamma-ray emission with respect to the pulsar. The distance to P2 have not yet been constrained. Ng et al. (2005) has argued a PWN distance $d=2\text{ to }5$ kpc, while Wang (2011) has claimed a much smaller distance $d=1.4\text{ to }1.9$ kpc. Thus, by looking at the gas distribution in various velocity ranges, we aim to highlight the distance which would support the PWN scenario.

From our CO(1–0) observations, we have detected several molecular complexes along the line of sight. We mostly focus at the velocity ranges $v_{\text{lsr}}=-65\text{ to }-55$ km/s, ($d\sim 5.6$ kpc), $\sim 50$ km/s (Scutum Crux arm, distance $d\sim 3.5$ kpc), and $\sim-5$ km/s (local arm, $d\sim 0.1$ kpc) respectively, as shown in Figure 6. In Figure 7, we also note a CO component at $v_{\text{lsr}}\sim-25$ km/s. As opposed to the other velocity ranges, no prominent CS detections have been found at $v_{\text{lsr}}\sim-25$ km/s in this region, and thus we suggest that the molecular cloud is located on the far distance.

At $v_{\text{lsr}}=-65\text{ to }-55$ km/s (see Figure 6a.), we note the bulk of the CO(1–0) emission is located towards the west and north side of HESS J1418$-$609. We also remark that the molecular emission shows little overlap with any of the TeV sources at this velocity range. We observed extended CS(1–0) emission, labelled ‘1’ to ‘5’, north of HESS J1418$-$609. Assuming a distance $d=5.6$ kpc, the combined mass of these clumps attains $M_{\text{H}_{2}}\left(\text{CS}\right)=1.1\times 10^{4}\textrm{M}_{\odot}$. Notably, the bulk of the CS(1–0) emission appears to wrap around the TeV emission (as expected for leptonic IC emission).

At $v_{\text{lsr}}=-52\text{ to }-42$ km/s, we note that the CO(1–0) emission overlaps the two TeV sources, and peaks west of HESS J1018$-$609 and south east of HESS J1420$-$607. However, due to the small velocity separation between the components and the components at $v_{\text{lsr}}=-65\text{ to }-55$ km/s, it is somewhat difficult to accurately describe the morphology of the diffuse molecular gas at these velocities. We particularly observe that the prominent CO emission in region ‘A’, with mass $M_{\text{H}_{2}}\left(\text{CO}\right)=3.5\times 10^{4}\textrm{M}_{\odot}$ and averaged density $n_{\text{H}_{2}}\left(\text{CO}\right)=5.6\times 10^{2}$ cm*-3*, anti-corresponds with the TeV source HESS J1418$-$609. Our CS(1–0) results highlight a filamentary structure (see regions ‘6’ to ‘10’ in Figure 6d.), with a combined CS mass $M_{\text{H}_{2}}\left(\text{CS}\right)=5.2\times 10^{3}\textrm{M}_{\odot}$ (CO mass $M_{\text{H}_{2}}\left(\text{CO}\right)=3.1\times 10^{4}\textrm{M}_{\odot}$) which crosses HESS J1420$-$607. Interestingly, the dense gas in regions ‘6’ to ‘8’ also appears in a shell-like arrangement centred towards the south-east of HESS J1420$-$607. CS(1–0) clumps, labelled ‘11’ to ‘13’ are also found towards HESS J1418$-$609. The regions ‘11’ and ‘12’ are in fact coincident with the TeV peak emission.

At $v_{\text{lsr}}=-10\text{ to }0\,$km/s (bottom panels), we found prominent CO emission south of HESS J1418$-$609, east of HESS J1420$-$607, and particularly between the two pulsars. Notably, the latter molecular region (labelled ‘B’) with mass $M_{\text{H}_{2}}\left(\text{CO}\right)=6.8\textrm{M}_{\odot}$ nests an extended dense clumps (see region ‘14’) with mass $M_{\text{H}_{2}}\left(\text{CS}\right)=2.7\textrm{M}_{\odot}$ and with averaged density reaching $n_{\text{H}_{2}}\left(\text{CS}\right)=3.8\times 10^{4}$ cm*-3*.

In Figure E.2, we compare the column density from CO and Hi measurements with the column density $N_{\text{H}}=2.7\pm 0.2\times 10^{22}$ cm*-2* towards the pulsar PSR J1418$-$6058 derived from X-ray measurements (Kishishita et al. 2012). Assuming that all components, except that at $v_{\text{lsr}}\sim 25$ km/s, are in the near distance, we observe in Figure E.2 that the column density value from X-ray measurements is matched at v${}_{\text{lsr}}\leq-50$ km/s, which infers a distance range $d\sim 3.5$ kpc. As a result, this column density X-ray study appears to favour a PWN location at $d\sim 3.5\,$kpc ($v_{\text{lsr}}=-52\text{ to }-42$ km/s), while our morphological study of the molecular ISM favours a PWN distance at $d\sim$5.6 kpc ($v_{\text{lsr}}=-65\text{ to }-55$ km/s), as it anti-corresponds with the TeV gamma-ray emission (see Figure 6).

3.5 HESS J1303-631

HESS J1303$-$631 was first classified as a ‘dark source’ due to its lack of any counterparts at other wavelengths (Aharonian et al. 2005; Acero et al. 2013). However, energy dependent morphology of the TeV source (see Abramowski et al. 2012b) unambiguously highlighted its association with the pulsar PSR J1301$-$6305 (P1 in Figure 8) with spin down energy $\dot{E}_{\text{SD}}=2.6\times 10^{36}$ erg s*-1*, a rotation period $P=184$ ms, and a characteristic age $\tau_{c}=11$ kyr. Follow-up observations with XMM-Newton revealed diffuse X-ray emission towards the pulsar with a power law spectral index $\Gamma_{X}=2.0^{+0.6}_{-0.7}$(Abramowski et al. 2012b). Acero et al. (2013) detected a GeV counterpart with a gamma-ray spectral index $\Gamma_{\gamma}=1.7$. Finally, from the 1.384 GHz ATCA observations, Sushch et al. (2017) recently announced the presence of a plausible SNR radio shell with radius $\sim 12^{\prime}$ next to the pulsar PSR J1301$-$6305, although this association appears unlikely.

Based on the dispersion measure of this pulsar, Cordes et al. (2002) suggested a distance $d\sim 6.6$ kpc, much closer than the previous distance $d\sim 12.6$ kpc (Taylor & Cordes 1993). From the Nanten CO(1–0) components identified (see Figure 9), we focus on several molecular complexes in the line of sight at $v_{\text{lsr}}=-35\text{ to }-25$ km/s (distance $d\sim 6.6$ kpc, Scutum Crux arm ), $v_{\text{lsr}}=-25\text{ to }-15$ km/s (distance $d\sim 1.5$ kpc, near Sagittarius-Carina arm), $v_{\text{lsr}}=-5\text{ to }5$ km/s ($d\sim 0.1$ kpc, local arm) and $v_{\text{lsr}}=25\text{ to }35$ km/s (distance $d\sim 12.6$ kpc, far Sagittarius-Carina arm) shown in Figure 8. From our 7mm CS observations, which cover the north-west part of the SNR towards the TeV source, we have identified several molecular clumps, that we have labelled ‘1 to 7’ (see Tables C.8 and D.4 for their physical parameters) but no extended CS(1–0) emission has been detected.

At all the aforementioned velocity ranges, it appears that some CO(1–0) emission always overlaps the TeV emission. It should be noted that the CO emission peaks south of HESS J1303$-$631 inside region ‘A’ at $v_{\text{lsr}}=-35\text{ to }-25$ km/s. At the local arm ($v_{\text{lsr}}=-5\text{ to }5\,$km/s), most of the CO(1–0) emission is distributed north east of the TeV source. At $v_{\text{lsr}}=15\text{ to }35$ km/s, prominent CO emission is found overlapping HESS J1303$-$631 (see regions ‘C’ and ‘D’). We also notice that the CO(1–0) emission does overlap the position of the SNR candidate represented in blue circle in Figure 8 at $v_{\text{lsr}}=-35\text{ to }-25$ km/s and $v_{\text{lsr}}=15\text{ to }35$ km/s.

Interestingly, we have found a CO(1–0) emission dip at $v_{\text{lsr}}=-25\text{ to }-15$ km/s localised towards the SNR candidate. From the position-velocity plots shown in Figure 10, we observe prominent CO emission which in fact appears to surround the SNR position (whose boundaries are shown in red dashed lines in Figure 10 bottom panels) between $v_{\text{lsr}}=-22\text{ to }-15$ km/s. From the CO(1–0) integrated intensity region at this velocity range, we also observe little spatial overlap between the molecular gas and the SNR candidate. Consequently, it may highlight the presence of a putative molecular shell at $d\sim 1.5$ kpc coincident with the recently observed SNR candidate. The green cross and ellipse in Fig. 10 indicate the position and the expansion speed ($v_{\text{exp}}\sim 4$ km/s) of the putative molecular shell surrounding the SNR candidate. From the dashed magenta circular region in Figure 10 centred towards the SNR candidate position, the maximum mass swept by the SNR candidate reaches $M_{\text{H}_{2}}=3.3\times 10^{4}\textrm{M}_{\odot}$. We obtain as a upper limit a required kinetic energy $E_{\text{kin}}=1/2M_{\text{H}_{2}}v_{\text{exp}}^{2}=5.1\times 10^{47}$ erg which represents $\sim 0.05\%$ of the total kinetic energy from powerful O stars stellar winds over 1 Myr time-scale (see Weaver et al. 1977 for detailed study on interstellar bubbles). Consequently, the molecular shell may have been produced by the SNR progenitor star. We thus suggest this new SNR candidate is at distance at $d\sim$1.5 kpc, supporting its non-physical association with PSR J1301$-$6305 claimed by Sushch et al. (2017).

3.6 HESS J1018–589

HESS J1018$-$589 was first reported by Abramowski et al. (2012a) and actually consists of two distinct sources. The gamma-ray binary 1FGL J1018.6$-$5856 appears to be responsible for the HESS J1018$-$589A TeV emission (Abramowski et al. 2012a) while HESS J1018$-$589B (shown in dashed black circle in Figure 11) is thought to be a PWN powered by the pulsar PSR J1016$-$587 (shown by a cyan diamond), with a rotation period $P=107$ ms, a spin down energy $\dot{E}_{\text{SD}}=2.6\times 10^{36}$ erg, and a characteristic age $\tau_{c}=21$ kyr (Camilo et al. 2001). It has been suggested that the pulsar, with dispersion measure distance $d$=8 kpc, is not associated with the nearby SNR G292$-$1.8 located at $d\sim 2.9$ kpc (Ruiz & May 1986).

From the Nanten CO(1–0) observations shown in Figure 11, we have identified CO emission at $v_{\text{lsr}}=-23\text{ to }-10$ km/s, inferring a near distance $d\sim 2.8$ kpc, matching the SNR distance. The CO(1–0) emission appears filamentary north of HESS J1018$-$589B , and partially overlaps the northern rim of SNR G292$-$1.8. The molecular gas in region ‘A’, with mass attaining $M_{\text{H}_{2}}\left(\text{CO}\right)=2.9\times 10^{3}\textrm{M}_{\odot}$, shows quite broad emission ($\Delta v\sim 12$ km/s, see Table C.9 and Figure 11). Although no CO emission was revealed at $v_{\text{lsr}}\sim 30$ km/s, we note from Figure E.3 significant Hi emission, which may suggest that the ISM surrounding HESS J1018$-$589B mostly consists of atomic gas.

4 Discussion of gamma-ray emission

In this section, we use results from our ISM studies to discuss whether the CRs and/or high energy electrons interacting with the ISM can contribute to the TeV emission or at least affect their morphology. In this section, we first introduce the method used to check whether hadronic cosmic-rays could contribute to the observed TeV emission. Then, we will briefly indicate how leptonic emission can also be affected by the ISM.

4.1 TeV emission from CRs

Based on the mass estimates towards molecular regions overlapping the TeV sources, we use eq. 10 from Aharonian (1991) to derive the cosmic-ray enhancement factor $k_{\text{CR}}=w_{\text{CR}}/w_{\odot}$ which represents the ratio between $w_{\text{CR}}$ the energy density of CRs towards a molecular cloud next to a TeV source, and $w_{\odot}\sim 1.0$ eV cm*-3* the energy density found in the solar neighbourhood. We use the combined atomic and molecular mass to obtain the total amount of target material available. Table 2 indicates the $k_{\text{CR}}$ values required for each molecular regions (partially) overlapping the TeV sources to account for the observed TeV fluxes. The $k_{\text{CR}}$ value can then be compared to typical cosmic-ray enhancement factors predicted in the vicinity of SNRs (and potentially towards PWNe) to check the plausibility of hadronic contribution.

Contribution from nearby SNRs?

Nearby SNRs are the most likely viable candidates to produce CRs energy densities up to $\sim 10^{3}$ eV cm*-3* (see Reynolds 2008 and references therein). CRs propagate along the magnetic field lines and they scatter from their interaction with magnetic field perturbations (provided the scale of the perturbation roughly equals the CR gyroradius). As the magnetic field turbulence is thought to be enhanced in molecular clouds, we here assume an isotropic diffusion of CRs and electrons as a first order approximation. If we assume an isotropic diffusion of CRs and neglect energy losses, we can estimate the energy density distribution of CRs at a distance $R$ from the SNR (here assumed as an impulsive source of CRs, see Aharonian & Atoyan 1996):

[TABLE]

with $m_{\text{p}}$ being the proton mass, $D_{10}=10^{28}$ cm${}^{2}\,$s*-1* being the diffusion coefficient of 10 GeV CRs, $\alpha$ the spectral index of the proton distribution, $\eta_{\text{pp}}$ the ratio of the total SNR energy $E_{\text{SNR}}$ transferred to CRs. $R_{\text{d}}\left(t\right)$ represents the diffusion radius travelled by CRs after a time $t$. The diffusion suppression factor $\chi$ accounts for slower diffusion of particles which can be caused, for instance, by streaming instabilities or perturbations caused by shocks (see Nava et al. 2016, Malkov et al. 2013 and references therein). Here we use, $\chi=0.01\text{ to }1$ which matches the slow and fast diffusion coefficient regimes defined by Aharonian & Atoyan (1996). Various studies of the cosmic-rays interacting with nearby molecular clouds W 28 (Giuliani et al. 2010; Li & Chen 2010; Gabici et al. 2010) have, in fact, suggested a suppression factor between $\chi=0.01\text{ and }0.1$. However, the value of $\chi$ remains poorly constrained. From Eq. 1, we can then obtain the total energy density of CRs $w_{\text{CR}}\left(R,t\right)$ using the following equation:

[TABLE]

with $\epsilon_{0}=m_{\text{p}}c^{2}$ is the proton energy at rest. Figure 12 illustrates the range of $k_{\text{CR}}$ produced by SNRs with initial energy $E_{\text{SNR}}=10^{51}$ erg as a function of the SNR age, using a standard proton spectral index $\alpha=2.2$ and $\eta_{\text{pp}}=0.1$. From the distance between the SNRs and the surrounding ISM regions (see Table 2), and the age of the SNRs, we can then check whether the required enhancement factors from Table 2 fall within the predicted range from nearby SNRs.

CR contribution from PWNe?

Additionally, a few authors (see Amato 2014 and references therein) have also argued that high energy hadrons could also be produced inside the pulsar environment, and be responsible for several features inside the PWN (e.g wisps in the Crab PWN, Gallant & Arons 1994). Providing high energy CRs have indeed been produced inside PWNe and not suffered heavy adiabatic losses, we will thus discuss whether the pulsars considered in this work could also generate the required cosmic-ray enhancement factors shown in Table 2. In order to model the high cosmic-ray energy density potentially produced by those pulsars, we account for the evolution of the spin down power $\dot{E}_{\text{SD}}$ as a function of time $t$, which can be described as :

[TABLE]

with $t_{\text{age}}$ being the current age of the PWN, $P\left(t\right)$ and $\dot{P}\left(t\right)$ being the current period and period derivative of a pulsar respectively at time $t$, $\dot{E}_{\text{SD}}\left(t\right)$ being the pulsar spin down power at time $t$, and $n_{\text{b}}$ being the pulsar braking index. In order to obtain the density of CRs at a given radius $R$ from the pulsar, we rewrite Eq. 1 using the source term $S\left(t\right)=\eta_{\text{pp}}\dot{E}_{\text{SD}}\left(t\right)/\left(m_{\text{p}}c^{2}\right)^{2-\alpha}$, with $\eta_{\text{pp}}$ here being the fraction of the spin down power transferred to CRs. We thus numerically solve the following equation :

[TABLE]

with $\xi=t_{\text{age}}-t$. As per the SNR scenario, we obtain the CR energy density at a given distance $R$ from pulsar using Eq. 3.

Figure 13 illustrates the CR energy density produced by the various pulsars as a function of the diffusion coefficient suppression factor $\chi$. Here, we assumed $n_{\text{b}}=3$ for all pulsars, except for PSR J1119$-$6127 in which we used the measured braking index $n_{\text{b}}=2.68$. From the spectral modelling towards several PWNe, Bucciantini et al. (2011) suggested that the maximum energy rate available to be transferred to CRs must be at most 0.20 of the total spin down power. We thus have used $\eta_{\text{pp}}=0.20$ and the $k_{\text{CR}}$ values shown in Figure 13 are consequently upper-limits. In the following subsections, we compare these predictions with the various required $k_{\text{CR}}$ (see Table 2) derived inside the molecular regions.

4.2 TeV emission from high energy electrons

As opposed to CRs, high energy electrons suffer heavy radiation losses as they diffuse inside the dense ISM because of the potentially enhanced magnetic field strength (see Crutcher et al. 2010). In the lack of intense radiation fields, synchrotron losses, with time-scale $\tau_{\text{sync}}=2.45\times 10^{7}\left(B^{2}_{\text{mg}}E/m_{\text{e}}c^{2}\right)^{-1}$ yr ($B_{\text{mg}}=B/1\text{mG}$ being the magnetic field strength and $E$ being the total energy of the electron) is likely to dominate over inverse-Compton losses inside molecular clouds. Although Bremsstrahlung emission, with time-scale $\tau_{\text{brem}}\sim 4.6\times 10^{7}\left(n/1\text{ cm}^{-3}\right)^{-1}$ yr, may contribute inside molecular clouds with densities $n_{\text{H}}>\text{ a few}\times 10^{3}$ cm*-3* in the 0.1 to 1 TeV band, we expect the TeV emission above 1 TeV to anti-correspond with the molecular ISM as the IC radiation is likely to dominate above 1 TeV. Table 3 indicates the synchrotron energy loss time-scale $\tau_{\text{sync}}$ for electrons with energy $E=10$ TeV, required to produce photons with energy $E_{\gamma}\sim 1$ TeV. To obtain $B$ values inside molecular clouds, we use the Crutcher et al. (2010) relation :

[TABLE]

We compare these time-scales with the age of nearby high energy sources and the diffusion time-scale $\tau_{\text{diff}}$ required for the particles with energy above 10 TeV to fully cross a molecular region. We here define $\tau_{\text{diff}}=\left(\langle r\rangle+d_{\text{source}}\right)^{2}/6D\left(E,B,\chi\right)$ with $\langle r\rangle$ being the mean radius of the molecular region, and $d_{\text{source}}$ the distance separating the centre of the high energy source to the centre of the molecular region. To check the leptonic scenario, we first select dense molecular regions which anti-correspond with each TeV source and/or its X-ray counterpart. Here, we preferably use the extended CS regions as most of the observed gas is confined. However, in the cases where no CS(1–0) emission was detected (i.e. HESS J1018$-$189B, HESS J1303$-$631 and HESS J1119$-$164), we use the CO regions. We again note that the densities $n_{\text{H}_{2}}$ from our CO analyses are upper-limits. Thus, these synchrotron time-scale values should be, in these cases, used as upper-limits. Then, we derive the synchrotron time-scale and compare it to $\tau_{\text{diff}}$ required for the high energy electrons to cross the cloud. In the case where the synchrotron time-scale is much smaller than the diffusion time-scale, we expect the high energy electrons to radiate most of their energy and thus become unable to produce X-rays and TeV emission inside the molecular cloud. An anti-correlation between the ISM and X-rays would support the physical connection between the TeV source and the ISM as well as the leptonic nature of the TeV source (see Sano et al. 2013).

4.3 Discussion of Individual Sources

HESS J1809$-$193

Castelletti et al. (2016), and more recently Araya (2018), have discussed the hadronic scenario in which they assume the SNR G011.0$-$0.0 to be located at a distance $d\sim$3.0 kpc. We focus here on how molecular clouds at $d\sim 3.7$ kpc can affect the TeV emission produced by the PWN powered by PSR J1809$-$1917 and SNR G011.0$-$0.0, in which we argue they are physically connected.

We remark in Table 3 that the synchrotron energy loss time-scale $\tau_{\text{sync}}$ for electrons above 10 TeV is significantly smaller than the diffusion time-scale through the various ISM regions. Consequently, most electrons producing gamma-rays above 1 TeV would not penetrate the dense molecular regions.

We also check whether the CRs in molecular clouds may produce significant TeV emission. To match the observed ISM regions, the TeV flux was scaled simplistically by area ratio, that is 16%, 16% and 12% towards regions ‘1’, ‘2’ and ‘3’ respectively. We conclude that a cosmic-ray enhancement factor $k_{\text{CR}}=10-30$ (see Table 2) is required to attain the TeV flux towards region ‘1’, ‘2’ and ‘3’ respectively. As the distances between these SNRs and the molecular regions range between 10 and 20 pc (see Table 3), we find, based on the red and green regions in Figure 12 that, at $t_{\text{age}}=51$ kyr, the predicted $k_{\text{CR}}$ ranges to be between $k_{\text{CR}}\sim 1\text{ and }300$. Consequently, the SNR G011.0$-$0.0 may be able to produce the required cosmic-ray enhancement and produce hadronic TeV emission. Interestingly, we find that the PWN powered by PSR J1809$-$1917 (see black line in Figure 13) could also significantly increase the predicted cosmic-ray enhancement factor towards nearby molecular clouds with $k_{\text{CR}}$ values peaking at $\sim 270$ (for diffusion suppression $\chi=0.02$) and thus generally exceeds the aforementioned $k_{\text{CR}}$ towards the molecular regions ‘1’ and ‘2’. Although these CR energy densities should be considered as upper-limits, the PWN powered by the PSR J1809$-$5158 may consequently be a viable laboratory to probe possible hadronic components originating from the PWN.

HESS J1026$-$583

It has already been stated that, based on energetics alone, the spin down energy of the pulsar PSR J1028$-$5819 could produce the HESS J1026$-$583 TeV emission (Abramowski et al. 2011).

From our CO(1–0) and CS(1–0) data, we have however identified a dense molecular region spatially coincident with the TeV emission. If the TeV emission is indeed of hadronic origin, we also require a proton spectrum $J\left(E\right)\propto E^{-1.94}$ to generate the photon spectrum towards HESS J1026$-$583 (see Abramowski et al. 2011).

Taking into account this hard gamma-ray spectra, the eq. 10 from Aharonian (1991) must be adjusted by a factor of $1.6/0.94=1.7$. We thus obtain the cosmic-ray enhancement factor $k_{\text{cr}}=60$. Using the combined mass of the regions ‘5’ to ‘7’ $M_{\text{H}_{2}}\left(\text{CS}\right)=5.0\times 10^{3}\textrm{M}_{\odot}$, the required cosmic-ray enhancement factor only increases to $k_{\text{CR}}=103$. We find that this value can only be attained if there is a SNR located art $d<10$ pc to the molecular region ‘B’ (see red region in Figure 12). However, no SNRs nor any cosmic-ray accelerators have been detected close to this molecular cloud so far.

HESS J1119$-$164

The centroid of the TeV emission appears to be located north-west of PSR J1119$-$6127, towards regions ‘X1’ and ‘X2’. The non-thermal X-ray emission, thought to be produced by high energy electrons inside the front shock, is unlikely to have also produced the TeV emission. The possible cause of the TeV gamma-ray emission may either come from high energy electrons from the PWN, or the cosmic-rays produced by the SNR front shock as they interact with the ISM in region ‘D‘. By modelling the GeV and TeV spectral energy distribution using Fermi-LAT and HESS data respectively, Acero et al. (2013) showed that electrons originating from the PWN could reproduce the observed spectrum. The offset of the TeV emission with respect to the pulsar position could suggest that the PWN may have been crushed by the reverse shock of the progenitor SNR. The Nanten CO(1–0) survey has however not revealed any significant emission east of PSR J1119$-$6127 between $v_{\text{lsr}}=20\text{ to }40\,$km/s. A higher resolution gas survey could resolve clumps which could be responsible for the asymmetry.

The hadronic scenario may also be plausible as the TeV gamma-ray emission corresponds with the CO(1–0) emission west of SNR G292.2$-$0.5. The thermal X-ray emission towards region ‘D’ suggest that some of the gas in region ‘D’ may have been heated as it interacted with the SNR. From Table 2, the required cosmic-ray enhancement factor $k_{\text{CR}}$ required for the cosmic-rays to produce the TeV gamma-ray flux is $\sim 143$. A SNR age $t=1.9\,$kyr suggests that the bulk of CRs are remained confined inside the SNR front shock. However, the steep photon spectral index $\Gamma=2.64$ above 1 TeV (Abdalla et al. 2018a) also indicates that CRs above 10 TeV have already escaped the SNR. With the increased sensitivity of CTA, one could look at the spatial and spectral evolution of the TeV gamma-ray emission inside the ISM towards region ‘D’. A hardening of the TeV gamma-ray spectrum with increasing distance of the SNR could favour the hadronic scenario.

Finally, as shown by the orange line in Figure 13, we also find that the pulsar is unlikely to have produced enough CRs to produce the TeV emission towards HESS J1119$-$164.

HESS J1418$-$609 and HESS J1420$-$607

From Table 2, we observe $k_{\text{CR}}$ values ranging between 50 and 90 towards the clumps ‘6’ to ‘10’, whose projected distance ranges between 22 pc and 35 pc (see Table 3). Assuming the progenitor SNR of PSR J1418$-$6058 is 1.6 kyr old, we note that escaping CRs are unlikely to produce significant TeV gamma-ray emission towards HESS J1420$-$607 (see blue region in Figure 12). As per HESS J1119$-$164, it is expected that CRs may still remain confined inside the SNR front shock. Similarly for the leptonic scenario (see Table 3), a fast diffusion coefficient ($\chi>0.1$) is required for electrons to reach the clumps ‘6’ to ‘10’ and the gas in region ‘A’. Kishishita et al. (2012) indicated that a slow diffusion of electrons inside the PWNe is required to produce the diffuse X-ray emission observed with Suzaku. Consequently, it is unlikely that high energy electrons escaping HESS J1418$-$609 can contribute to the TeV emission towards HESS J1420$-$607.

The various clumps north of HESS J1418$-$609 could affect its TeV emission. From Table 3, we note that, although some clumps could impede the propagation of multi-TeV electrons ($\tau_{\text{sync}}\ll\tau_{\text{diff}}$), a large diffusion coefficient would also be necessary for the high energy electrons to reach these clumps. Consequently, at this distance, the effect of the ISM on the morphology of HESS J1418$-$609 is also not significant.

HESS J1303$-$631

HESS J1303$-$631 has clearly been identified as a PWN based on energy dependent morphology at TeV energies (Abramowski et al. 2012b). We here check the potential contribution of the molecular regions ‘A’ to ‘D’ to the observed TeV gamma-ray emission.

From Table 3, we find that, assuming a suppression factor $\chi=0.1$, high energy electrons would only reach the molecular regions ‘A’ (located at $d\sim 6$ kpc) and ‘C’,‘D’ (both located at $d\sim 13.6$ kpc). In all cases, the shorter synchrotron time-scale indicates that any high energy electrons would lose most of their energy while diffusing into these molecular regions. Consequently, if associated with the PWN, less TeV emission spatially overlapping with these molecular regions should be expected.

Now, we discuss whether the selected molecular regions can also affect the morphology of the TeV emission from hadronic interactions. As per the case of HESS J1809$-$193, we scaled down the photon flux of HESS J1303$-$631 to 44%, 7% and 7% for regions ‘A’, ‘C’ and ‘D’ respectively, based on the size ratio between the molecular regions and the TeV emission. We then obtain $k_{\text{CR}}\sim 200$ (see Table 2). Based on the distances between the molecular clouds to the pulsar or SNR candidate (see Table 3), we find from Figures 12 and 13 that neither the pulsar (see purple line) nor the SNR candidate can provide sufficient CR energy density to produce hadronic TeV emission.

HESS J1018$-$589

We now discuss how the ISM in region ‘A’ may affect the TeV emission towards HESS J1018$-$589a. We find that the enhancement factor $k_{\text{CR}}=91$ (see Table 2) is required for the cosmic-ray to produce observed TeV flux. A 10 kyr old SNR is capable of producing such cosmic-ray enhancement if the distance between the SNR and the cloud is $\sim 10$ pc (see red region in Figure 12). A lack of TeV emission coincident with the region ‘A’ may suggest that the molecular cloud may be foreground/background to the SNR and thus may not be physically connected. Here we discuss only the possible contribution from CRs originating from the SNR G284.3$-$1.8. From our mass estimate in region ‘A’ (see Table C.9), we find that a slow diffusion of CRs escaping the SNR might lead to a significant contribution of the hadronic TeV emission towards HESS J1018$-$583A.

5 Conclusion

In this paper, we have mapped the molecular and atomic ISM towards pulsar wind nebulae (PWNe) and PWNe candidates, combining our new 7mm Mopra survey with the Nanten CO(1–0) data and the SGPS/GASS Hi survey. Except for the cases of HESS J1303$-$631 and HESS J1018$-$589B, our ISM studies provided additional information about the distance of astrophysical sources potentially related to the TeV source. We have also identified molecular clouds which could explain the morphology of the TeV gamma-ray emission as observed by CTA (Acharya et al. 2013), and also highlight cosmic-rays (CRs) originating from progenitor supernova remnants (SNRs). From the hypothesis that CRs could also be produced within the pulsar environment, and that at most 20% of the spin down power could be transferred to CRs (Bucciantini et al. 2011), we find that, among the studied PWNe, only CRs from PSR J1809$-$1917 could contribute to the TeV emission.

Towards HESS J1809$-$193, we have found several dense molecular regions at $v_{\text{lsr}}=10\text{ to }22$ km/s and $v_{\text{lsr}}=25\text{ to }38$ km/s adjacent to the pulsar PSR J1809$-$1918 and the SNRs G011.0$-$0.0 and G011.1$+$0.1. Notably, we have detected SiO(1–0) emission towards the dense molecular cloud south of HESS J1809$-$193 with no infra-red counterparts, suggesting a possible SNR-MC interaction at the pulsar dispersion measure distance $d\sim 3.7$ kpc. We found that the ambient density required to reconcile the projected radius of SNR G011.0$-$0.0 with the characteristic age of the pulsar PSR J1809$-$1917 is consistent with the averaged density $n_{\text{H}_{2}}$ obtained from our CO analysis. We argue that the PWN electrons could lose most of their energy before propagating deep inside the nearby molecular clouds, and would therefore become unable to produce coincident TeV gamma-ray emission. We also note that the CRs that have escaped the pulsar’s progenitor SNR may produce significant TeV emission towards these molecular clouds.

Towards HESS J1026$-$582, the molecular clouds east of the pulsar PSR J1028$-$5819 may explain its offset position with respect to the peak of the TeV emission at $d\sim 2.7$ kpc. However, we have also found a partial shell structure spatially coincident with the TeV emission at $d\sim 5$ kpc and highlight a potential hadronic origin, powered by an unknown CR source at $d<10$ pc from the molecular cloud.

Towards HESS J1119$-$164, the molecular clouds at the kinematic distances $d\sim 8.6\text{ to }9.7$ kpc found on the west side of SNR 292.2$-$0.5 show good correspondence with the thermal X-ray detections with Chandra and XMM-Newton, suggesting a possible SNR-MC interaction. Our comparative study of our column density estimates with the column density derived from X-ray measurements also supports the TeV source distance $d>8.6$ kpc.

Combining our column density study as well as our ISM morphology study towards HESS J1418$-$609, we argue the TeV source may be located at $d\sim 3.5-5.6$ kpc. We also claim that the various molecular clumps are too far from the pulsar PSR J1418$-$6058 to affect or contribute to the TeV gamma-ray emission.

Our ISM analysis has not constrained the distance of HESS J1303$-$631 as the CO(1–0) emission overlaps the TeV emission at all velocity ranges. Notably, at $v_{\text{lsr}}=-22\text{ to }-15$ km/s, however, we did find a CO dip spatially coincident with the SNR candidate found by Sushch et al. (2017), inferring a distance $d\sim 1.5$ kpc.

Although CO(1–0) emission north of the SNR G284.3$-$1.8 towards HESS J1018$-$589 has been found at $d\sim 2.8$ kpc, no CS clumps have been found at the pulsar PSR J1019$-$5857 dispersion measure distance $d\sim 8.6$ kpc, suggesting an atomic-dominated ISM surrounding HESS J1018$-$589B. An extension of high resolution Hi survey would shed more light concerning this TeV source.

As a conclusion, the arc-minute structure of the ISM plays an important role in understanding the morphology of the TeV gamma-ray emission produced by PWNe and progenitor SNR. The angular resolution and sensitivity of the Mopra CO(1–0) and 13CO surveys (Braiding et al. 2015) will refine the structure and dynamics of the diffuse molecular gas.

As future work, we aim to use the various properties of the identified clumps; to model the effect of the diffusion of CRs, and high energy electrons, escaping the PWNe and their progenitor SNRs; ; to model the interaction of these particles with the identified molecular cloud clumps; and to identify spectral and/or morphological signatures of the TeV gamma-ray emission. These could then be compared to upcoming observations with CTA (Acharya et al. 2013), which could further constrain the nature of the TeV sources.

Acknowledgments

The Mopra radio telescope is part of the Australia Telescope National Facility which is funded by the Australian Government for operation as a National Facility managed by CSIRO. Operations support was provided by the University of New South Wales and the University of Adelaide. The University of New South Wales Digital Filter Bank used for the observations with the Mopra Telescope (the UNSW–MOPS) was provided with support from the Australian Research Council LE160.100094 (ARC).

Appendix A Cleaning and integrated intensity mapping methods

A flow chart which outlines the various processes used to clean the various data-cubes and produce improved integrated intensity maps is shown in Figure A.1. The procedure to produce clean data-cubes is explained in Section 2.

A.1 Integrated intensity maps

From the unsmoothed cleaned data-cubes, we first apply a Gaussian smooth with FWHM$\sim 1^{\prime}$ using the Miriad task Convol. Then, from our cleaned data-cubes, we choose twelve background regions, each consisting of nine pixels, across the map, with no CS(1–0) emission. From the spectra obtained from these background regions, a binned distribution of the number of channels $N_{\text{ch}}\left(T_{\text{A}}^{*}\right)$ within the antenna temperature $T_{\text{A}}^{*}$ and $T_{\text{A}}^{*}+\Delta T_{\text{A}}^{*}$ is then plotted ($\Delta T_{\text{A}}^{*}=0.01$ K in this work). We then record the median $N_{\text{ch}}$ value for each bin to obtain our background distribution (see purple histogram in Figure A.2). Note that we also record the third quartile of the background distribution (see error bars in Figure A.2) to look at the fluctuation of the background distribution.

Now that we have estimated the background level distribution at the studied velocity range, the integrated intensity of each pixel is derived using the following method:

Obtain the binned distribution of the region, consisting of 9 pixels, centred on the current pixel (see blue histogram in Figure A.2). 2. 2.

Determine the threshold antenna temperature $T_{\text{A}}^{s}$ so that channels with antenna temperature $|T_{\text{A}}^{*}|<T_{\text{A}}^{s}$ are automatically discarded (see grey shaded area in Figure A.2). We here define $T_{\text{A}}^{s}$ as the lowest antenna temperature of the ’ON’ whose $N_{\text{ch}}\left(T_{\text{A}}^{*}\right)$ and $N_{\text{ch}}\left(T_{\text{A}}^{*}+\Delta T_{\text{A}}^{*}\right)$ exceed the third quartile of the background distribution (to account for fluctuation of the background distribution). It should be noted that the antenna temperature $T_{\text{A}}^{*}<-T_{\text{A}}^{s}$ are not discarded in order to mitigate potential excess from a region with increased $T_{\text{rms}}$. 3. 3.

For each bin, we subtract the background distribution from $N_{\text{ch}}^{\text{ON}}$ to obtain the excess number of channels $N_{\text{excess}}\left(T_{\text{A}}^{*}\right)$. If the value is negative, then we define $N_{\text{excess}}\left(T_{\text{A}}^{*}\right)=0$. 4. 4.

We multiply $N_{\text{excess}}\left(T_{\text{A}}^{*}\right)$ with $T_{\text{A}}^{*}$ and the velocity channels spacing ($\Delta v_{\text{lsr}}\sim 0.2$ km/s for 7mm Mopra observation) to obtain the excess integrated intensity $W_{\text{excess}}\left(T_{\text{A}}^{*}\right)$ at each bin. 5. 5.

We finally derive the integrated intensity at a given pixel, by first summing $W_{\text{excess}}\left(T_{\text{A}}^{*}\right)$ and then dividing by the number of pixels.

To account for possible spatial fluctuations, we then smooth our integrated intensity map with a Gaussian with FWHM$\sim 45^{\prime\prime}$, resulting in a total smoothing of $1.25^{\prime}$ (as per the data cubes).

Now that the integrated intensity maps are made, we need to establish the background level for each maps. Indeed, our new integrated intensity maps are not devoid of noise despite our aforementioned cuts. To do so, we simulate data consisting of noise with the same $T_{\text{rms}}$ as per our data-cubes. Using the same aforementioned steps, we produce integrated intensity maps within the same velocity ranges. We then define our noise level for our integrated intensity maps as the 99% distribution value of the ‘noise maps’ distribution.

A.2 comparison of methods

Figure A.3 compares the CS(1–0) integrated intensity maps using uncleaned data-cubes with the standard intensity maps after removing baseline ripples, and with the method used in this work, to obtain integrated intensity maps. We remark the striking difference between the uncleaned and cleaned data. We finally notice the removal of additional noise using our method.

Appendix B Galactic model

We used the Galactic model from Vallée (2013) (see Figure B.1) to help identify the Galactic arm in which each molecular cloud is located.

Appendix C Fit parameters

Tables C.4 to C.9 show the (RA, Dec) position, the semi-major axes and the Gaussian fitting parameters of the CS(1–0) , C34S(1–0), HC3N(5–4,F=4–3), SiO(1–0, v=0) and CO(1–0) detections found towards the studied TeV sources.

Appendix D Physical parameters

D.1 CS(1-0)

The averaged optical depth $\tau_{\text{CS(1--0)}}$ can be derived using the integrated intensity ratio $W_{\text{CS(1--0)}}$ and $W_{\text{C}^{34}\text{S(1--0)}}$ via :

[TABLE]

where $\alpha$=[32S]/[34S]=22 represents the abundance ratio between the two isotopologues based on terrestrial measurements (Fink 1981). We obtain the column density of the upper state $N_{\text{CS}_{1}}$:

[TABLE]

Assuming the gas to be in local thermal equilibrium at temperature $T_{\text{kin}}$, we can thus obtain the total CS column density $N_{\text{CS}}$:

[TABLE]

D.2 Mass and density

[TABLE]

with $\mu=2.8$ being the averaged atomic weight (accounting for 20% Helium), $m_{\text{H}}$ the hydrogen mass and $a$,$b$ the semi-minor and semi-major axis of the selected elliptic regions.

[TABLE]

D.3 Physical parameters of individual sources

Tables D.1 to D.6 show the column density $N_{\text{H}_{2}}$ the atomic and molecular mass $M_{\text{H}_{\textsc{i}}}$ and $M_{\text{H}_{2}}$ and, molecular density $n_{\text{H}_{2}}$ estimates based on the CS, CO, and H${}_{\textsc{i}}$ emission found towards the regions listed in Section C.

Appendix E Column densities and distance studies

Appendix F HESS J1809$-$193 additional figures