Free Access
Issue
A&A
Volume 652, August 2021
Article Number A33
Number of page(s) 10
Section Numerical methods and codes
DOI https://doi.org/10.1051/0004-6361/202140285
Published online 05 August 2021

© ESO 2021

1. Introduction

Even though more than 4000 exoplanets have been discovered, only less than one percent of these have been directly imaged. Nevertheless, high-contrast imaging surveys like SHINE (Vigan et al. 2021), the GPI exoplanet survey (Nielsen et al. 2013, 2019), and the MMT-Clio survey (Heinze et al. 2010) have assessed the demographics of giant exoplanets and brown dwarfs. To observe exoplanets with high-contrast imaging, state-of-the-art technology is required to overcome the two main challenges: high angular resolution and high contrast.

To achieve high angular resolution, large telescopes (such as the 8.2-m UT telescopes of the Very Large Telescope (VLT) observatory) are required because of the Rayleigh criterion. To achieve high contrasts, the light of the star must be suppressed, as the diffraction halos of the target stars are typically several orders of magnitude brighter than an accompanying exoplanet. At the VLT, the Spectro-Polarimetric High-contrast Exoplanet Research (SPHERE; Beuzit et al. 2019) instrument is equipped with coronagraphs. Its extreme adaptive optics (AO) system helps to reduce speckle noise introduced by atmospheric turbulence and to further increase the final contrast. To obtain spectra, SPHERE is equipped with an integral field spectrograph (IFS) which can measure low-resolution spectra (Claudi et al. 2008). The wavelength range can be set to 0.95−1.35 μm with a resolution of R ≈ 50, or it can be set to 0.95−1.68 μm with a spectral resolution of R ≈ 30. In both cases, there are 39 wavelength channels split over the wavelength range (Wahhaj et al. 2019). The field of view (FOV) of SPHERE/IFS is approximately 1.73″ × 1.73″.

State-of-the-art technology alone is not sufficient – advanced data-reduction algorithms are necessary to increase the achievable contrast performance in data post-processing. Angular differential imaging (ADI) as described by Marois et al. (2006) (see Sect. 2.1) lays the foundation for most advanced techniques used today. Examples of such algorithms are LOCI (Lafreniere et al. 2007), TLOCI (Marois et al. 2013, 2014), MLOCI (Wahhaj et al. 2015), ANDROMEDA Mugnier et al. (2009), Cantalloube et al. (2015), PACO (Flasseur et al. 2018), and TRAP Samland et al. (2021). Lastly, there are attempts to use machine-learning techniques to improve exoplanet detections in high-contrast imaging (Gonzalez et al. 2018; Gebhard et al. 2020).

Amara & Quanz (2012) and Stolker et al. (2019) developed the Python package PynPoint for processing and analysis of high-contrast imaging data taken by, for example, SPHERE or NACO (NACO; Lenzen et al. 2003; Rousset et al. 2003). PynPoint uses principal component analysis (PCA; Amara & Quanz 2012; Soummer et al. 2012) and contains many functionalities to further analyse and process data (e.g. calculate contrast limits); a full list can be found in the official PynPoint documentation1. However, the basic ADI assumption of a quasi-static residual stellar point spread function (PSF) only holds partially because of changes in the AO performance, mostly due to non-common path aberrations (NCPA). This leads to an imperfect PSF model and subtraction can create residuals that resemble a planet. Spectral differential imaging (SDI; Racine et al. 1999) on the other hand uses images taken simultaneously and therefore with the same NCPAs (see Sect. 2.1). Recently, the PACO-ASDI algorithm (Flasseur et al. 2020) was developed, which uses IFS data to add SDI (in combination with ADI) to PACO. These latter authors found that adding SDI to ADI results in improved sensitivity of their detection maps.

Other pipelines similar to PynPoint have been developed for an efficient analysis of high-contrast imaging data. Examples of such pipelines are the Vortex Image Processing Package (VIP; Gonzalez et al. 2017), the SPHERE speckle calibration tool (SpeCal; Galicher et al. 2018), and the automated data-processing architecture for the Gemini Planet Imager Exoplanet Survey (Wang et al. 2018).

Spectral differential imaging and combinations of SDI and ADI for exoplanet imaging surveys are already in use. For example, Biller et al. (2013), Nielsen et al. (2013), and Wahhaj et al. (2013) conducted a survey using SDI in combination with ADI. For their SDI, they took two images simultaneously at two different wavelengths. They set the first wavelength channel at a methane absorption line (1.652 μm) so that T-type substellar objects would not appear in the image. The second wavelength was chosen outside the methane band (1.578 μm). The first image was then used as a PSF model to correct the second one. Samland et al. (2017) and Zurlo et al. (2016) applied different combinations of ADI and SDI to detect and characterise the exoplanets of 51 Eridani and HR 8799, respectively, using 39 wavelength channels to achieve a more accurate PSF model. Biller et al. (2013), Nielsen et al. (2013), Wahhaj et al. (2013), Samland et al. (2017), and Zurlo et al. (2016) all used ADI and SDI one after the other. In contrast, the algorithm of Flasseur et al. (2020) combines ADI simultaneously with SDI. This is possible, as both techniques exploit the diversity of sky rotation of a potential exoplanet companion (see Sect. 2.1). Christiaens et al. (2019) compared the performance of applying ADI followed by SDI with the output when alternating the order of both processing methods (i.e. SDI followed by ADI). These authors used different post-processing techniques on PDS 70, a K7-type star that hosts two gas giant planets (Keppler et al. 2018; Haffert et al. 2019) which are embedded in its transitional disc. Injecting fake planets allowed them to study which differential imaging technique best recovers specific features.

In this paper, we looked into different combinations of ADI and SDI to determine which yields the highest signal-to-noise ratio (S/N) detection of directly imaged planets. Even though various ADI/SDI combinations have been developed (e.g. Zurlo et al. 2016; Samland et al. 2017; Flasseur et al. 2018), they have not yet been extensively compared to each other. To do so, we added IFS support to PynPoint allowing wavelength-dependent differential imaging techniques as described in Sect. 2. In Sect. 3, we use this implementation to reanalyse data on β Pictoris, 51 Eridani, and HR 8799 and to calculate the S/N of their exoplanets. We then discuss the results in Sect. 4 in order to decipher which combination of ADI and SDI achieves the highest S/N. In Sects. 3 and 4, we determine the dependence of the differential imaging techniques on field rotation, wavelength range, and total integration time. In Sect. 5, we summarise our findings.

2. Differential imaging

2.1. ADI and SDI

The primary goal is to understand how ADI and SDI can be effectively combined to create a PSF model that provides the best possible contrast performance. To use ADI, raw images are taken with the instrument derotator turned off. This causes the FOV to rotate around the target star during observations. The stellar PSF is, to first order, quasi-static but a potential exoplanet rotates around the stellar position during the observation. The images can be used to build a model of the stellar PSF with negligible contribution from the companion, which can then be subtracted from each individual image. After de-rotating and stacking the images, a possible companion should become visible. In this paper, ADI is used in combination with PCA to optimise planet detection. From the input images, an orthogonal basis is calculated where the components are sorted with decreasing variance. The elements of this basis set are called principal components (PCs). A companion should not influence the variance significantly. Thus, the first few PCs can be used to create a PSF model while minimising the influence of planetary signals. For more information on PCA the reader is referred to Jee et al. (2007). To perform ADI together with PCA, we used PynPoint as described by Amara & Quanz (2012) and Stolker et al. (2019).

To perform SDI, PynPoint was extended to support SPHERE/IFS data. Like ADI, our SDI implementation uses displacements of a possible companion to model the stellar PSF while minimising the contribution of a potential exoplanet.

Because of the diffraction caused by the telescope, the flux from a distant point source will be spread out into its PSF. The shape of the PSF is approximately the same for all wavelengths, but with different full width at half maximum (FWHM) according to the Rayleigh criterion:

(1)

where 𝒟 is the telescope diameter and λ the observed wavelength. Scaling images that are obtained at different wavelengths with a factor S ∝ 1/λ corrects for the wavelength dependence of the PSF size. At the same time, the scaling causes a radial shift of a possible companion. The conditions are equivalent to ADI which uses a quasi-constant stellar PSF over time and azimuthal displacements of a possible companion. Therefore, a similar processing framework can be used for ADI and SDI.

To scale the images and thus align the stellar PSF, the inverse of β(λ) can be used. Because only the relative scaling between the images matters, the scaling factor can be defined as:

(2)

We choose the reference wavelength λref as the longest observed wavelength to prevent under-sampling of the residual star light caused by downscaling images.

(3)

Using this reference wavelength results in an up-scaling of all images with a shorter wavelength (S(λ) ≥ 1, ∀λ). The scaling of the images leads to a reduced FOV as illustrated in Fig. 1. SPHERE/IFS has a FOV of 1.73 × 1.73 arcsec2 and a spectral range of 0.95−1.35 μm (or 0.95−1.68 μm; Wahhaj et al. 2019). Scaling the images with the scaling factor S(λ) leads to a final FOV of 1.22 arcsec (or 1.00 arcsec). Everything outside this range is not covered by all wavelengths and therefore the stellar PSF cannot be modelled well. This FOV reduction does not occur if only ADI is applied.

thumbnail Fig. 1.

Schematic of SDI (including the FOV loss through scaling). First, the images (Aλ) are scaled inversely proportional to their wavelength (Bλ). A stellar PSF model (C) is then created using PCA analysis. Afterwards the stellar PSF model (C) is rescaled for each wavelength channel proportional to its wavelength (Dλ) to have the same scaling as the original images (Aλ). Therefore, the outer region of the small wavelength images is not covered by the scaled stellar PSF models (Dλ). This reduces the FOV by a factor R = λmin/λmax. The reduced FOV is marked by a red border, the companion is marked as a red dot.

We did not select specific wavelength channels for SDI and used all the available images. This confers the advantage that no assumptions about a possible companion have to be made and the results can be compared between different observations. Other papers removed certain wavelengths to either reduce self-subtraction or exploit known absorption bands (e.g. Thatte et al. 2007; Samland et al. 2017). If using the differential imaging techniques discussed in this paper to search for exoplanets, observation specific wavelength selection for SDI could increase the performance.

2.2. Displacements

The radial displacements Δr of a companion at separation r are caused by the scaling of the images according to Eq. (2). The smaller the displacements, the more planet signal is present in the PSF model. Large radial displacements are therefore desirable. The total radial displacements Δr between the shortest (λmin) and longest (λmax) observed wavelength at the separation of a companion rpl can be calculated with the following equation:

(4)

where Δλ = λmax − λmin is the observed wavelength range.

The azimuthal displacements Δrot of a companion used in ADI are due to the field rotation during the observation. Because the companion shifts its position over the circular segment with an angle of Δα and radius rpl, Δrot can be calculated with the following equation:

(5)

In our analysis, we did not set a required minimum displacement. We analysed a range of different field rotations and wavelength ranges (see Sect. 2.4) to find their influence on the performance of differential imaging techniques.

2.3. Differential imaging techniques

Both ADI and SDI use displacements of possible companions to create a PSF model excluding most of the signal from the companion. Therefore, these two techniques can be combined in three different ways: using first ADI then SDI; using first SDI then ADI; or using both techniques simultaneously. The implementation of all five techniques is described below.

ADI. First, all images Iλ, t within the dataset ℋ = {Iλ, t ∀λ, t} are split into wavelength subsets . Afterwards, using only the images of the subset , a PSF model is created using PCA with n ∈ ℕ PCs2. Each image is then subtracted by the PSF model . Lastly, all images Iλ, t are derotated and stacked using the median in order to create the final image IADI.

SDI. First, all images Iλ, t are scaled according to λ using a scaling factor S(λ) = λmax/λ. Afterwards, the dataset ℋ = {Iλ, t ∀λ, t} is split into time subsets . Using PCA with n PCs2, a PSF model is created using only images from the subset . Each image is then subtracted by the PSF model . Lastly, all images Iλ, t are scaled back (using a factor 1/S(λ)) before they are derotated and stacked using the median in order to create the final image ISDI.

CODI (combined differential imaging). First, all images Iλ, t are scaled according to λ using a scaling factor S(λ) = λmax/λ. Afterwards, a PSF model 𝒫λ, t is created using PCA with n PCs3 The PSF model is then subtracted from each image . Lastly, all images Iλ, t are scaled back (using a factor 1/S(λ)) before they are derotated and stacked using the median in order to create the final image ICODI.

SADI. First SDI is performed as described above but without derotating and stacking after the reduction. Afterwards it performs ADI as described above to create the final image ISADI. Both steps use n/2 PCs where n is an even number.

ASDI. First ADI is performed as described above but without derotating and stacking after the reduction. Afterwards it performs SDI as described above to create the final image IASDI. Both steps use n/2 PCs where n is an even number.

The differential imaging techniques and the framework to process SPHERE/IFS data are available within the PynPoint library. Furthermore, all images were pre-processed using EsoReflex (Freudling et al. 2013). This includes flat-field corrections, dark-frame subtraction, and wavelength calibration. For the PCA, we use the entire subset to create the PCA basis for the PSF model. The first few elements of this basis are than fitted to each image of the subset individually to create the PSF models.

Theoretically, every SDI step can be followed by an ADI step, or vice versa, therefore creating many more possible combinations. We focused on techniques using a maximum of one ADI step and one SDI step, one after the other (namely SADI and ASDI). Furthermore, each combination also allows different numbers of PCs to be used in each ADI or SDI step. We focused on using the same number of PCs in all ADI and SDI steps to be able to test a large range of PCs.

As ADI splits the dataset into its wavelength channels, it will conduct a PCA for each wavelength. For example, if ADI is performed on SPHERE/IFS data, a total of 39 PCAs will be conducted. As SDI splits the dataset into time components, it will conduct a PCA for each time instance of the dataset. If SDI is applied to a dataset with 100 time instances (each containing 39 wavelength frames), this leads to 100 PCAs. CODI does not split the dataset at all, and therefore only uses 1 PCA. To compare the results between different differential imaging techniques, we calculated the model completeness Ω by dividing the number of PCs used by the total number of PCs available:

(6)

Therefore, the residuals are zero if a model completeness of 1 is used because the image is fitted perfectly.

For SDI, SADI, and ASDI, the number of wavelength channels is the limiting factor and therefore PCmax = 39. For ADI, PCmax = Ntot, λ, where Ntot, λ is the number of images per wavelength taken during an observation. For CODI, PCmax = 39 × Ntot, λ because it uses all images at once for a single PCA. A summary can be found in Table 2.

2.4. Differential imaging tests

We test for dependence on total integration time, observed rotation, and/or observed wavelength range with the following:

Integration time. We split each dataset into subsets with different total integration times. Each dataset starts at the first image and contains a number of subsequent images until the given integration time of the subset is reached. Afterwards, each differential imaging technique is applied separately to each subset. The same PC number was chosen for all subsets. We compare the results to a slope proportional to . If the S/N values follow this curve, they behave in a photon-noise-limited manner.

Field rotation. We split each dataset into subsets with an equal number of images M. In the first subset, we select the first M images (1, 2, …, M); in the second, every second image is selected (1, 3, …, 2M); in the third set every third image is selected (1, 4, …, 3M), and so forth. Therefore, all subsets have the same integration time but a different total field rotation Δα. Because small rotations lead to a PSF overlap of a potential exoplanet, it is expected that small rotations lead to lower S/N. To achieve comparable results, we related the angular rotation to displacements expressed in FWHM using Eq. (5) with Δrot = FWHM.

Wavelength range. We split each dataset into subsets with an equal number of wavelength frames from all images. In the first subset, we select the first F frames from each image (1, 2, …, F); in the second, every second frame of each image is selected (1, 3, …, 2F); in the third set every third frame of each image is selected (1, 4, …, 3F), and so forth. Similar to the rotational overlap, radial overlap (after scaling the images) can lead to a lower S/N. To achieve comparable results, we related the radial offsets to displacements expressed in FWHM using Eq. (4) with Δr = FWHM.

3. Results

3.1. Datasets

We used three stars with known exoplanets to compare the five different differential imaging techniques described in Sect. 2.3: β Pictoris b (Lagrange et al. 2010), 51 Eridani b (Macintosh et al. 2015), and HR 8799 e (Marois et al. 2010). The information on their observation can be found in Table 1.

Table 1.

All observations were performed with SPHERE/IFS (PI: J. Beuzit).

Among the three targets, β Pictoris b has the most favourable contrast. The exoplanet has a contrast of 10.4 ± 0.5 mag in the J-band (Bonnefoy et al. 2013). In the selected dataset, β Pictoris b is at a separation of 0.242 ± 0.0025 arcsec (Lagrange et al. 2019). As a second target, we chose 51 Eridani. Its exoplanet 51 Eridani b has a contrast of 15 ± 1.6 mag in the J-band (Samland et al. 2017) and is at a separation of 0.455 ± 0.006 arcsec (Rosa et al. 2015). The last test system we selected is HR 8799. From the four planets of the system, only HR 8799 ‘e’ and ‘d’ are within the FOV of SPHERE/IFS. After the scaling necessary for SDI (see Sect. 2.1 and Fig. 1), only HR 8799 e remains within the reduced FOV. This exoplanet has a contrast of 13.01 ± 0.21 mag in the J-band (Zurlo et al. 2016) and is at a separation of 0.381 ± 0.007 arcsec (Apai et al. 2016).

3.2. Calculating the S/N

Our goal was to quantify detections of exoplanets to compare the differential imaging techniques. We did not aim to assess the detection limitations of our techniques, nor did we aim to directly compare our detections with those of other studies. Therefore, we chose S/N calculations to compare the different differential imaging techniques.

The signal used in the S/N was calculated using aperture photometry at the position of a companion candidate. All reference apertures for the noise calculation were chosen at the same angular separation as the signal aperture (Mawet et al. 2014). The size of the apertures was chosen to be 1 FWHM of the PSF at the median wavelength observed (λ = 1.15 for β Pictoris and λ = 1.31 for HR 8799 and 51 Eridani). To determine the position, a 2D Gaussian is fitted to the PSF of the companion. To calculate the noise, we first calculate the sum of the pixel values within each reference aperture. The noise is then defined as the root mean square of the reference aperture values. The signal is calculated as the sum of the pixel values within the signal aperture. To reduce the influence of the planetary signal within the noise calculation, the apertures adjacent to the planetary signal were not used. The noise calculation is adjusted for small sample statistics (Mawet et al. 2014).

3.3. Differential imaging tests

For comparability between differential imaging techniques, the number of PCs used is given in model completeness Ω (see Eq. (6)). The values for PCmax can be seen in Table 2.

Table 2.

Maximum number of PCs, PCmax, that can be selected for PCA.

The model completeness Ω was limited to 0.6 (0.25 for 51 Eridani) because all differential imaging techniques showed a clearly reduced S/N for larger model-completion values. The resulting S/N for β Pictoris b, 51 Eridani b, and HR 8799 e depending on model completeness can be seen in Fig. 2. In all three targets, CODI achieves the highest S/N. SADI and ASDI achieve comparable results with the exception of SADI applied to β Pictoris b which results in a lower S/N than the other three techniques. Furthermore, the differential imaging techniques SADI, ASDI, and CODI always achieve higher S/N than ADI. SDI alone on the other hand always performs similarly to or worse than ADI. Regarding model completeness, SADI and ASDI show the highest S/N for Ω < 0.2 and CODI for Ω ⪅ 0.1.

thumbnail Fig. 2.

S/N as a function of model completeness using different differential imaging techniques. Top: β Pictoris b; middle: 51 Eridani b; bottom: HR 8799 e.

To quantify the residual speckle noise over the whole FOV, we plotted the noise at different angular separations (see Fig. 3). The noise is given as the root mean square of the pixel values within the reference aperture after PCA, averaged over all reference apertures at the same angular separation (see Sect. 3.2). In all three targets, CODI achieves the lowest noise over the entire angular separation range. The noise values of SADI and ASDI are similar for all three targets and over the entire range of tested angular separations. ADI achieves lower noise than SADI and ASDI at small angular separations (< 0.2 arcsec) in β Pictoris, at large angular separations (> 0.2 arcsec) in HR 8799, and over the entire range in 51 Eridani. Even though it has a considerably lower S/N for all three targets, SDI has a noise level comparable to that of CODI in 51 Eridani and that of SADI/ASDI in HR 8799.

thumbnail Fig. 3.

Noise as a function of angular separation using different differential imaging techniques. Top: β Pictoris b; middle: 51 Eridani b; bottom: HR 8799 e. The colours and symbols are the same as in Fig. 2. The noise calculation is described in Sect. 3.2.

For the tests on integration time, field rotation, and wavelength range, various subsets with fewer images and/or wavelength frames than the full dataset were used. Therefore, the optimal number of PCs could not be determined from Fig. 2. For all our tests, we chose the PC number that achieves the highest S/N over the whole range of tested integration times. The corresponding PC numbers can be seen in Table 3.

Table 3.

Number of principal components n that overall achieve the highest S/N during the integration time test.

The S/N for different integration times of β Pictoris b and 51 Eridani b can be seen in Fig. 4. Only S/N values above 0.5 are displayed. The tests showed that for integration times t > 1500 s, CODI, SADI, ASDI, and ADI approximately follow the photon noise curve with no significant deviation measured. In this case, the tests do not show significant differences in S/N for SADI, ASDI, and CODI.

thumbnail Fig. 4.

Dependence of S/N on time after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The grey line indicates the slope corresponding to photon-noise-limited behaviour.

The S/N for different field rotation of β Pictoris b and 51 Eridani b can be seen in Fig. 5. The same PC numbers as for the integration time tests were used (see Table 3). Each subset of β Pictoris contained 76 images and each subset of 51 Eridani contained 15 images. In both targets, a strong increase in S/N with increasing field rotation can be seen for Δα < 1 FWHM. After that, the increase begins to stagnate and after Δα ⪆ 1.5 FWHM no significant increase was measured.

thumbnail Fig. 5.

Dependence of S/N on field rotation after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The vertical grey lines indicate displacements of 1 and 2 FWHM according to Eq. (5).

The S/N for different wavelength ranges of β Pictoris b can be seen in Fig. 6. The same PC numbers as for the integration time tests were used (see Table 3).

thumbnail Fig. 6.

Dependence of S/N on wavelength range of β Pictoris b after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The vertical grey lines indicate displacements of 2 and 4 FWHMs according to Eq. (4).

Because the datasets were split into smaller subsets, none of the three dependence tests (integration time, field rotation, and wavelength range) for HR 8799 e achieved S/N of over 5σ. Therefore, no clear relation could be determined and the dependence analysis of HR 8799 e was not included in this paper.

4. Discussion

4.1. Differential imaging techniques

SADI, ASDI, and CODI all achieve higher S/N compared to using only ADI for the three test targets β Pictoris b, 51 Eridani b, and HR 8799 e. SDI alone does not result in S/N over 5σ in any of the three targets.

We compare our results with the work of Christiaens et al. (2019) on PDS 70. Of the three targets we tested, we focus on β Pictoris b because its angular separation is closest to PDS 70 b. Christiaens et al. (2019) subsequently used ADI (PCA on concentric annuli of 2 FWHM) and SDI (full-frame PCA). Additionally, these authors only used frames between which the companion would have a rotational offset of at least 1 FWHM for their PSF library. In their analysis, only PCA–ASDI (similar to our ASDI) was able to recover a point-like feature; PCA–SADI (similar to our SADI) was not. The authors did not apply ADI and SDI simultaneously. In our analysis of β Pictoris b, SADI achieved a higher S/N compared to ASDI. In contrast to our analysis, Christiaens et al. (2019) used the VLT instrument SINFONI (Eisenhauer et al. 2003; Bonnet et al. 2004) which observes in the near-infrared (1.1−2.45 μm; similar to SPHERE/IFS) but has a higher resolution (approximately 2000 in the J-band). Furthermore, PDS 70 b has a smaller angular separation (0.193 ± 0.0049 arcsec) than β Pictoris b (0.249 ± 0.001 arcsec). Because the quality of SDI strongly depends on the radial displacement (see Eq. (4)), the different separations can be a dominant factor in the image quality. Furthermore, SINFONI only achieves Strehl ratios of approximately 30% (whereas SPHERE can go to above 80% in the J band). This results in a larger FWHM in the SINFONI data with respect to the SPHERE data and therefore to smaller S/N at small separations. The observed angle variation strongly influences ADI. Their observations of PDS 70 b have a field rotation of 99.8° while the observations of β Pictoris have 41.1°. The comparison indicates that the first processing step (ADI in ASDI; SDI in SADI) has a stronger influence on the final images, but this may be limited to this specific case and may not be a general result.

For 51 Eridani b and HR 8799 e, no clear difference between the differential imaging techniques SADI, ASDI, and CODI could be found. Samland et al. (2017) analysed 51 Eridani b using classical SDI followed by ADI (applied using various algorithms). These authors selected reference wavelength channels to improve their SDI performance and achieved a S/N of close to 20, while our SADI achieves a S/N of roughly 10. The improvement could be due to the longer integration time used (≈16384 s) compared to our dataset (≈4928 s). Assuming photon-noise-limit behaviour, our dataset would achieve an S/N of roughly 18 and thus be comparable to the results of Samland et al. (2017) without requiring target-specific adjustments. On the other hand, the difference could also be influenced by their in-depth analysis of the different wavelength channels and subsequent improvement of their SDI step.

For all three test objects, CODI achieves slightly higher S/N than SADI or ASDI and requires the lowest model completeness. This can be explained by the higher resolution of the model completeness. The number of PCs used in SADI and ASDI can only be optimised up to ΔΩ = 1/PCmax = 0.026. The resolution of the model completeness is limited by the maximum number of PCs (see Table 2). The PC of CODI can be adjusted up to ΔΩ = 1/(39 × Ntot, λ) and in the case of 51 Eridani, this leads to ΔΩ = 10−4.

Similar to our CODI, Flasseur et al. (2020) developed PACO–ASDI which exploits the displacement caused by field rotation and wavelength-dependent scaling. In contrast to our work, PACO–ASDI does not use PCA but the PACO algorithm (Flasseur et al. 2018). These latter authors found that using angular and radial displacements simultaneously leads to a significant increase in S/N. This indicates that our results could be generalised to other ADI/SDI-based techniques. Flasseur et al. (2020) did not investigate other combinations of ADI and SDI (like SADI and ASDI) and did not test the dependence of their technique on integration time and wavelength range. They did test the dependence on separation, which can be related to field rotation. Even though field rotation and angular separation both influence self-subtraction linearly (see Eq. (5)), extending the angular separation decreases the influence of the stellar PSF which can additionally increase the S/N.

Our results were obtained using PCA, but a similar investigation could be done for other techniques like PACO, LOCI, or machine-learning algorithms. All of these algorithms rely on the angular displacement of a planetary companion within the quasi-static PSF halo of the primary star throughout the observing sequence. An additional radial displacement, as used in SDI-based processing methods, can also be used by all the aforementioned algorithms to enhance the diversity of planetary signal and stellar residuals. Therefore, an investigation into the quality of combining ADI and SDI could help determine the best observation strategy. The same holds true for our analysis on integration time, field rotation, and wavelength range.

4.2. Noise curves

The noise curves in Fig. 3 show that CODI more effectively reduces the residual speckle noise than the other differential imaging techniques. At separations larger than 0.5 arcsec, a steep drop in noise is observed. CODI is therefore expected to perform well at angular separations larger than 0.5 arcsec because self-subtraction effects also weaken with increased angular separation due to larger azimuthal and radial displacements (see Eqs. (4) and (5)).

Comparing the noise curves for ADI to SADI/ASDI of the three test targets shows no general trend. ADI can achieve noise values orders of magnitude lower than SADI/ASDI and vise versa. Nevertheless, SADI and ASDI achieve higher S/N in all three targets which can be explained by self-subtraction effects. Ideally, the signal from the companion is only present in the very last PCs. However, some of its signal will appear in lower order PCs if the companion has insufficient displacements between the images. Therefore, if more PCs are used (which leads to lower noise), more of the signal from the companion is reduced. Selecting the optimal PC number means balancing these two effects against each other. SDI never achieves high S/N even though the noise level can be comparable with the other techniques, indicating that SDI cannot efficiently reduce noise before self-subtraction effects start to become dominant.

In most cases, the noise becomes smaller with increasing angular separation, which is expected because the stellar flux decreases. For ADI and CODI in the dataset of β Pictoris, and for SADI and ASDI in the dataset of 51 Eridani, the noise increases for larger angular separations. This increase can be explained by wind-driven halos, which are caused by imperfect corrections of the AO-system due to atmospheric changes attributable to jet streams (Cantalloube et al. 2018, 2020). Because these effects are asymmetric, ADI might not correct them well, leading to increased noise at larger angular separation.

4.3. Dependence tests

4.3.1. Integration time

In both β Pictoris b and 51 Eridani b (Fig. 4), the S/N increase rapidly before following the photon noise curve. We find no significant difference between SADI, ASDI, and CODI. Compared to ADI, all three differential imaging techniques combining ADI and SDI (SADI, ASDI and CODI) require less integration time than ADI to achieve higher S/N.

Similar to our test, Janson et al. (2007) analysed the integration-time dependence of ϵ Eridani b using a combination of ADI and SDI similar to ASDI. In agreement with our results, these authors found no significant deviation from the photon noise curve with increasing integration time. This dependence is important because some processing algorithms only use a subset of the whole dataset to model the PSF to reduce self-subtraction effects (e.g. Zurlo et al. 2016; Samland et al. 2017). Our results indicate that the strategic selection of frames can be used as long as the total integration time of the selected frames still exceeds 103 s.

4.3.2. Field rotation

The results for β Pictoris b and 51 Eridani b (Fig. 5) show an increased S/N for larger angles if the field rotation is smaller than 1 FWHM. Beyond that, the S/N increase for larger rotations is significantly smaller. This suggests that observations should have a field rotation of at least 1 FWHM at the angular separation of the potential companion. Similar to the integration time analysis, no significant difference between SADI, ASDI, and CODI could be found. Because each angle set can contain different images, the S/N varies slightly depending on the overall quality of the images in the set. This is to be expected and is not caused by the field rotation.

Meshkat et al. (2013) analysed the S/N dependence of ADI on angular separation. As shown in Eq. (5), the rotational offset depends linearly on both field rotation and angular separation. In accordance with our results, their analysis showed a strong increase in S/N for larger angular separations. For larger angular separations (which correspond to larger angular displacements), these latter authors do not find a significant stagnation of the S/N. Even though field rotation and angular separation both influence self-subtraction linearly (see Eq. (5)), varying the angular separation additionally leads to a significant change in the influence of the stellar PSF.

4.3.3. Wavelength range

For 51 Eridani b (Fig. 6), ADI is approximately constant while SADI, ASDI, and CODI improve with increasing wavelength range. This confirms that an increase in observed wavelength range can lead to higher S/N if SDI is used as part of differential imaging techniques.

For β Pictoris b (Fig. 6), the S/N of ADI increases with larger wavelength ranges even though the technique itself does not depend on the wavelength range. This increase in S/N is explained with the wavelength-dependent brightness of β Pictoris b. We compared the S/N increase of ADI with the spectra of β Pictoris b obtained by Chilcote et al. (2017). These latter authors measured an approximately linear increase for the wavelength range 0.95−1.35 μm, which explains the S/N increase that we measure for ADI. SADI and ASDI only achieved reasonable detections (S/N > 5) if the whole wavelength range was used, and their performance is worse than ADI. Analysing the processed images shows strong self-subtraction features at the position of β Pictoris b which explains the reduced S/N. Images processed with CODI did not show signs of significant self-subtraction and achieve a reasonably high S/N in this test. This result suggests that CODI can be more robust in cases where self-subtraction is dominant. A further investigation of the effects of injecting a fake companion could be used to study this phenomenon in more detail.

5. Conclusion

We tested the performance of ADI, SDI, SADI (using SDI, then ADI), ASDI (using ADI, then SDI), and CODI (using ADI and SDI simultaneously) on three different targets: β Pictoris b, 51 Eridani b, and HR 8799 e. β Pictoris b and 51 Eridani b were used to test the dependence of the techniques on integration time, field rotation, and observed wavelength range. The differential imaging techniques used in this paper are available within the PynPoint library.

For these three targets, using a combination of ADI and SDI consistently achieved better results than using ADI alone. SADI and ASDI performed similarly but target-specific differences can occur. In all tests, CODI achieved the highest S/N and constantly performed slightly better than SADI or ASDI.

We analysed the dependence of differential imaging techniques on integration time by splitting the datasets of our test targets into subsets of different total integration time. In this test, SADI, ASDI, and CODI all showed a photon-noise-limited behaviour for integration times of > 1500 s if all images were considered. Choosing frames selectively to reduce self-subtraction might decrease this time limit.

We analysed the dependence of differential imaging techniques on field rotation by splitting the datasets of our test targets into subsets with different observed rotation but with the same integration time. In this test, SADI, ASDI, and CODI all showed a strong increase of S/N up to displacements of 1 FWHM. After that, increasing the field rotation only led to a small increase in S/N.

We analysed the dependence of differential imaging techniques on the observed wavelength range by splitting the datasets of our test targets into subsets with different wavelength ranges. Our tests show a general increase in S/N with increasing wavelength range, but no clear relation could be found. Furthermore, the spectra of a given exoplanet influence the results because its brightness can vary substantially in different wavelength channels.

To further investigate these techniques, SADI and ASDI should be compared to CODI while using different numbers of PCs in SDI and ADI. In combination with analysing more known exoplanets and brown dwarfs, this would help to decipher whether or not our findings can be generalised to a wider selection of targets. Unfortunately, both tasks are immensely computationally intensive. Computationally optimising the techniques could help to reduce the intensity and would make a large-scale search for exoplanets more feasible. Furthermore, the extraction of spectra using CODI should be investigated further. This requires careful treatment of the different wavelengths because of cross talk between adjacent spectral channels (Greco & Brandt 2016).

Our analysis of the differential imaging techniques ADI, SDI, CODI, SADI, and ASDI shows that using SDI in combination with ADI leads to higher S/N for exoplanet detections than using only ADI or SDI. Overall, CODI achieved the highest S/N in all three targets. This can be explained with the higher model completeness resolution which enables this technique to more precisely model the stellar PSF. We suggest that future observations should first be evaluated with CODI if it is not feasible to apply all of them. The observed wavelength range should be maximised (e.g. for SPHERE/IFS use 0.95−1.65 μm) and the field rotation should cover at least 2 FWHM at the angular separation where exoplanets are expected. An integration time of 103 s is sufficient to reach photon noise behaviour.


2

n has the same value for all subsets.

3

We note that in contrast to all other differential imaging techniques, CODI uses only one PCA using n PCs.

Acknowledgments

To achieve the scientific results presented in this article we made use of the Python, especially the SciPy (Virtanen et al. 2020), NumPy (Harris et al. 2020), Matplotlib (Hunter 2007), emcee (Foreman-Mackey et al. 2013), scikit-image (Walt et al. 2014), scikit-learn (Pedregosa et al. 2018), photutils (Bradley et al. 2016), and astropy (Astropy Collaboration 2013, 2018) packages. Part of this work has been carried out within the framework of the National Centre of Competence in Research PlanetS supported by the Swiss National Science Foundation. SPQ acknowledges the financial support of the SNSF.

References

  1. Amara, A., & Quanz, S. P. 2012, MNRAS, 427, 948 [NASA ADS] [CrossRef] [Google Scholar]
  2. Apai, D., Kasper, M., Skemer, A., et al. 2016, ApJ, 820, 40 [NASA ADS] [CrossRef] [Google Scholar]
  3. Astropy Collaboration (Robitaille, T. P., et al.) 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  4. Astropy Collaboration (Price-Whelan, A. M., et al.) 2018, AJ, 156, 123 [Google Scholar]
  5. Beuzit, J.-L., Vigan, A., Mouillet, D., et al. 2019, A&A, 631, A155 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  6. Biller, B. A., Liu, M. C., Wahhaj, Z., et al. 2013, ApJ, 777, 160 [Google Scholar]
  7. Bonnefoy, M., Boccaletti, A., Lagrange, A.-M., et al. 2013, A&A, 555, A107 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  8. Bonnet, H., Abuter, R., Baker, A., et al. 2004, The Messenger, 117, 17 [NASA ADS] [Google Scholar]
  9. Bradley, L., Sipocz, B., Robitaille, T., et al. 2016, Astrophysics Source Code Library [record ascl:1609.011] [Google Scholar]
  10. Cantalloube, F., Mouillet, D., Mugnier, L. M., et al. 2015, A&A, 582, A89 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Cantalloube, F., Por, E. H., Dohlen, K., et al. 2018, A&A, 620, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Cantalloube, F., Farley, O. J. D., Milli, J., et al. 2020, A&A, 638, A98 [CrossRef] [EDP Sciences] [Google Scholar]
  13. Chilcote, J., Pueyo, L., Rosa, R. J. D., et al. 2017, AJ, 153, 182 [NASA ADS] [CrossRef] [Google Scholar]
  14. Christiaens, V., Casassus, S., Absil, O., et al. 2019, MNRAS, 486, 5819 [NASA ADS] [CrossRef] [Google Scholar]
  15. Claudi, R. U., Turatto, M., Gratton, R. G., et al. 2008, Int. Soc. Opt. Photon., 7014, 70143E [Google Scholar]
  16. Eisenhauer, F., Abuter, R., Bickert, K., et al. 2003, Int. Soc. Opt. Photon., 4841, 1548 [Google Scholar]
  17. Flasseur, O., Denis, L., Thiébaut, E., & Langlois, M. 2018, A&A, 618, A138 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Flasseur, O., Denis, L., Thiébaut, E., & Langlois, M. 2020, A&A, 637, A9 [CrossRef] [EDP Sciences] [Google Scholar]
  19. Foreman-Mackey, D., Hogg, D. W., Lang, D., & Goodman, J. 2013, PASP, 125, 306 [Google Scholar]
  20. Freudling, W., Romaniello, M., Bramich, D. M., et al. 2013, A&A, 559, A96 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Galicher, R., Boccaletti, A., Mesa, D., et al. 2018, A&A, 615, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  22. Gebhard, T. D., Bonse, M. J., Quanz, S. P., & Schölkopf, B. 2020, ArXiv e-prints [arXiv:2010.05591] [Google Scholar]
  23. Gonzalez, C. A. G., Wertz, O., Absil, O., et al. 2017, AJ, 154, 7 [Google Scholar]
  24. Gonzalez, C. A. G., Absil, O., & van Droogenbroeck, M. 2018, A&A, 613, A71 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Greco, J. P., & Brandt, T. D. 2016, ApJ, 833, 134 [NASA ADS] [CrossRef] [Google Scholar]
  26. Haffert, S. Y., Bohn, A. J., de Boer, J., et al. 2019, Nat. Astron., 3, 749 [Google Scholar]
  27. Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [CrossRef] [PubMed] [Google Scholar]
  28. Heinze, A. N., Hinz, P. M., Kenworthy, M., et al. 2010, ApJ, 714, 1570 [Google Scholar]
  29. Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
  30. Janson, M., Brandner, W., Henning, T., et al. 2007, AJ, 133, 2442 [NASA ADS] [CrossRef] [Google Scholar]
  31. Jee, M. J., Blakeslee, J. P., Sirianni, M., et al. 2007, PASP, 119, 1403 [Google Scholar]
  32. Keppler, M., Benisty, M., Müller, A., et al. 2018, A&A, 617, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Lafreniere, D., Marois, C., Doyon, R., Nadeau, D., & Artigau, E. 2007, ApJ, 660, 770 [NASA ADS] [CrossRef] [Google Scholar]
  34. Lagrange, A.-M., Bonnefoy, M., Chauvin, G., et al. 2010, Science, 329, 57 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  35. Lagrange, A.-M., Boccaletti, A., Langlois, M., et al. 2019, A&A, 621, L8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Lenzen, R., Hartung, M., Brandner, W., et al. 2003, Int. Soc. Opt. Photon., 4841, 944 [Google Scholar]
  37. Macintosh, B., Graham, J. R., Barman, T., et al. 2015, Science, 350, 64 [NASA ADS] [CrossRef] [Google Scholar]
  38. Marois, C., Lafreniere, D., Doyon, R., Macintosh, B., & Nadeau, D. 2006, ApJ, 641, 556 [Google Scholar]
  39. Marois, C., Zuckerman, B., Konopacky, Q. M., Macintosh, B., & Barman, T. 2010, Nature, 468, 1080 [Google Scholar]
  40. Marois, C., Correia, C., Véran, J.-P., & Currie, T. 2013, Proc. Int. Astron. Union, 8, 48 [CrossRef] [Google Scholar]
  41. Marois, C., Correia, C., Galicher, R., et al. 2014, Int. Soc. Opt. Photon., 9148, 91480U [Google Scholar]
  42. Mawet, D., Milli, J., Wahhaj, Z., et al. 2014, ApJ, 792, 97 [NASA ADS] [CrossRef] [Google Scholar]
  43. Meshkat, T., Kenworthy, M., Quanz, S. P., & Amara, A. 2013, ApJ, 780, 17 [NASA ADS] [CrossRef] [Google Scholar]
  44. Mugnier, L. M., Cornia, A., Sauvage, J.-F., et al. 2009, J. Opt. Soc. Am. A, 26, 1326 [NASA ADS] [CrossRef] [Google Scholar]
  45. Nielsen, E. L., Liu, M. C., Wahhaj, Z., et al. 2013, ApJ, 776, 4 [NASA ADS] [CrossRef] [Google Scholar]
  46. Nielsen, E. L., Rosa, R. J. D., Macintosh, B., et al. 2019, AJ, 158, 13 [NASA ADS] [CrossRef] [Google Scholar]
  47. Pedregosa, F., Varoquaux, G., Gramfort, A., et al. 2018, ArXiv e-prints [arXiv:1201.0490] [Google Scholar]
  48. Racine, R., Walker, G. A. H., Nadeau, D., Doyon, R., & Marois, C. 1999, PASP, 111, 587 [Google Scholar]
  49. Rosa, R. J. D., Nielsen, E. L., Blunt, S. C., et al. 2015, ApJ, 814, L3 [NASA ADS] [CrossRef] [Google Scholar]
  50. Rousset, G., Lacombe, F., Puget, P., et al. 2003, Int. Soc. Opt. Photon., 4839, 140 [Google Scholar]
  51. Samland, M., Mollière, P., Bonnefoy, M., et al. 2017, A&A, 603, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  52. Samland, M., Bouwman, J., Hogg, D. W., et al. 2021, A&A, 646, A24 [CrossRef] [EDP Sciences] [Google Scholar]
  53. Soummer, R., Pueyo, L., & Larkin, J. 2012, ApJ, 755, L28 [NASA ADS] [CrossRef] [Google Scholar]
  54. Stolker, T., Bonse, M. J., Quanz, S. P., et al. 2019, A&A, 621, A59 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  55. Thatte, N., Abuter, R., Tecza, M., et al. 2007, MNRAS, 378, 1229 [NASA ADS] [CrossRef] [Google Scholar]
  56. Vigan, A., Fontanive, C., Meyer, M., et al. 2021, A&A, 651, A72 [EDP Sciences] [Google Scholar]
  57. Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nat. Meth., 17, 261 [Google Scholar]
  58. Wahhaj, Z., Liu, M. C., Nielsen, E. L., et al. 2013, ApJ, 773, 179 [NASA ADS] [CrossRef] [Google Scholar]
  59. Wahhaj, Z., Cieza, L. A., Mawet, D., et al. 2015, A&A, 581, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  60. Wahhaj, Z., Milli, J., Rodler, F., Girard, J., & Vigan, A. 2019, SPHERE User Manual [Google Scholar]
  61. Walt, S. V. D., Schönberger, J. L., Nunez-Iglesias, J., et al. 2014, PeerJ, 2, e453 [CrossRef] [PubMed] [Google Scholar]
  62. Wang, J. J., Perrin, M. D., Savransky, D., et al. 2018, J. Astron. Telesc. Instrum. Syst., 4, 018002 [NASA ADS] [Google Scholar]
  63. Zurlo, A., Vigan, A., Galicher, R., et al. 2016, A&A, 587, A57 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

All Tables

Table 1.

All observations were performed with SPHERE/IFS (PI: J. Beuzit).

Table 2.

Maximum number of PCs, PCmax, that can be selected for PCA.

Table 3.

Number of principal components n that overall achieve the highest S/N during the integration time test.

All Figures

thumbnail Fig. 1.

Schematic of SDI (including the FOV loss through scaling). First, the images (Aλ) are scaled inversely proportional to their wavelength (Bλ). A stellar PSF model (C) is then created using PCA analysis. Afterwards the stellar PSF model (C) is rescaled for each wavelength channel proportional to its wavelength (Dλ) to have the same scaling as the original images (Aλ). Therefore, the outer region of the small wavelength images is not covered by the scaled stellar PSF models (Dλ). This reduces the FOV by a factor R = λmin/λmax. The reduced FOV is marked by a red border, the companion is marked as a red dot.

In the text
thumbnail Fig. 2.

S/N as a function of model completeness using different differential imaging techniques. Top: β Pictoris b; middle: 51 Eridani b; bottom: HR 8799 e.

In the text
thumbnail Fig. 3.

Noise as a function of angular separation using different differential imaging techniques. Top: β Pictoris b; middle: 51 Eridani b; bottom: HR 8799 e. The colours and symbols are the same as in Fig. 2. The noise calculation is described in Sect. 3.2.

In the text
thumbnail Fig. 4.

Dependence of S/N on time after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The grey line indicates the slope corresponding to photon-noise-limited behaviour.

In the text
thumbnail Fig. 5.

Dependence of S/N on field rotation after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The vertical grey lines indicate displacements of 1 and 2 FWHM according to Eq. (5).

In the text
thumbnail Fig. 6.

Dependence of S/N on wavelength range of β Pictoris b after processing all images using different differential imaging techniques. Top: β Pictoris b; bottom: 51 Eridani b. The colours and symbols are the same as in Fig. 2. The vertical grey lines indicate displacements of 2 and 4 FWHMs according to Eq. (4).

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.