The SPHERE infrared survey for exoplanets (SHINE) -- II. Observations, Data reduction and analysis Detection performances and early-results

Over the past decades, direct imaging has confirmed the existence of substellar companions (exoplanets or brown dwarfs) on wide orbits (>10 au) from their host stars. To understand their formation and evolution mechanisms, we have initiated in 2015 the SPHERE infrared survey for exoplanets (SHINE), a systematic direct imaging survey of young, nearby stars to explore their demographics.} {We aim to detect and characterize the population of giant planets and brown dwarfs beyond the snow line around young, nearby stars. Combined with the survey completeness, our observations offer the opportunity to constrain the statistical properties (occurrence, mass and orbital distributions, dependency on the stellar mass) of these young giant planets.} {In this study, we present the observing and data analysis strategy, the ranking process of the detected candidates, and the survey performances for a subsample of 150 stars, which are representative of the full SHINE sample. The observations were conducted in an homogeneous way from February 2015 to February 2017 with the dedicated ground-based VLT/SPHERE instrument equipped with the IFS integral field spectrograph and the IRDIS dual-band imager covering a spectral range between 0.9 and 2.3 $\mu$m. We used coronographic, angular and spectral differential imaging techniques to reach the best detection performances for this study down to the planetary mass regime.}


Introduction and context
The discovery of the first brown dwarf companion Gl 229 B benefited from the combined technological innovation of infrared detectors and high contrast techniques (Nakajima et al. 1995).Following this discovery, the first generation of dedicated planet imagers on 10-m class telescopes (NaCo at VLT, NIRC2 at Keck, NICI at Gemini) conducted systematic surveys of young and nearby stars.They led the first direct detections of planetary mass companions in the early 2000's.These companions were detected at distances larger than several hundreds astronomical units (au) and/or with a mass ratio not much smaller than a tenth with respect to their host primaries (except for brown dwarf primaries), giving hints of formation via gravo-turbulent fragmentation (Hennebelle & Chabrier 2011) or gravitational disk instability (Boss 1997).Thanks to the improvement of direct imaging observation and data analysis techniques with ground-based adaptive optics systems (AO) or space telescopes, a few planetary mass objects and low-mass brown dwarfs have been detected since the first detection by Chauvin et al. (2004).Moreover, these developments, enabled the discoveries of giant planets within 100 au around young, nearby, and dusty early type stars like HR8799 b,c,d,e (Marois et al. 2008(Marois et al. , 2010)), β Pictoris b (Lagrange et al. 2009), and more recently HD 95086 b (Rameau et al. 2013a), and GJ 504 b (Kuzuhara et al. 2013a).
Direct imaging is the only viable technique to probe for planets at large separations with single epoch observations, but detecting them requires to overcome the difficulties caused by the angular proximity and the high contrast involved.With improved instruments and data reduction techniques, we are currently initiating the characterization of the giant planet population at wide orbits, between typically 10-100 au.More than a decade of direct imaging surveys targeting several hundred young, nearby stars have lead to the discovery of approximately a dozen sub-stellar companions excluding the companions to brown dwarfs located in the star vicinity (within 100 au).Despite the relatively small number of discoveries compared with other techniques, such as radial velocity and transit, each new imaged giant planet has provided unique clues on the formation, evolution and physics of young massive Jupiters.
Early surveys performed using the first generation of planet imagers enabled to conduct systematic surveys relatively modest in size sampling each less than hundreds of young nearby stars (Chauvin et al. 2018).Various strategies were followed for the target selection of these surveys: i/ complete census of given associations (Chauvin et al. 2003;Masciadri et al. 2005;Kasper et al. 2007;Chauvin et al. 2010), ii/ selection of young, intermediate-mass stars (Janson et al. 2011;Rameau et al. 2013b;Nielsen et al. 2013), iii/ or very low-mass stars (Delorme et al. 2012;Chauvin et al. 2012;Bowler et al. 2015;Lannier et al. 2016), iv/ application of figures of merit considering detection rate with toy models of planet population (Biller et al. 2013).Numerous direct imaging surveys to detect giant planet companions have reported null-detections (Masciadri et al. 2005;Kasper et al. 2007;Chauvin et al. 2010;Biller et al. 2007;Lafrenière et al. 2007;Ehrenreich et al. 2010;Chauvin et al. 2010;Janson et al. 2011;Delorme et al. 2012), but this allowed to set upper limits to the occurrence of giant planets.
Ongoing surveys target several hundreds stars, with the largest surveys to date due to be completed in the next years.SpHere INfrared survey for Exoplanets (SHINE, Chauvin et al. 2017b), Gemini GPIES (Macintosh et al. 2015) and SCExAO (Jovanovic et al. 2016) now combine dedicated extreme adaptive optics systems with coronagraphic and both angular and spectral differential image processing.With enhanced detection performances, the objective is to significantly increase the number of imaged planets to characterize, but also to provide better statistical constraints on the occurrence and the characteristics of exoplanets at wide orbits (>10 au).This should give us a more global picture of planetary systems architecture at all orbits to improve our understanding of planetary formation and evolution mechanisms.Early-results already confirm the gain in terms of detection performances compared with the first generation of planet imagers (Nielsen et al. 2019;Vigan et al. 2020).
The near-infrared wavelengths (H and K-bands) have been used intensively.They are a good compromise between lowbackground noise, high angular resolution and good Strehl correction.However, thermal imaging at L or M-bands has been very competitive in terms of detection performances as the planet-star contrast and the Strehl correction are more favorable in those wavelengths despite an increased thermal background and larger inner working angle.For instance, SPHERE could be in some cases less sensitive than NaCo for the detection of giant planets around young, nearby M dwarfs at typical separations larger than 500-1000 mas.
Both SHINE and GPIES surveys have recently discovered three new exoplanets (Macintosh et al. 2015;Chauvin et al. 2017b;Keppler et al. 2018), and a few additional higher mass brown dwarfs (Konopacky et al. 2016;Cheetham et al. 2018).Smaller surveys using SPHERE and GPI also discovered several substellar companions (Milli et al. 2017a;Wagner et al. 2020;Bohn et al. 2020).They offer unprecedented detection, astrometric and spectrophotometric capabilities which allow us to characterize fainter and closer giant planets, such as the recent discovery of 51 Eri b (2 MJup at 14 au, T5-type, of age 20 Myr; Macintosh et al. 2015;Samland et al. 2017), HIP 65426 b, a young, warm, and dusty L5-L7 massive jovian planet located at about 92 au from its host star (Chauvin et al. 2017b), and the young solar analogue PDS 70 now known to actually host two planets PDS 70 b discovered during the SHINE campaign (Keppler et al. 2018;Müller et al. 2018) and PDS 70 c by MUSE (Haffert et al. 2018).Such surveys also provide key spectral and orbital characterisation data for known exoplanets (e.g.De Rosa et al. 2016;Samland et al. 2017;Chauvin et al. 2018;Wang et al. 2018;Cheetham et al. 2019;Lagrange et al. 2019;Maire et al. 2019).Despite these new discoveries, SHINE and GPIES have yielded significantly fewer exoplanet detections than predicted using extrapolations of radial velocity planet populations to larger semi-major axes (e.g.Cumming et al. 2008).This results in setting strong statistical constraints on the distribution of giant exoplanets at separations >10 au from their stars, as well as sub-stellar companions to young stars (Nielsen et al. 2019;Vigan et al. 2020).As these systems are young (<100 Myr), and thus closer to their epoch of formation than for instance radial velocity planets (typically 1-10 Gyr), statistical analysis of large direct-imaging surveys can provide hints into the potential formation mechanism responsible for producing giant exoplanets orbiting at wide orbits.
This paper is part of a series of three papers describing early results obtained from analysis of a subset of the SHINE data.This paper describes the observations, data quality, analysis methods, and presents results in terms of detections and upper limits for this subset of SHINE data.Two associated papers describe the general characteristics of the survey sample (Desidera et al. submitted) and the statistical analysis and a discussion of the implication for formation scenarios (Vigan et al. 2020).
Following this introduction, we present in section 2 the observations that were performed to underpin this work.In Sect.3, we describe the reduction and calibration of the data sets.Section 4 is dedicated to the main results of the survey, including the detection limits of the survey.Finally, we summarize, in Sect 5, the characterization of newly discovered and known substellar companions along with a description of the circumstellar disks detected within the survey data.The conclusions and prospects are drawn in section 7. The relative astrometry and photometry of companion candidates from the sample is listed in Annex A.

The SHINE survey
The SpHere INfrared survey for Exoplanets (SHINE: Chauvin et al. 2017a) has been designed for 200 telescope nights, allocated in visitor mode, using the SPHERE consortium guaranteed time.SHINE has been designed by the SPHERE consortium to: (i) identify and characterized new planetary and brown dwarf companions; (ii) study the architecture of planetary systems (multiplicity and dynamical interactions); (iii) investigate the link between the presence of planets and disks (in synergy with the GTO program aimed at disk characterization); (iv) determine the frequency of giant planets beyond 10 au; (v) investigate the impact of stellar mass (and even age if possible) on the frequency and characteristics of planetary companions over the range 0.5 to 3.0 M .
The SHINE survey started in February 2015 and will be completed in mid-2021, with observations of 500 targets out of a larger sample of 800 nearby young stars aiming at searching for new sub-stellar companions.The sample is oversized with respect to the available telescope time by a factor of approximately two, on the basis of the adopted observing strategy, which relies on observations surrounding meridian passage in order to achieve the maximum field-of-view (FoV) rotation for optimal angular differential imaging.This requires some flexibility in the target list in order to optimize the scheduling.
The general design of the survey, the sample selection and the simulations performed for building it, the parameters of the individual targets used in this early statistical analysis, and the general properties of the sample used in this series of papers are described in detail in Desidera et al. (submitted).We describe in this paper early results obtained from the analysis of about a third of the large SHINE survey sample, considering only those targets whose first observations were done before Feb 2017.Several second epochs observations were also performed between 2016 and 2019 to discriminate companions from background stars on the selected sample described in this paper.This selected sample dataset including 150 targets, is already large enough for reviewing the survey efficiency, to discuss the incidence of massive planets at a separation > 10 au, and to have new indications about the formation scenarios for giant planets.

Observations setup
All the observations were performed with the Spectro-Polarimeter High-contrast Exoplanet REsearch (SPHERE; Beuzit et al. 2019) combining its SAXO extreme adaptive optics system (Fusco et al. 2006;Sauvage et al. 2016;Fusco et al. 2014) and its apodized pupil Lyot coronagraphs (Boccaletti et al. 2008;Carbillet et al. 2011;Guerri et al. 2011).Observations were acquired in either IRDIFS or IRDIFS-EXT mode, i.e. with both NIR sub-systems, IFS (Claudi et al. 2008) andIRDIS (Dohlen et al. 2008), observing in parallel (Zurlo et al. 2014;Mesa et al. 2015).The IFS covers a 1.7 ×1.7 FoV and IRDIS covers a more or less circular, unvignetted FoV of diameter ∼9 .The APLC_YJHs and APLC_Ks coronagraphic configurations were used for the IRDIFS and IRDIFS-EXT observations respectively, and in IRDIS all first epoch observations were performed with the DB-H23 and DB-K12 dual-band filter pairs (Vigan et al. 2010).Due to the presence of known companions and/or the detection of (new) candidate companions, some targets were observed multiple times for astrometric monitoring.In addition to the initial selected filters, follow-up observations were performed in different configurations with IRDIS, for example with the broad-band BB-H and dual-band DB-J23 filters, as listed in the filter column of Table 1.This resulted in a varying number of observations for each target.

Observations Optimized planning: SPOT
Given the large number of targets, each associated with a priority and an urgency (how soon the observation has to be made), and the various observing constrains, including those connected to Angular Differential Imaging (ADI: Marois et al. 2006) observations, we built a dedicated tool, SPOT, to deliver an optimised scheduling of the observations, both on long and short terms.SPOT is based on simulated annealing.It is described in in details in Lagrange et al. (2016).
In brief, SPOT uses as inputs the calendar of allocated observing nights, the list of targets available for e.g. the whole semester, together with their associated instrumental set up (that is associated with specific overheads) and any specific scheduling constrains (needed for second epoch observations), the targets coordinates and magnitude, the minimal coronagraphic exposition duration, the maximum air mass, the maximum proportion (in time) of the observation that is allowed either before or after the meridian crossing, the minimum exposure time, the amount of FoV rotation during the observation.It must be noted that the amount of FoV rotation depends on the target coordinates and on the actual time and duration of observations.As it also drives the exposure time, we also set a maximum exposition time, to avoid extending too much the duration of the exposure.
Classical, additional scheduling constraints are applied on all targets: minimal angular distance to the Moon, avoidance of zenith observations.Poor atmospheric conditions may also be taken into account through pointing restrictions that usually depend on the wind direction and speed and through magnitude restrictions in case of non photometric conditions.Finally, SPOT also takes into account the various overheads and the need for astrometric and spectro-photometric calibrations.The astrometric calibrations were to be observed as much as possible in the first night of each run, and the spectro-photometric in the following night.After optimisation SPOT returns a schedule, and produces the Observing Blocks (OBs) that can be automatically transferred to P2 (the ESO observing preparation tool).
In practice, we generally request the coronagraphic data to be obtained while the target is crossing the meridian at least 15 min before and after the meridian passage, and with at least 30 degrees of FoV rotation, unless it required more than 7200 s.The minimum exposure time was set to 3600 s.The scheduling of the targets was always optimal; the efficiency of the night could be as high as 100% when typically 2-3 time more targets were available than actually those scheduled.Yet, when not enough targets were available as inputs, short (typically 30 min) holes could be present in the schedules during the nights.These holes were used to observe fillers requiring short exposure times (binaries, stars with RV trends, astrometric and spectrophotometric standards).By comparing with other surveys executed in service such as BEAST (Janson et al., in preparation), we conclude that use of SPOT allowed an increase of more than 30% in the field rotation angles for identical observing time with respect to the service scheduling routinely used at ESO.

Observing conditions and Data quality
The SPHERE Instruments are fed by an extreme adaptive optics (AO) system called SAXO.It delivers a very high Strehl ratio, which reaches above 90% in H-band for very good observing conditions by correcting both perturbations induced by the atmospheric turbulence and from the internal aberrations of the instrument itself.A comprehensive description of the SAXO design can be found in Fusco et al. (2016) and Beuzit et al. (2019).We derived here the overall statistical AO telemetry data from SAXO and from the ESO MASS-DIMM measurements for the survey observations and relate these parameters to the Strehl ratio, raw contrast and processed contrast in order to evaluate the performance constraints from these observations.The telemetry data points are spread over 130 different nights and cover more than 200 different observations.The AO telemetry data (available for a large number of our survey observations) includes estimates from the real-time computer (called SPARTA) several quantities of interest that could be related to the final performances (hereafter SPARTA data) such as: the Strehl ratio and additional atmospheric parameters including the seeing and the coherence time.The Strehl is defined at 1.6 µm, while the seeing and coherence time are defined at 500 nm.
These quantities are also connected here with the brightness of the target (V and H magnitudes), retrieved from the Simbad database, and used as a proxy in V-Band for the number of photons received by the WFS1 .The distribution in magnitude for the sample considered here (for targets with telemetry measurements) is shown in Fig. 1 (Bottom Right).It illustrates the fact that our target selection as described in Desidera et al. (submitted) was based on setting a magnitude limit (V<12.5) in order to guarantee good AO performances.Our sample selection criteria also favors red targets leading to the faintest H-band observations to be below 8 magnitude.
The overall good AO performance for the selected range of target brightness is confirmed by the Strehl measurements shown in Fig. 1 (bottom left).This Figure shows the distribution of the Strehl ratio as a function of the star magnitude in V-band where the median Strehl is greater than 70% averaging around 80% in moderate seeing conditions.This Figure also shows a decrease in the Strehl beyond magnitude 7. The blue data points, representing coherence time smaller than 3 ms, highlight the impact of short coherence time on the performances with Strehl ratio below 80%.It is suggested in Milli et al. (2017b) that the Strehl or its estimation could be affected by the small number of photons reaching the WFS (for 7.5 -9 magnitude guide stars) while the AO loop is operating at full speed, wich could explain the strehl decrease around magnitude 8.For fainter stars, a slower AO loop frequency of 600 Hz is used in order to still collect a sufficient number of photons per AO cycle for high-order wavefront determination and correction.
This hypothesis is supported further by the greater decrease in Strehl for the low coherence time data points which explain most of the performance inflection around V=8 magnitude.The large scatter in the Strehl distributions indicates that the seeing is also, as expected, a parameter influencing the Strehl among others such as τ 0 and the star magnitude.Despite this large scatter the (Top Left) plot from the same Figure shows a linear decrease trend of the Strehl with the seeing.This decrease is on average 1.8% Strehl for an increase in seeing of 0.1" which is similar to the Strehl behavior versus seeing obtained by Milli et al. (2017b) when correction is applied according to the SPARTA versus DIMM measurements dependency.While the SPARTA seeing provides an estimate which is closer to the image quality in the science frames, it has been shown that the SPARTA seeing estimations are smaller than the ESO DIMM measurements, a fact that may be accounted for the turbulence outer scale.The overall impact of the performance limitation from the seeing on our survey is clearly important as illustrated by the large number of observations occurring for seeing greater than 0.8" (Fig. 1 -Top-Left).
The Top-Right part of Fig. 1 shows the distribution of the Strehl as a function of the coherence time for our full sample.There is a steep rise in Strehl ratio with the coherence time (for τ 0 < 3 ms) followed by plateau for larger coherence times (for τ 0 < 5 ms).This shows that for low coherence times, the AO performances are limited by the temporal bandwidth error as mentioned in Cantalloube et al. (2020) which also derives a 3 ms threshold from the system point of view.This is in agreement with the laboratory and first on-sky measurements described in details in Petit et al. (2012).The overall impact of this performance limitation on our survey is likely important as illustrated by the large number of observations occurring when τ 0 < 4 ms (i.e 70% of the observations).Further analysis of the impact of this effect should include the seeing contribution in the Strehl error budget to disentangle the seeing and coherence time correlations.
Summarizing, we have shown in this section that the AO performances are clearly related to the observing conditions for our sample.We will further describe the impact of these observational parameters on the high contrast performances in section 4.

Data Reduction and Analysis
Large surveys such as SHINE (especially because of the large IRDIS FoV) unveiled a large number of point sources (mainly background sources) from which true sub-stellar mass companions need to be distinguished.For this reason, both astrometry and photometry of the detected points sources need to be precisely and homogeneously calibrated and extracted to test their companionship and charactezerize their nature.We detail in the following sub-sections the SHINE strategy: i/ to calibrate both IRDIS and IFS on various important aspects (distorsion, plate scale, True North, parallactic angle determination, central start position) in sub-section 3.1, ii/ to pre-process the scientific observations and apply basic reduction steps (cosmetics, flat, frame registration and recentering, cross-talk and wavelength calibration) as described in sub-section 3.2, iii/ to apply advanced ADI/ASDI algorithms summarized in sub-section 3.3, and ultimately extract the companion relative astrometry and photometry as detailed in sub-section.The final error budget for both astrometry and photometry is described in this last subsection.

Astrometry
We present in this section the astrometric calibration of the SHINE survey which depends on several factors described in detail in the following paragraphs: i/ platescale, True North and distorsion correction, ii/ angular offset (between field and pupiltracking and parallactic angle correction), iii/ star-centering.This section also presents as a conclusion, iv/ a sanity check of the calibration startegy based on Gaia-DR2 results.
Platescale, True North and distorsion An extensive description of astrometry with SPHERE can be found in Maire et al. (2016b) and Beuzit et al. (2019).For what concerns in particular the SHINE survey, the astrometric calibration consisted in the correction for the instrument anamorphism (0.60±0.02% between the horizontal and vertical directions of the detector, i.e. 6 mas at 1 arcsec), correction for constant offset angles (between the IFS and IRDIS fields of view, between pupil-tracking mode and field-tracking mode), and determination of the values for the pixel scale and the correction to the true north.We find small but non negligible variations of the last quantities with time requiring dedicated calibration for each observing runs.Appropriate values were estimated for each run using several reference fields of view in clusters with a large number of stars and having accu-rate (sub-mas) astrometry from Hubble Space Telescope (HST) or ground-based diffraction limited observations: 47 Tuc (Bellini et al. 2014;Soto et al. 2017;Bellini et al. 2017); NGC 6380 (Bellini et al. 2011;Soto et al. 2017) and Noyola (private communication); NGC 3603 (Harayama et al. 2008;Khorrami et al. 2016;Rochau et al. 2010); and Trapezium B1-B4 (Close et al. 2012(Close et al. , 2013) ) (see Figure 2).
IRDIS was used for all astrometric calibrations, because its field of view allows observation of a large number of stars, in most cases between 50 and 100.Only seven stars were available for the Trapezium B1-B4 field, leading to a less accurate calibration of the field orientation (by 1%) and on the plate scale (by 0.1%), and was used only when the other calibration fields were not accessible from February 2015 to March 2016.Due to the less good accuracy of the catalog positions for NGC 6380 and the smaller number of calibrating objects, we used this field only twice in May-June 2015 and once in June 2017.Using several calibrators enable to perform calibrations throughout the year and to do cross-calibration.Because the coronagraph has a small effect on the pixel scale, the astrometric fields were observed with a coronagraphic plate (including an offset for the Trapezium B1-B4 field observation, to shift the B1 star out of the coronagraphic mask).
The pixel scale is slightly different for the H2, H3, K1, and K2 filters, with mean values of 12.255, 12.250, 12.267, and 12.263 mas/pixel, respectively.For IFS, we used a constant value of the pixel scale of 7.46 ± 0.02 mas/pixel.Variations in this case have less impact on the results due to the small field of view.The typical measurement accuracy of the pixel scale is ±0.012 mas/pixel using our two best fields 47 Tuc and NGC 3603, while the true north correction (weighted value -1.77 degree) typically has an error of ±0.07 degree.For IFS, an additional offset of 100.48 degree in the clockwise direction is applied to account for the orientation of the instrument FOV.This leads to uncertainties in the position of 3-4 mas at the edge of the IRDIS field of view and less than 1 mas for IFS.For comparison, typical accuracy of GPI using their best calibrations are ±0.021mas/pixel in the pixel scale and ±0.12 degree in the true north correction (De Rosa et al. 2019).The better calibration accuracy obtained for SPHERE is due to the wider field of view of IRDIS that allows the use of stellar cluster fields as calibrators with more accurate catalog positions.
Field/Pupil tracking and parallactic angle calibration While the SHINE science observations are taken in pupil-tracking mode for ADI purposes, the astrometric fields were observed in field-tracking mode.The offset between the pupil-tracking and field-tracking modes was measured at the beginning of the survey to be equal to -135.99±0.11degree.The parallactic angle is computed from the data FITS header: timestamp, RA/DEC of the derotator (INS4.DROT2.RA/DEC) which are more accurate than the "RA" and "DEC" keywords (J2000 coordinates of the target).In case these parameters are used it is possible to correct for the precession of coordinates between 2000 and the date of observation.In order to derive the precise parallactic angle of each DIT from the timestamp we also included a parametric model of the overheads.We observe a systematic error in the parallactic angle estimation due to backlash in the derotator mechanism of ∼0.05 • , as demonstrated in Beuzit et al. (2019).In pupil-stabilized mode, this leads to a ∼0.4 pixel difference in the position of an object located at the edge of the IRDIS FoV on either side of the meridian.

Star centering
The astrometry we consider in this paper is relative to the star.Since the star point spread function peak is hidden by the coronagraphic mask and to avoid concerns due to the non uniform intensity distribution of the coronagraphic leakage, its position was determined using a special calibration (STAR-CENTER) where four faint replicas of the star image are created by giving a bi-dimensional sinusoidal profile to the deformable mirror (see Beuzit et al. 2019;Makidon et al. 2005).The STAR-CENTER calibration was repeated before and after each science observation, and results center estimations were averaged.While this calibration greatly reduces uncertainties in the exact position of the star and the center of the field rotation, experience shows that small drifts of a few mas of the star center during long sequences of ∼1-2h can occur.For this reason in order to increase the star centering accuracy, we generally used, for the second epochs, when doing (candidate) companions follow up, the STAR-CENTER setup solely for the science exposure.
We performed specific measurements to estimate the accuracy of the central star position when the STAR-CENTER setup is not used for the science exposure (i.e most cases presented in this paper).To do so we measured the position of the central diffraction peak on the IFS datacubes collapsed in wavelengths, and making the mean over all Detector Integration Time (DIT); as mentioned above, errors are likely independent of the errors in the STAR-CENTER procedure.We found that the mean position of the center is offset with respect to the nominal position along the Left-Right direction in the pupil reference frame by a small but significant amount 0.52 ± 0.10 mas for Y-H observations, while there is no offset for the Y-J mode (0.03 ± 0.05 mas) or along the Top-Bottom direction in the pupil reference frame in both modes.The root mean square (rms) scatter of the residuals after a 3-σ clipping are 1.35 mas (1.23 mas) in RA and 1.24 mas (1.43 mas) in declination, for the Y-H (Y-J) mode.Most outliers are found in data sets that were not validated and a few of them are binaries.There are 5% of the validated observations that have much larger dispersion of the position of the peak than usual; these anomalous cases make ∼10% of the observations acquired before February 2016, while the fraction reduces to ∼2% after that epoch.This is likely due to improvements in the AO calibration, that resulted in less distorted diffraction peaks.A few of these residual cases may be unresolved binaries.We conclude that a reasonable estimate for the accuracy of the absolute central star position is ±1.5 mas for both IFS and IRDIS.
Cross-check with Gaia DR2 An external check of the accuracy of our astrometry is provided by wide companions in the Gaia Data Release 2 catalogue (Gaia Collaboration et al. 2018).The number of sources in the IRDIS field of view with a contrast adequate to be detected with Gaia is limited because stars with known bright companions were not included in our sample to avoid problems for the AO and heavy saturation of the detector.By comparison with the SPHERE results, we find that Gaia DR2 limiting contrast (in the visual G-band) corresponds to contrasts in the near infrared of ∆H=2.5 mag for separation sep< 2.5 arcsec and 6 mag for sep> 4 arcsec2 .For this comparison, we considered 34 IRDIS close companion candidates (some of them not included in the sample described in this paper) with ∆H2<6 mag; out of which twelve are in Gaia DR2.We did not considered one of the Gaia DR2 data because it has very large error bars.In addition we removed some SHINE objects because they are saturated, or located at the edge of the IRDIS FOV (sep>6.5 arcsec), both leading to inaccurate astrometry.Gaia measurements are most likely unreliable for very high contrasts.This is the case of HD1160C, for which the Gaia contrast is 8.24 mag.Binary periods for these stars are so long that they should not affect the result.However, a couple of the objects are not binaries but rather background stars (HIP82430 and PDS 70).We took into account the relative proper motion between the SHINE observation and the Gaia epoch (2015.5).At the end we have seven good comparative measurements.The separation measured by SHINE is slightly larger than that measured by Gaia DR2.The mean offset is −2.8 ± 1.5 mas, with a root mean square (rms) of 3.9 mas.The position angle measured by SHINE is similar to that measured by Gaia DR2: the mean offset is 0.06 ± 0.04 degree, with an rms= 0.11 degree.The rms agrees well with the expected uncertainties in these quantities in the SPHERE data.
We conclude that the minimum accuracy of the astrometric calibration of SHINE data is ±2 mas at separations < 1 arcsec, and ±3 mas at larger separations.These values are similar to the scatter typically observed for GPI astrometric calibrators (De Rosa et al. 2019).In practice, the astrometric accuracy on the position of faint substellar companions detected in the SHINE data is limited by the measurement uncertainties from the image post-processing (Sect.3.3) and by the companion magnitude.

Photometry
We describe below the strategy to derive the relative photometry of the SHINE candidates considering the unsaturated and coronographic observations of a scientific target, and two tests done by using IFS to validate our strategy for both IRDIS and IFS which have very similar photometric biases.

Strategy
In the SHINE survey, photometry of candidate companions and limiting magnitudes for non detections are relative in contrast to the central star.Since the star is behind the coronagraphic mask, simultaneous photometry is not possible.We thus include a flux calibration (STAR-FLUX) for both IRDIS and IFS that is acquired just before and after the science exposure by offsetting the Differential Tip-Tilt Stage (DTTS) by about 0.5 arcsec with respect the coronagraphic mask using the SPHERE tip/tilt mirror (Beuzit et al. 2019).When performing this calibration, suitable neutral density filters are inserted to avoid detector saturation.The transmissions of these neutral density filters were carefully calibrated from 0.9 to 2.3 microns (see Beuzit et al. 2019 and ESO website3 ) and are taken into account for the contrast estimation.This procedure works very well in stable conditions, but it may be affected by variations of the Strehl ratio from evolving observing conditions.Higher photometric accuracy can be achieved by using the waffle pattern continuously (for companion candidates follow up) during the observing sequence to monitor the Strehl variations.
Cross-check with catalogues and binary companions A first test of the photometric calibration accuracy is provided by comparing the peak counts of the diffraction image of the flux calibration (corrected for the integration time and for the neutral density filter transmission) to the apparent magnitude of the stars from the 2MASS catalogue (Skrutskie et al. 2006), corrected for atmospheric extinction.This procedure neglects stellar variabil-Fig.3: Spectra of the white dwarfs used as spectrophotometric standards compared with predictions from model atmospheres (solid lines).Different symbols are results obtained from different epochs.We normalized the observed spectra at the median flux in a narrow range around 1.25 micron.
ity and the impact of variability on the Strehl ratio.Using J-band IFS data, we obtain an rms of residuals of 0.15 mag for observation in the Y-H mode, and 0.19 mag for those in the Y-J mode.This has been computed after eliminating the expected variations due the Strehl ratio and after clipping outliers deviating greater than 2.5 standard deviation from the mean with an iterative procedure (about 12% of the data).A similar value is obtained by considering the rms scatter of the estimate of the coronagraph central transmission -measured as the ratio between the maximum counts of the diffraction peak in the coronagraphic image and in the flux calibration (corrected for the length of the integration and neutral density filter transmission).The average coronagraphic transmission measured in this way are 0.00168±0.00005for the Y-H mode, and 0.00208 ± 0.00001 for the Y-J mode, respectively, for a median SHINE observation.
A second test of the photometric calibration accuracy is provided by the photometry of three standard systems (HD2133, HD114174, and REJ1925-563) composed of a main sequence star and a white dwarf with separations in the range 0.2-0.7 arcsec, we have observed during the survey.These systems were observed several times, with typically at least one observation per run.For schedule optimization, these observations were acquired in ADI with the star quite far from meridian and in poorer observing conditions, than typical survey targets: we consider here only those observations that were obtained in fair to good atmospheric conditions.Photometry of the faint companions was obtained by inserting a negative scaled point spread function at the position of the companion in the raw datacube, i.e before combining the various DITs.We then measured the rms residuals of the differential image obtained through monochromatic principal component analysis (PCA: Soummer et al. 2012) around the companion.The procedure then minimizes these residuals by simultaneously adjusting the contrast and the position.The values we consider are the mean of the results obtained using 2 to 6 PCA modes.The average contrasts in the J-band obtained with this procedure using the main sequence star and the white dwarfs are 7.68, 10.04, and 6.10 mag for HD2133, HD114174, and REJ1925-563, respectively.The rms value for the contrasts in the J-band/H-band are 0.11/0.05,0.18/0.08,and 0.05/0.07mag for HD2133 (8 observations), HD114174 (14 observations) and REJ1925-563 (6 observations).The photometric errors as expected increase with magnitude and as a consequence are smaller in H-band than in J-band.This may be also attributed to a higher impact of the variation of the SR and of the speckle noise for fainter targets and shorter wavelengths.A comparison of the spectra determined by our method with atmospheric models is given in Figure 3.
We conclude from these measurements that a reasonable photometric errors estimate for both the limiting contrast and for companions characterisation (not including the ADI post processing error contribution), is around ±0.2 mag.It is dominated by the STAR-FLUX variation which is also the case for the survey data photometric error which is on average equal to ±0.25 mag.

IRDIS-only steps
The SHINE survey was reduced by the SPHERE Data Centre (hereafter SPHERE-DC)4 .For IRDIS data, the first reductions steps (dark/background, flat, and bad pixel correction) rely on the SPHERE Data Reduction and Handling (hereafter DRH) pipeline Pavlov et al. (2008) provided by ESO.We used the on sky background because there is a significant difference between day-time background calibrations and the on sky background recorded in the science frames.This is the case mostly in Kband with a systematic, spatially variable, offset (typically of the order of 100 ADU for 64s exposure time) which is due to the sky background contribution itself.A similar effect is also visible in H-band, on a smaller level (typically of only few ADU).
Since most of the SHINE observations use the pupil-tracking observing mode, a very accurate determination of the star center is needed in order to successfully use both the angular differential imaging and spectral differential imaging methods.Also since IRDIS is used in dual-band imaging, the star center in both IRDIS channels is used to combine them.By default, we use the DRH sph-ird-star-center routine to determine the star center position, using the waffle images acquired for this purpose just before and after the science observations.This very fast routine provides an accurate centering in many cases, but could fail for weak waffle spots, especially in the K-band where these spots can be hidden beneath a much stronger thermal background noise but also when the deformable mirror offsets to create this pattern are set too low.We therefore designed an automated way to check the quality of the DRH centering by comparing the two center positions derived out of the four waffles.When the distance between these two possible center positions was found to be greater than 0.9 pixel, we used a SPHERE-DC made IDL routine that is more robust to identify weak waffle spots.This routine is able to detect weaker waffles by locating them within small circular apertures located at the expected (wavelength-dependent) position and by using in combination a high-pass spatial filtering, sky-background subtraction and median stacks of all waffle images available in order to increase their SNR.
A small fraction of the SHINE datasets use continuous waffle mode observations, meaning that the waffle spots are activated during the entire observation, for better astrometric monitoring (with some localized loss in the limiting contrast because of the secondary spots).In this specific cases, we perform individual re-centering of each frame in the sequence using the SPHERE-DC dedicated star-centering routine described above.This improves the quality by correcting any drift or jitter of the targeted star behind the coronagraph, and as a result improves the quality of the astrometric measurements by removing these sources of error.

IFS-steps
For IFS data we use the Data Reduction and Handling (DRH Pavlov et al. 2008) pipeline but complement it with additional steps implemented at the SPHERE Data Center (Mesa et al. 2015;Delorme et al. 2017a) to improve the wavelength calibration, apply a correction for cross-talk, and improve the handling of bad pixels.The improvements of the wavelength calibration are obtained by using a cubic fit whose coefficients are estimated from the wavelength shifts of the spots generated by the STAR-CENTER calibration -that are known to scale linearly with wavelength.The cross-talk correction is performed using an iterative procedure that corrects for the spectrograph PSF, using coefficients derived using appropriate tests performed in the laboratory during the instrument assembly.Bad-pixels are corrected using a dedicated sky observation acquired at the end of each science exposure.This procedure yields more accurate results than the one based on the flat field calibration within the DRH.

Final steps common to IRDIS and IFS
After these pre-reduction steps, both IRDIS and IFS datasets are corrected for the instrument anamorphism and the astrometric solution (pixel scale and True north) estimated from the calibration described in Sec.3.1.1 is applied to the dataset.The output data is composed of a pre-reduced master cube combining all frames obtained during a given observation sequence, that is used as input for all Angular Differential Imaging (ADI Marois et al. 2006) algorithms.We also associate to this master cube a vector of de-rotation angles for each frames, using the accurate timing and overheads for each frame, to produce a frame-toframe determination of the parallactic angle.

ADI and ASDI postprocessing
For both IRDIS and IFS, we obtain a good-quality re-centered images gathered in a single master cube associated with their parallactic angle values.Subsequent steps follow to estimate and subtract the stellar halo from each images, followed by derotation and stacking of the residuals.The most critical step is the estimation of the stellar halo that drives the level of the noise residuals.We applied different Angular Differential Imaging (ADI) algorithms to optimize the detection performances and to identify associated biases.We rely mainly on the SpeCal (Galicher et al. (2018) software which offers various ADI options.We selected for the homogeneous reduction TLOCI (Marois et al. 2014;Galicher et al. 2018) and Principal Component Analysis (PCA: Soummer et al. 2012) for IRDIS and PCA Angular and Spectral Differential Imaging (ASDI: Mesa et al. 2015) for IFS.
In the TLOCI approach (Lafrenière et al. 2007), the PSFreference is estimated for each frame and each location.Linear combinations of all data are computed to minimize the residuals into an optimization zone, which is much bigger than the subtraction zone to avoid the self-removal of point-like sources.The SpeCal version of the TLOCI algorithm is derived from the one described in Galicher et al. (2018) assuming a flat planet spectrum in contrast.Adjustable parameters are used to select the frames and to describe the regions of interest.In SpeCal, the gap between this region and the region of interest is set to 0.5 Full Width Half-Maximum (FWHM).Hence, the optimizing region is far enough from the region of interest so that the flux of a source in the latter does not significantly bias the linear combination.Finally, an additional parameter sets the radial width of the optimizing region.We considered here a radial width of the subtraction zone of 1 FWHM in radius; a radial-to-azimuthal width ratio of 1.5; a standard surface of the optimization zone was N = 20 PSF FWHM and 10% is used for the minimum residual flux ratio due to self subtraction compared to the flux of a putative candidate.
For historical reasons, two PCA algorithms are implemented in SpeCal.The first version can be applied on IRDIS or IFS data using ADI or ASDI.This algorithm follows the equation of Soummer et al. (2012).In the ADI case, which is the option selected to reduce IRDIS data in this paper, the principal components are calculated for each spectral channel independently.Each frame is then projected onto a limited number of modes.In the ASDI case, the algorithm is the same but it works simultaneously on the spatial and spectral frames.The second version of PCA used to reduce IFS data in this paper, is very similar to the first version we have described but it was applied only on IFS data using the ASDI option (Mesa et al. 2015).In addition, it is worth mentioning that both PCA algorithms we use have no frame selection to minimize the self-subtraction of point-like sources when deriving the principal components.
For PCA and TLOCI algorithms that bias the photometry of off -axis point sources, SpeCal estimates the throughput at each position in the field by generating a datacube of fake planets for which the ratio of the flux in the resulting image to the flux of the fake planet is calculated to obtain the centro-symmetrical 1Dthroughput as a function of the angular separation which is then applied to the images.For PCA, Specal estimates the throughput by inserting fake planets injections which are only 10 times brighter than the local stellar residuals level after PCA.This makes it possible to be close to the level of interest while being sufficiently above the residual to minimize the bias in the throughput estimation.For TLOCI, the throughput is calculated from an analytical formula (see section 2.10.2 of Galicher et al. (2018).The contrast curves, we derived in the following, for each spectral channel are based on the azimuthal standard deviation calculated in annuli of 0.5 FWHM width.Finally, the 5σ detection limits are derived by taking into account the flux loss from ADI self-subtractions, the transmission of the coronagraph at short separations (close to the inner working angle (IWA: 0.1") and the transmission of the neutral-density filter if used when registering the PSF.These detection limits are thus normalized by the unsaturated PSF flux.Both 1D and 2D contrast maps are estimated following these steps for each star for each reduction technique.
All target stars were processed for each instrument in a homogeneous way using at least TLOCI or PCA with similar sets of parameters.We considered 50, 100 and 150 PCA modes for IFS and 5 for IRDIS.We inspected by eye at least three residual maps for each star and for both IFS and IRDIS to look for candidate companions (CC).We also included small number statistic correction.In addition to the standard SpeCal reductions, we also used, on several cases, ANDROMEDA (Cantalloube et al. 2015) and PACO (Flasseur et al. 2018(Flasseur et al. , 2020b) ) algorithms to search for points sources.Given their statistical robustness to derive detection limits and to better identify false detections these algorithms will become the main algorithms for the final analysis of the SHINE survey.

Astrometric and Photometric extraction for point sources
PCA and TLOCI algorithms are known to distort the images of any off-axis point-like source.To retrieve accurately the relative photometry and astrometry of the detected candidate companions (CC) with their uncertainties, SpeCal fits a model of an offaxis point source image to the detected point source and adjusts its position and flux to locally minimize the flux.After building a model of the point source using the technique described in (Galicher & Marois 2011), the flux and the position of this synthetic image are adjusted within a disk of diameter 3 FWHM so that it includes the positive and the negative parts of the point source image.To optimize the computation time, instead of calculating the synthetic image each time we test a new planet position, we shift the synthetic planet image to its rough position.Once the optimization is completed, we measure the excursion of each parameters that increases the minimum residual level by a factor of 1.15.These excursions correspond to the 1 σ accuracy due to the fitting errors in the SpeCal outputs.The spectrum extraction from IFS data is performed similarly by processing the wavelength channels separately.For the channels with no detection above 1 σ we provide an estimation of the upper limit only.
In conclusion, we summarize the astrometric and photometric error budget given by Specal considering both calibrations and ADI/ASDI extraction errors.
For astometry, this budget includes: 1. calibration with uncertainties of the detector distorsion, plate scale and True North, 2. determination of the correction for constant offset angles (between the IFS and IRDIS fields of view, between pupiltracking and field-tracking modes); 3. calculation of the parallactic angle variation (correction of precession and timestamp), 4. central star position, 5. determination of the companion/candidate relative position with Specal (Galicher et al. 2018).
For photometry, it includes: 1. Contrast estimation with proper calibration of neutral density, exposure time, and associated error related to the variation between Start/End STAR-FLUX calibrations but not including the strehl variation during the observing sequence (psf), 2. Temporal variations of the stellar coronographic flux between 30-50 pixels during the observing sequence (Seq), 3. Companion flux determination considering ASDI signature, taking into account flux cancellation from the algorithm used, and coronograph attenuation correction as described in (Galicher et al. 2018). .

Contrast curves
We derived the 5-σ IRDIS and IFS contrast curves of each observation for all the targets in the sample as presented on Figures 4 and 5.These detection limits are derived based on the noise in the speckle-subtracted image, compensated for the throughput of the algorithm (calibrated with fake planet injections), the transmission of the coronagraph (calibrated from measurements in SPHERE), and the small sample statistics (Mawet et al. 2014).
More details are provided in Galicher et al. (2018).Two types of contrast are discussed in detail in this section: the raw contrast, computed on the median coronagraphic image of each observation and the contrast after post-processing described above.As This dependency remains clearly visible on post-processed data especially at 500 mas, despite other factors also coming into considerations such as the field rotation and the stability of the conditions (see Fig B .2).At small separation, these scatter plots shows that one can easily gain one magnitude in raw contrast by increasing the Strehl by 10% or the decreasing the seeing by 30%.As a consequence it is clear that conducting the survey in visitor mode which does not offer the best seeing condition had some impact on the survey ultimate contrast performances.Outside the AO correction radius, at separation greater than 900 mas, there is still a smaller dependency of the contrast with the seeing due to the residual light scattering outside the AO control radius and the lower height of the diffraction peak.As part of this study, the dependency of the post-ADI contrast on parameters tracing the quality and stability of the conditions were also investigated as illustrated on Fig B.1,B.4 and B.2.We considered the dispersion in the seeing, coherence time, Strehl, during the duration of the pupil-stabilised sequence.No significant correlation could be drawn from this sample.It is worth noticing that at small separation (100 mas) the processed contrast is lower than the raw contrast because of the very small angular rotation, the small throughput of the algorithm (from self subtraction) and the stronger coronagraphic residuals.
The gain from raw contrast to post-processed contrast is clearly visible on Fig 6 which shows a typical improvement greater than 5 magnitudes at short separations for both IRDIS and IFS.It is worth noticing that the contrast gain in the IFS field of view is consistent for all separations and reach at least 7 magnitude.For the shorter separation located at the edge of the coronagraph (100 mas) there is marginal improvement for IFS from ASDI and no improvement for IRDIS in contrast due to the very small angular rotation, the small throughput of the algorithm and the stronger coronagraphic residuals for both instruments.At larger separation than the AO cutoff, the improvement from the post-processing range from 2 to 4 magnitudes.
We also plot on Fig B .3 the processed contrast as function of the star magnitude for both IFS and IRDIS at various separations.From these figures it is clear that there seem to be small correlation between these parameters when the post-processing if performed using TLOCI.On the contrary the PCA method for both IRDIS and IFS seems to be affected by the magnitude of the star at least for the largest separations which is likely related to the noise from the instruments background.We also highlight that one clear cause of contrast degradation in the post-ADI contrast is the presence of a smooth halo within the AO-corrected region which results from either bad seeing or from high-altitude wind related to the jet stream as discussed in (Cantalloube et al. 2018) leading to a non symmetrical halo in the direction of the wind (called the wind driven halo, WDH).This halo is rotating in the pupil-stabilised data set because it is fixed on the sky as described in (Cantalloube et al. 2020).In post-ADI frames, it therefore appears as a brighter elongation along the wind direction, with negative counterparts at 90 degrees.The impact of this effect does not appear to be very strong on the azimuthal averaged contrasts we have plotted but is noticeable on the 2D contrast maps.Increasing the temporal bandwidth of SAXO is the foreseen solution to help mitigate this effect.This is considered as one option in a forthcoming upgrade of the instrument, SPHERE+ (Boccaletti et al. 2020).

Mass detection limits
The contrast curves are converted into mass limits using massluminosity relationships.Whereas for old ( 1 Gyr) systems this relationship is essentially unique for gas giants at large separation (Burrows et al. 1997;Baraffe et al. 2003), at young ages the value of the post-formation luminosity still remains uncertain (Marley et al. 2007;Spiegel & Burrows 2012;Marleau & Cumming 2014).We present here a mass conversion of the contrast curves for a few specific targets using the canonical predictions of the COND-2003 evolutionary models (Baraffe et al. 2003).The impact of using other mass-luminosity relationships on the sensitivity of the SHINE survey are explored in more details in Vigan et al. (2020).For IRDIS contrast curves, we convert the contrast curves using the evolutionary models computed in the appropriate dual-band filter, while for the IFS contrast curves we use the predictions in the J-band filter for YJ data and in the Hband filter for YJH data.As demonstrated in Vigan et al. (2015), this approach for the IFS detection limits provides an accurate estimation of the detection limit.
For a selection of companions presented in 5.1, we show the detection of planet-like limits achieved by the survey in Fig. 10.Several targets appear to be relatively easy (HR 8799 bcd, HIP 65426 and HIP 64892) for this survey, allowing for exquisite characterisation.Other targets appear to be much more challenging mainly because they are at separations closer than 20 AU.The detection limits for these typical objects reach few Jupiter masses at such small separation while it can reach around 1 Jupiter mass at separations greater than 50 AU.(Chauvin et al. 2017b) and the brown dwarf companion to HIP 64892 B (Cheetham et al. 2018).PDS 70 was not originally part of the SHINE sam-Fig.5: Contrast curves at 5 σ obtained for a small subsample of targets dedicated to follow up for IRDIS mode observations.The solid line gives the median value of the contrast.The dash line gives the mean value of the contrast.Red color is for IRDIS data reduced in PCA ADI, Blue color is for IRDIS data reduced in TLOCI ADI Fig. 6: IRDIS and IFS processed contrast computed using respectively TLOCI and PCA-ASDI as function of raw contrast at various separations in irdifs mode.Blue is for 100 mas separation, orange for 200 mas separation and pink for 400 mas separations, black for 500 mas separation, , brown for 700 mas separation, red for 1000 mas separation , green for 2000 mas separation and yellow for 4000 mas separations.The black and brown line represent respectively a gain of 1 and 6.25 assuming 39 spectral channels.In fact So there are 15 independent spectral IFS channels leading to a maximum gain of 3.9 when neglecting the possible combination of the two IRDIS channels.

Exoplanet, Candidate Companion Detections
ple and was specifically targeted to explore the properties of the transition disk that led to the discovery of the first planet PDS 70 b in this system (Keppler et al. 2018;Müller et al. 2018), and the confirmation of the second planet PDS 70 c (Mesa et al. 2019).The other substellar companions, known from previous direct imaging campaigns, were originally blindly selected based on their host star properties, but then carefully characterized to study their orbital, spectral and physical properties in connection sometimes with their environment (the presence of additional planets or disk structures in the system).The SNR maps for a selection of candidates are shown in Fig. 7 for IFS data and Fig. 8 for IRDIS data.The extracted spectro-photometric data are reported in Fig. 11.The measured position, contrast, semi-major axis, predicted mass for each companion are listed in Table 1.The detection limits for a few of these systems are shown on Fig. 10 as illustration of the SPHERE performances.The Sco-Cen members HIP 71724 and HIP 73990 with stellar to substellar-mass companion previously discovered in sparseaperture masking by Hinkley et al. (2015) were observed, but the close companions were not resolved given their probable small separations at the epoch of our observations.In addition to the discovery and characterization papers on HIP 65426 b and HIP 64892 B summarized below, in-depths studies of various of these companions based on from the SHINE sample have been published in these early years of the survey: PZ Tel B (Maire et al. 2016a), HIP 107412 B (Milli et al. 2017a;Delorme et al. 2017b) for the brown dwarf companion, and HR 8799 b,c,d,e (Zurlo et al. 2016;Bonnefoy et al. 2016), 51 Eri b (Samland et al. 2017;Maire et al. 2019), HD 95086 b HIP 65426 b was the first planet discovered with SPHERE as part of the SHINE Survey (Chauvin et al. (2017b)).This warm, dusty giant planet is orbiting at a relatively large projected angular distance of of 830 mas (92 au projected) from its intermediate-mass primary (Fig. 7 and 8).Multi-epoch observations confirm that it shares common proper motion with HIP 65426, a young A2 member of the 17 Myr old Lower Centaurus-Crux association.Spectro-photometric measurements extracted with IFS and IRDIS, as shown in Fig. 11, between 0.95 and 2.2 µm indicate a warm, dusty atmosphere characteristics of young low-surface gravity L5-L7 dwarfs.Hot-start evolutionary models predict a luminosity consistent with a 6 − 12 M Jup , T eff = 1300 − 1600 K and R =1.5 R Jup for the planet.These results were later confirmed by Cheetham et al. (2019) combining SPHERE and NaCo observations.
Given its physical and spectral properties, HIP 65426 b represents a particularly interesting case to study the presence of clouds as a function of particle size, composition, and location in the atmosphere, to search for signatures of non-equilibrium chemistry, and finally to test models of planet formation and evolution (Cheetham et al. 2018).The planet location would not fa-vor a formation by core accretion unless HIP 65426 b formed significantly closer to the star followed by a planet-planet scattering event.Dedicated simulations by Marleau et al. (2019) showed that core formation at small separations from the star followed by outward scattering and runaway accretion at a few hundred astronomical units succeeds in reproducing the mass and separation of HIP 65426 b.In such a scenario, the planet is predicted to have a high eccentricity (≥ 0.5) and to be accompanied by one or several roughly Jovian-mass planets at smaller semimajor axes, which also could have a high eccentricity.SHINE detection limits setting relatively good constraints at close separation (2 M Jup beyond 20 au) as shown in Fig. 10 do not exclude the presence of unseen inner massive planets in that system that could have scattered out HIP 65426 b.

HIP 64892 B
A bright, brown dwarf companion was also discovered around the star HIP 64892, a B9.5V member of the Lower Centaurus-Crux association.The measured angular separation of the companion (about 1.27") corresponds to a projected distance of 159 au.We observed this target with IRDIS dual-band Imaging and in long-slit spectroscopy to estimate its spectral energy distribution and astrometry.More details can be found in the discovery paper by Cheetham et al. (2018).Luminosity and spectrum are consistent with a young (16 Myr), low-surface gravity brown dwarf with a spectral type of M9.From comparison with the BT-Settl atmospheric models its effective temperature was estimated to T eff = 2600 K, and comparison of the companion photometry to the COND evolutionary models yields a mass between 29 and 37 M Jup .Despite this, HIP 64892 is one of the highest mass stars around which a substellar companion has been detected, due to the challenges associated with observing such stars and the tendency for large surveys to focus on solar-like or low-mass stars.HIP 64892 is a rare example of high or intermediate mass stars with extreme mass ratio (q = 0.01) companions at large separations and will be useful for testing models relating to the formation and evolution of such low-mass objects.Notes.(Samland et al. 2017), while PZ Tel B is an M-type object (Maire et al. 2016a).These very red colors were quite unexpected for such faint objects, and indicate that young substellar objects have a spectrum quite different from those of older ones of the same luminosity, which can be interpreted as an indication for the presence of dust in their atmosphere (Chauvin et al. 2004;Bonnefoy et al. 2014).From the point of view of planet detection -the focus of this paper -the implication is that young substellar companions may be below detection limit in the Y and J-band and even faint in the H-band, while they may be detected quite easily in the K-band (extension of the datacubes to short wavelengths is however useful to model and better subtract speckles).An extreme case of this behaviour is represented by HD 95086 b (Chauvin et al. 2018).This has a significant impact on the detection limits of our survey and points toward the need of model atmospheres that include dust to properly derive masses from the observed magnitudes.

Ranking and identification of SHINE candidates
We identified 1483 unique candidates in the IRDIS images (all epochs), 1176 of which were re-detected at several epochs or using archival data and 307 of which were only detected at a single epoch.Most of the candidates have been observed with the H23 filters (> 95%).Sixty-nine candidates have K1-K2 photometry extracted from the IRDIFS_EXT observations.They were identified in follow-up observations of known systems (HD 95086, HIP107412, HIP 65426, HIP 64892, CD-35 2722) or for the most part in first epoch observations of stars belonging to the Sco-Cen association for which the detection of late-L dusty companions with red near-infrared colors is favored at K-band (e.g., HD95086b; Chauvin et al. 2018).
We used color-magnitude diagrams (CMDs) to identify the most promising candidates and perform an efficient follow-up campaign.The candidate absolute magnitudes were computed assuming they are bound and at the same distance as the observed star.Most contaminants are expected to be M-and K-type stars (e.g., Parravano et al. 2011).Their characteristic colors as well as the range of contrasts probed by the observations creates a locus which can fall in a distinctive color-luminosity space from the sequence of substellar objects in the CMDs.The diagrams are described in details in Bonnefoy et al. (2018) and are enriched with the recent photometry and distances of reference objects measured by Müller et al. (2018), Mesa et al. (2019), andDupuy et al. (2020).
The large number of candidates identified in the H23 observations allows for building a well defined locus of contaminants in the corresponding CMD (see Figure 9).The locus shows an elongation along a vector which may be related to the distribution of spectral type of the contaminants combined to interstellar extinction (the fainter the object with respect to the host star, the higher the optical extinction and color reddening).The locus intersects the sequence of substellar objects at the L/T transition and therefore can not discriminate contaminants from bound companions falling in that luminosity range.We could however define an exclusion zone (grey shaded region) where substellar companions are not expected to produce such colors and abso-lute magnitudes.The zone is defined to avoid missing i/ young mid-and late-T-type planets such as 51 Eri b (Macintosh et al. 2015) which can have reddened colors with respect to mature field T-dwarfs and ii/ young planets at the L/T transition such as HR8799 which can appear under-luminous with respect to older objects in the same temperature range.This criterion was used to identify 44% of the candidates with H23 photometry as likely contaminants for which a re-observation was not warranted.Consequently, we decided to prioritize the observations of candidates falling right on the sequence of late-M to mid-L dwarfs and for which the locus spans a narrow range of bluer colors in the same luminosity range.This allowed for identifying HIP 64892B and HIP 65426b as promising candidates before they could be confirmed as co-moving objects.We could also identify blindly known companions (e.g., HIP 78530B, η Tel B, PZ Tel B, AB Pic b, CD-35 2722B, HIP 73990B) as promising candidates.
Most of the remaining ambiguous candidates were observed at a second epoch to determine whether they are co-moving companions.Two crowded fields were observed with the J23 filter set.As shown in Bonnefoy et al. (2018), the locus of contaminants falls at a distinct location from the substellar objects in these CMDs and allowed to classify most of these candidates as contaminants.The locus in the K1K2 diagrams is more dispersed and the error bars on the photometry tend to be wider owing to the higher thermal background at these wavelengths.Therefore, detections in K1K2 requires follow-up observations to clarify the nature of the detections.
The final status of the candidates is given in Table A.1.The companions with a confirmed common-proper motion are flagged as 'C'.The contaminants either discriminated using the CMD or the comon-proper motion test are flagged as 'B'.Ambiguous cases are flagged as 'A'.We report in addition several not clear detections ('NC') in that table which correspond to the identification of candidates with unreliable astrometry and/or photometry and for which our standard classification scheme could not be applied.

Disks detected in SHINE
While not specifically optimized for this purpose, our observations allowed detection of circumstellar disks, thanks to the exquisite sensitivity and large FoV of IRDIS, twelve of SHINE targets, out of which three were new detections.By construction of our survey, they are debris disk -in fact, we did not consider stars with gas-rich disks to reduce the concern related to the impact of disks in planet detection.Table 2 collects the main data for the disks we detected; a gallery of their appearance from the IFS data is presented in Figure 12.Detailed analysis of each one of these disks is subject of specific studies, most of them already published (see Table 2).
Not surprisingly, stars with detected debris disks are typically young; the median age is 24 Myr (data from Paper I) and there is only one target with an age older than 50 Myr (HIP682=HD377).This result is similar to that obtained by the survey by (Esposito et al. 2020).We notice that most of the detected disks are seen at high inclination i, with a median value of i = 82 degree.The only disk that has a low inclination in our survey (TWA 7) was detected only in polarimetry (Olofsson et al. 2018).The strong dependence on inclination is because of the combination of a deeper optical depth (debris disks being optically thin), of higher efficiency of forward scattering, and of the use of angular differential imaging to reach deeper contrast (see Milli et al. 2012;Sissa et al. 2018).This selection effect would be less severe using polarimetric observations, as often done when studying disks (see the case of TWA 7; see also Esposito et al. 2020).The strong selection bias towards edge-on disks implies that debris disks at the typical age of our targets should be much more common than indicated by eleven detections over a sample of 150 stars.If we indeed consider only the 52 stars with age <30 Myr in our sample, we detected disk around 8 of them.In five cases the inclination is i > 80 degrees, an occurrence that has 17% of the probability to happen.This suggests that debris disks similar to the ones we detected surround a large fraction of young stars (see Sissa et al. 2018 for a similar argument about the M-stars in our survey or Meshkat et al. 2017).
Debris disks have a limited impact on detection limits of point sources.Low inclination debris disks have negligible effects because they are optically thin and well below detection limit.For what concerns high inclination disks, the effect is more relevant (see Milli et al. 2012 for a discussion).They typically cover a small fraction of the sky area; however, planet orbit is possibly coplanar with the disk, so that the sky region covered by the disk is where most likely the planets are.As discussed by Milli et al. (2012), the flux loss depends on the technique used, and may be large for some of them.On the other hand, in addition to flux losses, Milli et al. (2012) showed that disk features may lead to false positive point sources detections.Furthermore, the multi-wavelength approach allowed by simultaneously using both IFS and IRDIS provides more diagnostics for an appropriate classification of detected features and more particullarly point sources.The good sensitivity of ASDI-PCA to structures within disks is clearly shown by the case of AU Mic, where this technique allows to detect very faint point-like sources with contrasts of ∼ 14.5 mag at SNR> 10 along the edge-on disk (Boccaletti et al. 2018).The SNR of these detections agree with the expectation based on our estimate of the contrast limit at this separation based on the whole images.Note that these structures are however not planets: the spectrum is flat, inconsistent with that of young planets, and the objects are by far too bright to be reflecting planets.We carefully considered a similar possibility for all our targets, and we are confident that there is no false positive among our detections.
Concerning the detection of the disk structure themselves, the standard ADI methods used in this paper lead to several known causes of artifacts altering both the morphology and the photometry of the disks.In particular they suffer from positive / negative replicas and from the self-subtraction phenomenon which can only be alleviated by using a physical model of instrumental effects as recently proposed by (Pairet et al. 2020;Flasseur et al. submitted) but which are beyond the scope of this paper.

Conclusions and prospects
We process in a uniform manner more than 300 datasets from the SPHERE/SHINE Survey obtained at the VLT/ESO in visitor mode and assess SHINE survey typical sensitivity as a function of the host star and observing conditions.From these first 150 stars observed out of the planned 400-star survey, we reach typical contrasts of 10 6 at 1" separation and reach up to 10 7 contrast at 2" separation for the best observations.SPHERE/SHINE delivers one of the few largest and deepest direct imaging surveys for exoplanets conducted to date.Compared to GPIES, SHINE has been operating in priority visitor mode, the observing conditions for SHINE in regular visitor mode are highly variable and on average not as good leading to a wider spreading of the contrast performances.Sixteen substellar companions around twelve host stars have been detected within the first half of the SHINE survey, eight brown dwarfs and eight planetary-mass companions.The paper summarizes the measured position, contrast, and derived mass for each companion as well as it produces astro-photometric data for all the candidates detections.The two new discoveries in the first half of the SHINE survey are the planetary mass companion to HIP 65426 b (Chauvin et al. 2017b) and the brown dwarf companion to HIP 64892 B (Cheetham et al. 2018).All other substellar companions characterized during SHINE were discovered in previous direct imaging campaigns.
The SHINE survey is due to be completed in 2021, but will certainly extend over a few more years to become complete in terms of follow-up for all candidates within at least a 300 au and possibly even farther away.The final sample will include over 500 stars, which will make SHINE the largest high contrast imaging survey to date, covering from B to M stars in the solar neighborhood.The reanalysis of the complete SHINE data with advanced processing techniques (Cantalloube et al. 2015;Ruffio et al. 2017;Flasseur et al. 2018Flasseur et al. , 2020b,a;,a;Berdeu et al. 2020) will be performed in order to provide improved detection limits and reach the full power of this survey using the full sample.
There are a limited number of bright, young, nearby stars available to high-contrast direct imagers like SPHERE and GPI.As a result, at the current achievable contrasts we cannot expect a significant number of new detections of imaged giant planets following the completion of both GPIES and SPHERE/SHINE (Desidera et al. (submitted)) campaigns.An improvement in performance from instrument upgrades (Boccaletti et al. (2020), Chilcote et al. (2018) ), however, could unlock mass/separation phase space around these target stars that are currently inaccessible, which are in addition likely to host more planets Wagner et al. (2019).Additional planets discovered with upgraded GPI and SPHERE would allow us to test the robustness of the trends with stellar mass and planetary mass uncovered in this paper.The Gaia mission will identify astrometric signatures of planets that may be confirmed with direct imaging, which can also boost the number of companions in this separation range.Finally, cross checks with results from other indirect techniques, e.g.transits or high precision radial velocity, will allow a much better picture of the architecture of planetary systems.The spectra shown are obtained using ASDI PCA, and corrected for the attenuation factor using fake planets at similar separation.This method is reasonable when the contrast is very large, so that the planet cannot be seen well at individual wavelengths.Red symbols: contrast values obtained from IRDIS data.Note that there is some offset between IFS and IRDIS contrasts.This is at least in part due to the non optimal extraction of IFS spectra shown in this Figure .More accurate spectral extractions may be found in the papers about the individual targets.Fig. 12: Gallery of images of debris disks in the F100 sample detected with IFS.Some of detection is very marginal (e.g.HIP86598) .The disk is not detectable for TWA7 and HIP682 (it was however detected with SPHERE using other approaches).Whenever possible, we combined several images to have a cleaner detection.This is the reason the planet is not visible in the beta Pic image.Generally, we used ASDI PCA with 25 modes, that is a quite conservative method, and smoothed the images to have a better view.In each panel, the green cross mark the star position Fig. 13: Gallery of images of debris disks in the F100 sample detected with IRDIS obatained using PCA with 5 modes in H 2 .Some of detections are marginal they are however better detected with other signal extraction approaches.

IFS RAW CONTRASTS
IFS PROCESSED CONTRASTS

Fig. 1 :
Fig. 1: Observing conditions and AO performances in IRDIFS observing mode: (Top-Left) Strehl ratio as a function of the seeing measured by SPARTA.(Top-Right) Strehl Ratio as function of the atmospheric coherence time estimated by SPARTA.A linear fit is overlaid in red on the Top figures.(Bottom-Left) Strehl ratio as function of the star V magnitude.The blue data points represent the data taken with coherence time smaller than 3 ms.(Bottom-Right)Histogram illustrating the repartition of the observations (including multiple observations of the same target) as function of the star magnitude in V (dark blue) and H-bands (light blue).The apparent large number of observations of R=3.5 stars is an artificial effect due to the large number of observations targeting Beta Pictoris.

Fig. 2 :
Fig. 2: IRDIS H2 images of the fields used for the astrometric calibration.Upper left panel: 47 Tuc; upper right panel: NGC 3603; lower left panel: NGC 6380; lower right panel: Orion B1B4.A logarithmic intensity scale is used to show also faint stars.
illustrated on Fig B.1 and B.4, the raw coronagraphic contrast at various separations shows a strong dependency with the Strehl and seeing both estimated by the average of the SPARTA values during the coronagraphic sequence.The two smaller separations (from 100 to 700 mas) are within the AO control radius located at 840 mas separation radius in the H-band.

Fig. 4 :
Fig. 4: Contrast curves at 5 σ obtained for the full sample for irdifs-ext (Left) and irdifs (right) modes observations.The solid line gives the median value of the contrast.The dash line gives the mean value of the contrast.Red color is for IRDIS data reduced in PCA ADI, Blue color is for IRDIS data reduced in TLOCI ADI, and green color is for IFS with PCA ASDI reduction

5. 1 .
Brown dwarfs and exoplanets A total of sixteen substellar companions were imaged in the course of this part of the SHINE survey, including seven brown dwarf companions (PZ Tel B, η Tel B, CD -35 2722 B, HIP 78530 B, HIP 107412 B, GSC 8047-0232 B and HIP 64892 B), and ten planetary-mass companions (51 Eri b, β Pictoris b, HD 95086 b, HR8799 bcde, GJ 504 b, AB Pic b and HIP 65426 b).Two new companions have been discovered in this sample: the exoplanet HIP 65426 b

Fig. 7 :
Fig. 7: SNR maps for the sample data where at least a substellar close companion was detected with IFS.In each panel, the green cross mark the star position

Fig. 8 :Fig. 9 :
Fig. 8: SNR maps for the sample data where at least a substellar close companion were detected with IRDIS.In each panel, the white cross marks the star position Fig. 11 shows the spectra of the brown dwarf and exoplanet companions detected within the IFS field of view: 51 Eri b, β Pic b, HD 95086 b, HIP 65426 b, HR 8799 d and e, HIP 107412 B and PZ Tel B. Most of these spectra are very red and can be classified as of L-or early T-type.The only middle T-type object detected in the IFS field of view is 51 Eri b(Samland et al. 2017), while PZ Tel B is an M-type object(Maire et al. 2016a).These very red colors were quite unexpected for such faint objects, and indicate that young substellar objects have a spectrum quite different from those of older ones of the same luminosity, which can be interpreted as an indication for the presence of dust in their atmosphere(Chauvin et al. 2004;Bonnefoy et al. 2014).From the point of view of planet detection -the focus of this paper -the implication is that young substellar companions may be below detection limit in the Y and J-band and even faint in the H-band, while they may be detected quite easily in the K-band (extension of the datacubes to short wavelengths is however useful to model and better subtract speckles).An extreme case of this behaviour is represented by HD 95086 b(Chauvin et al. 2018).This has a significant impact on the detection limits of our survey and points toward the need of model atmospheres that include dust to properly derive masses from the observed magnitudes.

Fig. 10 :
Fig. 10: SPHERE IRDIS and IFS mass detection limits.The IFS H data are represented in blue, the IFS YJ data are represented in violet, while the IRDIS H23 and IRDIS K12 data are represented respectively in green and red.The limits have been converted into mass using the COND-2003 evolutionary models (Baraffe et al. 2003).

Fig. 11 :
Fig. 11: Black symbols: SPHERE contrast spectra of substellar companions detected within the sample of stars considered in this paper.Upper row: left: 51 Eri b; center: β Pic b; right: HD 95086 b; Middle row: left: HIP 65426 b; center: HR 8799 d; right: HR 8799 e. Lower row: left: HIP 107412 B; center: PZ Tel B. Results shown are obtained by making a median of results at different epochs with good quality results.Spectra obtained with the YJ and YH mode were kept separate.The error bar is either the result of the scatter of results at different epochs, or if only one was available in that epoch, the error bar obtained by considering other positions at the same separation from the star.The spectra shown are obtained using ASDI PCA, and corrected for the attenuation factor using fake planets at similar separation.This method is reasonable when the contrast is very large, so that the planet cannot be seen well at individual wavelengths.Red symbols: contrast values obtained from IRDIS data.Note that there is some offset between IFS and IRDIS contrasts.This is at least in part due to the non optimal extraction of IFS spectra shown in this Figure.More accurate spectral extractions may be found in the papers about the individual targets.

Fig. B. 4 :
Fig. B.4: IFS raw (Left) and processed (Right) contrasts in YJ bands at various separations computed using PCA ASDI estimated as function of the observing conditions in irdifs mode.Blue is for 100 mas separation, orange for 200 mas separation and pink for 400 mas separations, black for 500 mas separation, , brown for 700 mas separation.

Table 1 :
Exoplanets and brown dwarf companions observed during the SHINE campaign

Table 2 :
Debris disk detected from the SHINE survey (F150 sample)