A search for accreting young companions embedded in circumstellar disks: High-contrast H$\alpha$ imaging with VLT/SPHERE

TL;DR

This study utilized high-contrast Hα imaging with VLT/SPHERE to detect and analyze accreting young companions in circumstellar disks, successfully re-detecting known companions and setting limits on accretion activity around several young stars.

Contribution

It provides the first high-sensitivity Hα observations of known and candidate young companions, improving detection limits and estimating accretion rates in circumstellar environments.

Findings

Re-detected HD142527 B with highest S/N to date.

Estimated accretion rate of HD142527 B as 1-2×10^{-10} M_sun/yr.

Set upper limits on accretion fluxes around other young stars.

Abstract

Aims: We want to detect and quantify observables related to accretion processes occurring locally in circumstellar disks, which could be attributed to young forming planets. We focus on objects known to host protoplanet candidates and/or disk structures thought to be the result of interactions with planets. Methods: We analyzed observations of 6 young stars (age $3.5-10$ Myr) and their surrounding environments with the SPHERE/ZIMPOL instrument on the VLT in the H$\alpha$ filter (656 nm) and a nearby continuum filter (644.9 nm). Results: We re-detect the known accreting M-star companion HD142527 B with the highest published signal to noise to date in both H$\alpha$ and the continuum. We derive new astrometry ($r = 62.8^{+2.1}_{-2.7}$ mas and $\text{PA} = (98.7\,\pm1.8)^\circ$) and photometry ($\Delta$N_Ha=$6.3^{+0.2}_{-0.3}$ mag, $\Delta$B_Ha=$6.7\pm0.2$ mag and $\Delta$Cnt_Ha=$7.3^{+0.3}_{-0.2}$ mag) for the companion in agreement with previous studies, and estimate its mass accretion rate ($\dot{M}\approx1-2\,\times10^{-10}\,M_\odot\text{ yr}^{-1}$). A faint point-like source around HD135344 B (SAO206462) is also investigated, but a second deeper observation is required to reveal its nature. No other companions are detected. In the framework of our assumptions we estimate detection limits at the locations of companion candidates around HD100546, HD169142 and MWC758 and calculate that processes involving H$\alpha$ fluxes larger than $\sim8\times10^{-14}-10^{-15}\,\text{erg/s/cm}^2$ ($\dot{M}>10^{-10}-10^{-12}\,M_\odot\text{ yr}^{-1}$) can be excluded. Furthermore, flux upper limits of $\sim10^{-14}-10^{-15}\,\text{erg/s/cm}^2$ ($\dot{M}<10^{-11}-10^{-12}\,M_\odot \text{ yr}^{-1}$) are estimated within the gaps identified in the disks surrounding HD135344B and TW Hya.

Tables

Table 1. Table 1: Summary of observations.

Object	H$\alpha$	Obs. date	Prog. ID	DITb𝑏bb𝑏b$b$DIT = Detector integration time, i.e., exposure time per image frame.	# of	Field	Mean	${\tau }_0$c𝑐cc𝑐c$c$Coherence time	Mean
	Filtera𝑎aa𝑎a$a$Each dataset consists of data obtained in one of the two H$\alpha$ filters and simultaneous data taken with the continuum filter inserted in the other ZIMPOL camera.	[dd.mm.yyyy]		[s]	DITs	rotation [∘]	airmass	[ms]	seeingd𝑑dd𝑑d$d$Mean DIMM seeing measured during the observation. [as]
HD142527	B_Ha	31.03.2016	096.C-0248(B)	30	70	47.8	1.06	$2.7\pm 0.2$	$0.71\pm 0.06$
	N_Ha	31.03.2016	096.C-0248(B)	30	70	48.6	1.05	$2.7\pm 0.3$	$0.69\pm 0.07$
HD135344 B	N_Ha	31.03.2016	096.C-0248(B)	50	107	71.7	1.04	$4.4\pm 1.2$	$0.47\pm 0.17$
TW Hya	B_Ha	23.03.2016	096.C-0267(B)	80	131	134.1	1.16	$1.4\pm 0.4$	$1.33\pm 0.53$
HD100546	B_Ha	23.04.2015	095.C-0273(A)	10	1104e𝑒ee𝑒e$e$As we explain in Section 4.4 and Appendix E, for this dataset a frame selection was applied, which reduced the number of frames to 366 and the field rotation to ${20.7}^{\circ }$.	68.3e𝑒ee𝑒e$e$As we explain in Section 4.4 and Appendix E, for this dataset a frame selection was applied, which reduced the number of frames to 366 and the field rotation to ${20.7}^{\circ }$.	1.46	$1.7\pm 0.2$	$0.98\pm 0.28$
HD169142	B_Ha	09.05.2015	095.C-0298(A)	50	90	123.2	1.01	$1.4\pm 0.1$	$1.24\pm 0.04$
MWC 758	B_Ha	30.12.2015	096.C-0267(A)	60	194	54.8	1.63	$3.2\pm 0.8$	$1.39\pm 0.24$

Table 2. Table 2: Target sample.

Object	RA	DEC	Spec. type	$m_R$ [mag]	Distance [pc]	Age [Myr]
HD142527	15h56m41.89s	-42∘19′23$\stackrel{\prime \prime }{.}$27	F6III	7.91	$157.3\pm 1.2$	${8.1}_{-1.6}^{+1.9}$
HD135344 B	15h15m48.44s	-37∘09′16$\stackrel{\prime \prime }{.}$03	F8V	8.45	$135.9\pm 1.4$	$9\pm 2$
TW Hya	11h01m51.90s	-34∘42′17$\stackrel{\prime \prime }{.}$03	K6Ve	$10.43\pm 0.1$	$60.1\pm 0.1$	$\sim 10$
HD100546	11h33m25.44s	-70∘11′41$\stackrel{\prime \prime }{.}$24	B9Vne	8.78	$110.0\pm 0.6$	$7\pm 1.5$
HD169142	18h24m29.78s	-29∘46′49$\stackrel{\prime \prime }{.}$32	B9V	8.0	$114.0\pm 0.8$	$\sim 6$
MWC 758	05h30m27.53s	-25∘19′57$\stackrel{\prime \prime }{.}$08	A8Ve	$9.20\pm 0.01$	$160.3\pm 1.7$	$3.5\pm 2$

Table 3. Table 3: Summary of the stellar fluxes measured in the different filters in our ZIMPOL data and the derived H α 𝛼 \alpha line fluxes for our targets (last column). The extinction values A H ​ α subscript 𝐴 𝐻 𝛼 A_{H\alpha} were estimated as described in Section 4.1.4 from A V subscript 𝐴 𝑉 A_{V} .

Object	$A_V$ [mag]	$A_{H\alpha }$ [mag]	$F_{\text{F\_H}\alpha }^{\ast }$ [erg/s/cm2]	$F_{\text{Cnt\_H}\alpha }^{\ast }$ [erg/s/cm2]	$F_{H\alpha }^{\ast }$ [erg/s/cm2]
HD142527 (N_Ha)	¡0.05a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.04	$3.0\pm 0.8\times {10}^{-11}$	$6.1\pm 0.2\times {10}^{-11}$	$1.7\pm 0.8\times {10}^{-11}$
HD142527 (B_Ha)	¡0.05a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.04	$9.7\pm 0.8\times {10}^{-11}$	$6.1\pm 0.2\times {10}^{-11}$	$1.0\pm 0.5\times {10}^{-11}$
HD142527 B (N_Ha)	¡0.05a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.04	${9.1}_{-2.9}^{+3.5}\times {10}^{-14}$	${7.4}_{-2.1}^{+1.4}\times {10}^{-14}$	${7.6}_{-2.9}^{+3.5}\times {10}^{-14}$
HD142527 B (B_Ha)	¡0.05a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.04	$2.0\pm 0.4\times {10}^{-13}$	${7.4}_{-2.1}^{+1.4}\times {10}^{-14}$	${1.0}_{-0.4}^{+0.5}\times {10}^{-13}$
HD135344 B	0.23 a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.19	$3.1\pm 1.0\times {10}^{-11}$	$4.9\pm 0.6\times {10}^{-11}$	$1.8\pm 0.8\times {10}^{-11}$
TW Hya	0.0b𝑏bb𝑏b$b$Uyama et al. (2017).	0.0	$9.9\pm 0.4\times {10}^{-11}$	$1.5\pm 0.05\times {10}^{-11}$	$7.8\pm 0.3\times {10}^{-11}$
HD100546	¡0.05a𝑎aa𝑎a$a$Fairlamb et al. (2015).	0.04	$4.2\pm 0.2\times {10}^{-10}$	$1.6\pm 0.1\times {10}^{-10}$	$1.7\pm 0.2\times {10}^{-10}$
HD169142	0.43c𝑐cc𝑐c$c$Fedele et al. (2017).	0.35	$1.1\pm 0.1\times {10}^{-10}$	$7.4\pm 0.2\times {10}^{-11}$	$3.2\pm 4.4\times {10}^{-12}$
MWC758	0.22d𝑑dd𝑑d$d$van den Ancker et al. (1998).	0.18	$8.1\pm 0.7\times {10}^{-11}$	$5.3\pm 0.2\times {10}^{-11}$	$6.3\pm 3.7\times {10}^{-12}$

Table 4. Table 4: Summary of our detection limits for each target. While for HD100546, HD169142, and MWC 758 we consider the specific locations (separation and PA) of previously claimed companion candidates, we focused our analyses for HD135344B and TW Hya on separations related to disk gaps (hence no specific PA). Columns 5 and 6 give the mass and radius assumed for the accretion rate calculations, column 7 gives the contrast magnitude at the specific location and columns 8–11 report the values for the H α 𝛼 \alpha line flux, H α 𝛼 \alpha line luminosity, accretion luminosity, and mass accretion rate ignoring any possible dust around the companion.

Target	Sep.	PA	Ref.	Mass	Radius	$\Delta$Hα	$F_{H\alpha }^p$	$L_{H\alpha }$ [$L_{\odot }$]	$L_{\mathrm{acc}}$ [$L_{\odot }$]	$\dot{M}$ [$M_{\odot }{\text{ yr}}^{-1}$]
	[mas]	[∘]		[$M_J$]	[$R_J$]	[mag]	[erg/s/cm2]
HD135344B	$180$		(a)	10.2(h)ℎ(h)(h)ℎ(h)footnotemark: $\left(h\right)$	1.6(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>9.8$
TW Hya	$390$		(b)	2(k)𝑘(k)(k)𝑘(k)footnotemark: $\left(k\right)$	1.3(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>9.3$
HD100546	$480\pm 4$	$8.9\pm 0.9$	(c)	15(c)𝑐(c)(c)𝑐(c)footnotemark: $\left(c\right)$	2(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>11.4$
	$\sim 140$	$\sim 133$	(d)	15(l)𝑙(l)(l)𝑙(l)footnotemark: $\left(l\right)$	2(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>9.3$
HD169142	$\sim 340$	$\sim 175$	(e)	0.6(e)𝑒(e)(e)𝑒(e)footnotemark: $\left(e\right)$	1.4(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>10.7$
	$156\pm 32$	$7.4\pm 11.3$	(f)	10(f)𝑓(f)(f)𝑓(f)footnotemark: $\left(f\right)$	1.7(j)𝑗(j)(j)𝑗(j)footnotemark: $\left(j\right)$	$>9.9$
MWC 758	$111\pm 4$	$162\pm 5$	(g)	5.5(m)𝑚(m)(m)𝑚(m)footnotemark: $\left(m\right)$	1.7(n)𝑛(n)(n)𝑛(n)footnotemark: $\left(n\right)$	$>9.4$

Table 5. Table 5: Resulting signal flux and FPF for different beamsplitter behaviors.

Target	Signal flux	FPF
	($\times 1000$, arb. unit)	($\times {10}^{12}$)
Cnt flux 5% decreased	$5.45\pm 0.06$	$9.97\pm 4.72$
Cnt flux not changed	$3.62\pm 0.02$	$1.69\pm 1.47$
Cnt flux 5% enhanced	$2.84\pm 0.03$	$11.17\pm 2.86$

Equations

$$F_{\lambda }^{\ast }(Cnt\_Ha)=cts\cdot 10^{0.4(am\cdot k_1+m_{mode})}\cdot c_{zp}^{cont}(Cnt\_Ha),$$

$$cts(F\_Ha)=\frac{F_{\lambda }^{\ast }(Cnt\_Ha)}{c_{zp}^{cont}(F\_Ha)}\times 10^{-0.4(am\cdot k_1+m_{mode})},$$

$$F_{Cnt\_Ha}^p=7.4_{-2.1}^{+1.4}\times 10^{-14}{\text{ erg/s/cm}}^2,$$

$$F_{B\_Ha}^p=2.0\pm 0.4\times 10^{-13}{\text{ erg/s/cm}}^2,$$

$$F_{N\_Ha}^p=9.1_{-2.9}^{+3.5}\times 10^{-14}{\text{ erg/s/cm}}^2.$$

$$\log (L_{\mathrm{acc}})=b+a\log (L_{H\alpha }),$$

$${\dot{M}}_{\mathrm{acc}}={\left(1-\frac{R_c}{R_{in}}\right)}^{-1}\frac{L_{\mathrm{acc}}{R}_c}{G{M}_c}\sim 1.25\frac{L_{\mathrm{acc}}{R}_c}{G{M}_c}$$

$$L_{H\alpha }=4.7\times 10^{-6}{L}_{\odot }{\left(\frac{R_T}{R_J}\right)}^2\left(\frac{v_s}{59{\text{km s}}^{-1}}\right).$$

Full text

Methods.

Results.

*Aims. *We seek to detect and quantify observables related to accretion processes occurring locally in circumstellar disks, which could be attributed to young forming planets. We focus on objects known to host protoplanet candidates and/or disk structures thought to be the result of interactions with planets.

*Methods. *We analyzed observations of six young stars (age $3.5-10$ Myr) and their surrounding environments with the SPHERE/ZIMPOL instrument on the Very Large Telescope (VLT) in the H$\alpha$ filter (656 nm) and a nearby continuum filter (644.9 nm). We applied several point spread function (PSF) subtraction techniques to reach the highest possible contrast near the primary star, specifically investigating regions where forming companions were claimed or have been suggested based on observed disk morphology.

*Results. *We redetect the known accreting M-star companion HD142527 B with the highest published signal to noise to date in both H$\alpha$ and the continuum. We derive new astrometry ($r=62.8^{+2.1}_{-2.7}$ mas and $\text{PA}=(98.7\,\pm 1.8)^{\circ}$) and photometry ($\Delta$N_Ha=$6.3^{+0.2}_{-0.3}$ mag, $\Delta$B_Ha=$6.7\pm 0.2$ mag and $\Delta$Cnt_Ha=$7.3^{+0.3}_{-0.2}$ mag) for the companion in agreement with previous studies, and estimate its mass accretion rate ($\dot{M}\approx 1-2\,\times 10^{-10}\,M_{\odot}\text{ yr}^{-1}$). A faint point-like source around HD135344 B (SAO206462) is also investigated, but a second deeper observation is required to reveal its nature. No other companions are detected. In the framework of our assumptions we estimate detection limits at the locations of companion candidates around HD100546, HD169142, and MWC 758 and calculate that processes involving H$\alpha$ fluxes larger than $\sim 8\times 10^{-14}-10^{-15}\,\text{erg/s/cm}^{2}$ ($\dot{M}>10^{-10}-10^{-12}\,M_{\odot}\text{ yr}^{-1}$) can be excluded. Furthermore, flux upper limits of $\sim 10^{-14}-10^{-15}\,\text{erg/s/cm}^{2}$ ($\dot{M}<10^{-11}-10^{-12}\,M_{\odot}\text{ yr}^{-1}$) are estimated within the gaps identified in the disks surrounding HD135344 B and TW Hya. The derived luminosity limits exclude H$\alpha$ signatures at levels similar to those previously detected for the accreting planet candidate LkCa15 b.

1 Introduction

Providing an empirical basis for gas giant planet formation models and theories requires the detection of young objects in their natal environment, i.e., when they are still embedded in the gas and dust-rich circumstellar disk surrounding their host star. The primary scientific goals of studying planet formation are as follows: To understand where gas giant planet formation takes place, for example, at what separations from the host star and under which physical and chemical conditions in the disk; how formation occurs, i.e., via the classical core accretion process (Pollack et al., 1996) or a modified version of that process (e.g., pebble accretion, Lambrechts & Johansen, 2012) or direct gravitational collapse (Boss, 1997)); and the properties of the suspected circumplanetary disks (CPDs).

While in recent years high-contrast, high spatial resolution imaging observations of circumstellar disks have revealed an impressive diversity in circumstellar disk structure and morphology, the number of directly detected planet candidates embedded in those disks is still small (LkCa15 b, HD100546 b, HD169142 b, MWC 758 b, PDS 70 b; Kraus & Ireland, 2012; Quanz et al., 2013a; Reggiani et al., 2014; Biller et al., 2014; Reggiani et al., 2018; Keppler et al., 2018). To identify these objects, high-contrast exoplanet imaging can be used. These observations are typically performed at near- to mid-infrared wavelengths using an adaptive optics-assisted high-resolution camera. In addition to the intrinsic luminosity of the still contracting young gas giant planet, the surrounding CPD, if treated as a classical accretion disk, contributes significantly to fluxes beyond 3$\,\mu$m wavelength (Zhu, 2015; Eisner, 2015), potentially easing the detection of young forming gas giants at these wavelengths. While the majority of the forming planet candidates mentioned above were detected in this way, it has also been realized that the signature from a circumstellar disk itself can sometimes mimic that of a point source after PSF subtraction and image post-processing (e.g., Follette et al., 2017; Ligi et al., 2018). As a consequence, it is possible that some of the aforementioned candidates are false positives.

Another approach is to look for direct signatures of the suspected CPDs, such as their dust continuum emission or their kinematic imprint in high-resolution molecular line data (Perez et al., 2015; Szulágyi et al., 2018). In one case, spectro-astrometry using CO line emission was used to constrain the existence and orbit of a young planet candidate (Brittain et al., 2013, 2014). Moreover, Pinte et al. (2018) and Teague et al. (2018) suggested the presence of embedded planets orbiting HD163296 from local deviations from Keplerian rotation in the protoplanetary disk. A further indirect way to infer the existence of a young, forming planet is to search for localized differences in the gas chemistry of the circumstellar disk, as the planet provides extra energy to the chemical network in its vicinity (Cleeves et al., 2015).

Finally, it is possible to look for accretion signatures from gas falling onto the planet and its CPD. Accretion shocks are able to excite or ionize the hydrogen atoms, which then radiate recombination emission lines, such as H$\alpha$, when returning to lower energy states (e.g., Calvet & Gullbring, 1998; Szulágyi & Mordasini, 2017; Marleau et al., 2017). High-contrast imaging using H$\alpha$ filters was already successfully applied in three cases. Using angular spectral differential imaging (ASDI) with the Magellan Adaptive Optics System (MagAO), Close et al. (2014a) detected H$\alpha$ excess emission from the M-star companion orbiting the Herbig Ae/Be star HD142527, and Sallum et al. (2015) also used MagAO to identify at least one accreting companion candidate located in the gap of the transition disk around LkCa15. The accretion signature was found at a position very similar to the predicted orbital position of one of the faint point sources detected by Kraus & Ireland (2012), attributed to a forming planetary system. Most recently, Wagner et al. (2018) have claimed the detection of H$\alpha$ emission from the young planet PDS70 b using MagAO, albeit with comparatively low statistical significance (3.9$\sigma$).

In this paper we present a set of H$\alpha$ high-contrast imaging data for six young stars, aiming at the detection of potential accretion signatures from the (suspected) young planets embedded in the circumstellar disks of the stars. The paper is structured as follows: In Section 2 we discuss the observations and target stars. We explain the data reduction in Section 3 and present our analyses in Section 4. In Section 5 we discuss our results in a broader context and conclude in Section 6.

2 Observations and target sample

2.1 Observations

The data were all obtained with the ZIMPOL sub-instrument of the adaptive optics (AO) assisted high-contrast imager SPHERE (Beuzit et al., 2008; Petit et al., 2008; Fusco et al., 2016), which is installed at the Very Large Telescope (VLT) of the European Southern Observatory (ESO) on Paranal in Chile. A detailed description of ZIMPOL can be found in Schmid et al. (2018). Some of the data were collected within the context of the Guaranteed Time Observations (GTO) program of the SPHERE consortium; others were obtained in other programs and downloaded from the ESO data archive (program IDs are listed in Table 1). We focused on objects that are known from other observations to host forming planet candidates that still need to be confirmed (HD100546, HD169142, and MWC 758)222In the discussion (Section 5) we also include the analysis of a dataset of LkCa15 (PI: Huelamo) to set our results in context, but the data were poor in quality and hence not included in the main part of the paper., objects known to host accreting stellar companions (HD142527), and objects that have well-studied circumstellar disks with spatially resolved substructures (gaps, cavities, or spiral arms), possibly suggesting planet formation activities (HD135344 B and TW Hya). All data were taken in the noncoronagraphic imaging mode of ZIMPOL using an H$\alpha$ filter in one camera arm and a nearby continuum filter simultaneously in the other arm (Cont_Ha; $\lambda_{c}=644.9$ nm, $\Delta\lambda=3.83$ nm). As the data were observed in different programs, we sometimes used the narrow H$\alpha$ filter (N_Ha; $\lambda_{c}=656.53$ nm, $\Delta\lambda=0.75$ nm) and sometimes the broad H$\alpha$ filter (B_Ha; $\lambda_{c}=655.6$ nm, $\Delta\lambda=5.35$ nm). A more complete description of these filters can be found in Schmid et al. (2017). To establish which filter allows for the highest contrast performance, we used HD142527 and its accreting companion (Close et al., 2014a) as a test target and switched between the N_Ha and the B_Ha filter every ten frames within the same observing sequence. All datasets were observed in pupil-stabilized mode to enable angular differential imaging (ADI; Marois et al., 2006). The fundamental properties of the target stars are given in Table 2, while a summary of the datasets is given in Table 1.

We note that because of the intrinsic properties of the polarization beam splitter used by ZIMPOL, polarized light might preferentially end up in one of the two arms, causing a systematic uncertainty in the relative photometry between the continuum and H$\alpha$ frames. The inclined mirrors in the telescope and the instrument introduce di-attenuation (e.g., higher reflectivity for $I_{\perp}$ than $I_{\parallel}$) and polarization cross talks, so that the transmissions in imaging mode to the $I_{\perp}$ and $I_{\parallel}$ arm depend on the telescope pointing direction. This effect is at the level of a few percent (about $\pm 5\leavevmode\nobreak\ \%$), but unfortunately the dependence on the instrument configuration has not been determined yet. We discuss its potential impact on our analyses in Appendix A, even though we did not take this effect into account since it is small and could not be precisely quantified.

2.2 Target sample

*HD142527

*HD142527 is known to have a prominent circumstellar disk (e.g., Fukagawa et al., 2006; Canovas et al., 2013; Avenhaus et al., 2014b) and a close-in M star companion (HD142527 B; Biller et al., 2012; Rodigas et al., 2014; Lacour et al., 2016; Christiaens et al., 2018; Claudi2018) that shows signatures of ongoing accretion in H$\alpha$ emission (Close et al., 2014a). This companion orbits in a large, optically thin cavity within the circumstellar disk stretching from $\sim 0\aas@@fstack{\prime\prime}07$ to $\sim 1\aas@@fstack{\prime\prime}0$ (e.g., Fukagawa et al., 2013; Avenhaus et al., 2014b), and it is likely that this companion is at least partially responsible for clearing the gap by accretion of disk material (Biller et al., 2012; Price et al., 2018). Avenhaus et al. (2017) obtained polarimetric differential imaging data with SPHERE/ZIMPOL in the very broad band (VBB, as defined in Schmid et al., 2018) optical filter, revealing new substructures, and resolving the innermost regions of the disk (down to $0\aas@@fstack{\prime\prime}025$). In addition, extended polarized emission was detected at the position of HD142527 B, possibly due to dust in a circumsecondary disk. Christiaens et al. (2018) extracted a medium-resolution spectrum of the companion and suggested a mass of $0.34\pm 0.06\,M_{\odot}$. This value is a factor of $\sim 3$ larger than that estimated by spectral energy distribution (SED) fitting (Lacour et al., 2016, $M=0.13\pm 0.03\,M_{\odot}$). Thanks to the accreting close-in companion, this system is the ideal target to optimize the H$\alpha$ observing strategy with SPHERE/ZIMPOL and also the data reduction.

*HD135344 B

HD135344 B (SAO206462) is surrounded by a transition disk that was spatially resolved at various wavelengths. Continuum (sub-)millimeter images presented by Andrews et al. (2011) and van der Marel et al. (2016) revealed a disk cavity with an outer radius of $0\aas@@fstack{\prime\prime}32$. In polarimetric differential imaging (PDI) observations in the near-infrared (NIR), the outer radius of the cavity appears to be at $0\aas@@fstack{\prime\prime}18$, and the difference in apparent size was interpreted as a potential indication for a companion orbiting in the cavity (Garufi et al., 2013). Data obtained in PDI mode also revealed two prominent, symmetric spiral arms (Muto et al., 2012; Garufi et al., 2013; Stolker et al., 2016). Vicente et al. (2011) and Maire et al. (2017) searched for planets in the system using NIR NACO and SPHERE high-contrast imaging data, but did not find any. Using hot start evolutionary models these authors derived upper limits for the mass of potential giant planets around HD135344 B (3 MJ* beyond $0\aas@@fstack{\prime\prime}7$).

*TW Hya

*TW Hya is the nearest T Tauri star to Earth. Its almost face-on transitional disk ($i\sim 7\pm 1^{\circ}$; Qi et al., 2004) shows multiple rings and gaps in both dust continuum and scattered light data. Hubble Space Telescope (HST) scattered light images from Debes et al. (2013) first allowed the identification of a gap at $\sim 1\aas@@fstack{\prime\prime}48$. Later, Akiyama et al. (2015) observed in $H$-band polarized images a gap at a separation of $\sim 0\aas@@fstack{\prime\prime}41$. Using Atacama Large Millimeter Array (ALMA), Andrews et al. (2016) identified gaps from the radial profile of the 870 $\mu$m continuum emission at $0\aas@@fstack{\prime\prime}41$, $0\aas@@fstack{\prime\prime}68$ and $0\aas@@fstack{\prime\prime}80$. Finally, van Boekel et al. (2017) obtained SPHERE images in PDI and ADI modes at optical and NIR wavelengths, and identified three gaps at $0\aas@@fstack{\prime\prime}11$, $0\aas@@fstack{\prime\prime}39,$ and $1\aas@@fstack{\prime\prime}57$ from the central star. A clear gap was also identified by Rapson et al. (2015) at a separation of $0\aas@@fstack{\prime\prime}43$ in Gemini/GPI polarimetric images and the largest gap at $r\simeq 1\aas@@fstack{\prime\prime}52$ has also been observed in CO emission with ALMA (Huang et al., 2018).

*HD100546

*The disk around HD100546 was also spatially resolved in scattered light and dust continuum emission in different bands (e.g., Augereau et al., 2001; Quanz et al., 2011; Avenhaus et al., 2014a; Walsh et al., 2014; Pineda et al., 2014). The disk appears to be almost, but not completely, devoid of dusty material at radii between a few and 13 AU. This gap could be due to the interaction with a young forming planet, and Brittain et al. (2013, 2014) suggested the presence of a companion orbiting the star at $0\aas@@fstack{\prime\prime}13$, based on high-resolution NIR spectro-astrometry of CO emission lines. Another protoplanet candidate was claimed by Quanz et al. (2013a) using $L^{\prime}$ band high-contrast imaging data. The object was found at $0\aas@@fstack{\prime\prime}48\pm 0\aas@@fstack{\prime\prime}04$ from the central star, at a position angle (PA) of $(8.9\,\pm 0.9)^{\circ}$, with an apparent magnitude of $L^{\prime}$=$13.2\pm 0.4$ mag. Quanz et al. (2015) reobserved HD100546 in different bands ($L^{\prime}$, $M^{\prime}$, $K_{s}$) and detected the object again in the first two filters. Based on the colors and observed morphology these authors suggested that the data are best explained by a forming planet surrounded by a circumplanetary disk. Later, Currie et al. (2015) recovered HD100546 b from $H$-band integral field spectroscopy (IFS) with the Gemini Planet Imager (GPI; Macintosh et al., 2006) and identified a second putative point source c closer to the star ($r_{\mathrm{proj}}\sim 0\aas@@fstack{\prime\prime}14$) potentially related to the candidate identified by Brittain et al. (2013, 2014). More recently, Rameau et al. (2017) demonstrated that the emission related to HD100546 b appears to be stationary and its spectrum is inconsistent with any type of low temperature objects. Furthermore, they obtained H$\alpha$ images with the MagAO instrument to search for accretion signatures, but no point source was detected at either the b or c position, and they placed upper limits on the accretion luminosity ($L_{\mathrm{acc}}<1.7\times 10^{-4}\;L_{\odot}$). The same data were analyzed by Follette et al. (2017), together with other H$\alpha$ images (MagAO), $H$ band spectra (GPI), and Y band polarimetric images (GPI). Their data exclude that HD100546 c is emitting in H$\alpha$ with $L_{H\alpha}>1.57\times 10^{-4}L_{\odot}$.

*HD169142

*HD169142 is surrounded by a nearly face-on pre-transitional disk. Using PDI images, Quanz et al. (2013b) found an unresolved disk rim at $0\aas@@fstack{\prime\prime}17$ and an annular gap between $0\aas@@fstack{\prime\prime}28$ and $0\aas@@fstack{\prime\prime}49$. These results were confirmed by Osorio et al. (2014), who investigated the thermal emission ($\lambda=7$ mm) of large dust grains in the HD169142 disk, identifying two annular cavities ($\sim 0\aas@@fstack{\prime\prime}16-0\aas@@fstack{\prime\prime}21$ and $\sim 0\aas@@fstack{\prime\prime}28-0\aas@@fstack{\prime\prime}48$). The latter authors also identified a point source candidate in the middle of the outer cavity at a distance of $0\aas@@fstack{\prime\prime}34$ and PA $\sim 175^{\circ}$. Biller et al. (2014) and Reggiani et al. (2014) observed a point-like feature in NaCo $L^{\prime}$ data at the outer edge of the inner cavity (separation = $0\aas@@fstack{\prime\prime}11-0\aas@@fstack{\prime\prime}16$ and PA=$0^{\circ}-7.4^{\circ}$). Observations in other bands ($H$, $K_{S}$, $z_{p}$) with the Magellan Clay Telescope (MagAO/MCT) and with GPI in the $J$ band failed to confirm the detection (Biller et al., 2014; Reggiani et al., 2014), but revealed another candidate point source albeit with low signal-to-noise ratio (S/N; Biller et al., 2014). In a recent paper, Ligi et al. (2018) explained the latter Biller et al. (2014) detection with a bright spot in the ring of scattered light from the disk rim, potentially following Keplerian motion. Pohl et al. (2017) and Bertrang et al. (2018) compared different disk and dust evolutionary models to SPHERE $J$-band and VBB PDI observations. Both works tried to reproduce and explain the complex morphological structures observed in the disk and conclude that planet-disk interaction is occurring in the system, even though there is no clearly confirmed protoplanet identified to date.

*MWC 758

*MWC 758 is surrounded by a pre-transitional disk (e.g., Grady et al., 2013). Andrews et al. (2011) found an inner cavity of $\sim$55 AU based on dust continuum observations, which was, however, not observed in scattered light (Grady et al., 2013; Benisty et al., 2015). Nevertheless, PDI and direct imaging from the latter studies revealed two large spiral arms. A third spiral arm has been suggested based on VLT/NaCo $L^{\prime}$ data by Reggiani et al. (2018), together with the claim of the detection of a point-like source embedded in the disk at $(111\pm 4)$ mas. This object was observed in two separate datasets from 2015 and 2016 at comparable separations from the star, but different PAs, which was possibly due to orbital motion. The contrast of this object relative to the central star in the $L^{\prime}$ band is $\sim 7$ mag, which, according to the BT-Settl atmospheric models (Allard et al., 2012), corresponds to the photospheric emission of a 41-64 $M_{J}$ object for the age of the star. More recently, ALMA observations from Boehler et al. (2018) traced the large dust continuum emission from the disk. Two rings at $0\aas@@fstack{\prime\prime}37$ and $0\aas@@fstack{\prime\prime}53$ were discovered that are probably related to two clumps with large surface density of millimeter dust and a large cavity of $\sim 0\aas@@fstack{\prime\prime}26$ in radius. Finally, Huélamo et al. (2018) observed MWC 758 in H$\alpha$ with SPHERE/ZIMPOL, reaching an upper limit for the line luminosity of $L_{\mathrm{H_{\alpha}}}\lesssim 5\times 10^{-5}L_{\odot}$ (corresponding to a contrast of 7.6 mag) at the separation of the protoplanet candidate. No other point-like features were detected.

3 Data reduction

The basic data reduction steps were carried out with the ZIMPOL pipeline developed and maintained at ETH Zürich. The pipeline remapped the original 7.2 mas / pixel $\times$ 3.6 mas / pixel onto a square grid with an effective pixel scale of 3.6 mas / pixel $\times$ 3.6 mas / pixel (1024 $\times$ 1024 pixels). Afterward, the bias was subtracted and a flat-field correction was applied. We then aligned the individual images by fitting a Moffat profile to the stellar point spread functions (PSFs) and shifting the images using bilinear interpolation. The pipeline also calculated the parallactic angle for each individual frame and added the information to the image header. Finally, we split up the image stacks into individual frames and grouped them together according to their filter, resulting in two image stacks for each object: one for an H$\alpha$ filter and one for the continuum filter444For HD142527 we have four image stacks as we used both the N_Ha and the B_Ha filter during the observing sequence.. In general, all images were included in the analysis if not specifically mentioned in the individual subsections. The images in these stacks were cropped to a size of $1\aas@@fstack{\prime\prime}08\,\times 1\,\aas@@fstack{\prime\prime}08$ centered on the star. This allowed us to focus our PSF subtraction efforts on the contrast dominated regime of the images. The removal of the stellar PSF was performed in three different ways: ADI, spectral differential imaging (SDI), and ASDI (a two-step combination of SDI and ADI).

To perform ADI, we fed the stacks into our PynPoint pipeline (Amara & Quanz, 2012; Amara et al., 2015; Stolker et al., 2018). The PynPoint package uses principal component analysis (PCA) to model and subtract the stellar PSF in all individual images before they are derotated to a common field orientation and mean-combined. To investigate the impact on the final contrast performance for all objects, we varied the number of principal components (PCs) used to fit the stellar PSF and the size of the inner mask that is used to cover the central core of the stellar PSF prior to the PCA. No frame selection based on the field rotation was applied, meaning that all the images were considered for the analysis, regardless of the difference in parallactic angle. The SDI approach aims at reducing the stellar PSF using the fact that all features arising from the parent star (such as Airy pattern and speckles) scale spatially with wavelength $\lambda$, while the position of a physical object on the detector is independent of $\lambda$. The underlying assumption is that, given that $\lambda_{c}$ is similar in all filters, the continuum flux density is the same at all wavelengths. To this end, modified versions of the continuum images were created. First, they were multiplied with the ratio of the effective filter widths to normalize the throughput of the continuum filter relative to the H$\alpha$ filter555This approach ignores any potential color effects between the filters, which, given their narrow band widths, should, however, not cause any significant systematic offsets.. Then, they were spatially stretched using spline interpolation in radial direction, going out from the image center, by the ratio of the central wavelengths of the filters to align the speckle patterns. Because of the possibly different SED shapes of our objects with respect to the standard calibration star used in Schmid et al. (2017) to determine the central wavelengths $\lambda_{c}$ of the filters, it is possible that $\lambda_{c}$ is slightly shifted for each object. This effect, however, is expected to alter the upscaling factor by at most 0.4% for B_Ha (assuming the unrealistic case in which $\lambda_{c}$ is at the edge of the filter), which is the broadest filter we used. This is negligible at very small separations from the star, where speckles dominate the noise. Values for filter central wavelengths and filter equivalent widths can be found in Table 5 of Schmid et al. (2017). The modified continuum images were then subtracted from the images taken simultaneously with the H$\alpha$ filter, leaving only H$\alpha$ line flux emitted from the primary star and potential companions. As a final step, the images resulting from the subtraction are derotated to a common field orientation and mean-combined. It is worth noting that if, as a result of the stretching, a potential point-source emitting a significant amount of continuum flux moves by more than $\lambda/D$, the signal strength in the H$\alpha$ image is only marginally changed in the SDI subtraction step, and only the speckle noise is reduced. If this is not the case, this subtraction step yields a significant reduction of the source signal in addition to the reduction of the speckle noise. For SPHERE/ZIMPOL H$\alpha$ imaging, a conservative SDI subtraction without substantial signal removal is achieved for angular separations $\gtrsim 0\aas@@fstack{\prime\prime}90$ ($\sim 250$ pixels). Nevertheless, this technique is expected to enhance the S/N of accreting planetary companions even at smaller separations, since young planets are not expected to emit a considerable amount of optical radiation in the continuum. In this case, the absence of a continuum signal guarantees that the image subtraction leaves the H$\alpha$ signal of the companion unchanged and only reduces the speckle residuals. Therefore, for this science case, there is no penalty for using SDI.

To perform ASDI, the SDI (H$\alpha$-Cnt_H$\alpha$) subtracted images are fed into the PCA pipeline to subtract any remaining residuals. During the analysis we varied the same parameters as described for simple ADI. The HD142527 dataset was used to compare the different sensitivities achieved when applying ADI, SDI, and ASDI. The results are discussed in Section 4.1.1 and Appendix B.

With ZIMPOL in imaging mode, there is a constant offset of $(135.99\pm 0.11)^{\circ}$ between the parallactic angle and the PA of the camera in sky coordinates (Maire et al., 2016). A preliminary astrometric calibration showed, however, that this reference frame has to be rotated by $(-2.0\pm 0.5)^{\circ}$ to align images with north pointing to the top (Ginski et al., in preparation). This means that overall, for every PSF subtraction technique, the final images have to be rotated by $(134\pm 0.5)^{\circ}$ in the counterclockwise direction.

4 Analysis and results

4.1 HD142527 B: The accreting M-star companion

4.1.1 Comparing the performance of multiple observational setups

In this section, we quantitatively compare the detection performance for multiple filter combinations and PSF subtraction techniques and establish the best strategy for future high-contrast H$\alpha$ observations with SPHERE/ZIMPOL. For the analysis, the HD142527 dataset was used; during the data reduction, no further frame selection was applied. The final images of HD142527 clearly show the presence of the M-star companion east of the central star. The signal is detected in all filters with ADI (B_Ha, N_Ha, and Cnt_Ha) and ASDI (in both continuum-subtracted B_Ha and N_Ha images) over a broad range of PCs and also for different image and inner mask sizes (see Figure 1).

We used the prescription from Mawet et al. (2014) to compute the false positive fraction (FPF) as a metric to quantify the confidence in the detection. The flux is measured in apertures of diameter $\lambda/D$ (16.5 mas) at the position of the signal and in equally spaced reference apertures placed at the same separation but with different PAs, so that there is no overlap between these angles and the remaining azimuthal space is filled. These apertures sample the noise at the separation of the companion. Since the apertures closest to the signal are dominated by negative wings from the PSF subtraction process, they were ignored. Then, we used Equation 9 and Equation 10 from Mawet et al. (2014) to calculate S/N and FPF from these apertures. This calculation takes into account the small number of apertures that sample the noise and uses the Student t-distribution to calculate the confidence of a detection. The wider wings of the t-distribution enable a better match to a non-Gaussian residual speckle noise than the normal distribution. However, the true FPF values could be higher if the wings of the true noise distribution are higher than those of the t-distribution666As an example, Figure 7 of Mawet et al. (2014) shows how the t-distribution produces lower FPF values than the case where speckle noise follows more closely a modified Rician distribution. Nevertheless, it has been shown that applying ADI removes the correlated component of the noise leaving quasi-Gaussian residuals (Marois et al., 2008)..

The narrow N_Ha filter delivers a significantly lower FPF than the broader B_Ha filter over a wide range of PCs (see Figure 15 in Appendix B). Figure 15 also shows that the combination of SDI and ADI yields lower FPF values than only ADI for both filters. Applying ASDI on N_Ha images is hence the preferred choice for future high-contrast imaging programs with SPHERE/ZIMPOL in the speckle-limited regime close to the star. Furthermore, as shown in Figure 16 and explained in Appendix C, it is crucial to plan observations maximizing the field rotation to best modulate and subtract the stellar PSF and to achieve higher sensitivities.

In Figure 2 we show the resulting contrast curves for the three filters for a confidence level (CL) of 99.99995%. For each dataset (B_Ha, N_Ha, and Cnt_Ha) and technique (ADI and ASDI), we calculated the contrast curves for different numbers of PCs (between 10 and 30 in steps of 5) after removing the companion (see Section 4.1.2). From each set of curves, we only considered the best achievable contrast at each separation from the central star. The presence of H$\alpha$ line emission from the central star made SDI an inefficient technique to search for faint objects at small angular separations.

To derive the contrast curves, artificial companions with varying contrast were inserted at six different PAs (separated by 60∘) and in steps of $0\aas@@fstack{\prime\prime}03$ in the radial direction. As the stellar PSF was unsaturated in all individual frames, the artificial companions were obtained by shifting and flux-scaling the stellar PSFs and then adding these companions to the original frames. Also for the calculation of the ASDI contrast curves, the original H$\alpha$ filter images, containing underlying continuum and H$\alpha$ line emission, were used to create artificial secondary signals. For each reduction run only one artificial companion was inserted at a time to keep the PCs as similar as possible to the original reduction. The brightness of the artificial signals was reduced/increased until their FPF corresponded to a detection with a CL of 99.99995% (i.e., a FPF of 2.5$\times 10^{-7}$), corresponding to $\approx$5$\sigma$ whether Gaussian noise was assumed. An inner mask with a radius of $0\aas@@fstack{\prime\prime}02$ was used to exclude the central parts dominated by the stellar signal. The colored shaded regions around each curve represent the standard deviation of the contrast achieved at that specific separation within the six PAs.

It is important to note that, while in Figure 15 the N_Ha filter provides the lowest FPF for the companion, Figure 2 seems to suggest that the B_Ha filter provides a better contrast performance. However, this is an effect from the way the contrast analysis is performed. As described above, the stellar PSF was used as a template for the artificial planets, as it is usually done in high-contrast imaging data analysis. The flux distribution within a given filter can vary significantly depending on the object. In this specific case, HD142527 B is known to have H$\alpha$ excess emission, hence the flux within either H$\alpha$ filter is strongly dominated by line emission ($\sim$50% in B_Ha and $\sim$83% in N_Ha filter) and a contribution from the optical continuum can be neglected. The primary shows, however, strong and non-negligible optical continuum emission that contributes to the flux observed in the H$\alpha$ filters. Indeed, for the primary, only 10% and 56% of the flux in the B_Ha and N_Ha filters are attributable to line emission. Hence, when using the stellar PSF as template for artificial planets, we obtain a better contrast performance for the B_Ha filter as it contains overall more flux. In reality, however, if the goal is to detect H$\alpha$ line emission from low-mass accreting companions, the N_Ha filter is to be preferred. Finally, as found by Sallum et al. (2015) for the planet candidate LkCa15 b, the fact that ASDI curves reach a deeper contrast confirms that this technique, in particular close to the star, is more effective and should be preferred to search for H$\alpha$ accretion signals.

4.1.2 Quantifying the H$\alpha$ detection

The clear detection of the M-star companion in our images allows us to determine its contrast in all the filters and its position relative to the primary at the epoch of observation. For this purpose, we applied the Hessian matrix approach (Quanz et al., 2015) and calculated the sum of the absolute values of the determinants of Hessian matrices in the vicinity of the companion’s signal. The Hessian matrix represents the second derivative of an n-dimensional function and its determinant is a measure for the curvature of the surface described by the function. This method allows for a simultaneous determination of the position and the flux contrast of the companion and we applied a Nelder-Mead (Nelder & Mead, 1965) simplex algorithm to minimize the curvature, i.e., the determinants of the Hessian matrices. We inserted negative, flux-rescaled stellar PSFs at different locations and with varying brightness in the input images and computed the resulting curvature within a region of interest (ROI) around the companion after PSF subtraction777For this analysis we used an image size of $0\aas@@fstack{\prime\prime}36\times 0\aas@@fstack{\prime\prime}36$ to speed up the computation and an inner mask of 10.8 mas (radius).. To reduce pixel-to-pixel variations after the PSF-subtraction step and allow for a more robust determination of the curvature, we convolved the images with a Gaussian kernel with a full width at half maximum (FWHM) of 8.3 mas ($\approx 0.35$ of the FWHM of the stellar PSF, which was calculated to be 23.7 mas on average). To fully include the companion’s signal, the ROI was chosen to be $(43.2\times 43.2)$ mas around the peak flux detected in the original set of PSF subtracted images. Within the ROI, the determinants of the Hessian matrices in 10,000 evenly spaced positions on a fixed grid (every 0.43 mas) were calculated and summed up.

For the optimization algorithm to converge, we need to provide a threshold criterion: if the change in the parameters (position and contrast) between two consecutive iterations is less than a given tolerance, the algorithm has converged and the optimization returns those values for contrast and position. The absolute tolerance for the convergence was set to be 0.1888This is an absolute value, meaning that if the sum of the determinants can be lowered only using steps in pixels and contrast lower than 0.1, then the algorithm stops., as this value is the precision to which artificial signals can be inserted into the image grid. This value applies for all the investigated parameters (position and contrast). Errors in the separation and PA measurements take into account the tolerance given for the converging algorithm and the finite grid. Errors in the contrast magnitude only consider the uncertainty due to the tolerance of the optimization. To account for systematic uncertainties in the companion’s location and contrast resulting from varying self-subtraction effects in reductions with different numbers of PCs, we ran the Hessian matrix algorithm for reductions with PCs in the range between 13 and 29 and considered the average of each parameter as final result. This range of PCs corresponds to FPF values below $2.5\times 10^{-7}$ (see Figure 15). To quantify the overall uncertainties in separation, PA, and contrast in a conservative way, we considered the maximum/minimum value (including measurement errors) among the set of results for the specific parameter and computed its difference from the mean. In Figure 3, we present the results from this approach for the N_Ha dataset and show the comparison between the original residual image and the image with the companion successfully removed.

4.1.3 Astrometry

The previously described algorithm was used to determine the best combination of separation, PA, and magnitude contrast for HD142527 B. In the N_Ha data the companion is located at $63.3^{+1.3}_{-1.0}$ mas from the primary star, in the B_Ha dataset at $62.3^{+1.7}_{-2.2}$ mas, and in the Cnt_Ha data at $62.8^{+2.1}_{-1.9}$ mas. The corresponding PAs are $(97.8\pm 0.9)^{\circ}$, $(99.4^{+1.1}_{-1.5})^{\circ}$ and $(99.0^{+1.5}_{-1.6})^{\circ}$, respectively. Errors in the PA measurements also take into account the above mentioned uncertainty in the astrometric calibration of the instrument, which was added in quadrature to the PA error bars.

As within the error bars all filters gave the same results, we combined them and found that HD142527 B is located at a projected separation of $62.8^{+2.1}_{-2.7}$ mas from the primary star ($9.9^{+0.3}_{-0.4}$ AU at $157.3\pm 1.2$ pc) and has a PA of $(98.7\pm 1.8)^{\circ}$. The final values result from calculating the arithmetic mean of all the values obtained from the three different datasets, while their errors are calculated identically to those for each single dataset.

In Figure 4 we compare the positions previously estimated (Close et al., 2014a; Rodigas et al., 2014; Lacour et al., 2016; Christiaens et al., 2018) and that resulting from our analysis. Lacour et al. (2016) used a Markov chain Monte Carlo analysis to infer the orbital parameters of HD142527 B. Because the past detections were distributed over a relatively small orbital arc ($\sim 15^{\circ}$), it was difficult to constrain the parameters precisely. The high precision measurement added by our SPHERE/ZIMPOL data extends the arc to a range of $\sim 30^{\circ}$. An updated orbital analysis is provided in (Claudi2018). Figure 4 shows that HD142527 B is clearly approaching the primary in the plane of the sky.

4.1.4 Photometry

The Hessian matrix approach yields the contrasts between HD142527 A and B in every filter: $\Delta$N_Ha $=6.3^{+0.2}_{-0.3}\text{ mag}$ in the narrow band, $\Delta$B_Ha $=6.7\pm 0.2\text{mag}$ in the broad band, and $\Delta$Cnt_Ha$=7.3^{+0.3}_{-0.2}\text{ mag}$ in the continuum filter. To quantify the brightness of the companion and not only its contrast with respect to the central star, we determined the flux of the primary in the multiple filters. We measured the count rate ($cts$) in the central circular region with radius $\sim 1\aas@@fstack{\prime\prime}5$ in all frames of each stack and computed the mean and its uncertainty $\sigma/\sqrt{n}$, where $\sigma$ is the standard deviation of the count rate within the dataset and $n$ is the number of frames. No aperture correction was required because the same aperture size was used by Schmid et al. (2017) to determine the zero points for the flux density for the three filters from photometric standard star calibrations. To estimate the continuum flux density we used their Equation 4

[TABLE]

where $c_{zp}^{cont}(Cnt\_Ha)$ is the zero point of the Cnt_Ha filter, $\textit{cts}=1.105\,(\pm 0.001)\times 10^{5}$ ct/s is the count rate measured from our data, $am=1.06$ is the average airmass, $k_{1}$ is the atmospheric extinction at Paranal ($k_{1}(\lambda)=0.085$ mag/airmass for Cnt_Ha, $k_{1}(\lambda)=0.082$ mag/airmass for B_Ha and N_Ha; cf. Patat et al., 2011), and $m_{mode}=-0.23$ mag is the mode dependent transmission offset, which takes into account the enhanced throughput of the R-band dichroic with respect to the standard gray beam splitter. The flux density of the primary star in the continuum filter $F^{*}_{\lambda}(Cnt\_Ha)$ was then used to estimate the fraction of counts in the line filters due to continuum emission via

[TABLE]

where $c_{zp}^{cont}(F\_Ha)$ is the continuum zero point of the H$\alpha$ filter used in the observations (cf. Schmid et al., 2017). During this step, we assumed that the continuum flux density was the same in the three filters. The continuum count rate was subtracted from the total count rate in B_Ha and N_Ha, $\text{cts}(B\_Ha)=1.631\,(\pm 0.001)\times 10^{5}$ ct/s and $\text{cts}(N\_Ha)=3.903\,(\pm 0.003)\times 10^{4}$ ct/s, leaving only the flux due to pure H$\alpha$ emission. These were used, together with Equation (1) with line zero points, to determine the pure H$\alpha$ line fluxes (see fifth column in Table 3). For each filter, the continuum flux density was multiplied by the filter equivalent width, and the flux contribution from line emission was added for the line filters. As in Sallum et al. (2015), we assumed the B object to have the same extinction as A, ignoring additional absorption from the disk. Indeed, we considered an extinction of $A_{V}=0.05$ mag (Fairlamb et al., 2015) and, interpolating the standard reddening law of Mathis (1990) for $R_{V}=3.1$, we estimated the extinction at $\sim 650$ nm to be A${}_{H\alpha}=0.04$ mag. The stellar flux was found to be $6.1\,\pm 0.2\times 10^{-11}$ erg/s/cm2 in the Cnt_Ha filter, $9.7\,\pm 0.8\times 10^{-11}$ erg/s/cm2 in the B_Ha filter and $3.0\,\pm 0.8\times 10^{-11}$ erg/s/cm2 in the N_Ha filter (see Table 3).

With the empirically estimated contrasts, we calculated the companion flux, i.e., line plus continuum emission or continuum only emission, in the three filters as follows:

[TABLE]

[TABLE]

[TABLE]

We note that the contrast we calculated in the continuum filter is very similar to that obtained by Close et al. (2014a) of $\Delta\text{mag}=7.5\pm 0.25$ mag. The direct estimation of the brightness of the primary in each individual ZIMPOL filter led to a larger difference when comparing the companion’s apparent magnitude in our work ($m^{B}_{Cnt\_Ha}=15.4\pm 0.2$ mag) with that from Close et al. (2014a) ($m^{B}_{\text{Close}}=15.8\pm 0.3$ mag). Such values are possibly consistent within the typical variability of accretion of the primary and secondary at these ages. However, given the different photometry sources and filters used for the estimation of the stellar flux densities in the two works, the results cannot be easily compared.

4.1.5 Accretion rate estimates

The difference between the flux in the line filters and the continuum filter (normalized to the H$\alpha$ filter widths) represents the pure H$\alpha$ line emission for which we find for HD142527 B $f^{line}_{B\_Ha}=1.0^{+0.5}_{-0.4}\times 10^{-13}$ erg/s/cm2 and $f^{line}_{N\_Ha}=7.6^{+3.5}_{-2.9}\times 10^{-14}$ erg/s/cm2, respectively. The line flux is then converted into a line luminosity multiplying it by the GAIA distance squared (see Table 2), yielding $L_{B\_Ha}=7.7^{+4.0}_{-3.6}\times 10^{-5}\,L_{\odot}$ and $L_{N\_Ha}=6.0^{+2.8}_{-2.4}\times 10^{-5}\,L_{\odot}$. We then estimated the accretion luminosity with the classical T Tauri stars (CTTS) relationship from Rigliaco et al. (2012), in which the logarithmic accretion luminosity grows linearly with the logarithmic H$\alpha$ luminosity

[TABLE]

and $a=1.49\pm 0.05$ and $b=2.99\pm 0.16$ are empirically determined. We calculated the accretion luminosity for both datasets, yielding $L^{\mathrm{acc}}_{B\_Ha}=7.3^{+6.8}_{-6.4}\times 10^{-4}L_{\odot}$ and $L^{\mathrm{acc}}_{N\_Ha}=5.0^{+4.4}_{-4.0}\times 10^{-4}L_{\odot}$.

Following Gullbring et al. (1998) we finally used

[TABLE]

to constrain the mass accretion rate. The quantity $G$ is the universal gravitational constant, and $R_{c}$ and $M_{c}$ are the radius and mass of the companion, respectively. Assuming that the truncation radius of the accretion disk $R_{in}$ is $\sim 5R_{c}$, we obtain $\left(1-\frac{R_{c}}{R_{in}}\right)^{-1}\sim 1.25$. For the companion mass and radius, two different sets of values were considered: Lacour et al. (2016) fitted the SED of HD142527 B with evolutionary models (Baraffe et al., 2003) and calculated $M_{c}=0.13\pm 0.03\;M_{\odot}$ and $R_{c}=0.9\pm 0.15\;R_{\odot}$, while Christiaens et al. (2018) estimated from H+K band VLT/SINFONI spectra $M_{c}=0.34\pm 0.06\;M_{\odot}$ and $R_{c}=1.37\pm 0.05\;R_{\odot}$, in the presence of a hot circumstellar environment999They considered two different cases in which the companion may or may not be surrounded by a hot environment contributing in $H$+$K$. Because of the presence of accreting material shown in this work, we decided to consider the first case.. The accretion rates obtained from the H$\alpha$ emission line are $\dot{M}_{B\_Ha}=2.0^{+2.0}_{-1.9}\times 10^{-10}M_{\odot}/\text{yr}$ and $\dot{M}_{N\_Ha}=1.4^{+1.3}_{-1.2}\times 10^{-10}M_{\odot}/\text{yr}$ in the first case and $\dot{M}_{B\_Ha}=1.2\pm 1.1\times 10^{-10}M_{\odot}/\text{yr}$ and $\dot{M}_{N\_Ha}=0.8\pm 0.7\times 10^{-10}M_{\odot}/\text{yr}$ in the second case. Some H$\alpha$ flux loss from the instrument when the N_Ha filter is used might explain the lower value of $\dot{M}_{N\_Ha}$ compared to $\dot{M}_{B\_Ha}$. Indeed, according to Figure 2 and Table 5 from Schmid et al. (2017), the N_Ha filter is not perfectly centered on the H$\alpha$ rest wavelength, implying that a fraction of the flux could be lost, in particular if the line profile is asymmetric. Moreover, high temperature and high velocities of infalling material cause H$\alpha$ emission profiles of CTTS to be broad (Hartmann et al., 1994; White & Basri, 2003). Also, line broadening from the rotation and line shift of the object due to possible radial motion might be important, even though it is not expected to justify the $\sim$40% H$\alpha$ flux difference of HD142527B. We argue, therefore, that with the available data it is very difficult to estimate the amount of line flux lost by the N_Ha filter, and that the value given by the B_Ha filter is expected to be more reliable, since all line emission from the accreting companion is included.

As shown in PDI images from Avenhaus et al. (2017), dust is present at the separation of the secondary possibly fully embedding the companion or in form of a circumsecondary disk. During our calculations, we neglected any local extinction effects due to disk material. It is therefore possible that on the one hand some of the intrinsic H$\alpha$ flux gets absorbed/scattered and the actual mass accretion rate is higher than that estimated in this work; on the other hand, the material may also scatter some H$\alpha$ (or continuum) emission from the central star, possibly contributing in very small amounts to the total detected flux.

Although the results obtained in this work are on the same order of magnitude as those obtained by Close et al. (2014a), who derived a rate of $6\times 10^{-10}\,M_{\odot}\text{ yr}^{-1}$, it is important to point out some differences in the applied methods. Specifically, Close et al. (2014a) used the flux estimated in the H$\alpha$ filter to calculate $L_{H\alpha}$, while we subtracted the continuum flux and considered only the H$\alpha$ line emission. Moreover, we combined the derived contrast with the stellar flux in the H$\alpha$ filters obtained from our data, while Close et al. (2014a) used the $R$-band magnitude of the star. As HD142527 A is also accreting and therefore emitting H$\alpha$ line emission, this leads to a systematic offset. Finally, Close et al. (2014a) used the relationship found by Fang et al. (2009) and not that from Rigliaco et al. (2012), leading to a difference in the $L_{\mathrm{H\alpha}}-L_{\mathrm{acc}}$ conversion.

4.2 HD135344 B

Visual inspection of the final PSF-subtracted ADI images of HD135344B showed a potential signal north to the star. Given the weakness of the signal and the low statistical significance, we analyze and discuss it further in Appendix D.

In Figure 5 we plot the contrast curves obtained as explained in section 4.1.1 using the N_Ha and the Cnt_Ha datasets and applying ASDI. In addition to the $1\aas@@fstack{\prime\prime}08\times 1\aas@@fstack{\prime\prime}08$ images we also examined $2\aas@@fstack{\prime\prime}88\times 2\aas@@fstack{\prime\prime}88$ images to search for accreting companions beyond the contrast limited region and beyond the spiral arms detected on the surface layer of the HD135344 B circumstellar disk. However, no signal was detected. We paid special attention to the separations related to the reported disk cavities (Andrews et al., 2011; Garufi et al., 2013). We chose to investigate specifically the cavity seen in scattered light at $0\aas@@fstack{\prime\prime}18$. The outer radius of the cavity seen in millimeter continuum is larger, but small dust grains are expected to be located inside of this radius increasing the opacity and making any companion detection more difficult. Neglecting the small inclination ($i\sim 11^{\circ}$, Lyo et al., 2011), the disk is assumed to be face-on and the contrast value given by the curve of Figure 5 at $0\aas@@fstack{\prime\prime}18$ is considered ($\Delta$N_Ha = 9.8 mag). We derived the H$\alpha$ flux from the star in the N_Ha filter as presented in section 4.1.4 using the stellar flux values for the different filters given in Table 3, and calculated the upper limits for the companion flux, accretion luminosity, and mass accretion rate following Section 4.1.4 and Section 4.1.5. The accretion rate is given by Equation 4, assuming a planet mass of $M_{c}=10.2\,M_{J}$, the maximum mass that is nondetectable at those separations according to the analysis of Maire et al. (2017). Being consistent with their approach, we then used AMES-Cond101010AMES-Cond and BT-Settl models used through the paper where downloaded on Feb. 06, 2018, from https://phoenix.ens-lyon.fr/Grids/AMES-Cond/ISOCHRONES/ and https://phoenix.ens-lyon.fr/Grids/BT-Settl/CIFIST2011_2015/ISOCHRONES/, respectively. evolutionary models (Allard et al., 2001; Baraffe et al., 2003) to estimate the radius of the object $R_{c}=1.6\,R_{J}$ based on the age of the system. All values, sources, and models used are summarized in Table 3 and in Table 4 together with all the information for the other objects. The final accretion rate upper limit has been calculated to be $<2.4\times 10^{-12}\,M_{\odot}\text{ yr}^{-1}$ at an angular separation of $0\aas@@fstack{\prime\prime}18$, i.e., the outer radius of the cavity seen in scattered light.

4.3 TW Hya

The TW Hya dataset does not show any point source either in the $1\aas@@fstack{\prime\prime}08\times 1\aas@@fstack{\prime\prime}08$ images (see Figure 6) or in the $2\aas@@fstack{\prime\prime}88\times 2\aas@@fstack{\prime\prime}88$ images, which are large enough to probe all the previously reported disk gaps. The final contrast curves are shown in Figure 7. We also looked specifically at detection limits within the gaps observed by van Boekel et al. (2017) and focused in particular on the dark annulus at 20 AU ($0\aas@@fstack{\prime\prime}39$) from the central star, which has a counterpart approximately at the same position in 870 $\mu$m dust continuum observations (Andrews et al., 2016)

Since the circumstellar disk has a very small inclination, we considered the disk to be face-on and assumed the gaps to be circular. At $0\aas@@fstack{\prime\prime}39$, planets with contrast lower than 9.3 mag with respect to TW Hya would have been detected with the ASDI technique (cf. Figure 7). This value was then combined with the stellar flux calculated as described in section 4.1.4, to obtain the upper limit of the companion flux in the B_Ha filter. This yielded $\dot{M}<1.0\times 10^{-11}\,M_{\odot}\text{ yr}^{-1}$ (see Table 4) as the upper limit for the mass accretion rate based on our SPHERE/ZIMPOL dataset.

4.4 HD100546

The HD100546 dataset suffered from rather unstable and varying observing conditions, which resulted in a large dispersion in the recorded flux (see Figure 18 in Appendix E). We hence selected only the last 33% of the observing sequence, which had relatively stable conditions, for our analysis (see Appendix E). The H$\alpha$ data did not confirm either of the two protoplanet candidates around HD100546 (see Figure 6) and we show the resulting detection limits in Figure 8.

In order to investigate the detection limits at the positions of the protoplanet candidates, we injected artificial planets with increasing contrast starting from $\Delta$B_Ha = 8.0 mag until the signal was no longer detected with a CL of at least 99.99995%, and we repeated the process subtracting different numbers of PCs (from 10 to 30). At the position where Quanz et al. (2015) claimed the presence of a protoplanetary companion, we would have been able to detect objects with a contrast lower than 11.4 mag (using PC=14 and the ADI reduction). Consequently, if existing, a 15 MJ companion (Quanz et al., 2015) located at the position of HD100546 b must be accreting at a rate $<6.4\times 10^{-12}\,M_{\odot}\text{ yr}^{-1}$ in the framework of our analysis and assuming no dust is surrounding the object. We note that, in comparison to the accretion luminosity $L_{\mathrm{acc}}$ estimated by Rameau et al. (2017), our upper limit is one order of magnitude lower (cf. Table 4).

For the position of HD100546 c, we used the orbit given in Brittain et al. (2014) to infer the separation and PA of the candidate companion at the epoch of our observations, i.e., $\rho\simeq 0\aas@@fstack{\prime\prime}14$ and $\text{PA}\simeq 133^{\circ}$. At this position our data reach a contrast of 9.3 mag (using PC=14 on the continuum-subtracted dataset), implying an upper limit for the companion flux in the H$\alpha$ filter of $7.9\times 10^{-14}$ erg/s/cm2 and a mass accretion rate $<1.1\times 10^{-10}\,M_{\odot}\text{ yr}^{-1}$. This puts $\sim 2$ orders of magnitude stronger constraints on the accretion rate of HD100546 c than the limits obtained from the polarimetric H$\alpha$ images presented in Mendigutía et al. (2017) for a $15\,M_{J}$ planet. We note that owing to its orbit, HD100546 c is expected to have just disappeared or to disappear quickly behind the inner edge of the disk (Brittain et al., 2014). Therefore, extinction could play a major role in future attempts to detect this source.

4.5 HD169142

We analyzed the data with ADI and ASDI reductions (see Figure 6 for the ADI image). The latter was particularly interesting in this case because the stellar flux density in the continuum and H$\alpha$ filter is very similar and the continuum subtraction almost annihilated the flux from the central PSF, indicating that the central star has limited to no H$\alpha$ line emission (cf. Table 3 and see Grady et al. (2007)). We calculated the detection limits as explained in section 4.1.1 for both filters for a confidence level of 99.99995%, as shown in Figure 9.

We investigated with particular interest the positions of the candidates mentioned in Section 2 and derived specific detection limits at their locations, independent from the azimuthally averaged contrast curve. At the position of the compact source found by Osorio et al. (2014) (we call this potential source HD169142 c), our data are sensitive to objects 10.7 mag fainter than the central star (obtained by subtracting 16 PCs with ASDI reduction). At the position of HD169246 b (Reggiani et al., 2014; Biller et al., 2014) an object with a contrast as large as 9.9 mag could have been detected (PC=19; ASDI). For the compact source from Osorio et al. (2014) we found $\dot{M}<4.4\times 10^{-11}\,M_{\odot}\text{ yr}^{-1}$. Similarly, for the object detected by Biller et al. (2014) and Reggiani et al. (2014)111111Within the uncertainties in the derived positions, these objects are indistinguishable and hence we assume it is the same candidate. we found an upper limit for the mass accretion rate of $\dot{M}<7.6\times 10^{-11}\,M_{\odot}\text{ yr}^{-1}$.

4.6 MWC 758

Our analysis of the SPHERE/ZIMPOL images did not show an H$\alpha$ counterpart to the MWC 758 companion candidate detected by Reggiani et al. (2018) as shown in Figure 6. This is consistent with the recently published results from Huélamo et al. (2018). Nonetheless, we provide a detailed analysis and discussion of the same MWC 758 data to allow a comparison with the other datasets.

In Figure 10 we show the detection limits obtained with ADI for the B_Ha and Cnt_Ha dataset, and the results of the ASDI approach. At separations larger than $0\aas@@fstack{\prime\prime}25$, companions with a contrast smaller than 10 mag could have been detected. At the specific position of the candidate companion121212For our analysis we considered the position obtained from the first dataset in Reggiani et al. (2018) because the observing date was close to the epoch of the H$\alpha$ observations. we can exclude objects with contrasts lower than 9.4 mag (obtained subtracting 15 PCs using ASDI).

To explain the presence of a gap in dust-continuum emission without a counterpart in scattered light, a steady replenishment of $\mu$m-sized particle is required, which implies that a companion in the inner disk should not exceed a mass of $M_{c}=5.5\,M_{J}$ (Pinilla et al., 2015; Reggiani et al., 2018). In line with the analysis of Reggiani et al. (2018), we used the BT-Settl model to estimate the radius of the companion and we derived an upper limit for the mass accretion rate of $\dot{M}<5.5\times 10^{-11}\,M_{\odot}\text{ yr}^{-1}$ (see Table 4). Our analysis puts slightly stronger constraints on the mass accretion rate in comparison to that in Huélamo et al. (2018).

5 Discussion

5.1 SPHERE/ZIMPOL as hunter for accreting planets

The SPHERE/ZIMPOL H$\alpha$ filters allow for higher angular resolution compared to filters in the infrared regime and can, in principle, search for companions closer to the star. For comparison, a resolution element is 5.8 times smaller in the H$\alpha$ filter than in the $L^{\prime}$ filter, meaning that the inner working angle (IWA) is smaller by the same amount so that closer-in objects could be observed, if bright enough131313We note that SPHERE does not operate at similarly high Strehl ratios in the optical regime as it is able to do in the infrared.. An instrument with similar capabilities is MagAO (Close et al., 2014b; Morzinski et al., 2016), but as the Magellan telescope has a primary mirror of 6.5 m diameter, it has a slightly larger IWA than SPHERE at the 8.2 m VLT/UT3 telescope. A direct comparison of the HD142527 B detection shows that ZIMPOL reaches a factor $\sim 2.5$ higher S/N in one-third of total integration time and field rotation of MagAO under similar seeing conditions, even if the companion is located $\gtrsim 20$ mas closer to the star. The VAMPIRES instrument combined with Subaru/SCExAO will soon be a third facility able to perform H$\alpha$ imaging in SDI mode (Norris et al., 2012)

In terms of detection performance using different filters and reduction techniques, we re-emphasize that the N_Ha filter is more efficient in detecting H$\alpha$ signals in the contrast limited regime. The smaller filter width reduces the contribution of the continuum flux, which often dominates the signal in the B_Ha filter, particularly for the central star. Hence, assuming the planetary companion emits only line radiation, the N_Ha filter reduces the contamination by the stellar signal in the remaining speckles. Moreover, the subtraction of the stellar continuum from H$\alpha$ images reduces the speckles in both B_Ha and N_Ha filters. Hence, ASDI enhances the signal of potential faint companions, in particular at separations $<0\aas@@fstack{\prime\prime}3$ (cf. Figures 7, 9, 10), where companions 0.7 mag fainter appear accessible in comparison to using simple ADI. ASDI should always be applied during the analysis of SPHERE/ZIMPOL H$\alpha$ data.

What remains to be quantified is how longer detector integration times (DITs) or the broad band filter could improve the detection limits in the background limited regime (i.e., $>0\aas@@fstack{\prime\prime}3$ where the contrast curves are typically flattening out) or for fainter natural guide stars. At these separations narrow band data can be detector read noise limited and the B_Ha filter might be more suitable because of its higher throughput. However, as we show in Figure 11, it seems that at least for our HD142527 dataset this does not seem to be the case. Future studies conducted in both filters and on several objects are required to derive a more comprehensive understanding. Finding the sweetspot between longer integration times and the smearing of the PSF because of field rotation is also warranted. At least for the object considered in Figure 11, at large separations (usually $>0\aas@@fstack{\prime\prime}3$, in the background limited region) it is even possible to ignore completely ADI and simply apply field stabilized observations.

5.2 Constraining planet accretion

For our mass accretion rate estimates of HD142527 B we assumed that 100% of the H$\alpha$ flux originates from accretion processes involving circumstellar material. We note, however, that the values may be overestimated if we consider that chromospheric activity of the M star (White & Basri, 2003; Fang et al., 2009) can also contribute to the measured line flux. Furthermore, as mentioned in Section 4.1.5, we warn that the narrow width of the N_Ha filter might be too narrow to fully encompass all H$\alpha$ line emission from fast-moving, accreting material, and therefore the results may be underestimated. Finally, given the presence of dusty material at the projected position of HD142527 B (Avenhaus et al., 2017), H$\alpha$ flux might have been partially absorbed. It is beyond the scope of this paper to properly estimate a value for intrinsic extinction due to disk material and consider this value in the $\dot{M}$ estimation. Nevertheless, in Figure 12 we show the fraction of H$\alpha$ flux that is potentially lost because of extinction as a function of $A_{V}$, converted into $A_{H\alpha}$ as explained in Section 4.1.4. Only 2% of the H$\alpha$ signal remains if the disk material causes an extinction of $A_{V}=5$ mag. This plot quantifies the impact of dust on the measured flux and the detectability of H$\alpha$ emission from embedded objects.

For the other five objects studied in this work we were not able to detect any clear accretion signature located in the disks. Therefore, our data were not able to support the scenario in which protoplanets are forming in those disks. We put upper limits on the accretion luminosity and mass accretion rate. Two notes have to be made: (1) the fundamental quantities directly derived from the data are $F_{\mathrm{H\alpha}}$ and $L_{\mathrm{H\alpha}}$; they should be used for future comparisons with other datasets or objects; (2) the presented upper limits on $\dot{M}$ are only valid for an object with the mass and radius given in Table 4, while the $L_{\mathrm{acc}}$ upper limits refer to objects of any mass. In particular, assuming lower mass objects implies larger $\dot{M}$, as shown in Figure 13: on the y-axis the mass accretion rate upper limits decrease as a function of the companion mass, for which the corresponding radius was calculated using the evolutionary models reported in Table 4 and assuming the age listed in Table 2. The plot highlights that the assumed mass of the companion may change the final $\dot{M}_{acc}$ by more than one order of magnitude. Moreover, we overplot in violet the mass accretion rates of the three objects presented in Zhou et al. (2014, see also Section 5.3) as well as LkCa15 b and PDS70 b (Sallum et al., 2015; Wagner et al., 2018), and in gray the range of mass accretion rates for HD142527.

We stress that, similar to HD142527 B, we always assumed that the flux limit is completely due to H$\alpha$ line emission without any contribution from continuum or chromospheric activity. Furthermore, for our analysis we always neglected intrinsic extinction effects from disk material, which likely weaken the signal. In particular, at locations where no gap in small dust grains has been identified the extinction $A_{H\alpha}$ can be significant (see Figure 12). Models and precise measurements of the dust content in the individual disks would be required to properly include local extinction into our analysis. Finally, investigating the H$\alpha$ luminosity upper limits for the specific positions as a function of the separation from the central star, it can be noticed that the constraints are stronger at larger separations. The only exception is HD100546, for which higher upper limits were achieved. The combination of suboptimal weather conditions, under which the dataset was taken, and the small field rotation of the subsample analyzed in this work made those limits worse. A more stable dataset with larger field rotation should provide more constraining limits.

5.3 Comparison with other objects

The accretion rate of HD142527 B is in good agreement with the mass accretion rates found in Rigliaco et al. (2012) for low-mass TTauri stars in the $\sigma$ Ori star-forming region ($5\times 10^{-11}\,M_{\odot}\text{ yr}^{-1}<\dot{M}_{CTTS}<10^{-9}\,M_{\odot}\text{ yr}^{-1}$). A slightly broader mass accretion rate range was found by Alcalá et al. (2014), with $2\times 10^{-12}\,M_{\odot}\text{ yr}^{-1}<\dot{M}_{CTTS}<4\times 10^{-8}\,M_{\odot}\text{ yr}^{-1}$ in the Lupus star-forming region.

Zhou et al. (2014) reported three very low-mass objects (GSC 06214-00210 b, GQ Lup b and DH Tau b), which exhibit H$\alpha$ emission from accretion. Those objects have separations of 100-350 AU from their parent stars and $\dot{M}\sim 10^{-9}-10^{-11}\,M_{\odot}\text{ yr}^{-1}$ (see violet stars in Figure 13). The accretion rates measured in the paper are of the same order as the limits we found in our work. At projected distances similar to those of the three objects mentioned above, ZIMPOL would have been able to observe and detect H$\alpha$ emitting companions. However, closer to the star in the contrast limited regime, our data would not have detected accretion processes occurring with $\dot{M}\lesssim 10^{-11}\,M_{\odot}\text{ yr}^{-1}$.

The mass accretion rate of PDS70 b was estimated by Wagner et al. (2018) without considering any extinction effects and it is slightly lower than the limits we achieve for our sample (see violet square in Figure 13 and black star in Figure 14). The flux was calculated from the contrast in Wagner et al. (2018) assuming $R_{\text{PDS70\,b}}=11.7$ mag and estimating the MagAO H$\alpha$ filter widths assuming a flat SED141414 https://visao.as.arizona.edu/software_files/visao/html/group__reduction__users__guide.html#visao_filters. In order to properly compare our limits and their H$\alpha$ detection, the same confidence levels should be considered. We therefore estimated the contrast limit for a CL corresponding to a $4\sigma$ detection for HD142527 at the separation of PDS70 b, which was 0.3 mag lower than the limits corresponding to a CL of 99.99995%. Hence, to bring all the contrast curves from Figure 14 to a 4$\sigma$ confidence level at $\sim 0\aas@@fstack{\prime\prime}19$, a multiplication by a factor 0.76 is required. We note, however, that this scaling is just an approximation to provide a more direct comparison between the two studies.

We also compared the H$\alpha$ line luminosity upper limits obtained from our ZIMPOL H$\alpha$ sample with that estimated by Sallum et al. (2015) for LkCa15 b ($L_{H\alpha}\sim 6\times 10^{-5}\,L_{\odot}$). Our specific limits for the candidates around HD169142, HD100546, and MWC 758 are slightly lower, but, except for HD100546 b and the compact source in HD169142 found by Osorio et al. (2014), of the same order of magnitude. LkCa15 itself was observed with SPHERE/ZIMPOL during the science verification phase in ESO period P96. We downloaded and analyzed the data, which were, however, poor in quality and also in terms of integration time and field rotation. Only $\sim 1$ hr of data is available with a field rotation of $\sim 16^{\circ}$, a coherence time of $2.6\pm 0.8$ ms, and a mean seeing of $1\aas@@fstack{\prime\prime}64\pm 0\aas@@fstack{\prime\prime}37$. As we show in Figure 14, with deeper observations including more field rotation, ZIMPOL can potentially detect the signal produced by LkCa15 b (Sallum et al., 2015) with a CL of 99.99995%. However, the higher airmass at the Paranal Observatory and the fact that LkCa15 is a fainter guide star may complicate the redetection of the companion candidate, and therefore exceptional atmospheric conditions are required.

In addition to H$\alpha$ also other spectral features like Pa$\beta$ and Br$\gamma$ lines may indicate ongoing accretion processes onto young objects. As an example, Daemgen et al. (2017) used the absence of those lines in the spectrum of the low-mass companion HD106906 b to infer its mass accretion rate upper limits ($\dot{M}<4.8\times 10^{-10}\,M_{J}/\text{yr}^{-1}$). Their constraint is stronger than the ones we were able to put with our ZIMPOL H$\alpha$ data. Several other studies also detected hydrogen emission lines like Pa$\beta$ from low-mass companions (e.g., Seifahrt et al., 2007; Bowler et al., 2011; Bonnefoy et al., 2014), but unfortunately they did not calculate mass accretion rates.

5.4 Comparison with existing models

Two models for planetary accretion are currently used to explain the accreting phase of planet formation: magnetospheric accretion (Zhu, 2015) and boundary layer accretion (Owen & Menou, 2016). During magnetospheric accretion, the magnetic field truncates the CPD and hot ionized hydrogen in the closest regions of the disk falls onto the planet following the magnetic field lines. Recombination on the planet surface then produces H$\alpha$ flux. For protoplanets, these models predict H$\alpha$ luminosities at least three orders of magnitudes lower than in CTTS, according to equation 22 in Zhu (2015),

[TABLE]

This is mainly owing to a one order of magnitude smaller infall velocity $v_{s}$ and a one order of magnitude smaller truncation radius $R_{T}$ (squared in the $L_{H\alpha}$ equation) due to weaker magnetic fields than in stars. We combined the magnetospheric accretion models (Zhu, 2015) with existing detections in the infrared and evolutionary models. As an example, we present the case of HD100546 b. According to models (Zhu, 2015), the observed $L^{\prime}$ brightness could be emitted by a CPD with inner radius of $1-4\,R_{J}$ and $M_{p}\dot{M}$ of $0.2-2.9\times 10^{-6}\,M_{J}^{2}\text{ yr}^{-1}$. The mass accretion constraints obtained from H$\alpha$ ZIMPOL data would therefore imply that $M_{p}\gtrsim 31\,M_{J}$. This result is in conflict with that obtained by Quanz et al. (2015) and the AMES-Cond evolutionary models, since the object $L^{\prime}$ brightness excludes masses larger than $\sim 15\,M_{J}$. This is the mass expected in the case in which the $L^{\prime}$ flux is only from photospheric emission. Moreover, a 30 $M_{J}$ object would have significantly shaped the disk morphology and would have been clearly visible in other bands, such as the $K_{s}$-band, where Quanz et al. (2015) could only put upper limits to the companion brightness.

Szulágyi & Mordasini (2017) found that only a minimal fraction of the hydrogen in CPDs might be thermally ionized if the planet is massive and hot enough. Consequently, the disk does not get truncated and ionized material does not get accreted through magnetospheric accretion along the field lines. Then, disk material falls directly onto the planet (boundary layer accretion). The same authors showed that material falling from the circumstellar disk onto the CPD and the protoplanet shocks, and eventually produces H$\alpha$ line emission both from the CPD and the planet. The contribution to the H$\alpha$ flux is larger from the CPD than from the planet (Szulágyi & Mordasini, 2017). These authors also showed that the majority of the accreted gas, however, remains neutral, especially for planets $<10\,M_{J}$. Hence, the H$\alpha$ flux can only estimate the ionized gas accretion rate and not the total accreted material. According to their simulations, a 10 $M_{J}$ planet would be accreting at a rate of $5.7\times 10^{-8}\,M_{J}\text{ yr}^{-1}$, producing $L_{H\alpha}\sim 7\times 10^{-6}\,L_{\odot}$. This value is on the same order of the limits our data allow us to put on the H$\alpha$ luminosity from known forming protoplanet candidates. Since considering lower planetary masses enhances the mass accretion rate (see equation 4) and higher masses should be visible in other infrared bands, we conclude that either extinction from disk material plays a major role in the nondetection of the existing candidates, or they are false positives resulting from image post-processing.

The comparison of $L_{H\alpha}$ limits from Table 4 with Figure 7 from Mordasini et al. (2017) indicates that, assuming completely cold accretion, the observed objects may be low-mass ($0.1-1M_{J}$) medium accreters ($\dot{M}\sim 10^{-10}-10^{-9}M_{\odot}/\text{yr}$) or higher mass objects ($1-15M_{J}$) showing very little accretion ($\dot{M}<10^{-10.5}M_{\odot}/\text{yr}$). Mordasini et al. (2017) also suggested another possible reason for some of the nondetections in H$\alpha$. If some of the planets, such as HD100546 b, have not yet completely detached from the disk, they would be cooler and would not be accreting at high accretion rates. In a later phase, they will possibly be able to open a gap and accrete a large amount of material.

Another aspect that we did not consider is the effect of the circumplanetary disk inclination on the flux that is emitted. Zhu (2015) considered the disk inclination including a factor $1/\cos(i)$, where $i$ is the CPD inclination. Detailed accretion models should investigate the consequences of a tilted protoplanetary disk on $L_{H\alpha}$.

6 Conclusions

Imaging in H$\alpha$ is one of the promising techniques to detect forming planets at very small separations. In this context, the SPHERE/ZIMPOL instrument will play a major role in investigating local accretion signatures in circumstellar disks. An important next step is to redetect the previous discoveries of MagAO of H$\alpha$ emission from LkCa15 b and PDS 70 b and to study potential accretion variability. None of the possible protoplanet candidates discovered in the infrared (HD169142 b, MWC758 b, and HD100546 b and c) could be confirmed in this study searching for accretion signatures, implying several possible scenarios. Their mass accretion rates could be lower than our limits and therefore they are currently not detectable. Other explanations are that protoplanetary accretion shows variability and some of the objects are currently going through a period of quiescence, or that extinction effects from disk material absorb a considerable fraction of the light. The study of NIR line diagnostics might reduce the effects of absorption and allow the detection of accretion processes. Furthermore, it is possible that the observed candidates are disk features that have been enhanced by image post-processing (Follette et al., 2017; Ligi et al., 2018), or our understanding of accretion processes during the formation of giant planets is not correct and, as an example, the use of the CTTS scaling relation is not correct. In order to investigate this, precise simulations of protoplanetary accretion, as well as of disk intrinsic effects (via full radiative transfer), have to be developed and combined with multiwavelength observations spanning from the optical to the (sub)millimeter.

The estimation of upper limits are of particular importance for the study of accretion variability of protoplanets in the future. Continuing surveys for accreting planets could possibly detect H$\alpha$ signatures and combine these with detection limits provided by this work to investigate variability in the accretion processes. Finally, we emphasize that although a lot of effort was put into the calculation of mass accretion rate upper limits, those values are model and parameter dependent. The H$\alpha$ flux upper limits are, however, the fundamental quantities that were measured from the data and can be directly compared with future observations.

Acknowledgements.

SPHERE is an instrument designed and built by a consortium consisting of IPAG (Grenoble, France), MPIA (Heidelberg, Germany), LAM (Marseille, France), LESIA (Paris, France), Laboratoire Lagrange (Nice, France), INAF - Osservatorio di Padova (Italy), Observatoire de Genéve (Switzerland), ETH Zurich (Switerland), NOVA (Netherlands), ONERA (France), and ASTRON (Netherlands), in collaboration with ESO. SPHERE also received funding from the European Commission Sixth and Seventh Framework Programmes as part of the Optical Infrared Coordination Network for Astronomy (OPTICON) under grant number RII3-Ct-2004-001566 for FP6 (2004-2008), grant number 226604 for FP7 (2009-2012), and grant number 312430 for FP7 (2013-2016). This work has been carried out within the frame of the National Center for Competence in Research PlanetS supported by the Swiss National Science Foundation. SPQ and HMS acknowledge the financial support of the SNSF. GC and SPQ thank the Swiss National Science Foundation for financial support under grant number 200021_169131. FMe and GvdP acknowledge fundings from ANR of France under contract number ANR-16-CE31-0013. This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. The authors thank Arianna Musso-Barcucci for the preliminary analysis on HD142527.

Appendix A Influence of the beamsplitter on flux measurements

We investigated the effects of the throughput uncertainties of the two ZIMPOL arms resulting from instrument polarization effects. It is currently not known how the overall throughput to the individual ZIMPOL arms depends on the telescope and instrument configurations. However, it is easy to estimate the overall effect because the Nasmyth mirror of the VLT introduces an instrument polarization of about 4 %. This is reduced by the first mirror in SPHERE to about 0 to 3 %, while the following mirrors in the instrument add further positive or negative polarization contributions of about 2 %, while polarization cross talks (linear $\rightarrow$ circular polarization) reduce the linear polarization. Thus, it is safe to adopt a maximum error of 5 % for the relative difference (e.g., $T(H_{\alpha})=(1\pm 0.05)T(cont)$) in throughput between the two channels. We therefore tested the impact of an enhancement/decrease by at most 5% in the continuum flux, analyzing the consequences on the detection of HD142527B and on the contrast performances of our pipeline. The signal flux is measured in an aperture of radius 8.3 mas, and the FPF was calculated as explained in section 4.1.1. The results averaged over a range of PCs (PC=10,15,20,25,30) on the ASDI-processed B_Ha dataset are shown in Table 5; a central mask of 21.6 mas was applied. As one may expect, the signal flux shows a strong variation of $20-50\%$ in the ASDI images, which is mainly due to the stronger/weaker subtraction of the continuum. The relative difference in this case is increased from the initial 5% by the ASDI processing, but it should be noted that together with the signal, the noise level also gets increased/decreased, causing the FPF to be less subject to variations. Indeed, regarding the FPF values, we argue that depending on the arm 1 to arm 2 transmission the confidence of the detection is lower by a factor of $\sim$10 in both extreme cases, which corresponds approximately to a maximum variation of $\sim$0.1 mag in $\Delta$mag. Therefore, we do not expect this effect to have a large impact on the detection limits estimated in this work. Nevertheless, it is important to keep in mind that the calculation of the mass accretion rate of HD142527 B does not consider this effect and a more accurate description of the instrument behavior is required to correct for it.

Appendix B FPF analysis of the HD142527 B dataset

In order to identify the best strategy for future SPHERE/ZIMPOL observations, we compared the FPF calculated after the subtraction of different numbers of PCs using different techniques and datasets. For the ADI technique we considered three datasets: B_Ha, N_Ha, and Cnt_Ha, of which the last set contained all the images taken with the continuum filter. For SDI and ASDI, we considered the subtraction of the Cnt_Ha images from the respective H$\alpha$ filter images. All the images had a size of $1\aas@@fstack{\prime\prime}08\times 1\aas@@fstack{\prime\prime}08$. For each case, we applied an inner mask of varying size (10.8 mas, 21.6 mas, 32.4 mas) and chose the smallest FPF value as representative value for the detection. The FPF calculation (see Section 4.1.1) followed the prescription suggested in Mawet et al. (2014). Because of the strong negative wings of the companion in the PSF subtracted images, we decided not to consider the two background apertures closest to the signal as they are not representative of the background and speckle noise.

In the four panels of Figure 15 we analyzed the FPFs of HD142527 B, obtained using different combinations of techniques and datasets. In the top panels we compare the detection from different filters using the same technique: ADI on the left and ASDI on the right. For the ADI analysis, the B_Ha and Cnt_Ha datasets show similar values with a stronger detection in B_Ha for fewer subtracted components, while the FPF values obtained with the N_Ha filter are, for a wide range of PCs (from 11 to 32), $\sim$ 5 orders of magnitude lower. The detections with the ASDI technique show a similar trend; there is a stronger detection in N_Ha, particularly between 10 and 27 PCs. The normal SDI technique, which is not presented in the image, was not efficient enough to properly subtract the stellar PSF and did not reveal the companion. This is probably for two reasons: (1) the central star is actively accreting material and emitting strong H$\alpha$ flux, which cannot be subtracted accurately with the Cnt_Ha images, impeding the detection of the companion, and (2) PSF shapes are slightly different for different filters due to nonmatching bandpasses.

In the lower panels, we consider the results from the B_Ha (left) and N_Ha (right) datasets for ADI and ASDI. In both cases ASDI seems to be more efficient in detecting signals. A larger gain is obtained for the B_Ha filter, while FPFs obtained with the N_Ha filter have more similar values, probably due to the minor impact of the continuum subtraction on images taken with the narrow filter with respect to the broad filter. We conclude that the best observing strategy to look for accreting companions in the contrast limited regime with SPHERE/ZIMPOL is to take images in the N_Ha filter and Cnt_Ha filter simultaneously and to perform ASDI. It is of particular interest that in the case of HD142527, ASDI also performs better than ADI. Indeed, we could expect that the presence of a clear signal in the continuum would have strongly compromised the detection with ASDI. On the contrary, the detection is even stronger, implying that the subtraction of the stellar pattern is much more important than the self-subtraction of the companion, boosting its S/N. We note, however, that observing fainter objects might cause the data to be readout noise limited. In this case, the B_Ha filter might be preferred to the N_Ha filter. This hypothesis, however, should be confirmed with a fainter source than the bright M-dwarf HD142527B.

Appendix C Impact of field rotation and total integration time

In addition to the best instrumental setup for the detection of accreting objects in H$\alpha$ imaging, we also investigated the effect of two observational parameters on the achieved upper limits: the field rotation and integration time on target. Two subsets were created from the HD142527 B_Ha data. The first was composed of every second frame of the dataset, while the second only included the first half of the dataset frames. In this way, the field rotation of the first subset is twice that of the second subset, while the integration time is the same for the two subsets. Figure 16 shows the resulting contrast curves, calculated in the same way as described in Section 4.1.1. The dashed blue line represents the subset composed of the first half of the dataset, which allows us to reach $\sim 9.4$ mag of contrast at $0\aas@@fstack{\prime\prime}2$. The green solid line shows the contrast limits estimated from the subset composed of every second frame. It is clear that at all separations, this subset allows us to detect fainter objects than the other subset and at $0\aas@@fstack{\prime\prime}2$ the difference reaches 1.1 mag. Finally, the entire dataset allows us to go, at the same distance, another 0.3 mag deeper. At least for this dataset, the field rotation seems to play a very important role, allowing a better modeling and subtraction of the stellar PSF.

Appendix D Is a companion candidate orbiting HD135344B?

We visually inspected the final PSF-subtracted ADI images of HD135344 B, which showed a potential signal north of the star. The feature is persistent in the N_Ha and Cnt_Ha datasets for different mask radii (e.g., $0\aas@@fstack{\prime\prime}02$, $0\aas@@fstack{\prime\prime}03$, $0\aas@@fstack{\prime\prime}04$, $0\aas@@fstack{\prime\prime}05$, and $0\aas@@fstack{\prime\prime}06$) and over a wide range of PCs (6-21). In particular, when using larger mask radii, the close-in speckles are removed and the signal appears to be stronger. We then investigated smaller images (101$\times$101 pixels) with the same technique and confirmed the signal for different mask radii and PCs. Next, we examined the ASDI images and found that the signal is present once again in different reductions, but appears fainter. If the signal is from a physical source, this is expected from an accreting object emitting H$\alpha$ line radiation. This signal is shown in Figure 17 in the N_Ha and in the Cnt_Ha filter for the parameter setups that yield the lowest FPFs, which are $5.9\times 10^{-5}$ for N_Ha and 0.0015 for Cnt_Ha.

A careful look at the bright signal in the Cnt_Ha images raises doubts on the nature of its source as it is very compact and does not have a PSF-like shape. Furthermore, the signal has high FPFs, with a minimal value of $\sim 0.0015$, which is not statistically significant enough to claim a detection.

The signal in the N_Ha frames has a morphology resembling that of a faint physical source. Varying the number of PCs seems to influence the apparent shape and location of the signal, as expected from faint close-in objects with low S/N when subtracting the stellar PSF. Even though for 9 PCs the FPF reaches a minimal value of $\sim 5.9\times 10^{-5}$, the FPFs for 7-17 PCs are in the range $10^{-3}-10^{-2}$, which does not give us sufficient confidence to claim a detection.

As a final check, we used the Hessian matrix approach as described in Section 4.1.2 to perform a signal characterization. We ran the algorithm for PCs between 7 and 17 (where the final images showed a clear signal) with a central mask of radius 57.6 mas and a ROI of $8\times 8$ pixels. The other parameters were kept identical to the analysis performed on HD142527 B. The signal appears to be located at a separation of $71.1^{+4.8}_{-4.2}$ mas with a PA of $(19.1^{+2.2}_{-2.8})^{\circ}$. The contrast was measured to be $8.1\pm 0.4$ mag. As visible in the error bars, the positions found are spread over a range of $\sim 9$ mas, which corresponds to $\sim 2.5$ pixels. Normally, a physical point source should be less affected by systematics introduced by the PSF subtraction process. However, a low S/N object at a separation of 71.1 mas is more difficult to measure properly and a larger spread in the recovered positions could be the result. A similar note can be made for the contrast values, which span over a range of $\sim 0.8$ mag. We conclude that to settle this issue and fully understand the origin of the signal in the H$\alpha$ filter, a dataset with higher S/N would be required.

Appendix E Frame selection for the HD100546 dataset

As briefly described in Section 4.4, the large HD100546 dataset (1104 frames, cf. Table 1) was taken in unstable conditions, which made a frame selection necessary. To determine a frame selection metric, we plotted the mean count value per image (image dimensions $1\aas@@fstack{\prime\prime}08\times 1\aas@@fstack{\prime\prime}08$ pixels, see Figure 18). It turned out that three phases could be identified within the observing run: a short initial phase of stability with some outliers (120 frames), a long period of 619 frames where the mean count values spanned a range between $\sim 0$ and $\sim 55$ counts per pixels, and, finally, a large amount of stable frames at the end of the observations (mean pixel value $\sim 55$ in B_Ha). We decided only to keep the images of the last stable period, composed of the frames $739-1104$ to perform our analysis. This subsample has a total on-target integration time of 61 minutes and its field rotation is $\sim 20.7^{\circ}$.