© ESO 2019

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, 2018; Biller et al. 2014; 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 μ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.

Anotherapproach 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 inferthe 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α, when returning to lower energy states (e.g., Calvet & Gullbring 1998; Szulágyi & Mordasini 2017; Marleau et al. 2017). High-contrast imaging using Hα 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α 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α emission from the young planet PDS70 b using MagAO, albeit with comparatively low statistical significance (3.9σ).

In this paper we present aset of Hα 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 Sect. 2 we discuss the observations and target stars. We explain the data reduction in Sect. 3 and present our analyses in Sect. 4. In Sect. 5 we discuss our results in a broader context and conclude in Sect. 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)¹, 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α filter in one camera arm and a nearby continuum filter simultaneously in the other arm (Cont_Ha; λ_c = 644.9 nm, Δλ = 3.83 nm). As the data were observed in different programs, we sometimes used the narrow Hα filter (N_Ha; λ_c = 656.53 nm, Δλ = 0.75 nm) and sometimes the broad Hα filter (B_Ha; λ_c = 655.6 nm, Δλ = 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 notethat 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α frames. The inclined mirrors in the telescope and the instrument introduce di-attenuation (e.g., higher reflectivity for I_⊥ than I_∥) and polarization cross talks, so that the transmissions in imaging mode to the I_⊥ and I_∥ arm depend on the telescope pointing direction. This effect is at the level of a few percent (about ± 5%), 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

2.2.1 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; Claudi et al. 2019) that shows signatures of ongoing accretion in Hα emission (Close et al. 2014a). This companion orbits in a large, optically thin cavity within the circumstellar disk stretching from ~ 0′′. 07 to ~ 1′′. 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′′.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 ± 0.06 M_⊙. This value is a factor of~3 larger than that estimated by spectral energy distribution (SED) fitting (Lacour et al. 2016, M = 0.13 ± 0.03 M_⊙). Thanks to the accreting close-in companion, this system is the ideal target to optimize the Hα observing strategy with SPHERE/ZIMPOL and also the data reduction.

2.2.2 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′′.32. In polarimetric differential imaging (PDI) observations in the near-infrared (NIR), the outer radius of the cavity appears to be at 0′′. 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 evolutionarymodels these authors derived upper limits for the mass of potential giant planets around HD135344 B (3 M_J beyond 0′′. 7).

2.2.3 TW Hya

TW Hya is the nearest T Tauri star to Earth. Its almost face-on transitional disk (i ~ 7 ± 1°; 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 ~ 1′′. 48. Later, Akiyama et al. (2015) observed in H-band polarized images a gap at a separation of ~0′′.41. Using Atacama Large Millimeter Array (ALMA), Andrews et al. (2016) identified gaps from the radial profile of the 870 μm continuum emission at 0′′.41, 0′′. 68 and 0′′. 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′′.11, 0′′. 39, and 1′′. 57 from the central star. A clear gap was also identified by Rapson et al. (2015) at a separation of 0′′. 43 in Gemini/GPI polarimetric images and the largest gap at r ≃ 1′′.52 has also been observed in CO emission with ALMA (Huang et al. 2018).

2.2.4 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′′. 13, based on high-resolution NIR spectro-astrometry of CO emission lines. Another protoplanet candidate was claimed by Quanz et al. (2013a) using L′ band high-contrast imaging data. The object was found at 0′′.48 ± 0′′.04 from the central star, at a position angle (PA) of ${(8.9\pm 0.9)}^{°}$, with an apparent magnitude of L′ = 13.2 ± 0.4 mag. Quanz et al. (2015) reobserved HD100546 in different bands (L′, M′, K_s) and detected theobject 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_proj ~ 0′′.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α images withthe 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_acc < 1.7 × 10⁻⁴ L_⊙). The same data were analyzed by Follette et al. (2017), together with other Hα images (MagAO), H band spectra (GPI), and Y band polarimetric images (GPI). Their data exclude that HD100546 c is emitting in Hα with L_Hα > 1.57 × 10⁻⁴L_⊙.

2.2.5 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′′. 17 and an annular gap between 0′′.28 and 0′′.49. These results were confirmed by Osorio et al. (2014), who investigated the thermal emission (λ = 7 mm) of large dust grains in the HD169142 disk, identifying two annular cavities (~ 0′′.16 − 0′′.21 and ~ 0′′.28 − 0′′.48). The latter authors also identified a point source candidate in the middle of the outer cavity at a distance of 0′′. 34 and PA ~ 175°. Biller et al. (2014) and Reggiani et al. (2014) observed a point-like feature in NaCo L′ data at the outer edge of the inner cavity (separation = 0′′.11–0′′.16 and PA = 0°−7.4°). 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.

2.2.6 MWC 758

MWC 758 issurrounded by a pre-transitional disk (e.g., Grady et al. 2013). Andrews et al. (2011) found an inner cavity of ~ 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′ data by Reggiani et al. (2018), together with the claim of the detection of a point-like source embedded in the disk at (111 ± 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′ band is ~ 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′′.37 and 0′′.53 were discovered that are probably related to two clumps with large surface density of millimeter dust and a large cavity of ~ 0′′.26 in radius. Finally, Huélamo et al. (2018) observed MWC 758 in Hα with SPHERE/ZIMPOL, reaching an upper limit for the line luminosity of $L_{{\text{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 × 3.6 mas/pixel onto a square grid with an effective pixel scale of 3.6 mas × 3.6 mas (1024 × 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α filter and one for the continuum filter². 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′′.08 × 1′′.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. 2019). 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 λ, while the position of a physical object on the detector is independentof λ. The underlying assumption is that, given that λ_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α filter³. 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 respectto the standard calibration star used in Schmid et al. (2017) to determine the central wavelengths λ_c of the filters, it is possible that λ_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 λ_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α filter, leaving only Hα 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 λ∕D, the signal strength in the Hα 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 asignificant reduction of the source signal in addition to the reduction of the speckle noise. For SPHERE/ZIMPOL Hα imaging, a conservative SDI subtraction without substantial signal removal is achieved for angular separations ≳ 0′′.90 (~ 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α 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α-Cnt_Hα) 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 Sect. 4.1.1 and Appendix B.

With ZIMPOL in imaging mode, there is a constant offset of ${(135.99\pm 0.11)}^{°}$ 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)}^{°}$ to align images with north pointing to the top (Ginski et al., in prep.). This means that overall, for every PSF subtraction technique, the final images have to be rotated by ${(134\pm 0.5)}^{°}$ 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α observationswith 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 Fig. 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 λ∕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 thatthere 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 Eqs. (9) and (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-distribution⁴.

The narrow N_Ha filter delivers a significantly lower FPF than the broader B_Ha filter over a wide range of PCs (see Fig. B.1). Figure B.1 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 Fig. C.1 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 Fig. 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 Sect. 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α 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′′. 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α filter images, containing underlying continuum and Hα 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 × 10⁻⁷), corresponding to ≈5σ whether Gaussian noise was assumed. An inner mask with a radius of 0′′. 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 Fig. B.1 the N_Ha filter provides the lowest FPF for the companion, Fig. 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 artificialplanets, 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α excess emission, hence the flux within either Hα filter is strongly dominated by line emission (~50% in B_Ha and ~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α 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α 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α accretion signals.

Fig. 1 Final ADI and ASDI reduced images of HD142527. Top row: B_Ha, Cnt_Ha, and N_Ha filter images resulting in the lowest FPFs (1.5 × 10−11, 2.2 × 10−9, and < 10−17, corresponding to S/Ns of 13.1, 9.8, and 26.6, respectively). Bottom row: final images after ASDI reduction for B_Ha-Cnt_Ha and N_Ha-Cnt_Ha frames (4.4 × 10−16 and < 10−17, corresponding to S/Ns of 22.7 and 27.6). We give the number of subtracted PCs and the radius of the central mask in milliarcseconds in the top left corner of each image. The color scales are different for the two rows. Because all images of the top row have the same color stretch, the detection appears weaker in the continuum band.

Fig. 2 Contrast curves for HD142527. The colored shaded regions around each curve represent the standard deviation of the achieved contrast at the 6 azimuthal positions considered at each separation. The markers (red diamond, orange circle, and violet star) represent the contrast of HD142527 B.

4.1.2 Quantifying the Hα 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 subtraction⁵. 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 (≈ 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 × 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.1⁶, 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 × 10⁻⁷ (see Fig. B.1). 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 Fig. 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.

Fig. 3 Image of HD142527 before (top panel) and after (bottom panel) the insertion of the negative companionresulting from the Hessian matrix algorithm. The image flux scale is the same in both images. In this case 14 PCs were subtracted and a mask of 10.8 mas (radius) was applied on the 101 × 101 pixels images of the N_Ha stack.

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.0}^{+1.3}$ mas from the primary star, in the B_Ha dataset at ${62.3}_{-2.2}^{+1.7}$ mas, and in the Cnt_Ha data at ${62.8}_{-1.9}^{+2.1}$ mas. The corresponding PAs are ${(97.8\pm 0.9)}^{°}$, ${({99.4}_{-1.5}^{+1.1})}^{°}$ and ${({99.0}_{-1.6}^{+1.5})}^{°}$, 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.7}^{+2.1}$ mas from the primary star(${9.9}_{-0.4}^{+0.3}$ AU at 157.3 ± 1.2 pc) and has a PA of ${(98.7\pm 1.8)}^{°}$. 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 Fig. 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 (~ 15°), 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 ~ 30°. An updated orbital analysis is provided in Claudi et al. (2019). Figure 4 shows that HD142527 B is clearly approaching the primary in the plane of the sky.

Fig. 4 Position of HD142527 B based on NaCo sparse aperture masking (red pentagons), MagAO (cyan triangles), GPI non-redundant masking (dark green diamonds) and VLT/SINFONI (blue circle) data from Rodigas et al. (2014), Close et al. (2014a), Lacour et al. (2016), and Christiaens et al. (2018), together with the SPHERE/ZIMPOL observation presented in this work (light green square). The position of HD142527 A is shown with the yellow star at coordinates (0,0).

4.1.4 Photometry

The Hessian matrix approach yields the contrasts between HD142527 A and B in every filter: ΔN_Ha$={6.3}_{-0.3}^{+0.2}\text{\hspace{0.17em}mag}$ in the narrow band, ΔB_Ha = 6.7 ± 0.2mag in the broad band, and ΔCnt_Ha$={7.3}_{-0.2}^{+0.3}\text{\hspace{0.17em}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 ~ 1′′.5 in all frames of each stack and computed the mean and its uncertainty $\sigma /\sqrt{n}$, where σ 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 Eq. (4) $$F_{\lambda }^{*}(\text{Cnt}\_\text{Ha})=\text{\hspace{0.17em}cts}\times {10}^{0.4(\text{am}\times k_1+m_{\text{mode}})}\cdot c_{\text{zp}}^{\text{cont}}(\text{Cnt}\_\text{Ha}),$$(1)

where $c_{\text{zp}}^{\text{cont}}(\text{Cnt}\_\text{Ha})$ is the zero point of the Cnt_Ha filter, cts = 1.105 (± 0.001) × 10⁵ ct s⁻¹ is the count rate measured from our data, am = 1.06 is the average airmass, k₁ is the atmospheric extinction at Paranal (k₁(λ) = 0.085 mag/airmass for Cnt_Ha, k₁(λ) = 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 }^{*}(\text{Cnt}\_\text{Ha})$ was then used to estimate the fraction of counts in the line filters due to continuum emission via $$\text{cts}(\text{F}\_\text{Ha})=\frac{F_{\lambda }^{*}(\text{Cnt}\_\text{Ha})}{c_{\text{zp}}^{\text{cont}}(F\_\text{Ha})}\times {10}^{-0.4(\text{am}\cdot k_1+m_{\text{mode}})},$$(2)

where $c_{\text{zp}}^{\text{cont}}(\text{F}\_\text{Ha})$ is the continuum zero point of the Hα 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, cts (B_Ha) = 1.631 (±0.001) × 10⁵ ct s⁻¹ and cts (N_Ha) = 3.903 (±0.003) × 10⁴ ct s⁻¹, leaving only the flux due to pure Hα emission. These were used, together with Eq. (1) with line zero points, to determine the pure Hα 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 ~650 nm to be A_Hα = 0.04 mag. The stellar flux was found to be 6.1 ± 0.2 × 10⁻¹¹ erg s⁻¹ cm⁻² in the Cnt_Ha filter, 9.7 ± 0.8 × 10⁻¹¹ erg s⁻¹ cm⁻² in the B_Ha filter and 3.0 ± 0.8 × 10⁻¹¹ erg s⁻¹ cm⁻² 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: $$F_{\text{Cnt}\_\text{Ha}}^p={7.4}_{-2.1}^{+1.4}\times {10}^{-14}\text{erg\hspace{0.17em}}{\text{s}}^{-1}{\text{cm}}^{-2},$$ $$F_{\text{B}\_\text{Ha}}^p=2.0\pm 0.4\times {10}^{-13}\text{erg\hspace{0.17em}}{\text{s}}^{-1}{\text{cm}}^{-2},$$ $$F_{\text{N}\_\text{Ha}}^p={9.1}_{-2.9}^{+3.5}\times {10}^{-14}\text{erg\hspace{0.17em}}{\text{s}}^{-1}{\text{cm}}^{-2}.$$

We notethat the contrast we calculated in the continuum filter is very similar to that obtained by Close et al. (2014a) of Δmag = 7.5 ± 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_{\text{Cnt}}^B\_\text{Ha\hspace{0.17em}}=15.4\pm 0.2$ mag) with that from Close et al. (2014a; $m_{\text{Close}}^B=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α filter widths) represents the pure Hα line emission for which we find for HD142527 B $f_{\text{B}\_\text{Ha}}^{\text{line}}={1.0}_{-0.4}^{+0.5}\times {10}^{-13}$ erg s⁻¹ cm⁻² and $f_{\text{N}\_\text{Ha}}^{\text{line}}={7.6}_{-2.9}^{+3.5}\times {10}^{-14}$ erg s⁻¹ cm⁻², respectively. The line flux is then converted into a line luminosity multiplying it by the Gaia distance squared (see Table 2), yielding $L_{\text{B}\_\text{Ha}}={7.7}_{-3.6}^{+4.0}\times {10}^{-5}{L}_{\odot }$ and $L_{\text{N}\_\text{Ha}}={6.0}_{-2.4}^{+2.8}\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α luminosity $$\log (L_{\text{acc}})=b+a\log (L_{\text{H}\alpha }),$$(3)

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

Following Gullbring et al. (1998) we finally used $${\dot{M}}_{\text{acc}}={\left(1-\frac{R_c}{R_{\text{in}}}\right)}^{-1}\frac{L_{\text{acc}}{R}_c}{G{M}_c}~1.25\frac{L_{\text{acc}}{R}_c}{G{M}_c}$$(4)

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 ~ 5R_c, we obtain ${\left(1-\frac{R_c}{R_{\text{in}}}\right)}^{-1}~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 ± 0.03 M_⊙ and R_c = 0.9 ± 0.15 R_⊙, while Christiaens et al. (2018) estimated from H + K band VLT/SINFONI spectra M_c = 0.34 ± 0.06 M_⊙ and R_c = 1.37 ± 0.05 R_⊙, in the presence of a hot circumstellar environment⁷. The accretion rates obtained from the Hα emission line are ${\dot{M}}_{\text{B}\_\text{Ha}}={2.0}_{-1.9}^{+2.0}\times {10}^{-10}{M}_{\odot }{\text{yr}}^{-1}$ and ${\dot{M}}_{\text{N}\_\text{Ha}}={1.4}_{-1.2}^{+1.3}\times {10}^{-10}{M}_{\odot }{\text{yr}}^{-1}$ in the first case and Ṁ_{B_Ha} = 1.2 ± 1.1 × 10⁻¹⁰ M_⊙ yr⁻¹ and Ṁ_{N_Ha} = 0.8 ± 0.7 × 10⁻¹⁰ M_⊙∕yr in the second case. Some Hα flux loss from the instrument when the N_Ha filter is used might explain the lower value of Ṁ_{N_Ha} compared to Ṁ_{B_Ha}. Indeed, according to Fig. 2 and Table 5 from Schmid et al. (2017), the N_Ha filter is not perfectly centered on the Hα 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α 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 ~ 40% Hα 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α 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 someHα (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 × 10⁻¹⁰ M_⊙yr⁻¹, it is important to point out some differences in the applied methods. Specifically, Close et al. (2014a) used the flux estimated in the Hα filter to calculate L_Hα, while we subtracted the continuum flux and considered only the Hα line emission. Moreover, we combined the derived contrast with the stellar density flux in the Hα 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α 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_Hα − L_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 Fig. 5 we plot the contrast curves obtained as explained in Sect. 4.1.1 using the N_Ha and the Cnt_Ha datasets and applying ASDI. In addition to the 1′′. 08 × 1′′.08 images we also examined 2′′.88 × 2′′.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′′.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 ~ 11°, Lyo et al. 2011), the disk is assumed to be face-on and the contrast value given by the curve of Fig. 5 at 0′′. 18 is considered (ΔN_Ha = 9.8 mag). We derived the Hα flux from the star in the N_Ha filter as presented in Sect. 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 Sects. 4.1.4 and 4.1.5. The accretion rate is given by Eq. (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-Cond⁸ 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 Tables 3 and 4 together with all the information for the other objects. The final accretion rate upper limit has been calculated to be < 2.4 × 10⁻¹² M_⊙yr⁻¹ at an angular separation of 0′′.18, i.e., the outer radius of the cavity seen in scattered light.

Fig. 5 Contrast curves for HD135344 B. The vertical lines indicate the outer radii of the cavities in small and large dust grains presented in Garufi et al. (2013) and Andrews et al. (2011), respectively.

4.3 TW Hya

The TW Hya dataset does not show any point source either in the 1′′. 08 × 1′′.08 images (see Fig. 6) or in the 2′′.88 × 2′′.88 images, which are large enough to probe all the previously reported disk gaps. The final contrast curves are shown in Fig. 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′′. 39) from the central star, which has a counterpart approximately at the same position in 870 μ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′′.39, planets with contrast lower than 9.3 mag with respect to TW Hya would have been detected with the ASDI technique (cf. Fig. 7). This value was then combined with the stellar flux calculated as described in Sect. 4.1.4, to obtain the upper limit of the companion flux in the B_Ha filter. This yielded Ṁ < 1.0 × 10⁻¹¹ M_⊙yr⁻¹ (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 Fig. E.1). We hence selected only the last 33% of the observing sequence, which had relatively stable conditions, for our analysis (see Appendix E). The Hα data did notconfirm either of the two protoplanet candidates around HD100546 (see Fig. 6) and we show the resulting detection limits in Fig. 8.

In order to investigate the detection limits at the positions of the protoplanet candidates, we injected artificial planets with increasing contrast starting from Δ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 todetect objects with a contrast lower than 11.4 mag (using PC = 14 and the ADI reduction). Consequently, if existing, a 15 M_J companion (Quanz et al. 2015) located at the position of HD100546 b must be accreting at a rate < 6.4 × 10⁻¹² M_⊙yr⁻¹ in the framework of our analysis and assuming no dust is surrounding the object. We note that, in comparison to the accretion luminosity L_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., ρ ≃ 0′′.14 and PA ≃ 133°. 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α filter of 7.9 × 10⁻¹⁴ erg s⁻¹ cm⁻² and a mass accretion rate < 1.1 × 10⁻¹⁰ M_⊙yr⁻¹. This puts ~2 orders of magnitude stronger constraints on the accretion rate of HD100546 c than the limits obtained from the polarimetric Hα 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.

Fig. 6 Final PSF subtracted ADI images of TW Hya, HD100546, HD169142, and MWC 758. We applied a central mask with radius 32.4 mas and 18 PCs were removed. No companion candidates were detected. All images have a linear, but slightly different, color scale.

Fig. 7 Contrast curves for TW Hya. The vertical line indicates the gap at 0′′. 39 detected in both scattered light (Akiyama et al. 2015; van Boekel et al. 2017) and submillimeter continuum (Andrews et al. 2016).

4.5 HD169142

We analyzed the data with ADI and ASDI reductions (see Fig. 6 for the ADI image). The latter was particularly interesting in this case because the stellar flux density in the continuum and Hα 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α line emission (cf. Table 3 and see Grady et al. 2007). We calculated the detection limits as explained in Sect. 4.1.1 for both filters for a confidence level of 99.99995%, as shown in Fig. 9.

We investigated with particular interest the positions of the candidates mentioned in Sect. 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 Ṁ < 4.4 × 10⁻¹¹ M_⊙yr⁻¹. Similarly, for the object detected by Biller et al. (2014) and Reggiani et al. (2014)⁹ we found an upper limit for the mass accretion rate of Ṁ < 7.6 × 10⁻¹¹ M_⊙yr⁻¹.

4.6 MWC 758

Our analysis of the SPHERE/ZIMPOL images did not show an Hα counterpart to the MWC 758 companion candidate detected by Reggiani et al. (2018) as shown in Fig. 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 Fig. 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′′. 25, companions with a contrast smaller than 10 mag could have been detected. At the specific position of the candidate companion¹⁰ 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 μ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 Ṁ < 5.5 × 10⁻¹¹ M_⊙yr⁻¹ (see Table 4). Our analysis puts slightly stronger constraints on the mass accretion rate in comparison to that in Huélamo et al. (2018).

Fig. 8 Contrast curves for HD100546. The gray dashed vertical line shows the separation of the outer gap edge cavity presented in Avenhaus et al. (2014a), while the solid blue lines indicate the separations of the forming planet candidates around HD100546 (Quanz et al. 2013a; Brittain et al. 2014).

5 Discussion

5.1 SPHERE/ZIMPOL as hunter for accreting planets

The SPHERE/ZIMPOL Hα 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α filter than in the L′ filter, meaning that the inner working angle (IWA) is smaller by the same amount so that closer-in objects could be observed, if bright enough¹¹. 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 ~ 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 ≳ 20 mas closer to the star. The VAMPIRES instrument combined with Subaru/SCExAO will soon be a third facility able to perform Hα 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α 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α 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′′.3 (cf. Figs. 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α 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′′.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 Fig. 11, it seems that at least for our HD142527 dataset this does not seem to be the case. Futurestudies 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 Fig. 11, at large separations (usually > 0′′.3, in the background limited region) it is even possible to ignore completely ADI and simply apply field stabilized observations.

Fig. 9 Contrast curves for HD169142. The shaded region represents the annular gap observed in scattered light (Quanz et al. 2013b) and in millimeter continuum (Osorio et al. 2014). The blue vertical lines represent the separation of the companion candidates (Reggiani et al. 2014; Biller et al. 2014; Osorio et al. 2014).

5.2 Constraining planet accretion

For our mass accretion rate estimates of HD142527 B we assumed that 100% of the Hα 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 Sect. 4.1.5, we warn that the narrow width of the N_Ha filter might be too narrow to fully encompass all Hα 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α 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 Ṁ estimation. Nevertheless, in Fig. 12 we show the fraction of Hα flux that is potentially lost because of extinction asa function of A_V, converted into A_Hα as explained in Sect. 4.1.4. Only 2% of the Hα 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α 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_Hα and L_Hα; they should be used for future comparisons with other datasets or objects; (2) the presented upper limits on Ṁ are only valid for an object with the mass and radius given in Table 4, while the L_acc upper limits refer to objects of any mass. In particular, assuming lower mass objects implies larger Ṁ, as shown in Fig. 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 Ṁ_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 Sect. 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 B.

We stress that, similar to HD142527 B, we always assumed that the flux limit is completely due to Hα 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α can be significant (see Fig. 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α 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 combinationof 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.

Fig. 10 Contrast curves for MWC 758. The gray dashed line shows the outer edge of the dust cavity observed by Andrews et al. (2011). The blue solid line indicates the separation at which Reggiani et al. (2018) found a candidate companion.

Fig. 11 Apparent flux detection limits as a function of the angular separation from HD142527 for both B_Ha and N_Ha filters.

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 T Tauri stars in the σ Ori star-formingregion (5 × 10⁻¹¹ M_⊙yr⁻¹ < Ṁ_CTTS < 10⁻⁹ M_⊙yr⁻¹). A slightly broader mass accretion rate range was found by Alcalá et al. (2014), with 2 × 10⁻¹² M_⊙yr⁻¹ < Ṁ_CTTS < 4 × 10⁻⁸ M_⊙yr⁻¹ in the Lupus star-formingregion. Zhou et al. (2014) reported three very low-mass objects (GSC 06214-00210 b, GQ Lup b and DH Tau b), which exhibit Hα emission from accretion. Those objects have separations of 100–350 AU from their parent stars and Ṁ ~ 10⁻⁹−10⁻¹¹ M_⊙yr⁻¹ (see violet stars in Fig. 13). The accretion rates measured in the paper are of the same order as thelimits 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α emitting companions. However, closer to the star in the contrast limited regime, our data would not have detected accretion processes occurring with Ṁ ≲ 10⁻¹¹ M_⊙yr⁻¹.

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 Fig. 13 and black star in Fig. 14). The flux was calculated from the contrast in Wagner et al. (2018) assuming R_{PDS70 b} = 11.7 mag and estimating the MagAO Hα filter widths assuming a flat SED¹². In order to properly compare our limits and their Hα detection, the same confidence levels should be considered. We therefore estimated the contrast limit for a CL corresponding to a 4σ 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 Fig. 14 to a 4σ confidence level at ~0′′.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α line luminosity upper limits obtained from our ZIMPOL Hα sample with that estimated by Sallum et al. (2015) for LkCa15 b (L_Hα ~ 6 × 10⁻⁵ L_⊙). 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 ~ 1 h of data is available with a field rotation of ~16°, a coherence time of 2.6 ± 0.8 ms, and a mean seeing of 1′′.64 ± 0′′.37. As we show in Fig. 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 Observatoryand 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α also other spectral features like Paβ and Brγ 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 (Ṁ < 4.8 × 10⁻¹⁰ M_J yr⁻¹). Their constraint is stronger than the ones we were able to put with our ZIMPOL Hα data. Several other studies also detected hydrogen emission lines like Paβ 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.

Fig. 12 Fraction of Hα flux absorbed as a function of the disk extinction AV assuming the extinction law of Mathis (1990) as explained in Sect. 4.1.4.

Fig. 13 Mass accretion rate upper limits as a function of the planetary mass for all the candidate forming planets investigated in this work. The violet stars represent the values reported in Zhou et al. (2014), while the violet squares indicate PDS70 b (Wagner et al. 2018) and LkCa15 b (Sallum et al. 2015). The gray shaded area represents the mass accretion rate of HD142527 B and is shown for mass accretion rate comparisonpurposes only. Indeed, the mass of the object is much larger than what is reported on the x-axis of the plot.

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α flux. For protoplanets, these models predict Hα luminosities at least three orders of magnitudes lower than in CTTS, according to Eq. (22) in Zhu (2015), $$L_{\text{H}\alpha }=4.7\times {10}^{-6}{L}_{\odot }{\left(\frac{R_T}{R_J}\right)}^2\left(\frac{v_s}{59{\text{kms}}^{-1}}\right).$$(5)

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α 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′ brightness could be emitted by a CPD with inner radius of 1–4 R_J and M_pṀ of $0.2\text{–}2.9\times {10}^{-6}{M}_J^2{\text{yr}}^{-1}$. The mass accretion constraints obtained from Hα ZIMPOL data would therefore imply that M_p ≳ 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′ brightness excludes masses larger than ~15 M_J. This is the mass expected in the case in which the L′ 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α line emission both from the CPD and the planet. The contribution to the Hα 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α 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 × 10⁻⁸ M_Jyr⁻¹, producing L_Hα ~ 7 × 10⁻⁶ L_⊙. This value is on the same order of the limits our data allow us to put on the Hα luminosity from known forming protoplanet candidates. Since considering lower planetary masses enhances the mass accretion rate (see Eq. (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α limits from Table 4 with Fig. 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 (Ṁ ~ 10⁻¹⁰–10⁻⁹M_⊙ yr⁻¹) or higher mass objects (1−15M_J) showing very little accretion (Ṁ < 10^−10.5M_⊙ yr⁻¹). Mordasini et al. (2017) also suggested another possible reason for some of the nondetections in Hα. 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α.

Fig. 14 Detection limits in apparent flux obtained for a 99.99995% CL in this work, together with limits achieved with the available ZIMPOL dataset for LkCa15 b (red dashed line) and the result presented in Sallum et al. (2015) and Wagner et al. (2018). A deeper dataset is required to redetect LkCa15 b with ZIMPOL, but this detection is feasible.

6 Conclusions

Imaging in Hα 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α 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α 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α 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. S.P.Q. and H.M.S. acknowledge the financial support of the SNSF. G.C. and S.P.Q. thank the Swiss National Science Foundation for financial support under grant number 200021_169131. F.M. and G.v.d.P. 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 → 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_α) = (1 ± 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 Sect. 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 A.1; 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 ~ 10 in both extreme cases, which corresponds approximately to a maximum variation of ~ 0.1 mag in Δ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 considerthis 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α filter images. All the images had a size of 1′′.08 × 1′′.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 Sect. 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 thetwo background apertures closest to the signal as they are not representative of the background and speckle noise.

In the four panels of Fig. B.1 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 forfewer subtracted components, while the FPF values obtained with the N_Ha filter are, for a wide range of PCs (from 11 to 32), ~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α 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 filterwith 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.

Fig. B.1 Performance comparison for Hα imaging with SPHERE/ZIMPOL using different filters (narrow and broad Hα) and reduction techniques (ADI and ASDI) based on the HD142527 dataset. In all panels the FPF obtained for HD142527 B is shown, as a function of the number of subtracted PCs used in PynPoint. On the right side of each panel we give the scale of the S/N to improve understanding of the plot and to compare different instrumental setups. We note again that this does not correspond to the classical σ notation. The gray regions indicates a confidence level for the detection of HD142527 B of at least 99.99995%, i.e., >5σ in case of Gaussian noise. Because of the applied corrections for small sample statistics, the border of the gray area does not correspond to an S/N of 5.

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α 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 thedataset, 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 C.1 shows the resulting contrast curves, calculated in the same way as described in Sect. 4.1.1. The dashed blue line represents the subset composed of the first half of the dataset, which allows us to reach ~ 9.4 mag of contrast at 0′′.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′′. 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.

Fig. C.1 Contrast curves calculated for the “first half” (blue dashed line) and the “very second frame” subsets (green solid line), and for the entire B_Ha dataset of HD142527 (orange dotted line).

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 thestar. The feature is persistent in the N_Ha and Cnt_Ha datasets for different mask radii (e.g., 0′′. 02, 0′′. 03, 0′′. 04, 0′′. 05, and 0′′. 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 × 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α line radiation. This signal is shown in Fig. D.1 in the N_Ha and in the Cnt_Ha filter for the parameter setups that yield the lowest FPFs, which are 5.9 × 10⁻⁵ 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 ~ 0.0015, which is notstatistically significant enough to claim a detection.

Fig. D.1 Lowest FPF images of HD135344 B (top panel: N_Ha filter; bottom panel: Cnt_Ha filter). The radius r of the inner mask and the number of subtracted PCs are given in each panel. The location of the tentative companion candidate is indicated by the arrow (see Sect. 4.1.2).

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 ~ 5.9 × 10⁻⁵, the FPFs for 7–17

PCs are in the range 10⁻³−10⁻², which does not give us sufficient confidence to claim a detection.

As a final check, we used the Hessian matrix approach as described in Sect. 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 × 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.2}^{+4.8}$ mas with a PA of ${({19.1}_{-2.8}^{+2.2})}^{°}$. The contrast was measured to be 8.1 ± 0.4 mag. As visible in the error bars, the positions found are spread over a range of ~ 9 mas, which corresponds to ~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 ~ 0.8 mag. We conclude that to settle this issue and fully understand the origin of the signal in the Hα filter, a dataset with higher S/N would be required.

Appendix E Frame selection for the HD100546 dataset

As briefly described in Sect. 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′′.08 × 1′′.08 pixels, see Fig. E.1). 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 ~ 0 and ~55 counts per pixels, and, finally, a large amount of stable frames at the end of the observations (mean pixel value ~ 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 min and its field rotation is ~ 20.7°.

Fig. E.1 B_Ha (blue circles) and Cnt_Ha (green crosses) mean count rates as a function of the image number in the observing sequence. The shaded region at the end represents the subset of frames that was chosen for the analysis.

References

All Tables

All Figures

Fig. 1 Final ADI and ASDI reduced images of HD142527. Top row: B_Ha, Cnt_Ha, and N_Ha filter images resulting in the lowest FPFs (1.5 × 10−11, 2.2 × 10−9, and < 10−17, corresponding to S/Ns of 13.1, 9.8, and 26.6, respectively). Bottom row: final images after ASDI reduction for B_Ha-Cnt_Ha and N_Ha-Cnt_Ha frames (4.4 × 10−16 and < 10−17, corresponding to S/Ns of 22.7 and 27.6). We give the number of subtracted PCs and the radius of the central mask in milliarcseconds in the top left corner of each image. The color scales are different for the two rows. Because all images of the top row have the same color stretch, the detection appears weaker in the continuum band.

In the text

	Fig. 2 Contrast curves for HD142527. The colored shaded regions around each curve represent the standard deviation of the achieved contrast at the 6 azimuthal positions considered at each separation. The markers (red diamond, orange circle, and violet star) represent the contrast of HD142527 B.
In the text

	Fig. 3 Image of HD142527 before (top panel) and after (bottom panel) the insertion of the negative companionresulting from the Hessian matrix algorithm. The image flux scale is the same in both images. In this case 14 PCs were subtracted and a mask of 10.8 mas (radius) was applied on the 101 × 101 pixels images of the N_Ha stack.
In the text

Fig. 4 Position of HD142527 B based on NaCo sparse aperture masking (red pentagons), MagAO (cyan triangles), GPI non-redundant masking (dark green diamonds) and VLT/SINFONI (blue circle) data from Rodigas et al. (2014), Close et al. (2014a), Lacour et al. (2016), and Christiaens et al. (2018), together with the SPHERE/ZIMPOL observation presented in this work (light green square). The position of HD142527 A is shown with the yellow star at coordinates (0,0).

In the text

	Fig. 5 Contrast curves for HD135344 B. The vertical lines indicate the outer radii of the cavities in small and large dust grains presented in Garufi et al. (2013) and Andrews et al. (2011), respectively.
In the text

	Fig. 6 Final PSF subtracted ADI images of TW Hya, HD100546, HD169142, and MWC 758. We applied a central mask with radius 32.4 mas and 18 PCs were removed. No companion candidates were detected. All images have a linear, but slightly different, color scale.
In the text

	Fig. 7 Contrast curves for TW Hya. The vertical line indicates the gap at 0′′. 39 detected in both scattered light (Akiyama et al. 2015; van Boekel et al. 2017) and submillimeter continuum (Andrews et al. 2016).
In the text

	Fig. 8 Contrast curves for HD100546. The gray dashed vertical line shows the separation of the outer gap edge cavity presented in Avenhaus et al. (2014a), while the solid blue lines indicate the separations of the forming planet candidates around HD100546 (Quanz et al. 2013a; Brittain et al. 2014).
In the text

	Fig. 9 Contrast curves for HD169142. The shaded region represents the annular gap observed in scattered light (Quanz et al. 2013b) and in millimeter continuum (Osorio et al. 2014). The blue vertical lines represent the separation of the companion candidates (Reggiani et al. 2014; Biller et al. 2014; Osorio et al. 2014).
In the text

	Fig. 10 Contrast curves for MWC 758. The gray dashed line shows the outer edge of the dust cavity observed by Andrews et al. (2011). The blue solid line indicates the separation at which Reggiani et al. (2018) found a candidate companion.
In the text

	Fig. 11 Apparent flux detection limits as a function of the angular separation from HD142527 for both B_Ha and N_Ha filters.
In the text

	Fig. 12 Fraction of Hα flux absorbed as a function of the disk extinction AV assuming the extinction law of Mathis (1990) as explained in Sect. 4.1.4.
In the text

Fig. 13 Mass accretion rate upper limits as a function of the planetary mass for all the candidate forming planets investigated in this work. The violet stars represent the values reported in Zhou et al. (2014), while the violet squares indicate PDS70 b (Wagner et al. 2018) and LkCa15 b (Sallum et al. 2015). The gray shaded area represents the mass accretion rate of HD142527 B and is shown for mass accretion rate comparisonpurposes only. Indeed, the mass of the object is much larger than what is reported on the x-axis of the plot.

In the text

	Fig. 14 Detection limits in apparent flux obtained for a 99.99995% CL in this work, together with limits achieved with the available ZIMPOL dataset for LkCa15 b (red dashed line) and the result presented in Sallum et al. (2015) and Wagner et al. (2018). A deeper dataset is required to redetect LkCa15 b with ZIMPOL, but this detection is feasible.
In the text

Fig. B.1 Performance comparison for Hα imaging with SPHERE/ZIMPOL using different filters (narrow and broad Hα) and reduction techniques (ADI and ASDI) based on the HD142527 dataset. In all panels the FPF obtained for HD142527 B is shown, as a function of the number of subtracted PCs used in PynPoint. On the right side of each panel we give the scale of the S/N to improve understanding of the plot and to compare different instrumental setups. We note again that this does not correspond to the classical σ notation. The gray regions indicates a confidence level for the detection of HD142527 B of at least 99.99995%, i.e., >5σ in case of Gaussian noise. Because of the applied corrections for small sample statistics, the border of the gray area does not correspond to an S/N of 5.

In the text

	Fig. C.1 Contrast curves calculated for the “first half” (blue dashed line) and the “very second frame” subsets (green solid line), and for the entire B_Ha dataset of HD142527 (orange dotted line).
In the text

	Fig. D.1 Lowest FPF images of HD135344 B (top panel: N_Ha filter; bottom panel: Cnt_Ha filter). The radius r of the inner mask and the number of subtracted PCs are given in each panel. The location of the tentative companion candidate is indicated by the arrow (see Sect. 4.1.2).
In the text

	Fig. E.1 B_Ha (blue circles) and Cnt_Ha (green crosses) mean count rates as a function of the image number in the observing sequence. The shaded region at the end represents the subset of frames that was chosen for the analysis.
In the text