# The Host Galaxy of the Recoiling Black Hole Candidate in 3C 186: An Old
Major Merger Remnant at the Center of a z=1 Cluster

T. Morishita  Space Telescope Science Institute, 3700 San Martin Drive,
Baltimore, MD 21218, USA  IPAC, California Institute of Technology, MC 314-6,
1200 E. California Boulevard, Pasadena, CA 91125, USA;
<EMAIL_ADDRESS>M. Chiaberge  Space Telescope Science Institute for
the European Space Agency (ESA), ESA Office, 3700 San Martin Drive, Baltimore,
MD 21218, USA  The William H. Miller III Department of Physics and Astronomy,
Johns Hopkins University, Baltimore, MD 21218, USA  B. Hilbert  Space
Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA
E. Lambrides  The William H. Miller III Department of Physics and Astronomy,
Johns Hopkins University, Baltimore, MD 21218, USA  L. Blecha  University of
Florida, Department of Physics, 2001 Museum Rd., Gainesville, FL 32611, USA
S. Baum  Department of Physics and Astronomy, University of Manitoba,
Manitoba, Canada  S. Bianchi  Dipartimento di Matematica e Fisica, Università
degli Studi Roma Tre, via della Vasca Navale 84, 00146 Roma, Italy  A. Capetti
INAF - Osservatorio Astrofisico di Torino, Pino Torinese, Italy  G. Castignani
Dipartimento di Fisica e Astronomia “Augusto Righi”, Alma Mater Studiorum
Università di Bologna, Via Gobetti 93/2, I-40129 Bologna, Italy  INAF -
Osservatorio di Astrofisica e Scienza dello Spazio di Bologna, via Gobetti
93/3, I-40129, Bologna, Italy  F. D. Macchetto  Space Telescope Science
Institute, 3700 San Martin Drive, Baltimore, MD 21218, USA  G. K. Miley
Leiden Observatory, University of Leiden, PO Box 9513, 2300 RA Leiden, The
Netherlands  C. P. O’Dea  Department of Physics and Astronomy, University of
Manitoba, Manitoba, Canada  C. A. Norman  The William H. Miller III Department
of Physics and Astronomy, Johns Hopkins University, Baltimore, MD 21218, USA
Space Telescope Science Institute, 3700 San Martin Drive, Baltimore, MD 21218,
USA

(Accepted April 23, 2022; Received )

###### Abstract

3C 186, a radio-loud quasar at $z=1.0685$, was previously reported to have
both velocity and spatial offsets from its host galaxy, and has been
considered as a promising candidate for a gravitational wave recoiling black
hole triggered by a black hole merger. Another possible scenario is that is in
an on-going galaxy merger, exhibiting a temporary displacement. In this study,
we present analyses of new deep HST/WFC3-IR and ACS images, aiming to
characterize the host galaxy and test this alternative scenario. We carefully
measure the light-weighted center of the host and reveal a significant spatial
offset from the quasar core ($11.1\pm 0.1$ kpc). The direction of the
confirmed offset aligns almost perpendicularly to the radio jet. We do not
find evidence of a recent merger, such as a young starburst in disturbed
outskirts, but only marginal light concentration in F160W at $\sim 30$ kpc.
The host consists of matured ($\,\hbox to0.0pt{\lower 3.5pt\hbox{$\mathchar
0\relax\sim$}\hss}\raise 1.0pt\hbox{$>$}\,200$ Myr) stellar populations and
one compact star-forming region. We compare with hydro-dynamical simulations
and find that those observed features are consistently seen in late-stage
merger remnants. Taken together, those pieces of evidence indicate that the
system is not an on-going/young merger remnant, suggesting that the recoiling
black hole scenario is still a plausible explanation for the puzzling nature
of

quasars: galaxies: active – galaxies: jets – gravitational waves - quasars:
individual (3C 186)

††journal: ApJ

##  1 Introduction

In the current paradigm of $\Lambda$-CDM cosmology, galaxy-galaxy mergers
carry many important roles in galaxy evolution (Volonteri et al., 2003;
Springel, 2005). Supermassive black hole (SMBH) merger is one of such expected
as a result of major galaxy mergers. To explain tight relationships such as
the $M$–$\sigma$ relation (Ferrarese & Merritt, 2000; Gebhardt et al., 2000),
an intimate relationship between galaxy growth and black hole growth during
mergers seems to be required (e.g., Peng, 2007).

Black hole mergers are thought to happen in three phases: the two black holes
are pulled towards the center of the gravitational potential of the merged
galaxy by dynamical friction; they become a binary system by losing angular
momentum via gravitational slingshot interaction with stars that have
appropriate angular momentum in the region of the parameter space (the so-
called loss cone; Frank & Rees, 1976; Begelman et al., 1980) and gas-driven
inspiral (e.g., Escala et al., 2005; Dotti et al., 2007); in the final phase,
the bound pair may lose the remaining angular momentum via the emission of
gravitational waves (GWs), and eventually, the two black holes coalesce. In
most cases, these GWs are emitted anisotropically. Depending on both the
relative orientation of the spins of the merging black holes and their mass
ratio, the merged black hole may receive a recoil kick (Centrella et al.,
2010; Blecha et al., 2011; Komossa, 2012) with velocities as high as $\sim
5000$ km s-1 (e.g., Campanelli et al., 2007; Tichy & Marronetti, 2007; Lousto
& Zlochower, 2011). As a result, the merged black hole may show displacement
from the center of the host galaxy both in position and velocity.

However, the details of the mechanisms that bring the two black holes to the
distance at which GW emission becomes substantial are still poorly understood.
This is called the final-parsec problem (e.g., Milosavljević & Merritt, 2003).
For example, if the loss cone is not replenished quickly enough for the SMBH
pair to get sufficiently close and lose energy via GWs, the pair may stall and
never merge; those may be observed as dual quasars or pairs of less active AGN
(e.g., Komossa et al., 2003; Hennawi et al., 2006; Goulding et al., 2019;
Silverman et al., 2020; O’Neill et al., 2021). If this is often the case and
prevents SMBH pairs to merge, our current understanding of the mechanisms for
black hole growth may need significant revisions. It is thus important to
find, either direct (GWs) or indirect (recoiling black holes), evidence of
SMBH mergers. While we still need to wait for next generation space-based GW
interferometers and pulsar timing arrays for the direct confirmation at the
black hole mass range of $M_{\rm BH}\,\hbox to0.0pt{\lower
3.5pt\hbox{$\mathchar 0\relax\sim$}\hss}\raise
1.0pt\hbox{$>$}\,10^{8}M_{\odot}$ (see, e.g., Moore et al., 2015), identifying
recoiling black holes helps us to advance our understanding of such extreme
events.

As of today, there are several candidates of recoiling black holes, with
velocity (Komossa et al., 2008; Steinhardt et al., 2012; Comerford & Greene,
2014; Pesce et al., 2018), spatial (Batcheldor et al., 2010; Koss et al.,
2014; Lena et al., 2014; Barrows et al., 2016; Skipper & Browne, 2018), or
both offsets (Pesce et al., 2021; Hogg et al., 2021). In particular, high-
velocity offsets ($>1000$ km s-1) are expected to be rare, but they are more
likely to be observed in combination with large spatial offsets (e.g., Blecha
et al., 2016) and thus can be used as a proxy for ideal followup targets.

For example, CID-42, a radio-quiet AGN at $z=0.36$, shows a velocity offset
between the narrow and broad components of the ${\rm H\beta}$ emission line of
$\sim 1300$ km/s, with two nuclei in HST images displaced by $\sim 2.5$ kpc in
projected distance (Civano et al., 2010, 2012; Blecha et al., 2013; Novak et
al., 2015). The presence of a second, obscured radio quiet AGN, however,
cannot be excluded (but see also Kim et al., 2017). Another intriguing case is
found in A2261-BCG, the brightest cluster galaxy (BCG) in Abell 2261 at
$z=0.2248$. This galaxy was found to have an unusually flat core profile, with
an offset of $\sim 0.7$ kpc from the photocenter of the surrounding envelope
(Postman et al., 2012). The followup study presented spectroscopic
measurements of three of the four knots identified in the central region, to
test their hypothesis that the massive black hole was ejected from the core
and harbour within one of those knots (Burke-Spolaor et al., 2017). Their
observations could not identify the precise location of the expected SMBH
($\sim 10^{10}\,M_{\odot}$), leaving the conclusion still pending.

is another promising candidate, but by far at higher redshift than others. The
system is a radio-loud quasar ($\sim 10^{5}$ yrs; Murgia et al., 1999;
Siemiginowska et al., 2005) of the compact-steep spectrum class (Fanti et al.,
1985; O’Dea, 1998; O’Dea & Saikia, 2021), located in the center of a cluster
of galaxies at $z=1.0685$, with a SMBH of $\sim 3$-$6\times 10^{9}\,M_{\odot}$
(Kuraszkiewicz et al., 2002; Chiaberge et al., 2017).

was recently observed as one of 22 3C radio galaxies/quasars in a snapshot
program of HST (Hilbert et al., 2016). Chiaberge et al. (2017) reported its
spatial offset by $1.3\pm 0.1$ arcsec, or $\sim 11$ kpc, with respect to the
center of the host galaxy seen in their medium-deep image of WFC3-IR F140W. In
the same study, was also revealed to have a significant velocity offset ($\sim
2140$ km/s) of broad lines associated with the quasar with respect to the
systemic redshift of the host measured by two narrow lines, [O ii] and [Ne
iii], suggestive of a binary black hole merger of a mass ratio
$m_{1}/m_{2}>0.25$ (Lousto et al., 2017). In addition, the followup
observations with Keck/OSIRIS integral field unit (IFU) revealed a tentative
but significant velocity shift ($>1800$ km/s) in the broad component of ${\rm
H\beta}$ line from narrow lines (Chiaberge et al., 2018, see also Sec. 4.4).

While all previous results support the interpretation of as a GW recoiling
black hole, there are two other possible scenarios that may explain the
observed properties; one is that we are seeing an on-going galaxy merger,
where its photocenter is temporarily displaced due to disturbances in the host
morphology; the other scenario is that the quasar is associated with an under-
massive galaxy that is superimposed to another brighter galaxy.

With the depth of the previous F140W image, those scenarios were not
completely excluded, as tidal features and recent star formation in the
outskirts can easily be missed (e.g., Koss et al., 2018). Deeper imaging with
HST is thus important to unveil any evidence of on-going merger and constrain
the current stage of the system.

The aim of this paper is to test these alternative hypotheses by analysing
deeper IR and optical images newly taken in the HST Cycle 25. In particular,
high photometric sensitivity of HST at $\sim 1$ kpc resolution provides
details of the underlining stellar population of the host galaxy and allows us
to further constrain the time scale since the last major merger and associated
star formation.

Throughout, we adopt the AB magnitude system (Oke & Gunn, 1983; Fukugita et
al., 1996), cosmological parameters of $\Omega_{m}=0.3$,
$\Omega_{\Lambda}=0.7$, $H_{0}=70\,{\rm km~{}s^{-1}}\,{\rm Mpc}^{-1}$, and the
Salpeter (1955) initial mass function. We refer to magnitude for the HST
filters used in this study, F606W, F110W, and F160W as $V_{606}$, $J_{110}$,
and $H_{160}$ respectively.

Figure 1: Pseudo color image of the field (F606W for blue, F110W green, F160W red). A zoomed image of ($\sim 27\times 27$ arcsec2) is shown in the two right panels, in the original (top) and PSF-subtracted (bottom) images. The compact star-forming blob, SF Blob1 (Sec. 3.3) is indicated in the right-bottom panel. The images have been rotated so that north is up and east to the left.  Table 1: Summary of HST observations in the field. Instrument | Filter | Total Exposure Time
---|---|---
|  | sec
ACS-WFC | F606W | 5567
WFC3-IR | F110W | 5042
WFC3-IR | F160W | 7563

##  2 Data

###  2.1 HST Observations and Data Reduction

was observed in a GO program (PID 15254; PI. M. Chiaberge) during the HST
Cycle 25, with WFC3-IR and ACS. Two filters, F110W and F160W, are used in the
WFC3 observations and one filter, F606W, in the ACS observations (Table 1).
The observing position angle was carefully designed during the phase II
process, so to minimize possible contamination on the host galaxy by stellar
spikes from the quasar and surrounding objects. To avoid saturation of a few
pixels in the quasar core, we started with one short exposure ($\sim 30$ sec)
in each of the 5 orbit ACS observations. Then, the following exposures were
taken at the same location and rotation angle of the first dither point in
each orbit, to ensure an identical PSF profile over the entire exposures.

To reduce imaging data, we follow the procedure presented in Hilbert et al.
(2016). For ACS data, we begin with files in the flc format. These images are
corrected for charge transfer efficiency (CTE) effects through the CTE
correction algorithm (Anderson & Bedin, 2010). For WFC3-IR data, we begin
files in the flt format, as the multiple non-destructive detector readouts
within each observation allowed for easy cosmic ray identification and removal
within the calwf3 data reduction pipeline. However, pixels subjected to high
flux levels in preceding observations are affected by persistence and need to
be corrected. To remove any persistence signal present in our data, we begin
by retrieving the persistence masks and persistence-corrected flt files of our
observations from MAST. These files have had the persistence signal modeled
and subtracted from them following the model described in Section 8.3 of
version 4.0 of the WFC3 Data Handbook, and are therefore different from the
standard flt files available in the archive. Before proceeding with the data
reduction, we also manually mask the satellite trail present in one of our
images.

Since our focus in this study is the host galaxy of we attempt to subtract the
light from the bright quasar nucleus in each distortion uncorrected image
(those in flt/flc formats) before the geometric correction and final stacking
steps (see Sec. 2.2).

The next step in our data reduction is to remove the geometric distortion from
all of the PSF-subtracted flc/flt files and to combine images into a final
image for each filter. We use Astrodrizzle to accomplish both of these steps.

For each filter, we combine the PSF-subtracted images into a final image with
a pixel scale of 0.045”/pixel. We find a value of ”5.5 5.0” for the
_driz_cr_snr_ parameter results in good cosmic ray rejection during the
stacking process. We also calculate an optimal _final_pixfrac_ value of 0.75
for our data and final pixel scale.

We then use Tweakreg (also part of the Drizzlepac software package) to align
this final image to the same world coordinate system (WCS) present in the
corresponding ACS drizzled image. At this point, it is possible to overlay the
ACS and WFC3 images for a particular object and compare the morphology and
brightness in the two observation bands.

The resulting final image is shown in Fig. 1. The PSF-subtracted image reveals
a clear spatial offset between the quasar and the host galaxy (Sec. 3.1). Our
deep images also reveal faint outskirts of the host, out to $\sim 60$ kpc.

###  2.2 PSF Subtraction

In order to subtract the PSF component without over-subtracting the host light
component, we fit the observed light profile with PSF+single Sérsic profile
using GALFIT (Peng et al., 2002a). In an effort to achieve the best model
fits, we experiment with the GALFIT inputs. We create a custom bad pixel file
for each exposure, where we remove all nearby astrophysical sources as well as
bad pixels and cosmic rays. By supplying this file to GALFIT, we reduce the
chances that the final fit is contaminated by nearby sources.

The PSFs used for our PSF subtractions are created using TinyTim software
(Krist et al., 2011). We use TinyTim PSFs due to a lack of unsaturated, high-
SNR, isolated PSFs in our data to use as models. In order to create realistic
PSFs, we used characteristics from our observations when creating the model
PSFs. We begin by finding the mean secondary mirror focus offset for the date
and time of our observation. This information is collected from the HST focus
model website111http://www.stsci.edu/hst/observatory/focus/FocusModel.

Other information used as input to TinyTim includes the HST instrument and
filter, the target’s pixel location on the detector, the amount of jitter (0
for this work), and the size in arcseconds of the output model PSF (7
arcseconds). Finally, we also provide a quasar spectrum, in order for Tiny Tim
to determine wavelength-dependent PSF characteristics. For simplicity, we used
the composite quasar near-IR spectrum from Glikman et al. (2007) redshifted to
$z\sim 1$. With this information in hand for each exposure, we run all three
stages of TinyTim. The final output is one geometrically distorted model PSF
corresponding to each of flc/flt images.

We also fit the data using a PSF at the nominal pixel scale, as well as a
supersampled PSF. We find that the fit with the supersampled PSF tends to
over-subtract the source flux, and therefore keep to the nominal pixel scale
for our final fittings. We also fit the data with TinyTim PSFs from a range of
focus values, in order to examine any differences. The best modeling and
subtraction came from values that match the focus values indicated by the
focus model website. Finally, as the data being fit are in units of electrons
per second, rather than the electrons that GALFIT assumes, we also experiment
with fits where we provide the ERR array of the flt/flc file to GALFIT as a
sigma image. We find in this case that the resulting fits again tend to over-
subtract the source’s flux. With our final set of GALFIT parameters and
inputs, we obtain fits with reduced-$\chi^{2}$ between 0.75 and 1.5. The best-
fit PSF component is then subtracted from each of the images before they are
processed for geometric correction and final stacking with Astrodrizzle.

Figure 2:  Left: PSF-subtracted image of in F160W, in a $500\times 500$ pixel
postage stamp. Regions that are masked and patched are indicated with purple
apertures. Middle: One of the patched images, where all masked regions are
filled by random values taken from nearby pixels (Sec. 3.1). The median light-
weighted center (red point), calculated over 100 such realizations, is clearly
off from the flux peak position of the quasar (cyan cross). Fitted ellipses
are shown (orange lines). Light concentration is seen in the South-West
direction from the quasar. Right: Zoom-in image of the middle panel with
different color stretch.

##  3 Analysis and Results

###  3.1 Light-weighted center of the host galaxy

We estimate the center of mass of the host galaxy by using F160W, which is the
deepest among our images here. It is noted that the light distribution in the
F160W image is considered to trace the stellar component of the host well, as
at the redshift of the filter corresponds to rest-frame $\sim 0.8\,\mu$m and
is clean from strong emission lines (see also Sec. 3.2).

While the extended light from the host galaxy is clearly visible in the PSF-
subtracted image (Fig. 1), it is still challenging to accurately estimate its
light-weighted center (LWC) due to the fact that fore/background objects are
superposed on the host. Especially, since is in a crowded field, the latter
remains as a critical issue because flux contribution from other sources can
easily affect the estimate.

To eliminate such flux contamination, we first manually mask regions that are
likely to belong to surrounding objects (including the compact star-forming
blob discussed in Sec. 3.3) and residuals at the position of quasar, by
placing circular apertures of arbitrary sizes in the PSF-subtracted F160W
image. By doing so, we ensure to mask all possible contaminants, whereas
automated software (e.g., SExtractor; Bertin & Arnouts 1996) may not
necessarily recognize all objects in crowded regions or mistakenly mask the
part of the host. The image with those aperture masks is shown in the left
panel of Fig. 2. To mask the residuals at the location of the quasar, we place
an aperture with radius of $r=0.\\!^{\prime\prime}9$.

The next step is then to fill those masked regions. This process is necessary,
because leaving those masked pixels could result in a biased estimate of LWC
depending on the choice of position and size used of each circular mask. For
each pixel of the masked regions, we use a Metropolis-Hastings-like algorithm
to randomly select one of the neighboring pixels that are not masked, to fill
in the target pixel. This filling process is as follows;

  * 1\. For the pixel of interest in a masked region, assign a radius ($r_{rand}$) randomly taken from a range of $[0,r_{max}]$, where we here set $r_{max}=50$ pixel.

  * 2\. Compare the probability of $r_{rand}$, assigned by $p(r)=\exp{(-r/r_{max})}$, with a random float number, $p_{r}$, drawn from a range of [0,1].

  * 3\. If $p(r_{rand})>p_{r}$, then select a random pixel from unmasked regions at $r_{rand}<r<r+\Delta r$, where we set $\Delta r=2$ pixel; otherwise return to 1 and repeat the process until it passes the criteria.

  * 4\. Repeat step 1-3 for all masked pixels.

  * 5\. Repeat step 1-4 for $n_{iter}$ times, to get $n_{iter}$ images with the masked regions filled respectively.

We here set $n_{iter}=100$. By adopting an exponential function for $p(r)$,
the filling flux is likely to be chosen from one of the nearby pixels, making
the reproduced light profile contiguous, which is reasonable given the primary
purpose here. One of the realizations is shown in Fig. 2 as an example. While
several pixels can occasionally have a distinct value from the average value
of surrounding pixels, contributions from such pixels to the determination of
LWC are smoothed out after repeating $n_{iter}$ times of realization.

At each realization of $n_{iter}$, we calculate the LWC by the following
equations;

$\displaystyle\overline{x}={{\sum xf_{\rm 160}(x,y)/e_{\rm
160}(x,y)^{2}}\over{\sum f_{\rm 160}(x,y)/e_{\rm 160}(x,y)^{2}}},$ (1)
$\displaystyle\overline{y}={{\sum yf_{\rm 160}(x,y)/e_{\rm
160}(x,y)^{2}}\over{\sum f_{\rm 160}(x,y)/e_{\rm 160}(x,y)^{2}}}$

where $f_{\rm 160}(x,y)$ and $e_{\rm 160}(x,y)$ represent flux value and its
associated RMS value at the position of [x,y], and the sum goes over pixels
with $f_{\rm 160}/e_{\rm 160}>{\rm SNR_{limit}}$ within $r_{\rm lim}<250$
pixel ($\sim 91$ kpc) from the target. We here set ${\rm SNR_{limit}}=0.7$,
which is equivalent to the detection threshold (DETECT_THRESH of SExtractor)
used for source detection in the F160W image. We empirically find that the
value in general provides reasonable detection out to faint outskirt of
sources without being contaminated by background pixels.

Figure 3: Central region of The RGB image is in the same configuration as in
Fig. 1. The position of the light-weighted center is shown (red circle). The
8.44 GHz radio continuum emission, retrieved from the NRAO VLA archive
(project ID AA0129; cyan contour), is overlaid. Two dashed lines are shown to
highlight the inferred axes of the offset (red) and the jet (cyan),
respectively. The image has been rotated so that north is up and east to the
left.

We calculate the median LWC over $n_{iter}$ realizations and find
$[x,y]=[-29.42_{-0.08}^{+0.07},7.44_{-0.05}^{+0.06}]$ (in pixel) with respect
to the position of in the image of Fig. 2, or 7:44:17.586,$+$37:53:18.165 in
the sky coordinate. The 16th/84th percentiles are associated as uncertainties.
The result reveals a significant offset $1.37\pm 0.01$ arcsec, or $11.1\pm
0.1$ kpc, in the projected distance, which is consistent with the estimate
with a shallower F140W image by Chiaberge et al. (2017). The angle from the
LWC toward the quasar core, measured from the west to the clockwise direction,
is $\sim+35^{\circ}$, which lies almost perpendicularly to the direction of
the radio jet (Fig. 3). We will discuss this in Sec. 4.1.

The choice of $r_{\rm lim}$ and the limiting SNR has only a minor effect on
the estimate as long as ${\rm SNR_{limit}}$ is set to $>0.7$; having lower
limiting SNR would significantly increase the estimated error due to noisy
pixels across the entire image. On the other hand, the inferred offset becomes
even larger when a higher value for ${\rm SNR_{limit}}$ is chosen. As in the
middle panel of Fig. 2, ellipse fit to the F160W light profile shows a shift
of the centroid toward the North-East direction when only the central part of
is used. This is due to the light concentration in the outskirts. We discuss
the light concentration in relation to merger-induced morphological
disturbances, or lack thereof, in Sec. 4.

Despite the depth that our new images reach, we might still be missing very
low-surface brightness components ($>25$ mag arcsec-2). However, we argue that
such light components only marginally affect the determined photocenter. Our
images here already reveal a substantial amount of stellar mass222Such an
amount of mass may be found in the diffuse stellar components but only in the
most massive clusters (e.g., Morishita et al., 2017). ($\log
M_{*}/M_{\odot}\sim 11.4$; see Sec. 3.2). The amount is $\times 10$ more than
what is found in the diffuse outskirts of low-$z$ elliptical galaxies (e.g.,
Huang et al., 2018).

###  3.2 Stellar populations of the host galaxy

We aim to investigate the underlying stellar population in the host. The
spatial distribution of stellar age enables us to estimate the timeline since
the last merger and associated star formation activity.

To get sufficient SNR for reliable estimates out to faint outskirts of the
host, we define sub-regions of the host by using a Voronoi tessellation method
(Cappellari & Copin, 2003). We set the minimum SNR to $15$, so that every
tessellated region has SNR greater than this value. As a result, at $\sim 6$
arcsec ($\sim 48$ kpc) from the system LWC we reach down to $\sim 25$ mag
arcsec-2 in F160W. We set the boundary of the host by using the segmentation
mask generated by SExtractor with the detection threshold of $0.7$, which is
consistent with the threshold used for the LWC calculation above. Note that in
this analysis, as it is not necessary for our purposes, we do not fill the
masked regions defined in Sec. 3.1 and leave those empty. We show the defined
mask and radial distribution of tessellated pixels in Fig. 4.

The boundary for tessellation defined in the F160W image are then applied to
both F110W and F606W images, to collect flux from the same sub-regions. Each
of the F110W and F606W images is convolved to the PSF size of the F160W image
beforehand. The convolution kernels were generated by providing median-stacked
PSFs to a python software, pypher (Boucaud et al., 2016), in the same way as
in Morishita (2021). Fluxes of each sub-region are summed and compiled into a
catalog for a spectral energy distribution (SED) analysis in the following
section.

Figure 4:  Left: Segmentation map, where each tessellated region is color-
coded differently, overlaid on the F160W image. IDs of each segment used in
the main text are shown (also, Table 2). The position of is marked by a cyan
cross symbol. Right: Central region of the same Voronoi segmentation map. The
region of the blue blob at [x,y] = [210,260] (Sec. 3.3) is excluded from the
host population analysis.

Figure 5:  Top: Schematic of our star formation history model, where star
formation initiates $T_{0}$ ago, constantly lasts for over time of $\tau$, and
then terminates at $t_{Q}$ ago. Bottom: Examples of SED fitting results of
three regions: Voronoi ID 4 (top), Voronoi ID 13 (middle), and SFB1 (bottom).
The best-fit SEDs (solid lines) and convolved data points (red squares) are
shown. Despite the number of observed data points (3; blue diamonds), the data
constrain the model by capturing the Balmer break ($V_{\rm 606}-J_{\rm 110}$;
sensitive to stellar age) and optical slope ($J_{\rm 110}-H_{\rm 160}$;
sensitive to dust). It is noted that at the redshift of F160W is clean from
strong emission lines such as ${\rm H\alpha}$, ${\rm H\beta}$, and [O iii],
making it a good tracer of the stellar component.

The flux catalog collected from the tessellated images is provided to an SED
fitting software, FAST (Kriek et al., 2009). Due to the number of broadband
filters available here (3), we are limited to a small number of fitting
parameters. Since our primary focus here is to estimate the epoch of the last
significant star formation (which could potentially be associated with major
merger), we simply assume a single burst started at $T_{0}$, with a length of
$\tau$. The schematic of this star formation model is shown in Fig. 5.

We use the Bruzual & Charlot (2003) stellar population model, with the
redshift fixed to the one estimated with the SDSS DR6, $z=1.07$ (Hewett &
Wild, 2010), and the Calzetti et al. (2000) dust attenuation model. It is
generally challenging to determine metallicity with only broad band
photometry, and thus we fix it to solar. This is reasonable given that
galaxies in this mass range are already metal-enriched, as seen at similar or
earlier epochs (e.g., Thomas et al., 2010; Onodera et al., 2015; Morishita et
al., 2019; Kriek et al., 2019).

In Fig. 5, we show examples of the SED results in three different regions; the
outer region ($\sim 0.9$ arcsec from the LWC; Voronoi ID 4), the inner ($\sim
0.2$ arcsec from the LWC; Voronoi ID 13), and the blue blob (see Sec. 3.3).
These examples demonstrate that three data points are sufficient to constrain
relative differences in age of different regions, by capturing both the Balmer
break by $V_{\rm 606}-J_{\rm 110}$ and the rest-optical slope by $J_{\rm
110}-H_{\rm 160}$. The best-fit parameters for all tessellated regions are
summarized in Table. 2.

A caveat is that with the number of data points used, degeneracy between age
and dust is still not completely resolved, which makes inference on absolute
values challenging. Our SED fitting results indicate small to moderate
attenuation ($A_{V}\,\hbox to0.0pt{\lower 3.5pt\hbox{$\mathchar
0\relax\sim$}\hss}\raise 1.0pt\hbox{$<$}\,2.0$ mag). However, in the case of
the host we do not expect significant dust attenuation, because no dust lane-
like feature is seen in our HST images — such feature would be clearly seen in
the presence of significant dust attenuation (see, e.g., Skipper & Browne,
2018; Hogg et al., 2021, for offset AGN candidates with clear dust lanes). In
Fig. 6, we show the two-dimensional distribution of $J_{\rm 110}-H_{\rm 160}$
color of the host. While the color does not necessarily represent the exact
dust attenuation due to the degeneracy, from the smooth morphology across the
host it is unlikely to expect significant dust structures.

Figure 6:  Two-dimensional $J_{\rm 110}-H_{\rm 160}$ color map of the host of
Note that the quasar component is subtracted from each of the images used
here. The image exhibits smooth morphology across the host, except two blue
regions near the quasar (mostly residual) and the blue blob (Sec. 3.3). The
outer regions of the host are excluded and set to blank.

In addition, Podigachoski et al. (2015) reported tentative ($<3\sigma$) or
non-detection of in FIR bands, at 160, 250, 350, and 500 $\mu$m by the
Herschel Space Observatory, implying no significant amount of dust is present
in the system. Fixing dust attenuation to zero would increase the age
estimates, to account for red components; therefore, our age estimates,
especially those with non-zero $A_{V}$, should be taken as lower limits.

In Fig. 7, we show the spatial distribution of each of the derived parameters.
Our interest here is $T_{0}$, the time at which primary star formation
started, and $t_{Q}$, which we introduce to refer to duration since star
formation was truncated ($\equiv T_{0}-\tau$). There is no clear radial trend
for $T_{0}$, with a median value of $160_{-40}^{+200}$ Myr, with the
associating errors capturing the 16 / 84 th percentiles. Similarly, no clear
spatial correlation nor radial trend is observed for $t_{Q}$. The parameter
ranges from $\sim 30$ Myr to $\sim 500$ Myr, with a median value of
$t_{Q}=90_{-30}^{+110}$ Myr, except a few sub-regions with on-going star
formation (i.e. $t_{Q}\leq 0$). The total mass, derived by summing the best-
fit estimate over the Voronoi segments, is found to be $\log
M_{*}/M_{\odot}=11.4$.

Figure 7:  Spatial and radial distributions of two parameters, $T_{0}$ (Top
rows) and $t_{Q}$ (bottom) are shown. Colors in the radial distribution plots
correspond to those in the spatial distribution plots. Four of the tessellated
pixels show on-going star formation (i.e. $t_{Q}=0$) and are shown at an
arbitrary value (blue triangle).

###  3.3 Compact blue star-forming region in the host galaxy

As we highlight in the color-composite image of Fig. 1, there is a compact
star-forming blob (SFB1 hereafter) located at $\sim 2$ arcsec from the quasar
to the East-North direction. This blob was seen in the early data and reported
in Chiaberge et al. (2017), but its properties remained unclear. Our SED
modeling revealed its age, $T_{0}\sim 160$ Myr, and stellar mass, $\sim
8\times 10^{9}\,M_{\odot}$ (Fig. 5). We estimate its star formation rate,
$65\pm 20\,M_{\odot}{\rm yr^{-1}}$, by averaging over the last 100 Myr of the
best-fit star formation history, which is consistent with the upper limit of
the whole system derived from non-detection in Herschel data,
$<80\,M_{\odot}\mathrm{yr^{-1}}$ (Podigachoski et al., 2015).

The best-fit SFR converts to star formation rate surface density of $\sim
1.2\,M_{\odot}{\rm yr^{-1}kpc^{-2}}$. Considering its stellar surface density
($\log\Sigma_{*}\sim 8.1\,M_{\odot}{\rm kpc^{-2}}$), this value is $\sim 0.8$
dex higher than the average value of resolved sub-regions in star-forming
galaxies at $z\sim 1$ (e.g., Wuyts et al., 2013), indicating intense star
formation activity may still be ongoing in the host despite relatively old
stellar populations elsewhere. Such intense star formation is actually seen in
a post-merged system (e.g., Barrows et al., 2018) but lasts only for a short
time scale as $\sim 100$ Myr. We revisit this in Sec. 4.

Since the template library used in the SED fitting in Section 3.2 does not
include emission lines, we check the fidelity by running another SED fitting
code, gsf (Morishita et al., 2019), with emission lines included. While an
exponential declining star formation history is used this time, due to the
limitation of the code, it derives a consistent value for star formation rate,
$70\pm 10\,M_{\odot}{\rm yr^{-1}}$.

While the IFU observations presented in Chiaberge et al. (2018) covered SFB1,
no line was detected, as the observations were originally designed for the
brightness of the quasar and exposure was not deep enough to detect any
emission lines in SFB1. Recent observations at NOEMA detected extended CO
emission near the flux peak position of SFB1 at the same redshift as the
optical narrow lines (Castignani et al., 2022).

##  4 Discussion

###  4.1 On the dynamics of the GW recoil

As we introduced in Sec 1, based on its spectroscopic and imaging properties,
was proposed as a clear example of a GW recoiling black hole candidate by
Chiaberge et al. (2017). With our deep HST observations, we both confirmed and
more accurately determined the presence of a spatial offset between the quasar
and the host photocenter. Based on our results (Sec. 3.1), we see that the
direction of the putative GW kick lies perpendicularly with respect to the
radio jet axis. In Fig. 3, we show the 8.44 GHz radio continuum emission
overlaid onto the central region of the host galaxy. The direction of the kick
is at the position angle of $\sim+33^{\circ}$, while the radio jet axis has a
position angle of $\sim-57^{\circ}$. As pointed out by, e.g., Lousto et al.
(2012), the maximum velocity of GW kick ($>2000$ km/s) may be more likely
achieved if the direction of the kick is aligned with the angular momentum of
the orbital plane of the merging black hole binary. Since the radio jet axis
likely indicates the orientation of the spin of the merged black hole, these
results may impose further constraints on the properties of the progenitor
binary, in particular with respect to the spin amplitude and relative
orientation. A detailed modeling of the binary black hole system and the GW
recoil kick lies beyond the scope of this work.

Figure 8:  A sequence of a merging system VELA04 (Camera05) at $2.0<z<2.7$,
taken from the VELA-Sunrise simulation suite. Redshift, cosmic time, and
specific star formation rate of each frame are shown in the inset. The RGB
composite images consist of F814W, F125W, and F160W, to compensate its higher
redshift than s. One of the star-forming regions (see Sec. 4) is indicated by
a blue arrow in the last frame.

###  4.2 Insight into the Current Evolutionary Stage

In this subsection, we discuss the main results of this work, which focuses on
the observed properties of the host. This allows us to speculate the current
evolutionary stage of this system. Chiaberge et al. (2017) reported a possible
shell or tidal feature in the host light profile and speculated its origin as
a major merger that occurred about $\sim 1$ Gyr ago from comparison with
numerical simulations. Our F160W image presented here revealed even fainter
stellar components out to $\sim 60$ kpc from the LWC of the host. While a
slight light concentration is seen at the outermost part in the south-east
direction, we did not see clear evidence of disturbed morphology or young star
clusters in its outskirts, which would be characteristic of on-going merger
(e.g., Mulia et al., 2015), implying that the system is a relatively old,
post-merged remnant. Furthermore, Podigachoski et al. (2015) derived high SFR
in other 3C galaxies at similar redshifts, most of which are associated with
on-going mergers, whereas an upper limit was estimated for 3C186 from non-
detection in their Herschel data ($<80\,M_{\odot}\mathrm{yr}^{-1}$). This
suggests that is relatively older than other 3C systems in terms of recent
star formation too.

To further investigate the current evolutionary phase of we compare the
observed properties with hydro-dynamical simulations. In Fig. 8, we show a
time sequence of a merging system from the VELA-Sunrise project (Simons et
al., 2019). The VELA-Sunrise project provides a set of mock images of
simulated galaxies that were originally taken from the VELA simulation suite
(Ceverino et al., 2014; Zolotov et al., 2015). The simulation suite consists
of a variety of galaxies at different evolutionary stages, including both
isolated and merging systems. While the original simulations do not include
specific recipes for black hole mergers, the datasets offer detailed
morphological information of, e.g., a galaxy during merger sequences, at the
resolution of HST. Here we select VELA04, one of the merging systems available
in the VELA suite. Specifically, the dataset of VELA04 captures the timeline
of pre- to post-merger and allows us to investigate the evolutionary stage of
by comparing their morphological properties.

The example shown in Fig. 8 is a sequence that captures the pre-merger to
post-merged phases of two galaxies at $z\sim 2$. The sequence exhibits several
features that resemble the observed properties of the host. The last frame of
Fig. 8 captures elongated (but not disturbed) morphology of the host galaxy
and compact star-forming regions in the main body of the host. We also see the
presence of young, compact blue regions across the entire sequence. Those blue
regions do not sustain more than one simulation frame, implying that their
life time is less than its time resolution, $\,\hbox to0.0pt{\lower
3.5pt\hbox{$\mathchar 0\relax\sim$}\hss}\raise 1.0pt\hbox{$<$}\,200$ Myr. In
the last frame of the sequence in Fig. 8, we highlight one such compact blue
region. Its specific star formation rate is $\sim 40$ Gyr-1 (cf. $\sim 8$
Gyr-1 of the blue blob in . The rest of the stellar populations of the galaxy
in the frame is relatively old, with the median age of $\sim 1$ Gyr (cf.
$\,\hbox to0.0pt{\lower 3.5pt\hbox{$\mathchar 0\relax\sim$}\hss}\raise
1.0pt\hbox{$>$}\,200$ Myr of . While we also notice morphological differences
between the two systems (e.g., more clumps and less relaxed), those physical
properties are broadly consistent with the observed properties of (Sec. 3.3).

The spatial offset found in this study is encouraging regarding a typical
timescale after a recoiling event. By using the derived spatial offset of the
quasar from the LWC (Sec 3.1) and the velocity offset derived in Chiaberge et
al. (2017), we estimate that the recoil event occurred $\,\hbox to0.0pt{\lower
3.5pt\hbox{$\mathchar 0\relax\sim$}\hss}\raise 1.0pt\hbox{$>$}\,5$ Myr ago,
though the exact number depends on the assumed kick angle. Theoretically, it
is possible for a recoiling BH to remain active for such amount of time. By
using the formula in Loeb (2007) with assumed radiative efficiency of
$\epsilon=0.1$ as in Chiaberge et al. (2017), the lifetime of an accretion
disk for the black hole mass of is estimated to be $\sim 100$ Myr. In fact,
Murgia et al. (1999) derived the age of the radio source in $\sim 0.1$ Myr,
from its synchrotron spectrum, which implies that the radio AGN turned on
after the recoiling event.

In addition, Blecha et al. (2011) presented several examples of simulated
recoiling black holes and estimated a typical lifetime of $\sim 10$ Myr (see
also Blecha et al., 2016). We note, however, that uncertainties in absolute
timescale in those simulations are not negligible. Such uncertainties
primarily originate in the resolution element, which is typically much larger
than the scale of black hole merger (several kpc), making the comparison with
observations challenging. This is still the case with more recent simulations
that record the exact time of the black hole merger (e.g., TNG; Hani et al.,
2020). Sub-grid models of black hole binary evolution (e.g., Kelley et al.,
2017) show a wide range of black hole merger timescales, from $\sim 1$ Gyr to
more than a Hubble time, depending on the model assumptions and the specific
conditions in the galaxy nuclei.

###  4.3 Alternative Scenario

As we mentioned in Sec. 1, there is an alternative scenario for the observed
spatial offset, where the quasar is hosted by another, fainter galaxy (“galaxy
B” for convenience) that is superposed on a galaxy at different redshift
(“galaxy A” i.e. the one we have analysed in this study as the host). However,
this scenario is unlikely because; (1) we did not confirm any stellar light
concentration near the position of the quasar in any filters and (2) there is
evidence of narrow absorption at the same redshift of the narrow emission
lines.

For (1), assuming the scaling relation between $M_{*}$ and $M_{\rm BH}$ (e.g.,
Kormendy & Ho, 2013), stellar mass $M_{*}\sim 10^{12}\,M_{\odot}$ is expected
for the black hole mass of $M_{\rm BH}\sim 3$-$6\times 10^{9}M_{\odot}$. Such
an amount of mass should be noticeably seen in our F160W image regardless of
its light profile. Chiaberge et al. (2017) tested this with various structural
parameters for unseen galaxy B by using galfit (Peng et al., 2002a) but
concluded that it is not likely not to see such an amount of mass even in
their F140W image. The recoiling black hole scenario, which assumes galaxy A
($\log M_{*}/M_{\odot}\sim 11.4$; Sec. 3.2) as the host, offers more
reasonable agreement in terms of the scaling relation for the stellar and
black hole masses observed. To further investigate this scenario, we repeat
the structural analyses by using our deep F160W image. We provide the original
image with the quasar component unsubtracted as the input image to galfit. We
fit the light profile of the system with one PSF component at the position of
the quasar, one Sérsic component fixed at the LWC position derived above, and
one additional Sérsic component (i.e. galaxy B) at the quasar position. We fix
the position angle and axis ratio of the host to 36 degree and 0.37 (Chiaberge
et al., 2017), respectively, and leave all other parameters as free. We find
that the second Sérsic component does not improve the fit. The best-fit
effective radius of this additional component is unphysically small, $54$ mas.
In addition, the total flux of this component remains to be $\sim 7\%$ of the
host. This converts to stellar mass of $\sim 1.8\times 10^{10}\,M_{\odot}$ by
adopting the same mass-to-light ratio of the host, which is significantly
under massive for the black hole mass of We therefore conclude that this
additional component, if not an artifact such as flux residuals of the PSF
component, is unlikely to be the unseen host galaxy of the quasar.

For (2), Chiaberge et al. (2017) found narrow absorption features both in
${\rm Ly\alpha}$ and [C iv] lines at the same redshift as the narrow emission
lines ([O ii] and [Ne iii]). This indicates the presence of gas associated
within galaxy A and thus secures its redshift (see also Castignani et al.,
2022, who identified CO emission at the same redshift of the narrow lines
within the host). Determining the systemic redshift of galaxy A directly from
stellar absorption lines will be a critical step to test the scenario.

###  4.4 Future Prospects

Lastly, we note that the observed blueshift measured in the permitted broad
lines should carefully be interpreted, as these have large uncertainties
primarily due to broad absorption components and blending with other lines
(Chiaberge et al., 2017). For example, such broad absorption is often
associated with extreme outflow and seen in a number of local AGN (Kaastra et
al., 2014; Ebrero et al., 2016), making it challenging to attribute the
observed blueshift solely to a recoiling event. An ideal line is ${\rm
H\beta}$, which is a good tracer of the kinematics of the accretion disk/inner
broad line region and is isolated from other emission line; thus it can
provide a clean measurement of the velocity shift of the accretion disk
(Chiaberge et al., 2018). However, this line is located at $\sim 1\,\mu$m for
the redshift (plus possible large blueshift) of and the access to the entire
line profile with the current spectrographs may be challenging due to
atmospheric absorption and the gap between CCD and IR-detector. The IFU mode
of JWST’s NIRSpec would be an ideal choice, allowing us to confirm the
velocity shift and to locate spatial distribution simultaneously.

High-resolution sub-mm/mm observations of dense gas in would also be a
critical step to confirm the physical origin of the blueshift. This is because
dense gas in the molecular torus is considered to be orbiting at large
distances and thus not gravitationally bounded to the recoiling black hole
itself. Indeed, Decarli et al. (2014) tested this on a candidate, QSO
J0927+2943 at $z=0.697$, and identified molecular gas traced by the CO(2-1)
line at the redshift of the quasar’s broad line, arguing the need for another
scenario to explain the properties of the quasar (see also Heckman et al.,
2009; Shields et al., 2009).

##  5 Summary

In this study, we investigated the properties of the host of a recoiling black
hole candidate at $z=1.0685$. By carefully analyzing newly taken deep images
by HST, we confirmed the previously reported spatial offset, $11.1\pm{0.1}$
kpc in projected distance, between the core of and the light-weighted center
of the host galaxy. We did not find evidence of a recent merger, such as a
young starburst in disturbed outskirts. We analyzed the underlining stellar
populations of the host and then compared the observed properties with
numerical simulations, finding that the observed properties of are consistent
with an old merger remnant. Based on those pieces of evidence, we concluded
that the recoiling black hole scenario is still a plausible explanation for
the puzzling nature of

Lastly, we wish to stress that is among the highest redshift candidates of a
recoiling black hole that have been identified as of today. Given the higher
merger rate at high redshift (e.g., Lotz et al., 2011; Snyder et al., 2017),
it is reasonable to expect that more candidates will be identified in
systematic high-resolution imaging surveys of radio quasars at similar
redshifts. Future observations with high-resolution imaging cameras of JWST,
even with a short amount of exposure time available in, e.g., Survey
programs333https://jwst-docs.stsci.edu/jwst-opportunities-and-policies/jwst-
call-for-proposals-for-cycle-1/jwst-cycle-1-proposal-categories/general-
observer-go-proposals/survey-go-proposals, will allow us to identify more
recoiling black hole candidates and investigate the ubiquity of this extreme
phenomenon.

We thank the anonymous referee for the constructive comments. We thank
Emanuele Berti and Davide Gerosa for fruitful discussion on the physical
interpretation of our findings, and Markos Georganopoulos, Andrea Marinucci,
Eileen Meyer, Eric S. Perlman, William B. Sparks, Grant R. Tremblay for their
contribution to the HST observing program GO-15254. MC thanks Marta Volonteri
for providing constructive comments on our study. Support for this study was
provided by NASA through a grant HST-GO-15254.005-A from the Space Telescope
Science Institute, which is operated by the Association of Universities for
Research in Astronomy, Inc., under NASA contract NAS 5-26555. This work is
based in part on observations taken by the VLA, operated by the U.S. National
Radio Astronomy Observatory, which is a facility of the National Science
Foundation, operated under cooperative agreement by Associated Universities,
Inc. LB acknowledges support from NSF award AST-1909933. SB and CO acknowledge
support from the Natural Sciences and Engineering Research Council (NSERC) of
Canada. GC acknowledges the support from the grant ASI n.2018-23-HH.0.

Table 2:  Positions of tessellated segments and the best-fit SED parameters.

Voronoi ID | SNR | $r^{\dagger}$ | $\log M_{*}$ | $\log T_{0}$ | $\log\tau$ | $\log t_{Q}$ | $A_{V}$
---|---|---|---|---|---|---|---
|  | arcsec | $\log M_{\odot}$ | yr | yr | yr | mag
0 | 19.0 | 0.59 | 9.2 | 8.1 | 7.8 | 7.8 | 1.7
1 | 22.8 | 0.80 | 9.6 | 8.1 | 7.0 | 8.1 | 1.5
2 | 26.7 | 1.91 | 10.0 | 8.0 | 7.2 | 7.9 | 1.2
3 | 23.1 | 0.86 | 9.6 | 8.5 | 7.8 | 8.4 | 1.2
4 | 23.9 | 4.40 | 10.7 | 9.4 | 9.8 | —∗ | 1.3
5 | 21.0 | 0.92 | 9.4 | 8.8 | 8.2 | 8.7 | 0.0
6 | 23.1 | 2.51 | 9.9 | 8.3 | 7.0 | 8.3 | 1.1
7 | 23.4 | 0.35 | 9.2 | 8.0 | 7.6 | 7.8 | 1.3
8 | 22.3 | 7.06 | 10.6 | 8.3 | 7.6 | 8.2 | 2.1
9 | 18.7 | 0.20 | 9.0 | 7.9 | 7.0 | 7.8 | 1.0
10 | 19.3 | 0.12 | 9.2 | 8.3 | 7.6 | 8.2 | 0.9
11 | 22.5 | 1.84 | 10.0 | 8.1 | 8.0 | 7.4 | 1.7
12 | 23.4 | 1.36 | 9.7 | 8.8 | 9.4 | —∗ | 1.6
13 | 15.7 | 0.18 | 8.9 | 8.0 | 7.0 | 8.0 | 0.9
14 | 20.3 | 0.55 | 9.3 | 8.6 | 8.2 | 8.4 | 0.0
15 | 18.2 | 0.39 | 9.0 | 7.9 | 7.2 | 7.8 | 1.0
16 | 16.4 | 0.36 | 9.1 | 8.5 | 8.0 | 8.3 | 0.4
17 | 22.1 | 2.60 | 10.3 | 8.1 | 7.6 | 7.9 | 2.1
18 | 17.9 | 0.49 | 9.2 | 8.0 | 7.0 | 8.0 | 1.2
19 | 21.5 | 0.81 | 9.5 | 7.9 | 7.0 | 7.8 | 1.6
20 | 26.0 | 1.52 | 9.7 | 8.0 | 7.8 | 7.6 | 1.8
21 | 20.8 | 1.21 | 9.6 | 7.9 | 7.4 | 7.7 | 1.8
22 | 26.1 | 3.87 | 10.4 | 8.1 | 7.6 | 7.9 | 1.6
23 | 22.3 | 2.51 | 10.2 | 8.3 | 7.8 | 8.1 | 1.9
24 | 22.4 | 0.86 | 9.5 | 8.5 | 7.0 | 8.5 | 0.0
25 | 23.8 | 1.52 | 9.8 | 8.4 | 9.6 | —∗ | 1.7
26 | 18.1 | 2.86 | 10.1 | 8.7 | 10.0 | —∗ | 1.8
27 | 15.4 | 4.21 | 10.3 | 8.2 | 7.6 | 8.1 | 2.2
28 | 24.9 | 2.12 | 10.1 | 8.3 | 8.0 | 8.0 | 0.6

Note. — $\dagger$: Distance from the light-weighted center.

$\ast$: Regions with on-going star formation (i.e. $t_{Q}=0$).

Software: Astropy (Astropy Collaboration et al., 2013, 2018), Astrodrizzle
(Hack et al., 2012), FAST (Kriek et al., 2009), GALFIT (Peng et al., 2002b,
2010), gsf (Morishita et al., 2019), LACOSMIC (van Dokkum et al., 2001), numpy
(Oliphant, 2006; Van Der Walt et al., 2011), Pypher (Boucaud et al., 2016),
SExtractor (Bertin & Arnouts, 1996).

## References

  * Anderson & Bedin (2010) Anderson, J., & Bedin, L. R. 2010, PASP, 122, 1035
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33
  * Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., SipHocz, B. M., et al. 2018, aj, 156, 123
  * Barrows et al. (2018) Barrows, R. S., Comerford, J. M., & Greene, J. E. 2018, ApJ, 869, 154
  * Barrows et al. (2016) Barrows, R. S., Comerford, J. M., Greene, J. E., & Pooley, D. 2016, ApJ, 829, 37
  * Batcheldor et al. (2010) Batcheldor, D., Robinson, A., Axon, D. J., Perlman, E. S., & Merritt, D. 2010, ApJ, 717, L6
  * Begelman et al. (1980) Begelman, M. C., Blandford, R. D., & Rees, M. J. 1980, Nature, 287, 307
  * Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393
  * Blecha et al. (2013) Blecha, L., Civano, F., Elvis, M., & Loeb, A. 2013, MNRAS, 428, 1341
  * Blecha et al. (2011) Blecha, L., Cox, T. J., Loeb, A., & Hernquist, L. 2011, MNRAS, 412, 2154
  * Blecha et al. (2016) Blecha, L., Sijacki, D., Kelley, L. Z., et al. 2016, MNRAS, 456, 961
  * Boucaud et al. (2016) Boucaud, A., Bocchio, M., Abergel, A., et al. 2016, A&A, 596, A63
  * Bruzual & Charlot (2003) Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000
  * Burke-Spolaor et al. (2017) Burke-Spolaor, S., Gültekin, K., Postman, M., et al. 2017, ApJ, 849, 59
  * Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682
  * Campanelli et al. (2007) Campanelli, M., Lousto, C. O., Zlochower, Y., & Merritt, D. 2007, Phys. Rev. Lett., 98, 231102
  * Cappellari & Copin (2003) Cappellari, M., & Copin, Y. 2003, MNRAS, 342, 345
  * Castignani et al. (2022) Castignani, G., Meyer, E., Chiaberge, M., et al. 2022, arXiv e-prints, arXiv:2204.05882
  * Centrella et al. (2010) Centrella, J., Baker, J. G., Kelly, B. J., & van Meter, J. R. 2010, Reviews of Modern Physics, 82, 3069
  * Ceverino et al. (2014) Ceverino, D., Klypin, A., Klimek, E. S., et al. 2014, MNRAS, 442, 1545
  * Chiaberge et al. (2018) Chiaberge, M., Tremblay, G. R., Capetti, A., & Norman, C. 2018, ApJ, 861, 56
  * Chiaberge et al. (2017) Chiaberge, M., Ely, J. C., Meyer, E. T., et al. 2017, A&A, 600, A57
  * Civano et al. (2010) Civano, F., Elvis, M., Lanzuisi, G., et al. 2010, ApJ, 717, 209
  * Civano et al. (2012) —. 2012, ApJ, 752, 49
  * Comerford & Greene (2014) Comerford, J. M., & Greene, J. E. 2014, ApJ, 789, 112
  * Decarli et al. (2014) Decarli, R., Dotti, M., Mazzucchelli, C., Montuori, C., & Volonteri, M. 2014, MNRAS, 445, 1558
  * Dotti et al. (2007) Dotti, M., Colpi, M., Haardt, F., & Mayer, L. 2007, MNRAS, 379, 956
  * Ebrero et al. (2016) Ebrero, J., Kriss, G. A., Kaastra, J. S., & Ely, J. C. 2016, A&A, 586, A72
  * Escala et al. (2005) Escala, A., Larson, R. B., Coppi, P. S., & Mardones, D. 2005, ApJ, 630, 152
  * Fanti et al. (1985) Fanti, C., Fanti, R., Parma, P., Schilizzi, R. T., & van Breugel, W. J. M. 1985, A&A, 143, 292
  * Ferrarese & Merritt (2000) Ferrarese, L., & Merritt, D. 2000, ApJ, 539, L9
  * Frank & Rees (1976) Frank, J., & Rees, M. J. 1976, MNRAS, 176, 633
  * Fukugita et al. (1996) Fukugita, M., Ichikawa, T., Gunn, J. E., et al. 1996, AJ, 111, 1748
  * Gebhardt et al. (2000) Gebhardt, K., Richstone, D., Kormendy, J., et al. 2000, AJ, 119, 1157
  * Glikman et al. (2007) Glikman, E., Djorgovski, S. G., Stern, D., Bogosavljević, M., & Mahabal, A. 2007, ApJ, 663, L73
  * Goulding et al. (2019) Goulding, A. D., Pardo, K., Greene, J. E., et al. 2019, ApJ, 879, L21
  * Hack et al. (2012) Hack, W. J., Dencheva, N., Fruchter, A. S., et al. 2012, in American Astronomical Society Meeting Abstracts, Vol. 220, American Astronomical Society Meeting Abstracts #220, 135.15
  * Hani et al. (2020) Hani, M. H., Gosain, H., Ellison, S. L., Patton, D. R., & Torrey, P. 2020, MNRAS, 493, 3716
  * Heckman et al. (2009) Heckman, T. M., Krolik, J. H., Moran, S. M., Schnittman, J., & Gezari, S. 2009, ApJ, 695, 363
  * Hennawi et al. (2006) Hennawi, J. F., Strauss, M. A., Oguri, M., et al. 2006, AJ, 131, 1
  * Hewett & Wild (2010) Hewett, P. C., & Wild, V. 2010, MNRAS, 405, 2302
  * Hilbert et al. (2016) Hilbert, B., Chiaberge, M., Kotyla, J. P., et al. 2016, ApJS, 225, 12
  * Hogg et al. (2021) Hogg, J. D., Blecha, L., Reynolds, C. S., Smith, K. L., & Winter, L. M. 2021, MNRAS, 503, 1688
  * Huang et al. (2018) Huang, S., Leauthaud, A., Greene, J. E., et al. 2018, MNRAS, 475, 3348
  * Kaastra et al. (2014) Kaastra, J. S., Kriss, G. A., Cappi, M., et al. 2014, Science, 345, 64
  * Kelley et al. (2017) Kelley, L. Z., Blecha, L., & Hernquist, L. 2017, MNRAS, 464, 3131
  * Kim et al. (2017) Kim, D. C., Yoon, I., Privon, G. C., et al. 2017, ApJ, 840, 71
  * Komossa (2012) Komossa, S. 2012, Advances in Astronomy, 2012, 364973
  * Komossa et al. (2003) Komossa, S., Burwitz, V., Hasinger, G., et al. 2003, ApJ, 582, L15
  * Komossa et al. (2008) Komossa, S., Zhou, H., & Lu, H. 2008, ApJ, 678, L81
  * Kormendy & Ho (2013) Kormendy, J., & Ho, L. C. 2013, ARA&A, 51, 511
  * Koss et al. (2014) Koss, M., Blecha, L., Mushotzky, R., et al. 2014, MNRAS, 445, 515
  * Koss et al. (2018) Koss, M. J., Blecha, L., Bernhard, P., et al. 2018, Nature, 563, 214
  * Kriek et al. (2009) Kriek, M., van Dokkum, P. G., Labbé, I., et al. 2009, ApJ, 700, 221
  * Kriek et al. (2019) Kriek, M., Price, S. H., Conroy, C., et al. 2019, ApJ, 880, L31
  * Krist et al. (2011) Krist, J. E., Hook, R. N., & Stoehr, F. 2011, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 8127, Optical Modeling and Performance Predictions V, ed. M. A. Kahan, 81270J
  * Kuraszkiewicz et al. (2002) Kuraszkiewicz, J. K., Green, P. J., Forster, K., et al. 2002, ApJS, 143, 257
  * Lena et al. (2014) Lena, D., Robinson, A., Marconi, A., et al. 2014, ApJ, 795, 146
  * Loeb (2007) Loeb, A. 2007, Phys. Rev. Lett., 99, 041103
  * Lotz et al. (2011) Lotz, J. M., Jonsson, P., Cox, T. J., et al. 2011, ApJ, 742, 103
  * Lousto & Zlochower (2011) Lousto, C. O., & Zlochower, Y. 2011, Phys. Rev. Lett., 107, 231102
  * Lousto et al. (2017) Lousto, C. O., Zlochower, Y., & Campanelli, M. 2017, ApJ, 841, L28
  * Lousto et al. (2012) Lousto, C. O., Zlochower, Y., Dotti, M., & Volonteri, M. 2012, Phys. Rev. D, 85, 084015
  * Milosavljević & Merritt (2003) Milosavljević, M., & Merritt, D. 2003, ApJ, 596, 860
  * Moore et al. (2015) Moore, C. J., Cole, R. H., & Berry, C. P. L. 2015, Classical and Quantum Gravity, 32, 015014
  * Morishita (2021) Morishita, T. 2021, ApJS, 253, 4
  * Morishita et al. (2017) Morishita, T., Abramson, L. E., Treu, T., et al. 2017, ApJ, 835, 254
  * Morishita et al. (2019) —. 2019, ApJ, 877, 141
  * Mulia et al. (2015) Mulia, A. J., Chandar, R., & Whitmore, B. C. 2015, ApJ, 805, 99
  * Murgia et al. (1999) Murgia, M., Fanti, C., Fanti, R., et al. 1999, A&A, 345, 769
  * Novak et al. (2015) Novak, M., Smolčić, V., Civano, F., et al. 2015, MNRAS, 447, 1282
  * O’Dea (1998) O’Dea, C. P. 1998, PASP, 110, 493
  * O’Dea & Saikia (2021) O’Dea, C. P., & Saikia, D. J. 2021, A&A Rev., 29, 3
  * Oke & Gunn (1983) Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713
  * Oliphant (2006) Oliphant, T. E. 2006, A guide to NumPy, Vol. 1 (Trelgol Publishing USA)
  * O’Neill et al. (2021) O’Neill, S., Kiehlmann, S., Readhead, A. C. S., et al. 2021, arXiv e-prints, arXiv:2111.02436
  * Onodera et al. (2015) Onodera, M., Carollo, C. M., Renzini, A., et al. 2015, ApJ, 808, 161
  * Peng (2007) Peng, C. Y. 2007, ApJ, 671, 1098
  * Peng et al. (2002a) Peng, C. Y., Ho, L. C., Impey, C. D., & Rix, H.-W. 2002a, AJ, 124, 266
  * Peng et al. (2002b) —. 2002b, AJ, 124, 266
  * Peng et al. (2010) Peng, Y.-j., Lilly, S. J., Kovač, K., et al. 2010, ApJ, 721, 193
  * Pesce et al. (2018) Pesce, D. W., Braatz, J. A., Condon, J. J., & Greene, J. E. 2018, ApJ, 863, 149
  * Pesce et al. (2021) Pesce, D. W., Seth, A. C., Greene, J. E., et al. 2021, ApJ, 909, 141
  * Podigachoski et al. (2015) Podigachoski, P., Barthel, P. D., Haas, M., et al. 2015, A&A, 575, A80
  * Postman et al. (2012) Postman, M., Lauer, T. R., Donahue, M., et al. 2012, ApJ, 756, 159
  * Salpeter (1955) Salpeter, E. E. 1955, ApJ, 121, 161
  * Shields et al. (2009) Shields, G. A., Rosario, D. J., Smith, K. L., et al. 2009, ApJ, 707, 936
  * Siemiginowska et al. (2005) Siemiginowska, A., Cheung, C. C., LaMassa, S., et al. 2005, ApJ, 632, 110
  * Silverman et al. (2020) Silverman, J. D., Tang, S., Lee, K.-G., et al. 2020, ApJ, 899, 154
  * Simons et al. (2019) Simons, R. C., Kassin, S. A., Snyder, G. F., et al. 2019, ApJ, 874, 59
  * Skipper & Browne (2018) Skipper, C. J., & Browne, I. W. A. 2018, MNRAS, 475, 5179
  * Snyder et al. (2017) Snyder, G. F., Lotz, J. M., Rodriguez-Gomez, V., et al. 2017, MNRAS, 468, 207
  * Springel (2005) Springel, V. 2005, MNRAS, 364, 1105
  * Steinhardt et al. (2012) Steinhardt, C. L., Schramm, M., Silverman, J. D., et al. 2012, ApJ, 759, 24
  * Thomas et al. (2010) Thomas, D., Maraston, C., Schawinski, K., Sarzi, M., & Silk, J. 2010, MNRAS, 404, 1775
  * Tichy & Marronetti (2007) Tichy, W., & Marronetti, P. 2007, Phys. Rev. D, 76, 061502
  * Van Der Walt et al. (2011) Van Der Walt, S., Colbert, S. C., & Varoquaux, G. 2011, Computing in Science & Engineering, 13, 22
  * van Dokkum et al. (2001) van Dokkum, P. G., Franx, M., Kelson, D. D., & Illingworth, G. D. 2001, ApJ, 553, L39
  * Volonteri et al. (2003) Volonteri, M., Haardt, F., & Madau, P. 2003, ApJ, 582, 559
  * Wuyts et al. (2013) Wuyts, S., Förster Schreiber, N. M., Nelson, E. J., et al. 2013, ApJ, 779, 135
  * Zolotov et al. (2015) Zolotov, A., Dekel, A., Mandelker, N., et al. 2015, MNRAS, 450, 2327