Subaru Hyper Suprime-Cam Survey for An Optical Counterpart of GW170817

We perform a $z$-band survey for an optical counterpart of a binary neutron star coalescence GW170817 with Subaru/Hyper Suprime-Cam. Our untargeted transient search covers $23.6$ deg$^2$ corresponding to the $56.6\%$ credible region of GW170817 and reaches the $50\%$ completeness magnitude of $20.6$ mag on average. As a result, we find 60 candidates of extragalactic transients, including J-GEM17btc (a.k.a. SSS17a/DLT17ck). While J-GEM17btc is associated with NGC 4993 that is firmly located inside the 3D skymap of GW170817, the other 59 candidates do not have distance information in the GLADE v2 catalog or NASA/IPAC Extragalactic Database (NED). Among 59 candidates, 58 are located at the center of extended objects in the Pan-STARRS1 catalog, while one candidate has an offset. We present location, $z$-band apparent magnitude, and time variability of the candidates and evaluate the probabilities that they are located inside of the 3D skymap of GW170817. The probability for J-GEM17btc is $64\%$ being much higher than those for the other 59 candidates ($9.3\times10^{-3}-2.1\times10^{-1}\%$). Furthermore, the possibility, that at least one of the other 59 candidates is located within the 3D skymap, is only $3.2\%$. Therefore, we conclude that J-GEM17btc is the most-likely and distinguished candidate as the optical counterpart of GW170817.


Introduction
The existence of gravitational waves (GWs) is predicted in the theory of general relativity. Although the existence is indirectly demonstrated by the energy loss of a binary pulsar system (Hulse & Taylor 1975;Taylor & Weisberg 1982), the direct observation of GWs had not been realized owing to its small amplitudes. The first direct detection is achieved with the Advanced Laser Interferometer Gravitational-Wave Observatory (LIGO) on Sep 14, 2015 (Abbott et al. 2016a). The first GW source originates from the coalescence of two black holes, each ∼ 30M⊙. The discovery is important not only for the direct probe of the strong field dynamics of general relativity, but also for the first evidence of a black hole binary. LIGO and Advanced Virgo subsequently detect three GW signals and one candidate signal, all from the coalescence of black-hole binaries (Abbott et al. 2016b;Abbott et al. 2016c;Abbott et al. 2017a;Abbott et al. 2017b). These discoveries open the era of "gravitational wave astronomy".
However, the conclusive identification of the GW sources on the sky remains challenging because of the poor sky localization with the gravitational wave observations. The sky localization areas of four GW sources are about 230 − 1160 deg 2 (90% credible region) with two detectors of LIGO and about 60 deg 2 (90% credible region) even with three detectors including the Advanced Virgo. Since there are many galaxies in the area, it is impossible to determine the host galaxy of a GW source only with the GW observations. Therefore, multiwavelength searches for electromagnetic (EM) counterparts are initiated after the alerts of GW detection from the LIGO-Virgo networks. So far, no firm EM counterparts have been found (e.g., Soares-Santos et al. 2016;Smartt et al. 2016;Kasliwal et al. 2016;Morokuma et al. 2016;Yoshida et al. 2017), except for a report of the putative detection for GW150914 with Fermi/GBM (Connaughton et al. 2016, but questioned by Savchenko et al. 2016;Greiner et al. 2016).
The non-detection of EM counterparts is not surprising because the four GWs are originated from mergers of black holes, although several theoretical studies try to explain the putative Fermi/GBM emission (e.g., Yamazaki et al. 2016). On the other hand, first-principle numerical simulations with general relativity demonstrate that binary coalescence including at least one neutron star (NS) can eject materials as dynamical ejecta (e.g., Rosswog et al. 1999;Goriely et al. 2011;Hotokezaka et al. 2013;Bauswein et al. 2013) and post-merger ejecta (e.g., Dessart et al. 2009;Fernández & Metzger 2013;Shibata et al. 2017). The ejecta dominantly consist of r-process elements (e.g., Lattimer & Schramm 1974;Eichler et al. 1989;Korobkin et al. 2012;Wanajo et al. 2014), and thus the decay of radioactive isotopes produced by the r-process nucleosynthesis heats up and brightens the ejecta. The EM-bright object is called "kilonova" or "macronova" (Li & Paczyński 1998;Kulkarni 2005;Metzger et al. 2010), and regarded as a promising EM counterpart of a GW (Metzger & Berger 2012;Kasen et al. 2013;Tanaka & Hotokezaka 2013;Metzger & Fernández 2014;Tanaka et al. 2014;Kasen et al. 2015;Metzger 2017). Also, the central engine of a short gamma-ray burst, which is believed to originate from a binary neutron star coalescence, is a possible energy source of EM counterparts through its jet and gamma/X-ray emission (e.g., Kisaka et al. 2016).
On Aug 17, 2017, 12:41:04 UTC, Advanced LIGO and Advanced Virgo detected a GW candidate from a binary NS coalescence, being coincident with a gamma-ray detection with Fermi/GBM (The LIGO Scientific Collaboration & the Virgo Collaboration 2017a; The LIGO Scientific Collaboration & the Virgo Collaboration 2017b). The sky localization with the three detectors is as narrow as 28 deg 2 for a 90% credible region centered at R.A. = 13 h 08 m , decl. = −22 • 30 ′ (J2000.0) (Abbott et al. 2017c). And the localization is overlapped with the error regions of gamma-ray detection with Fermi/GBM and INTEGRAL (Connaughton et al. 2017;Savchenko et al. 2017a;Savchenko et al. 2017b). The GW observation reveals the luminosity distance to the GW source, named GW170817, as 40 +8 −14 Mpc (90% probability) (Abbott et al. 2017c). Although GW170817 appeared at the position close to the Sun, the first significant alert of a binary NS coalescence and the narrow sky localization area initiate many EM follow-up observations (Abbott et al. 2017d).
Along with the EM follow-up observation campaign of GW170817, the Japanese collaboration for Gravitational wave ElectroMagnetic follow-up (J-GEM) performed a survey with Hyper Suprime-Cam (HSC, Miyazaki et al. 2012), which is a wide-field imager installed on the prime focus of the 8.2m Subaru telescope. Its FoV of 1.77 deg 2 is largest among the currently existing 8-10 m telescopes, and thus it is the most efficient instrument for the optical survey. In this paper, we summarize the observation with Subaru/HSC and properties of discovered candidates. Throughout the paper, we correct the Galactic reddening (Schlafly & Finkbeiner 2011) 1 , and all the magnitudes are given as AB magnitudes.

Observation and data analysis
We started HSC observation from Aug 18.23, 2017 (UTC), corresponding to 0.7 days after the GW detection, and also performed HSC observation on Aug 19, 25, and 27. All the observations were carried out in the z-band. The poor visibility of GW170817 from Maunakea compels us to conduct the survey during the astronomical twilight. The observations on Aug 25 and 27 concentrate on one field because the target fields set immediately after the sunset. The survey pointings are selected from HEALPix grid with resolution of NSIDE=64 by following criteria: higher probability of GW170817 sky localization and larger number of nearby galaxies in the GLADE catalog 2 (Table 1). We also choose the pointings located in footprints of Pan-STARRS1 (PS1, Chambers et al. 2016) and use the PS1 catalog and images for astrometric calibration and image subtraction, respectively. As the fields with smaller right ascension and declination set earlier, we conduct the observations in order of reaching the elevation limit earlier. The observed area is 28.9 deg 2 corresponding to the 66.0% credible region of GW170817 ( Figure 1). Exposures used in the following analysis are listed in Table 2.
The data are analyzed with hscPipe v4.0.5, which is a standard reduction pipeline of HSC (Bosch et al. 2017). It provides full packages for data analyses of images obtained with HSC, including bias subtraction, flat fielding, astrometry, flux calibration, mosaicing, warping, stacking, image subtraction, source detection, and source measurement. The astrometric and photometric calibration is made relative to the PS1 catalog with a 4.0 arcsec (24 pixel) aperture diameter. Further, in order to select variable sources, we perform image subtraction between the HSC and archival PS1 z-band images using a package in hscPipe based on an algorithm proposed by Alard & Lupton (1998). The PS1 images are adopted as the reference images and convolved to fit the point spread function (PSF) shape of the HSC images.
We measure the FWHM sizes of PSF in the stacked images with hscPipe. These scatters in a wide range from 0.7 to 1.8 arcsec depending on the pointings, especially on the elevation, and the median is ∼ 1.2 arcsec (Figure 2). The PSF size statistics is summarized in Table 3. The median FWHM size is slightly worse than that of the image quality of the PS1 3π survey (Magnier et al. 2016a), and the PSF convolution of the PS1 image for the image subtraction works well.
After the image subtraction, the 5σ limiting magnitudes in the difference images are estimated by measuring standard deviations of fluxes in randomly distributed apertures with a diameter of twice the FWHM of PSF, and scatter from 18.3 mag to 22.5 mag with a median of 21.3 mag ( Figure 3 and Table 4). The 5σ limiting magnitudes are mainly determined by the depths of HSC images which are typically shallower than those of the PS1 image. In particular, the depths in the pointings observed early on Aug 19 are quite shallow. We also evaluate completeness of detection by a random injection and detection of artificial point sources with various magnitude (dashed lines in Figure 4). The magnitude of artificial point sources are fixed in time. The large diversity in the depth of images taken on Aug 19 causes the shallow dependence of completeness on the PSF magnitude of artificial sources. The median of 5σ limiting magnitude is roughly comparable to the 70% completeness magnitude.
As the detected sources include many bogus detection, candidate selection is performed as done in Utsumi et al. (2017a). Criteria for the detection in a single difference image are (1) |(S/N )PSF| > 5, (2) (b/a)/(b/a)PSF > 0.65 where a and b are the lengths of the major and minor axes of a shape of a source, respectively, (3) 0.7 < FWHM/(FWHM) PSF < 1.3, and (4) PSF-subtracted residual < 3σ. These criteria confirm a high confidence level of detection and a stellar-like shape of a source. Further, we impose the sources to be detected in both of the difference images on Aug 18 and 19, and find 1551 sources. We also evaluate the completeness of this candidate selection with the artificial point sources (solid lines in Figure 4). The candidate selection makes the 50% completeness magnitudes shallower by 0.7 − 0.8 mag. The completeness of the two-epoch detection is comparable to that on Aug 19 because the observation on Aug 19 is shallower than that on Aug 18. The 50% completeness magnitude for two-epoch detection is 20.6 mag.
The two-epoch detection is only possible for the fields with the archival PS1 images and the HSC images on both of Aug 18 and 19. The resultant area for the transient search is 23.6 deg 2 corresponding to the 56.6% credible region of GW170817.

Source screening
Since the 1551 sources include sources unrelated to GW170817, we need to screen them in order to pick up candidates that may be related to GW170817. We adopt a procedure shown in the flowchart ( Figure 5).
First of all, the flux of optical counterpart of GW170817 needs not to be negative on Aug 18 and 19. We exclude sources having significantly negative fluxes (< −3σ) on Aug 18 or 19. We also rule out sources associated with stellar-like objects in the PS1 catalog (Magnier et al. 2016b;Flewelling et al. 2016) 3 with a separation of < 1.0 arcsec. Here we adopt the larger separation, similar to the typical seeing size, than the astrometric error in order to remove bogus detection that frequently appears around a bright star. According to the number density of stellarlike objects in the PS1 catalog, this exclusion reduces only 0.2% of the survey fields. After these screening, 384 sources remain. While 322 sources are located at the center of extended objects in the PS1 catalog, 62 sources have separations with > 1.0 arcsec to any objects in the PS1 catalog.    We further exclude sources associated with PS1 objects that is firmly located outside of the 3D skymap derived from the GW observations (LALInference v2.fits.gz; Abbott et al. 2017c), adopting the GLADE v2 catalog and NASA/IPAC Extragalactic Database (NED) 4 . While we primarily employ the distance in the GLADE catalog, we replace it with the redshiftindependent distance in NED if the associated PS1 objects or one of a galaxy pair containing the associated PS1 objects have information (Tully 1988;Willick et al. 1997;Freedman et al. 2001;Theureau et al. 2007;Sorce et al. 2014;Springob et al. 2014), and with the redshift-dependent distance in NED (Mould et al. 2000) if no distance information is available in the GLADE catalog. We search for possibly associated galaxies in the GLADE catalog or NED with a separation of < 2.0 arcsec for the 322 sources at the center of extended PS1 objects, which is smaller than the criteria to identify duplicate galaxies in the GLADE catalog (3.6 arcsec), and with a separation of < 15.0 arcsec for the 62 off-center sources, which corresponds to a separation of < 3 kpc at a distance of 40 Mpc. If the 3D probability of GW170817 occurrence at the location and distance of the associated PS1 object with a HEALPix 3D grid with resolution of NSIDE=1024 is less than 10 −3 of the maximum probability, the source is ruled out. This screening reduces the number of sources to 224 sources at the center of extended PS1 objects and 59 off-center sources.
There is only one source (J-GEM17btc) associated with a PS1 object that is located in the 3D skymap. The detail of J-4 https://ned.ipac.caltech.edu/ GEM17btc is described in the next subsection. On the other hand, the other 282 sources do not have distance measurement in the GLADE catalog or NED. After the catalog matching, four of the authors remove bogus and high proper motion stars by visual inspection (Figure 6). The number of final candidates, that may be related to GW170817, is 60 (Table 5). We note that 58 candidates are located at the center of extended PS1 objects and that some of them could be active galactic nuclei (AGN) or indistinguishable residuals resulting from different instrumental signatures between PS1 and HSC, but we conservatively hold them as candidates.

Properties of candidates
We investigate properties of remaining 60 candidates. Figure 7 shows the candidate with the associated PS1 object within the 3D skymap of GW170817. J-GEM17btc is located at The nearest object in the PS1 catalog is PSO J130947.744-232257.366 at R.A. = 13 h 09 m 47. s 74, decl. = −23 • 22 ′ 57. ′′ 4 (J2000.0) with a separation of 6.0 arcsec to J-GEM17btc, which is superposed on NGC 4993 and located at 4.6 arcsec north of the center of NGC 4993. According to an archival Hubble Space Telescope (HST) ACS image (Bellini et al. 2017), the PSF shape of PSO J130947.744-232257.366 is consistent with stellar-like sources surrounding it and PSO J130947.744-232257.366 is unlikely to be relevant to J-GEM17btc. Thus, we conclude that the second closest object NGC 4993 located 10.0 arcsec away from J-GEM17btc, well within the separation criterion of 15.0 arcsec, is associated with J-GEM17btc. NGC 4993 is an S0 galaxy at the distance of ∼ 40 Mpc (Freedman et al. 2001).
Among remaining 59 candidates, one candidate (J-GEM17bog) at R.A. = 13 h 04 m 44. s 11, decl. = −22 • 37 ′ 07. ′′ 2 (J2000.0) is registered as an off-center transient (Figure 8). However, we marginally find a persistent object overlapping with a galaxy in the archival PS1 i-band image, which is not registered in the PS1 catalog. The other 58 candidates are located at the center of extended PS1 objects. Among them, two candidates are associated with X-ray sources in the ROSAT catalog (Boller et al. 2016) with separations of 15.6 and 16.2 arcsec, and two other candidates are associated with radio sources in the NRAO VLA Sky Survey (NVSS, 1.4 GHz) catalog (Condon et al. 1998) with separations of 0.9 and 3.0 arcsec. These four candidates could be AGN showing optical variability. We also check the 3XMM-DR7 catalog (Rosen et al. 2016) but there are no associated sources in the 3XMM-DR7 catalog. Although some of them might have little possibility of the association with GW170817, we cannot exclude them from candidates of the optical counterpart of GW170817.
The apparent magnitude on Aug 18 and time variability between Aug 18 and 19 of the 60 candidates are shown in Figure 9. It is remarkable that J-GEM17btc is much brighter than the other candidates with ∼ 19 − 22 mag. Although this diagram illustrates the distinguished feature of J-GEM17btc, this fact alone is not conclusive evidence that J-GEM17btc is the most likely counterpart of GW170817.
For comparison, theoretical models of kilonovae at 1-5 days after the merger are also plotted in Figure 9 (Tanaka et al. 2014;Tanaka et al. 2017). Although the observation takes place at 1 day after the merger, we adopt a wide range of time after the merger in order to take into account theoretical uncertainties. Although the rapid time evolution is believed to be a clue for identification of a kilonova, the time variability at early epochs of the theoretical models can take any values between −1.0 and +1.0 mag day −1 , which is consistent with the properties of all the candidates. On the other hand, there is a discrepancy in the apparent magnitude. However, it can be explained by different ejecta masses, which could result from uncertainties of the equation of state and different efficiency of viscus heating (e.g., Shibata et al. 2017). Therefore, we cannot rule out the possibility of any candidates as an EM counterpart of GW170817 from the time variability and the brightness at 1 day after GW170817.
Since the distance to the candidates are unknown except for J-GEM17btc, we evaluate the probability P3D that the associated PS1 object is located inside of the 3D skymap of GW170817, with a luminosity function of galaxies at a rest wavelength λ, φ(λ,M ), derived from rest-frame U BV RI luminosity functions (Ilbert et al. 2005)  where Dmean and σD are the mean and standard deviation of the distance to GW170817 at the position, respectively, M (mj , D) is the absolute magnitude of a galaxy with observer-frame jband apparent magnitude mj at a distance of D, λ(λj, D) is the rest wavelength redshifted from the observed wavelength λj with a distance of D, and A(D) is the surface area of observed region at a distance of D.
We evaluated P3D for the PS1 objects associated with the 59 candidates using the rand/or i-band Kron magnitude in the PS1 catalog. We also estimate the probability of NGC 4993 (R = 12.09 mag, Lauberts & Valentijn 1989) associated with J-GEM17btc, which is P3D = 64%. On the other hand, the probability of J-GEM17bog is 9.3 × 10 −3 % and the probabilities of the other 58 candidates range from 1.2 × 10 −2 % to 2.1 × 10 −1 %. Furthermore, the possibility, that at least one of the 59 candidates including J-GEM17bog is located in the 3D skymap of GW170817, is only 3.2%. Therefore, we conclude that J-GEM17btc is more likely, by more than an order of magnitude, to be the optical counterpart of GW170817 than the other candidates. The large difference between J-GEM17btc and the other 59 candidates stems from the faintness of the associated objects of the other 59 candidates, which prevent them from being registered in the GLADE catalog or NED. Given the luminosity function of galaxy and the comoving volume, the faint objects are likely to be distant objects and thus P3D of them are small. We note that the integrand of the denominator in equation (1) is nearly zero at a redshift of z ≥ 0.7 for all of the 60 candidates and that these results are almost independent on the adopted filters.

Conclusions
We have performed the survey for the optical counterpart of GW170817 with Subaru/HSC. Our untargeted transient search covers 23.6 deg 2 corresponding to the 56.6% credible region of GW170817 and reaches the 50% completeness magnitude of 20.6 mag. We find 1551 sources with two-epoch detection, and screen them with the catalog matching and the visual inspection. The number of our final candidates is 60.
We find only one candidate J-GEM17btc with an associated object firmly located within the 3D skymap of GW170817. On the other hand, the other 59 candidates do not have distance information of associated objects. The candidates include one off-center candidate other than J-GEM17btc, but it is associated with the marginally-detected persistent object in the archival PS1 i-band image. The other 58 candidates are located at the center of extended PS1 objects and could be AGN. Four of them are actually associated with the ROSAT X-ray sources or NVSS radio sources. However, we can not rule out the other 59 candidates from our observations because the kilonova model can have any time variability of −1.0 to +1.0 mag day −1 at the early epochs.
Hence, we evaluate the probability that the PS1 object associated with the candidate is located inside of the 3D skymap of GW170817. The probability of NGC 4993 associated with J-GEM17btc is 64%, while the possibility, that at least one of the other 59 candidates is located in the 3D skymap, is only 3.2%. Therefore, we conclude that J-GEM17btc (a.k.a. SSS17a/DLT17ck) is the most-likely and distinguished candi-date as the optical counterpart of GW170817. The same conclusion is brought by the other untargeted wide-field survey with the Dark Energy Camera (DECam, Soares-Santos et al. 2017). We note that J-GEM17btc is intensively observed by many telescopes, satellites, and instruments (e.g., Abbott et al. 2017d;Utsumi et al. 2017b).