###### Abstract

We study the merger history of primordial black holes (PBHs) in a scenario
where they represent the dominant dark matter component of a typical dwarf
galaxies’ core. We investigate the possibility of a sequence of collisions
resulting in a hierarchical merger of black holes, and look at the final mass
spectrum in such clusters, which initially present a monochromatic (single-
mass) PBH population. Our study shows that the merging process results in the
transfer of about $40\%$ of the total mass of the core to the merger products
regardless of the initial mass of PBHs, with about $5\%$ of energy radiated
out in the form of gravitational waves. We find that, in the lighter mass
limit, black holes up to eight times more massive than the original population
can be formed within a Hubble time.

Hierarchical Merger of Primordial Black Holes in Dwarf Galaxies

Encieh Erfani${}^{1,\,*}$, Tadeo D. Gomez-Aguilar${}^{2,\,3,\,{\dagger}}$ and
Juan Carlos Hidalgo${}^{3,\,{\ddagger}}$

1Department of Physics, Institute for Advanced Studies in Basic Sciences
(IASBS), Zanjan 45137-66731, Iran
2Facultad de Ciencias en Física y Matemáticas, Universidad Autónoma de
Chiapas, 29050, Tuxtla Gutiérrez, Chiapas, México
3Instituto de Ciencias Físicas, Universidad Nacional Autónoma de México,
62210, Cuernavaca, Morelos, México

000∗erfani@iasbs.ac.ir000†tadeo.dga<EMAIL_ADDRESS>
##  1 Introduction

###

Observation of gravitational waves (GWs) within the last six years is
transforming the field of astronomy and our understanding of compact objects
[1, 2, 3, 4, 5, 6]. The LIGO and Virgo experiments have so far detected around
90 binary compact object mergers composed of black holes (BHs) and neutron
stars (NSs) [6], with a few events challenging the models of stellar
evolution. LIGO/Virgo’s third observing run [4] announced the most massive
binary BH (BBH) merger (GW190521) [7] with the primary component mass of $\sim
85\,M_{\odot}$. This observation challenges the formation and evolution
mechanisms of stellar BHs since models predict no objects larger than about
$65\,M_{\odot}$ due to the pulsational pair-instability (PI) process [8]. This
sets a mass gap for BHs between $70-150\,M_{\odot}$ [9]. Primordial black
holes (PBHs) formed from the collapse of large density fluctuations in the
early Universe [10, 11, 12, 13], could provide an alternative origin for such
events. Another challenge is to model detected pairs which present a large
asymmetry in the progenitors masses (e.g. GW190412), with mass ratios of
nearly 4–to–1 ($30\,M_{\odot}+8\,M_{\odot}$) [14]. A scenario where more than
one species of stars or BHs is confined to a small volume could be explained
by a sequence of mergers.

Several mechanisms for the origin of the progenitors in BBHs have been
considered in the literature; for example, 1) binary star evolution through
common envelope [15, 16, 17], 2) dynamical process in triples [18], 3) pairing
of PBHs [19, 20], or 4) as a product of hierarchical mergers [21, 22, 23, 24,
25, 26, 27]. In this article we propose a combination of two mechanisms (3 &
4); i.e. PBHs of a common mass pair to form a binary and the “product” of
their merger will pair with another (P)BH, and so on. We explore the
possibility of repetition of this process in a “hierarchical merger” sequence
throughout the age of the Universe, considering the radiation of GWs, as they
could potentially provide a detectable population of binaries with one or both
component masses in the mass gap of stellar BHs.

The astrophysical environments suitable for hierarchical mergers, where
multiple GW events will be produced, require a large escape velocity so that
merger remnants are efficiently retained. A compact (dense) region will
prevent merger remnants from being ejected by GW recoils. Scenarios previously
considered for the efficient production of hierarchical mergers include
globular clusters [28, 29, 30, 31], nuclear star clusters [32, 33, 34, 35],
and accretion disks surrounding active galactic nuclei [36, 37, 38]. Here we
propose a hierarchical merger of PBHs in Dwarf Galaxies (DGs) for two reasons:
First, DGs are dominated by DM [39], which we assume to be in the form of
PBHs111PBHs with a mass larger than $10^{15}$ g are candidates for Dark Matter
(DM) [40]., and second, these environments are dense enough for the
hierarchical merger of PBHs [39].

We explore the possible formation of systems like GW190412 and GW190521 as a
result of our model (hierarchical merger is also considered in [41]). We also
investigate the formation of BHs in the mass gap and also merged pairs with a
large mass ratio [42].

The paper is organized as follows: in section 2 we discuss DGs as hosts of
PBHs. We also summarize the basic equations of PBH as a candidate for DM [43].
In section 3, we explain our model and in section 4 we present our results.
Finally, section 5 is devoted to our conclusions.

##  2 Dwarf Galaxies and Primordial Black Holes

###

Dwarf Galaxies (DGs) are the oldest, and least chemically evolved stellar
systems known [39]. These galaxies resemble globular clusters in density
profiles, however, they are distinguished from star clusters since they
present a dynamical mass that is substantially larger than the mass inferred
from the luminous stellar population [44]. The stellar kinematics of DGs
confirm that they contain a dominant component of DM which is in contrast with
the star clusters. Altogether, these properties indicate that DGs are the
astrophysical systems where DM is most abundant.

Since these galaxies are highly DM dominated, with mass-to-light ratios $\sim
1000\,M_{\odot}/L_{\odot}$ [45], they are extremely valuable laboratories for
characterizing the nature of DM, mostly as probes of DM on the smaller end of
the large scale structure spectrum [46].

Dozens of DGs have recently been discovered by several surveys like the SDSS
[47, 48, 49] as satellites of the Milky Way, Andromeda, and as members of the
Local Group [50]222Recall that DGs have been a source of a tension known as
missing satellites problem since predictions for the abundance of massive
satellite galaxies in simulations present significantly larger numbers than
the observed objects [51, 52].. The ultra-faint Milky Way satellites have
masses ranging from just over $10^{6}M_{\odot}$ (Coma Berenices) up to
$2.8\times 10^{7}M_{\odot}$ (Canes Venatici I) [45]. Discovered DGs have
Plummer (half-light) radii as small as $\sim 20$ pc [53], with the total mass
inside this radius in the range of one to three orders of magnitude larger
than their stellar masses [54]. Consequently, the measured central densities,
and the density profiles are essentially those of DM halo. The central
densities of DGs range from $\sim 0.08\,M_{\odot}$/pc3 up to $\sim
2.1\,M_{\odot}$/pc3 [45].

The identity of the dominant DM in these environments, however, remains a
mystery [55]. One appealing possibility is that the DM consists of black holes
formed in the early Universe, known as Primordial Black Holes (PBHs) [19, 56].
Constraints on their abundance imposed in DG environments are reported in
Refs. [57, 58]333Recently, some studies examined whether DM candidates in the
form of PBHs can solve the cusp-core problem in low-mass DGs [59].. In the
following, we will briefly explain the properties of PBHs and their
abundances.

The most common mechanism for PBHs formation is the collapse of density
fluctuations larger than a threshold that reenters the horizon after inflation
in the early Universe [10, 11, 12]444See however [60] for an alternative, more
likely mechanism during a long-lasting period of reheating.. Thus the mass of
PBH is roughly the mass of the horizon, $M_{\rm H}$

$M_{\rm PBH}\approx M_{\rm H}\sim 10^{15}\left(\dfrac{t}{10^{-23}\,{\rm
s}}\right){\rm g}\,.$ (2.1)

PBHs radiate thermally due to the Hawking radiation [40], but massive enough
PBHs survive for more than a Hubble time. In consequence, their mass spectrum
spans many orders of magnitude, from $\sim 10^{15}$ g to well over $10^{50}$
g.

A parameter which represents the abundance (the energy density fraction) of
PBHs in the epoch of their formation is

$\beta\equiv\frac{\rho_{\rm PBH}(t)}{\rho(t)}\,.$ (2.2)

Assuming an adiabatic expansion of the Universe, the above can be related to
the present abundance mass fraction of PBHs, $f_{\rm PBH}\equiv\Omega_{\rm
PBH}/\Omega_{\rm DM}$, by [43]

$\beta\simeq 3.7\times
10^{-9}\,\left(\frac{g_{*}}{10.75}\right)^{1/4}\left(\frac{M_{\rm
PBH}}{M_{\odot}}\right)^{1/2}\,f_{\rm PBH}\,,$ (2.3)

where $g_{*}$ is a number of relativistic degree of freedom at the time of
PBHs’ formation.

A population of PBHs with a specific mass is subject to observational
constraints which in many cases reduces the contribution of $f_{\rm PBH}$ to
Dark Matter. However, DM can be completely formed of PBHs in the mass window
around $(5\times 10^{-16}-2\times 10^{-14})\,M_{\odot}$ [43, 61]. The recent
detection of GWs by LIGO/Virgo also imposes constraints on the abundance of
PBHs with masses of order $10$ to $50\,M_{\odot}$. Searches for compact
objects of sub-Solar masses also impose constraints on PBH abundance [62].
There are also constraints on sub-Solar masses by gravitational microlensing
observations [63]. These constraints are imposed without any assumption about
the distribution of PBHs at their formation, their clustering, or the
environment, which can alter the bounds significantly [64, 65]. In the next
section, in light of the observational constraints, we will consider DM PBHs
with four different masses.

##  3 Model for the Merger Rate

###

In what follows we describe the details of the model employed to produce a
succession of mergers of BHs in a DG. We consider a population of PBHs at the
core of the DG, with $N_{\rm PBH}$ elements, which we call a cluster. We
assume PBHs are the only component of DM, with a single initial mass
(monochromatic population). We test four different initial masses of PBHs,
where $N_{\rm PBH}$ is determined by the mass of DG. Throughout this work we
consider a DG with the total mass, $M_{\rm DG}=10^{9}\,M_{\odot}$, and radius,
$R_{\rm DG}\sim 10$ pc.

Figure 1: The hierarchical tree of (P)BH mergers from an initial population of
PBHs with mass $10\,M_{\odot}$. In the age of the Universe, collisions reach
up to a fourth generation of BHs. Different colors indicate different
generations, and approximate masses, in units of $M_{\odot}$, consider the
mass-loss due to GW radiation.

In order for PBHs to merge more than once, a dense environment is required. We
thus focus on core of DGs with mass, $M_{\rm c}=10^{5}\,M_{\odot}$ confined in
a core radius, $R_{\rm c}=0.9$ pc, which is admitted by observations [45]. The
number of PBHs is proportional to the density of DM, $\rho_{\rm DM}$,

$N_{\rm PBH}=\dfrac{4\pi}{3}\dfrac{\rho_{\rm DM}\,R_{\rm c}^{3}}{m_{\rm
PBH}}\,.$ (3.1)

In our analysis we consider four different initial masses for PBHs, which we
define as first generation (1G);
$10^{-14}\,M_{\odot},\,10^{-2}\,M_{\odot},\,1\,M_{\odot}$ and 10 $M_{\odot}$.
We choose $10^{-14}\,M_{\odot}$ since as mentioned there is no confirmed
observational constraints for this mass and PBHs can be whole DM. Our interest
in other mass values comes from the GWs observations. From Eq. (3.1), we get
the number of PBHs at the core (see Table 1 below).

According to our assumption, the components of the first merger will have
equal masses. Therefore, the merger of 1G+1G will produce the second
generation BHs (2G) with roughly double the mass. Subsequently, from the third
generation (3G) onward, BHs can form from 1G+2G and 2G+2G; two populations of
BHs will form, and therefore, a variety of masses is possible (see Figure 1).

For estimating the merger rate, we know that the cross-section of two objects
with masses $m_{i}$ and $m_{j}$ which merge after first being captured into
binaries by the emission of gravitational radiation is [66, 67]

$\sigma=\dfrac{\sigma(m_{i},m_{j})}{|v_{i}-v_{j}|^{18/7}}\,,$ (3.2)

with

$\sigma(m_{i},m_{j})=2\pi\left(\frac{85\pi}{6\sqrt{2}}\right)^{2/7}\frac{G^{2}(m_{i}+m_{j})^{10/7}\left(m_{i}m_{j}\right)^{2/7}}{c^{10/7}}\,.$
(3.3)

Here $v_{i}$ represents the velocity of the components in the binary, which is
given by the root mean-squared (rms) velocity,
$\overline{v_{i}}\equiv\sqrt{\overline{v_{i}^{2}}}$

$\overline{v_{i}^{2}}(r)=\frac{4\pi}{n_{i}(r)}\int_{0}^{\phi(r)}\,f_{i}(E)\left[2\left(\phi(r)-E\right)\right]^{3/2}\,dE\,,$
(3.4)

where $n_{i}(r)=\dfrac{\rho_{i}(r)}{m_{i}}$ is the number density, and the DM
density is modelled by a Plummer sphere [53] at all times. This prescribes the
following expressions for the density and the gravitational potential,
respectively [53]

$\displaystyle\rho_{i}(r)$ $\displaystyle=$ $\displaystyle{\rho}_{i}^{\rm
c}\left(1+\frac{r^{2}}{R_{\rm
c}^{2}}\right)^{-5/2}=\frac{3\,m_{i}\,N_{i}}{4\pi\,R_{\rm
c}^{3}}\left(1+\frac{r^{2}}{R_{\rm c}^{2}}\right)^{-5/2}\,,$ (3.5)
$\displaystyle\phi(r)$ $\displaystyle=$ $\displaystyle{\phi}^{\rm
c}\left(1+\dfrac{r^{2}}{R_{\rm c}^{2}}\right)^{-1/2}=\dfrac{GM_{\rm c}}{R_{\rm
c}}\left(1+\dfrac{r^{2}}{R_{\rm c}^{2}}\right)^{-1/2}\,,$ (3.6)

where we defined the core density of the $i$-th species, ${\rho}_{i}^{\rm
c}\equiv\frac{3\,m_{i}\,N_{i}}{4\pi\,R_{\rm c}^{3}}$, with $N_{i}$
representing the total number of the BHs with mass $m_{i}$, as well as a core
potential, ${\phi}^{\rm c}\equiv\frac{GM_{\rm c}}{R_{\rm c}}$. The
distribution function of each BH population, $f_{i}(E)$ is given by [66]

$f_{i}(E)=\frac{32\sqrt{2}}{7\pi^{2}}\,(\phi^{c})^{-5}\,n_{i}^{c}\,E_{i}^{7/2}\,,$
(3.7)

where $n_{i}^{c}\equiv n_{i}(r=R_{c})$, and according to the virial theorem,
the energy of each species of the cluster is given by

$E_{i}=\dfrac{G\,m_{i}^{2}\,N_{i}^{2}}{2\,R_{\rm c}}\,.$ (3.8)

An approximate analytic expression for the merger rate in this case is
estimated from the Fokker-Planck equation, in particular, from the terms that
account for the loss and gain of BHs due to mergers with other BHs [66] (see
also [68]). Once the velocities take the rms value, this can be expressed as

$\Gamma_{j}=\frac{14\pi}{3}\sum_{i}\sigma(m_{i},\,m_{j})\int\;dr\,r^{2}\,\frac{n_{i}}{\overline{v}_{i}}\frac{n_{j}}{\overline{v}_{j}}\left[\left(\overline{v}_{i}+\overline{v}_{j}\right)^{3/7}-|\overline{v}_{i}-\overline{v}_{j}|^{3/7}\right]\,.$
(3.9)

The rms velocity is the solution of Eq. (3.4), which in light of the above
definitions, can be expressed as

$\overline{v_{i}}^{2}(r)=\frac{1}{2}\phi(r)\left(\frac{\phi(r)}{\phi^{c}}\right)^{5}\left(\frac{n_{i}(r)}{n^{c}_{i}}\right)^{-1}\,.$
(3.10)

In order to define the merging components, we take as different elements the
species, $i$ or $j$, defined by their masses, and we also divide the galaxy
core in shells which components will collide with a distance-dependent merger
rate. From the expression in Eq. (3.9), the merger rate in region between
shells at radii $r_{A}$ and $r_{B}$ is given by

$\begin{split}\Gamma&(r_{A},\,r_{B})=\frac{14\pi}{3}\sum_{i}\sigma(m_{i},\,m_{j})\,\frac{n_{j}(r_{B})n_{i}^{c}}{\overline{v_{j}}(r_{B})}\,\left(\frac{\phi^{c}}{2}\right)^{-2/7}\int\,dr_{A}\,r^{2}_{A}\,(1+A)^{-9/4}\\\
&\times\left\\{\left[\left(\frac{1}{(1+B)^{1/4}}\right)+\left(\frac{1}{(1+A)^{1/4}}\right)\right]^{3/7}-\left|\left(\frac{1}{(1+B)^{1/4}}\right)-\left(\frac{1}{(1+A)^{1/4}}\right)\right|^{3/7}\right\\}\,,\end{split}$
(3.11)

where $A=(r_{A}/R_{\rm c})^{2}$ and $B=(r_{B}/R_{\rm c})^{2}$, label the shell
to which each of the merging components belong.

For each example, we have sliced the galaxy core in 10 shells and we proceed
to count the number of mergers by dividing the cosmic time in merger epochs of
constant $\Gamma(r_{A},\,r_{B})$. The lasting time for such epochs, $\tau_{e}$
is defined as one tenth times the inverse of the largest merging rate value
among the pairs of shells. Mathematically,

$\tau_{e}=\frac{1}{10}\mathrm{min}\left\\{\Gamma(r_{A},\,r_{B})^{-1}\right\\}\,.$
(3.12)

We count the number of mergers, starting with a single population of PBHs at
redshift $z=20$. This initial condition sets an early enough time to discard
the stellar origin for the initial BHs, and late enough for the formation of
DGs. For the first iteration, the initial single-mass population is divided
into two identical sets, taken as two different species as an input of our
merging algorithm. In every epoch, we evaluate the probability of encounters
from the merger rate in Eq. (3.11) for a given combination of two species at
two given radii ($r_{A}$ and $r_{B}$), verifying that the merger time is
shorter than the period each merger epoch lasts. We thus count the number of
mergers by considering the number of BHs available in each shell for both
species. For the next epoch, we reset the merging rates at time
$t_{i}+\tau_{e}$ for the updated population of each species. Thus the above
procedure is repeated to compute the mergers at each shell. In Figures 2 and
3, we report results for the total number of BHs formed in the whole DG core
(after integrating all shells).

In each merger, we also consider the mass loss due to the emission of
gravitational radiation. The binary is formed through gravitational capture
where the energy is released through the emission of GW. Such dissipative
binary formation is followed by the coalescence of the binary with a well-
known rate of energy loss. The time-averaged energy loss rate of the binary in
the Keplerian orbit is given by [69]

$\bigg{\langle}\dfrac{dE}{dt}\bigg{\rangle}=-\dfrac{32}{5}\dfrac{G^{4}\,(m_{i}m_{j})^{2}\,(m_{i}+m_{j})}{a^{5}}F(e)\,,$
(3.13)

where $e$ and $a$ are eccentricity and semi-major axis of the orbit of binary,
respectively, and the explicit dependence on the eccentricity is given by

$F(e)=\frac{1}{(1-e^{2})^{7/2}}\left(1+\dfrac{73}{24\,}\,e^{2}+\dfrac{37}{96}\,e^{4}\right)\,.$
(3.14)

Since the orbits follow initially a parabolic path, which is also required by
the cross-section expressed in Eq. (3.2), we take the eccentricity as $e=0.99$
as an approximation to the parabola and the semi-major axis as the initial
(periastron) separation. The energy emission rate can be expressed as a
function of the semi-major axis decrease [69]

$\frac{da}{dt}=-\frac{64}{5}\frac{G^{3}}{c^{5}}\frac{m_{i}m_{j}\,(m_{i}+m_{j})}{a^{3}}F(e)\,,$
(3.15)

Therefore the energy loss in GWs throughout the merger is given by

$\Delta E=\int_{a_{i}}^{a_{\rm
merge}}\frac{dE}{da}\,\frac{da}{dt}\,dt=\frac{1}{2}m_{i}m_{j}\left(\frac{1}{a_{\rm
merge}}-\frac{1}{a_{i}}\right)\,,$ (3.16)

where the semi-major axis is integrated from an initial separation, $a_{i}$ up
to the sum of the Schwarzschild radii of the progenitors, $a_{\rm merge}$. The
energy emitted by gravitational radiation is accounted for by the effective
mass loss. This results in the merger product with a mass smaller than the sum
of the progenitors’ masses, as illustrated in Figure 1. The accumulated effect
on the total mass of the cluster is also taken into account in our evaluation
of merger rates (see Table 1). Our prescription for the mass loss is in
agreement with events observed by the LIGO/Virgo collaboration555The computed
mass via Eq. (3.16) from the progenitor components deviates in average only 2%
from the final mass estimated for the events in the three runs of LIGO/Virgo
(well within the determined errors).. Finally, it is worth mentioning that we
do not consider the disruption of BHs binaries by encounters with other
compact objects (BHs and stars). We also neglect the ejection of BHs due to
the recoil kicks of GW radiation. Note that even if we take these effects into
account, they will only affect an insignificant fraction of the total number
of merger products.

##  4 Results

###

In this section, we present the results of a succession of merger epochs
according to the model described above. We consider DG cores dominated by DM
PBHs, with four cases of a single initial PBH mass given by $m_{\rm
PBH}/M_{\odot}=\left(10^{-14},\,10^{-2},\,1,\,10\right)$. We keep the initial
mass of the core as a constant, $M_{\rm c}(t_{i})=10^{5}M_{\odot}$, which
implies that the initial number of PBHs and number density is sensitive to the
initial mass of (P)BHs (see Table 1). For illustration purposes, in Figure 1
we show the merger tree of (P)BHs, for the specific case of initial mass,
$m_{\rm PBH}=10\,M_{\odot}$ (blue circles). Each merger gives way to a BH of
the $n$-th generation if at least one of its progenitors belongs to the
$n-1$-th generation. Thereby, the second generation is constituted by the
1G+1G progenitors (yellow circles), while mergers of 1G+2G and 2G+2G BHs
constitute the third generation (red circles), and so forth. In our approach,
only four epochs fit within the age of the Universe and thus, BHs merge up to
the fourth generation, because the $n$-th epoch will show mergers from the
first up to the $n$-th generation. As stated earlier, due to the GW radiation
in each merger, the mass of each merger product is smaller than the sum of the
components of the binary.

$m_{\rm PBH}$ | $N_{\rm PBH}(t_{i})$ | $M_{\rm c}(t_{0})/M_{\rm c}(t_{i})$ | $\sum m_{\rm PBH}(t_{0})/M_{\rm c}(t_{0})$
---|---|---|---
$10^{-14}\,M_{\odot}$ | $10^{19}$ | $0.953$ | $0.5897$
$10^{-2}\,M_{\odot}$ | $10^{7}$ | $0.952$ | $0.5897$
$1\,M_{\odot}$ | $10^{5}$ | $0.952$ | $0.5896$
$10\,M_{\odot}$ | $10^{4}$ | $0.953$ | $0.5891$
Table 1: The first and second columns show the mass and number of the original
population of PBHs (1G), respectively. The third column represents the
remaining mass of the core after four epochs. The final column shows the final
(present) percentage of the original PBH population in terms of mass.

Table 1 indicates the number of PBHs in the initial population (at time
$t_{i}$) and the remaining mass of the core after the four merger epochs at
the present time, $t_{0}$, ($M_{\rm c}(t_{0})/M_{\rm c}(t_{i})$)–third column.
The last column shows the mass fraction of the initial PBH population which
did not collide with another BH up to $t_{0}$. It is worth noting that these
results are practically independent of the initial mass of PBHs.

Figure 2 shows the number count of merger products at each epoch for the
different initial populations of PBHs considered. Note that for smaller PBH
masses, merger products of mass up to $8\,m_{\rm PBH}$ are possible. On the
other hand, for large initial PBH masses, only BHs of mass $4.3\,m_{\rm PBH}$
are significantly produced. This is due to the number of PBHs initially
present, as listed in Table 1. Note also that in Figure 2, the dominant
product populations in numbers are those of the first and second generation of
mergers. This is true for all investigated cases and is due to the fact that
the relative merger rates are independent of the PBH mass but depend on the
number of PBHs.

(a)

(b)

(c)

(d)

Figure 2: Number of BHs formed from mergers after each epoch, shown as a
function of the age of the Universe, starting from $t=0.18~{}\mathrm{Gyr}$
($z=20$). Each epoch is marked with a vertical dashed line. The panels
correspond to the following masses for the initial populations: (a)
$10\,M_{\odot}$, (b) $1\,M_{\odot}$, (c) $10^{-2}\,M_{\odot}$, (d)
$10^{-14}\,M_{\odot}$.

Focusing on the galaxy core mass fraction of the merger products, we plot
their contribution in Figure 3, where we show the relative mass abundance of
each population in a DG core initially formed of PBHs with mass $m_{\rm PBH}$.
We find that the resulting mass fractions are independent of the value of the
original PBH mass (see Table 1). The left plot of this figure shows that the
fraction of PBHs decreases steadily with time, increasing the mass fraction of
merger products after each merger epoch. After the four iterations, we find
that more than 40% of the PBHs have formed at least one binary and merged to
form larger BHs. The right panel of Figure 3 shows the mass fraction in merger
products. Note that the mass fraction is dominated by the more massive species
even when the number count of such species is subdominant (see the cyan
component in Figures 2 and 3).

(a)

(b)

Figure 3: Percentage of the total mass in each population (of a given mass)
for the most abundant populations. This is plotted for the case where $m_{\rm
PBH}=10\,M_{\odot}$ but the fractions are largely independent of the PBH mass.
The original population is included in Panel (a), where part of the missing
mass is in smaller species, and only $4.7\%$ is lost in the radiation of
gravitational waves. As shown in Panel (b), the percentage of mass in the
three most numerous populations generated is of order $10\%$ each, while the
most massive population constitutes more than $15\%$ of the final mass. Note
that these percentages are the same for a range of initial masses of PBHs in
the cluster. Finally, note that the plots account only for the mass
contribution of each population, and therefore we make no distinction among
the generation in which they are formed.

##  5 Summary and Discussion

###

Primordial black holes formed in the early Universe, before the formation of
stars, can exist as dark matter and also contribute to the black hole merger
events observed through gravitational waves. GW observations have demonstrated
that BHs mergers may be more frequent than expected. If merger products form
new binaries, they may subsequently merge as detectable GW sources. For
sequential mergers of BHs, a dense, DM-dominated environment is required.
Dwarf galaxies are ideal scenarios to host a hierarchical merger of BHs.

In this paper, we studied sequence mergers of BH starting with a monochromatic
cluster of PBHs at the core of DGs. Our study, featuring stages of constant
merger rate, represents a first step towards the full numerical analysis of
the Fokker-Planck equation. Since PBHs span several decades of mass, we
considered four different initial masses. Starting the evolution of the system
at redshift $z=20$, we find that in cores with massive PBHs ($\sim
10\,M_{\odot}$), BHs with up to four times the initial mass can be formed
significantly, while in the small PBH mass limit the masses of products are up
to eight times the initial PBH mass. Our results also show that the total mass
loss at the DG core from GW emission, and the mass fraction of BH that undergo
collisions are mostly independent of the initial PBH mass. Figure 3 in
particular shows that the original population of PBHs is reduced to little
more than 50% of the total mass. The proportion of larger mass BH populations
may be tested by future GW surveys.

Our results are not only relevant to the rate and spectrum of GW events. We
have shown that binary formation in the dense cores of DM-dominated systems
can give way to more than one population of PBHs. Our results are also in
agreement with studies investigating systems of multiple stars/BHs in the
central regions of clusters with an accretion disk [70, 71], and without it
[31, 41, 72]. There are important differences in the populations of BH
binaries that may distinguish our merger scenario (see e.g. [22]). Ultimately,
it is important to assess if the stochastic GWs of the proposed mechanism
yield a detectable signal (see e.g. [73]). There is, however, room for
important complements to our assumptions before producing accurate forecasts.
For example, since the DGs are faint with low stellar mass, we have neglected
the effect of mergers on stars, but a more detailed study contemplating such
collisions would include electromagnetic signals which will constrain the
parameters of the model. The recoil velocity of the product BH after a merger
event is also an ingredient to be included in future studies [74, 75, 76]. In
the meantime, our results indicate that a series of sequential mergers may
take place at the cores of DGs.

## Acknowledgments

We are grateful to Luis Padilla for his insightful comments on the manuscript.
EE thanks the TWAS 2021 Fellowship for Research and Advanced Training. This
work is sponsored by CONACyT grant CB-2016-282569 and by Program UNAM-PAPIIT
Grant IN107521 “Sector Oscuro y Agujeros Negros Primordiales”.

## References

  * [1] F. Acernese et al. [VIRGO], “Advanced Virgo: a second-generation interferometric gravitational wave detector,” Class. Quant. Grav. 32, no.2, 024001 (2015) [arXiv:1408.3978 [gr-qc]].
  * [2] B. P. Abbott et al. [LIGO Scientific and Virgo], “Observation of Gravitational Waves from a Binary Black Hole Merger,” Phys. Rev. Lett. 116, no.6, 061102 (2016) [arXiv:1602.03837 [gr-qc]].
  * [3] B. P. Abbott et al. [LIGO Scientific and Virgo], “Binary Black Hole Mergers in the first Advanced LIGO Observing Run,” Phys. Rev. X 6, no.4, 041015 (2016) [erratum: Phys. Rev. X 8, no.3, 039903 (2018)] [arXiv:1606.04856 [gr-qc]].
  * [4] R. Abbott et al. [LIGO Scientific and Virgo], “GWTC-2: Compact Binary Coalescences Observed by LIGO and Virgo During the First Half of the Third Observing Run,” Phys. Rev. X 11, 021053 (2021) [arXiv:2010.14527 [gr-qc]].
  * [5] R. Abbott et al. [LIGO Scientific and Virgo], “Population Properties of Compact Objects from the Second LIGO-Virgo Gravitational-Wave Transient Catalog,” Astrophys. J. Lett. 913, no.1, L7 (2021) [arXiv:2010.14533 [astro-ph.HE]].
  * [6] R. Abbott et al. [LIGO Scientific, VIRGO and KAGRA], “GWTC-3: Compact Binary Coalescences Observed by LIGO and Virgo During the Second Part of the Third Observing Run,” [arXiv:2111.03606 [gr-qc]].
  * [7] R. Abbott et al. [LIGO Scientific and Virgo], “GW190521: A Binary Black Hole Merger with a Total Mass of $150M_{\odot}$,” Phys. Rev. Lett. 125, no.10, 101102 (2020) [arXiv:2009.01075 [gr-qc]].
  * [8] A. Heger, C. L. Fryer, S. E. Woosley, N. Langer and D. H. Hartmann, “How massive single stars end their life,” Astrophys. J. 591 (2003), 288-300 [arXiv:astro-ph/0212469 [astro-ph]].
  * [9] K. Belczynski, A. Heger, W. Gladysz, A. J. Ruiter, S. Woosley, G. Wiktorowicz, H. Y. Chen, T. Bulik, R. O’Shaughnessy and D. E. Holz, et al. “The Effect of Pair-Instability Mass Loss on Black Hole Mergers,” Astron. Astrophys. 594, A97 (2016) [arXiv:1607.03116 [astro-ph.HE]].
  * [10] Y. B. Zel’dovich and I. D. Novikov, “The Hypothesis of Cores Retarded during Expansion and the Hot,” Cosmological Model, Soviet Astronomy 10 (1967), 602.
  * [11] S. Hawking, “Gravitationally collapsed objects of very low mass,” Mon. Not. Roy. Astron. Soc. 152 (1971), 75.
  * [12] B. J. Carr and S. Hawking, Mon. Not. Roy. Astron. Soc. 168 (1974), 399.
  * [13] B. J. Carr, “The primordial black hole mass spectrum,” Astrophys. J. 201 , 1 (1975).
  * [14] R. Abbott et al. [LIGO Scientific and Virgo], “GW190412: Observation of a Binary-Black-Hole Coalescence with Asymmetric Masses,” Phys. Rev. D 102, no.4, 043015 (2020) [arXiv:2004.08342 [astro-ph.HE]].
  * [15] A. Tutukov and L. Yungelson, “Evolution of massive close binaries,” Nauchnye Informatsii 27, 70 (1973).
  * [16] H. A. Bethe and G. E. Brown, “Evolution of binary compact objects which merge,” Astrophys. J. 506, 780-789 (1998) [arXiv:astro-ph/9802084 [astro-ph]].
  * [17] K. Belczynski, V. Kalogera and T. Bulik, “A Comprehensive study of binary compact objects as gravitational wave sources: Evolutionary channels, rates, and physical properties,” Astrophys. J. 572, 407-431 (2001) [arXiv:astro-ph/0111452 [astro-ph]].
  * [18] M. Arca-Sedda, G. Li and B. Kocsis, “Order in the chaos - Eccentric black hole binary mergers in triples formed via strong binary–binary scatterings,” Astron. Astrophys. 650, A189 (2021) [arXiv:1805.06458 [astro-ph.HE]].
  * [19] B. Carr, F. Kuhnel and M. Sandstad, “Primordial Black Holes as Dark Matter,” Phys. Rev. D 94, no.8, 083504 (2016) [arXiv:1607.06077 [astro-ph.CO]].
  * [20] Y. Ali-Haïmoud, E. D. Kovetz and M. Kamionkowski, “Merger rate of primordial black-hole binaries,” Phys. Rev. D 96, no.12, 123523 (2017) [arXiv:1709.06576 [astro-ph.CO]].
  * [21] M. Fishbach, D. E. Holz and B. Farr, “Are LIGO’s Black Holes Made From Smaller Black Holes?,” Astrophys. J. Lett. 840, no.2, L24 (2017) [arXiv:1703.06869 [astro-ph.HE]].
  * [22] D. Gerosa and E. Berti, “Are merging black holes born from stellar collapse or previous mergers?,” Phys. Rev. D 95, no.12, 124046 (2017) [arXiv:1703.06223 [gr-qc]].
  * [23] L. Liu, Z. K. Guo and R. G. Cai, “Effects of the merger history on the merger rate density of primordial black hole binaries,” Eur. Phys. J. C 79 (2019) no.8, 717 [arXiv:1901.07672 [astro-ph.CO]].
  * [24] Z. Doctor, D. Wysocki, R. O’Shaughnessy, D. E. Holz and B. Farr, “Black Hole Coagulation: Modeling Hierarchical Mergers in Black Hole Populations,” [arXiv:1911.04424 [astro-ph.HE]].
  * [25] Y. Wu, “Merger history of primordial black-hole binaries,” Phys. Rev. D 101, no.8, 083008 (2020) [arXiv:2001.03833 [astro-ph.CO]].
  * [26] C. Kimball, C. Talbot, C. P. L. Berry, M. Zevin, E. Thrane, V. Kalogera, R. Buscicchio, M. Carney, T. Dent and H. Middleton, et al. “Evidence for Hierarchical Black Hole Mergers in the Second LIGO–Virgo Gravitational Wave Catalog,” Astrophys. J. Lett. 915, no.2, L35 (2021) [arXiv:2011.05332 [astro-ph.HE]].
  * [27] D. Gerosa and M. Fishbach, “Hierarchical mergers of stellar-mass black holes and their gravitational-wave signatures,” Nature Astron. 5, no.8, 749-760 (2021) [arXiv:2105.03439 [astro-ph.HE]].
  * [28] C. L. Rodriguez, S. Chatterjee and F. A. Rasio, “Binary Black Hole Mergers from Globular Clusters: Masses, Merger Rates, and the Impact of Stellar Evolution,” Phys. Rev. D 93, no.8, 084029 (2016) [arXiv:1602.02444 [astro-ph.HE]].
  * [29] A. Askar, M. Szkudlarek, D. Gondek-Rosińska, M. Giersz and T. Bulik, “MOCCA-SURVEY Database – I. Coalescing binary black holes originating from globular clusters,” Mon. Not. Roy. Astron. Soc. 464, no.1, L36-L40 (2017) [arXiv:1608.02520 [astro-ph.HE]].
  * [30] J. Samsing, “Eccentric Black Hole Mergers Forming in Globular Clusters,” Phys. Rev. D 97, no.10, 103014 (2018) [arXiv:1711.07452 [astro-ph.HE]].
  * [31] M. Mapelli, M. Dall’Amico, Y. Bouffanais, N. Giacobbo, M. A. Sedda, M. C. Artale, A. Ballone, U. N. Di Carlo, G. Iorio and F. Santoliquido, et al. “Hierarchical black hole mergers in young, globular and nuclear star clusters: the effect of metallicity, spin and cluster properties,” Mon. Not. Roy. Astron. Soc. 505, no.1, 339-358 (2021) [arXiv:2103.05016 [astro-ph.HE]].
  * [32] C. Kimball, C. Talbot, C. P. L.Berry, M. Carney, M. Zevin, E. Thrane and V. Kalogera, “Black Hole Genealogy: Identifying Hierarchical Mergers with Gravitational Waves,” Astrophys. J. 900, no.2, 177 (2020) [arXiv:2005.00023 [astro-ph.HE]].
  * [33] G. Fragione and J. Silk, “Repeated mergers and ejection of black holes within nuclear star clusters,” Mon. Not. Roy. Astron. Soc. 498, no.4, 4591-4604 (2020) [arXiv:2006.01867 [astro-ph.GA]].
  * [34] K. Kritos, V. De Luca, G. Franciolini, A. Kehagias and A. Riotto, “The Astro-Primordial Black Hole Merger Rates: a Reappraisal,” JCAP 05, 039 (2021) [arXiv:2012.03585 [gr-qc]].
  * [35] G. Fragione, B. Kocsis, F. A. Rasio and J. Silk, “Repeated mergers, mass-gap black holes, and formation of intermediate-mass black holes in nuclear star clusters,” [arXiv:2107.04639 [astro-ph.GA]].
  * [36] Y. Yang, I. Bartos, V. Gayathri, S. Ford, Z. Haiman, S. Klimenko, B. Kocsis, S. Márka, Z. Márka and B. McKernan, et al. “Hierarchical Black Hole Mergers in Active Galactic Nuclei,” Phys. Rev. Lett. 123, no.18, 181101 (2019) [arXiv:1906.09281 [astro-ph.HE]].
  * [37] H. Tagawa, B. Kocsis, Z. Haiman, I. Bartos, K. Omukai and J. Samsing, “Eccentric Black Hole Mergers in Active Galactic Nuclei,” Astrophys. J. Lett. 907, no.1, L20 (2021) [arXiv:2010.10526 [astro-ph.HE]].
  * [38] H. Tagawa, B. Kocsis, Z. Haiman, I. Bartos, K. Omukai and J. Samsing, “Mass-gap Mergers in Active Galactic Nuclei,” Astrophys. J. 908, no.2, 194 (2021) [arXiv:2012.00011 [astro-ph.HE]].
  * [39] J. D. Simon, “The Faintest Dwarf Galaxies,” ARA&A 57 (2019) 375, [arXiv:1901.05465 [astro-ph.GA]].
  * [40] S. Hawking, “Black hole explosions?,” Nature 248 (1974) 30.
  * [41] C. L. Rodriguez, K. Kremer, M. Y. Grudic, Z. Hafen, S. Chatterjee, G. Fragione, A. Lamberts, M. A. S. Martinez, F. A. Rasio and N. Weatherford, et al. “GW190412 as a Third-Generation Black Hole Merger from a Super Star Cluster,” Astrophys. J. Lett. 896, no.1, L10 (2020) [arXiv:2005.04239 [astro-ph.HE]].
  * [42] D. Gerosa and E. Berti, “Escape speed of stellar clusters from multiple-generation black-hole mergers in the upper mass gap,” Phys. Rev. D 100, no.4, 041301 (2019) [arXiv:1906.05295 [astro-ph.HE]].
  * [43] B. Carr, K. Kohri, Y. Sendouda and J. Yokoyama, “Constraints on primordial black holes,” Rept. Prog. Phys. 84 (2021) no.11, 116902 [arXiv:2002.12778 [astro-ph.CO]].
  * [44] A. Drlica-Wagner et al. [DES], “Eight Ultra-faint Galaxy Candidates Discovered in Year Two of the Dark Energy Survey,” Astrophys. J. 813 (2015) no.2, 109 [arXiv:1508.03622 [astro-ph.GA]].
  * [45] J. D. Simon and M. Geha, “The Kinematics of the Ultra-Faint Milky Way Satellites: Solving the Missing Satellite Problem,” Astrophys. J. 670 (2007), 313-331 [arXiv:0706.0516 [astro-ph]].
  * [46] R. Errani, J. Peñarrubia, M. G. Walker, “Systematics in virial mass estimators for pressure-supported systems,” Mon. Not. Roy. Astron. Soc. 481, 5073-5090 (2018) [arXiv:1805.00484 [astro-ph.GA]].
  * [47] B. Willman, M. R. Blanton, A. A. West, J. J. Dalcanton, D. W. Hogg, D. P. Schneider, N. Wherry, B. Yanny and J. Brinkmann, “A New Milky Way companion: Unusual globular cluster or extreme dwarf satellite?,” Astron. J. 129, 2692-2700 (2005) [arXiv:astro-ph/0410416 [astro-ph]].
  * [48] N. F. Martin, R. A. Ibata, S. C. Chapman, M. Irwin and G. F. Lewis, “A Keck/DEIMOS spectroscopic survey of faint Galactic satellites: Searching for the least massive dwarf galaxies,” Mon. Not. Roy. Astron. Soc. 380, 281-300 (2007) [arXiv:0705.4622 [astro-ph]].
  * [49] D. Homma, M. Chiba, S. Okamoto, Y. Komiyama, M. Tanaka, M. Tanaka, M. N. Ishigaki, K. Hayashi, N. Arimoto and J. A. Garmilla, et al. “Searches for New Milky Way Satellites from the First Two Years of Data of the Subaru/Hyper Suprime-Cam Survey: Discovery of Cetus~III,” Publ. Astron. Soc. Jap. 70, 18 (2018) [arXiv:1704.05977 [astro-ph.GA]].
  * [50] B. P. M. Laevens, N. F. Martin, E. J. Bernard, E. F. Schlafly, B. Sesar, H. W. Rix, E. F. Bell, A. M. N. Ferguson, C. T. Slater and W. E. Sweeney, et al. “Sagittarius II, Draco II and Laevens 3: Three new Milky way Satellites Discovered in the Pan-starrs 1 3$\pi$ Survey,” Astrophys. J. 813, no.1, 44 (2015) [arXiv:1507.07564 [astro-ph.GA]].
  * [51] A. A. Klypin, A. V. Kravtsov, O. Valenzuela and F. Prada, “Where are the missing Galactic satellites?,” Astrophys. J. 522, 82-92 (1999) [arXiv:astro-ph/9901240 [astro-ph]].
  * [52] B. Moore, S. Ghigna, F. Governato, G. Lake, T. R. Quinn, J. Stadel and P. Tozzi, “Dark matter substructure within galactic halos,” Astrophys. J. Lett. 524, L19-L22 (1999) [arXiv:astro-ph/9907411 [astro-ph]].
  * [53] H. C. Plummer, “On the Problem of Distribution in Globular Star Clusters: (Plate 8.)” Mon. Not. Roy. Astron. Soc. 71 (1911) 460-470.
  * [54] A. W. McConnachie, “The observed properties of dwarf galaxies in and around the Local Group,” Astron. J. 144, 4 (2012) [arXiv:1204.1562 [astro-ph.CO]].
  * [55] N.Aghanim et al. [Planck], “Planck 2018 results. VI. Cosmological parameters,” Astron. Astrophys. 641, A6 (2020) [erratum: Astron. Astrophys. 652, C4 (2021)] [arXiv:1807.06209 [astro-ph.CO]].
  * [56] A. M. Green and B. J. Kavanagh, “Primordial Black Holes as a dark matter candidate,” J. Phys. G 48 (2021) no.4, 4 [arXiv:2007.10722 [astro-ph.CO]].
  * [57] Q. Zhu, E. Vasiliev, Y. Li and Y. Jing, “Primordial Black Holes as Dark Matter: Constraints From Compact Ultra-Faint Dwarfs,” Mon. Not. Roy. Astron. Soc. 476 (2018) no.1, 2-11 [arXiv:1710.05032 [astro-ph.CO]].
  * [58] J. Stegmann, P. R. Capelo, E. Bortolas and L. Mayer, “Improved constraints from ultra-faint dwarf galaxies on primordial black holes as dark matter,” Mon. Not. Roy. Astron. Soc. 492 (2020) no.4, 5247-5260 [arXiv:1910.04793 [astro-ph.GA]].
  * [59] P. Boldrini, Y. Miki, A. Y. Wagner, R. Mohayaee, J. Silk and A. Arbey, “Cusp-to-core transition in low-mass dwarf galaxies induced by dynamical heating of cold dark matter by primordial black holes,” Mon. Not. Roy. Astron. Soc. 492 (2020) no.4, 5218-5225 [arXiv:1909.07395 [astro-ph.CO]].
  * [60] L. E. Padilla, J. C. Hidalgo and K. A. Malik, “A new mechanism for primordial black hole formation during reheating,” [arXiv:2110.14584 [astro-ph.CO]].
  * [61] G. Franciolini, A. Maharana and F. Muia, “The Hunt for Light Primordial Black Hole Dark Matter with Ultra-High-Frequency Gravitational Waves,” [arXiv:2205.02153 [astro-ph.CO]].
  * [62] B. P. Abbott et al. [LIGO Scientific and Virgo], “Search for Subsolar Mass Ultracompact Binaries in Advanced LIGO?s Second Observing Run,” Phys. Rev. Lett. 123 (2019) no.16, 161102 [arXiv:1904.08976 [astro-ph.CO]].
  * [63] H. Niikura, M. Takada, N. Yasuda, R. H. Lupton, T. Sumi, S. More, T. Kurita, S. Sugiyama, A. More, M. Oguri and M. Chiba, “Microlensing constraints on primordial black holes with Subaru/HSC Andromeda observations,” Nature Astron. 3 (2019) no.6, 524-534 [arXiv:1701.02151 [astro-ph.CO]].
  * [64] S. Clesse and J. García-Bellido, “The clustering of massive Primordial Black Holes as Dark Matter: measuring their mass distribution with Advanced LIGO,” Phys. Dark Univ. 15 (2017), 142-147 [arXiv:1603.05234 [astro-ph.CO]].
  * [65] K. M. Belotsky, V. I. Dokuchaev, Y. N. Eroshenko, E. A. Esipova, M. Y. Khlopov, L. A. Khromykh, A. A. Kirillov, V. V. Nikulin, S. G. Rubin and I. V. Svadkovsky, “Clusters of primordial black holes,” Eur. Phys. J. C 79 (2019) no.3, 246 [arXiv:1807.06590 [astro-ph.CO]].
  * [66] G. D. Quinlan, S. .L. Shapiro, “Dynamical Evolution of Dense Clusters of Compact Stars” Astrophys. J. 343 (1989) 725.
  * [67] H. Mouri and Y. Taniguchi, “Runaway merging of black holes: analytical constraint on the timescale,” Astrophys. J. Lett. 566, L17-L20 (2002) [arXiv:astro-ph/0201102 [astro-ph]].
  * [68] V. D. Stasenko and A. A. Kirillov, “The Merger Rate of Black Holes in a Primordial Black Hole Cluster,” MDPI Physics 3, no.2, 372-378 (2021) [arXiv:2105.14523 [astro-ph.CO]].
  * [69] P. C. Peters, “Gravitational Radiation and the Motion of Two Point Masses,” Phys. Rev.136 (1964) B1224.
  * [70] B. McKernan, K. E. S. Ford, B. Kocsis, W. Lyra and L. M. Winter, “Intermediate-mass black holes in AGN discs – II. Model predictions and observational constraints,” Mon. Not. Roy. Astron. Soc. 441, no.1, 900-909 (2014) [arXiv:1403.6433 [astro-ph.GA]].
  * [71] J. Miralda-Escude and A. Gould, “A cluster of black holes at the galactic center,” Astrophys. J. 545, 847 (2000) [arXiv:astro-ph/0003269 [astro-ph]].
  * [72] M. Martinelli, F. Scarcella, N. B. Hogg, B. J. Kavanagh, D. Gaggero and P. Fleury, “Dancing in the dark: detecting a population of distant primordial black holes,” [arXiv:2205.02639 [astro-ph.CO]].
  * [73] J. García-Bellido, S. Jaraba and S. Kuroyanagi, “The stochastic gravitational wave background from close hyperbolic encounters of primordial black holes in dense clusters,” Phys. Dark Univ. 36 (2022), 101009 [arXiv:2109.11376 [gr-qc]].
  * [74] M. J. Fitchett, “The influence of gravitational wave momentum losses on the centre of mass motion of a Newtonian binary system,” Mon. Not. Roy. Astron. Soc. 203, (1983) 1049-1062.
  * [75] J. A. Gonzalez, U. Sperhake, B. Bruegmann, M. Hannam and S. Husa, “Total recoil: The Maximum kick from nonspinning black-hole binary inspiral,” Phys. Rev. Lett. 98 (2007), 091101 [arXiv:gr-qc/0610154 [gr-qc]].
  * [76] V. Varma, S. Biscoveanu, T. Islam, F. H. Shaik, C. J. Haster, M. Isi, W. M. Farr, S. E. Field and S. Vitale, “Evidence of Large Recoil Velocity from a Black Hole Merger Signal,” Phys. Rev. Lett. 128, no.19, 191102 (2022) [arXiv:2201.01302 [astro-ph.HE]].


# Characterisation of the MUSIC ASIC for large-area silicon photomultipliers
for gamma-ray astronomy

Nicolas De Angelis111Corresponding author.
<EMAIL_ADDRESS>DPNC, Université de Genève, 24 Quai Ernest
Ansermet, CH-1205 Geneva, Switzerland  David Gascón  Institute of Cosmos
Sciences - University of Barcelona (ICC-UB), Martí i Franquès, 1, 08028,
Barcelona, Spain  Sergio Gómez  Serra Húnter Fellow, Polytechnic University of
Catalonia (UPC), Eduard Maristany, 16, 08019, Barcelona, Spain  Institute of
Cosmos Sciences - University of Barcelona (ICC-UB), Martí i Franquès, 1,
08028, Barcelona, Spain  Matthieu Heller  DPNC, Université de Genève, 24 Quai
Ernest Ansermet, CH-1205 Geneva, Switzerland  Teresa Montaruli  DPNC,
Université de Genève, 24 Quai Ernest Ansermet, CH-1205 Geneva, Switzerland
Andrii Nagai  DPNC, Université de Genève, 24 Quai Ernest Ansermet, CH-1205
Geneva, Switzerland

(August 28, 2024)

## Abstract

Large-area silicon photomultipliers (SiPMs) are desired in many applications
where large surfaces have to be covered. For instance, a large area SiPM has
been developed by Hamamatsu Photonics in collaboration with the University of
Geneva, to equip gamma-ray cameras employed in imaging atmospheric Cherenkov
telescopes. Being the sensor about 1 cm2, a suitable preamplification
electronics has been investigated in this work, which can deal with long
pulses induced by the large capacitance of the sensor. The so-called Multiple
Use SiPM Integrated Circuit (MUSIC), developed by the ICCUB (University of
Barcelona), is investigated as a potential front-end ASIC, suitable to cover
large area photodetection planes of gamma-ray telescopes. The ASIC offers an
interesting pole-zero cancellation (PZC) that allows dealing with long SiPM
signals, the feature of active summation of up to 8 input channels into a
single differential output and it can offer a solution for reducing power
consumption compared to discrete solutions. Measurements and simulations of
MUSIC coupled to two SiPMs developed by Hamamatsu are considered and the ASIC
response is characterized. The 5th generation sensor of the Low Cross Talk
technology coupled to MUSIC turns out to be a good solution for gamma-ray
cameras.

Keywords: Photon detectors for UV, visible and IR photons (solid-state) (PIN
diodes, APDs, Si-PMTs, G-APDs, CCDs, EBCCDs, EMCCDs, CMOS imagers, etc);
Front-end electronics for detector readout; Cherenkov detectors

###### Contents

  1. 1 Introduction: Current challenges of SiPMs in gamma-ray astronomy
  2. 2 Experimental and simulation setups
  3. 3 LCT5 pulse shapes, photoelectron spectra, and linearity measurements with the MUSIC high gain channel
  4. 4 Comparison between simulations and measurements using LCT5 sensors coupled to the MUSIC ASIC
  5. 5 Response prediction of a hexagonal sensor with LCT5 technology
  6. 6 Conclusions and discussion

##  1 Introduction: Current challenges of SiPMs in gamma-ray astronomy

The next generation of cameras for gamma-ray astronomy are all moving toward
silicon photo-multipliers (SiPMs) [1, 2, 3, 4, 5, 6] for their stability,
robustness and higher sensitivity compared to photomultiplier tubes (PMTs).
Additionally, SiPMs are more and more evolving toward tailored technologies
for a given application. It is not possible to foresee that each new
technology is tried out on the field or even sometimes in the laboratory.
Therefore, it is of high relevance that each user is able to foresee by means
of simulations the response of the readout chain to a change in sensor
technology. This is even more true when the sensor has a custom size or shape
[7, 8] which requires most of the time a custom photo-mask and therefore
prohibitive cost when just a trial sample is needed.

Our work was triggered by the needs of the SST-1M camera, described in [3],
though it can be generalized. The SST-1M camera is a SiPM-based fully
digitizing camera, developed for Cherenkov imaging. It was originally
developed in the frame of the Cherenkov Telescope Array (CTA) project [9]. It
is composed of a photo-detection plane (PDP), responsible for detecting
Cherenkov light and pre-amplifying the signal, and of a digitizing system
(DigiCam) that applies a trigger selection and digitizes out the signal to the
camera server. The entrance of the camera is protected by a filtering window
to avoid dust from entering the PDP and to filter background-dominated
wavelengths222A typical Cherenkov photon has a wavelength from near-UV up to
$\sim$550 nm. Cutting off high wavelengths (above 540 nm for the entrance
window) allows for reducing the night-sky background (NSB), whose main
emission is above 600 nm. The PDP is based on hexagonal silicon photo-
multipliers (SiPMs) coupled to conically-shaped hollow light guides. The SiPMs
are read out by a two stages front-end electronics (FEE) [10]. The FEE
consists of a preamplification board and a slow control board that routes the
analog signals from the preamplification stage to DigiCam. An exploded view of
both the camera and a single PDP module is shown in Figure 1.

Figure 1: Exploded CAD view of the SST-1M camera and of a PDP module [3].

The overall optical efficiency of the camera [3] (conversion factor from
photons to photoelectrons) was found to be 18.7%. In order to increase such
value to more than 20%, as the light guide design and the window transmittance
are already optimized for high efficiency, we investigate an upgrade of the
sensor to a higher photo-detection efficiency (PDE). Two sensors from
Hamamatsu (HPK) were originally considered for this upgrade: the so-called
LCT5 sensor (S13360), based on the same low cross-talk technology as the
already adopted LCT2 sensor in the current camera, and the LVR3 (S14520),
based on the low voltage resistor technology. The results are compared with
the LCT2 sensor, with which the current camera is instrumented. The results
are presented in Figure 2. The measurements and analysis procedure are
described in [5, 11]. The LCT5 sensor, operated at the same optical crosstalk
value as the LCT2 of 8%, reaches 25.4% overall optical efficiency, while LVR3
reaches 28.0%. These two sensors not only have a higher PDE but also a higher
capacitance and quenching resistance, inducing long-tail signals. The LVR3 has
a 25% larger diode capacitance and 200% higher quenching capacitance which
translates into slower response time and a larger electronic noise. As a
result, the width of the signal integration window has to be larger for the
LVR3 than for the LCT5 resulting in higher contamination from NSB photons
hence a worse charge resolution. In addition, signals from NSB photons will
more easily pile up increasing the probability to trigger on noise for a given
threshold. The consequence would therefore be an increase in the energy
threshold. Therefore, we decided to select the LCT5 device for this study.
This decision is even more motivated when foreseeing their use in even larger
telescopes where pulses as short as 3 ns FWHM are targeted [12].

Figure 2: Comparison of SIPM peak $PDE$ (left) and crosstalk probability
$P_{XT}$ (right) as a function of over-voltage $\Delta V$.

Fine-tuning of the current analog front-end electronics of the SST-1M camera,
in order to get sharp enough signals, is not sufficient. Hence, in order to
profit from the new sensors we describe the redesign of the preamplification
stage of the front-end electronics that can deal with these higher-capacitance
SiPMs.

Given these constraints, the only ASIC at the time of this study that allows
coping with the large sensor capacitance is the Multiple Use SiPM Integrated
Circuit (MUSIC) ASIC [13, 14]. This 8-channel pre-amplifying ASIC has been
developed at the University of Barcelona. It contains several functionalities:
(1) summation of up to eight input channels and providing a differential
output; (2) individual readout of the time over threshold response, and (3)
individual analog response using a single-ended output driver. Both output
drivers are capable of driving 50 $\Omega$. Each operating mode, summation and
single input processing includes a configurable dual-gain option. Moreover,
the summation mode implements a dual-gain output providing four different
gains for this operating mode. Figure 3 shows the architecture and the main
functionalities of the MUSIC ASIC.

Figure 3: MUSIC ASIC architecture illustrating the different operating modes.
[14].

The ASIC also offers the possibility to apply a pole-zero cancellation in any
operating mode, reducing the pulse full width at half-maximum (FWMH). For the
particular case studied here, this allows reducing the signal tail and summing
from four channels of a pixel as it is done in the SST-1M front-end
electronics. This is particularly interesting for the current SiPM
technologies which exhibit large capacitance. Furthermore, this ASIC would
also reduce the power consumption of the FEE. The idea would be to use a MUSIC
board per pixel in high gain (HG) mode to readout, sum, and amplify the
signals of the 4 channels of a pixel and to send a single signal to the slow
control board. The integration of the MUSIC ASIC in the current camera
topology is shown in figure 4.

Figure 4: MUSIC ASIC as preamplification electronics for the SST-1M camera
[15].

The measurement setup used to characterize the MUSIC ASIC with the LCT5 sensor
is first described, followed by the resulting pulse shapes and photo-electron
spectra obtained by coupling the MUSIC ASIC with SiPMs of area 3$\times$3 mm2.
A study of the linearity of the HG channel is then presented. Finally, the
measured pulses are compared with the simulated MUSIC response. The study is
completed by the simulation of the LCT5 sensor in hexagonal shape enabling
performance extrapolation for future cameras.

##  2 Experimental and simulation setups

Figure 5: Schematic view of the experimental setup [15].

The measurement setup used to characterize the LCT5 sensor together with the
MUSIC is shown in Figure 5. Inside the dark box, an LCT5 sensor is soldered on
a small dedicated PCB connected to the MUSIC board. Two other boards are
connected to the MUSIC PCB, one to communicate with the computer and configure
the MUSIC ASIC, and the other to readout the single-ended and differential
outputs using SMA cables. An LED with a wavelength of 525 nm is placed in the
dark box and connected to a pulse generator placed outside the box in order to
control the light source. The single-ended and differential channels are read
out using a Teledyne Lecroy 620Zi oscilloscope. The oscilloscope is controlled
by the computer via a LabViewTM script in order to automate the data
acquisition. Two power supplies are connected to the MUSIC board: the High
Voltage power supply to bias the sensors ($\sim$40-60 V) and the Low Voltage
power supply in order to bias the MUSIC (6.5 V).

Figure 6: Picture of the experimental setup [15]

A picture of the measurement setup is shown in Figure 6. The light source
excitation pulse is set with the shortest possible width (20 ns) and a
frequency of 1 kHz. The MUSIC is configured via the command line tool of the
computer with a windows executable specifying the channels that have to be
enabled both in single-ended and differential modes. The internal gain of the
analog processing chain can also be tuned by configuring an internal register.
Finally, the PZC can be enabled and tuned using two parameters R (3-bits) and
C (5-bits) that adjust the values of the resistors and capacitor in the shaper
circuit [14].

Figure 7: SiPM model employed in the simulations.

Simulations have been performed using the Cadence Design environment with
Spectre simulator333Cadence Design Systems, Inc. https://www.cadence.com/ in
order to compare to experimental results. The simulations have been performed
with both LCT5 and LVR3 technologies, the former to compare with experimental
measurements and to extrapolate the behavior to larger sensors, and the latter
to address the concerns about the SNR expressed in Section 1. Each SiPM is
simulated following the circuit structure depicted in Figure 7.

$R_{q}$ and $C_{q}$ are the quenching resistor and capacitance respectively,
$R_{d}$ and $C_{d}$ are the diode resistance and capacitance, $C_{g}$ the
parasitic capacitance, $N_{tot}$ is the number of micro-cells, $N_{f}$ the
number of fired cells and $V_{bd}$ is the breakdown voltage of the SiPM.
$S_{A}$ is a switch implemented in verilog-A code that models the arrival time
of the photons. A complete description of the SiPM circuit model can be found
in [16]. In this simulation, the SiPM is connected to the bias voltage from
one side and is DC coupled to the input of the MUSIC chip from the other side.

It is important to highlight that the simulations have been done at the
transistor level by using the MUSIC circuit. The employment of the ASIC itself
and not of a model provides a more accurate response of the system compared to
the usage of commercial amplifiers, which contain simplified electrical
models. Note that the combination of both the MUSIC and the SiPM has been
simulated for different settings of the ASIC.

##  3 LCT5 pulse shapes, photoelectron spectra, and linearity measurements
with the MUSIC high gain channel

The LCT5 technology has been measured with all possible PZC configurations.
Figure 8 shows the measured pulse shape of a 3$\times$3 mm2 LCT5 sensor with
different R and C configurations.

Figure 8: Pulse shape for the 3$\times$3 mm2 LCT5 sensor at 4.4 V over-voltage
for different PZC configurations compared to the pulse without cancellation

The initial pulse has a length of 120 ns, while applying a PZC allows to
shorten it below 40 ns even reaching a recovery time smaller than 10 ns, thus
satisfying the required integration time (see Figure 9). As described in [14],
measurements performed using the 6$\times$6 mm2 LCT5 sensor show that it is
possible to count photons even with this larger area sensor. More measurements
with the LCT5 technology will be detailed in section 13. The LCT5 coupled to
MUSIC could therefore be a suitable solution for sensor/electronics upgrade.

Figure 9: Recovery time on the z-axis as a function of R and C for the LCT5
pulse and for each PZC configuration (left, adapted from [15]) and photo-
electron spectrum for the 3$\times$3 mm2 LCT5 sensor at 4.4 V over-voltage
(right, [15]).

The evolution of the pulse shape as a function of the light level has also
been measured in order to study the behavior in saturation of the LCT5
3$\times$3 mm2 sensor coupled to MUSIC. The differential outputs are used for
these measurements. The idea is to see if the HG channel has a predictive
behavior in saturation in order to use it to reconstruct signals up to a few
thousand of photons with a good resolution. Indeed, for cameras that would
already adopt a single gain channel, as it is for DigiCam [10], or for cameras
with a large amount of pixels with demanding power-consumption constraints,
doubling the number of channels, using both HG and LG outputs, may be
unfeasible. Figure 10 shows the evolution of the LCT5 pulse shapes with the
bias voltage of the LED in HG mode. The pulses are normalised to the maximum
amplitude in order to bring out the shape evolution. As expected, when
amplitude saturation is reached, the pulse duration, and therefore the charge,
are still increasing.

Figure 10: Pulse shapes of the LCT5 sensor as a function of light level for
the HG differential output [15] without PZC.

The amplitude and the charge are extracted from the previous pulses and are
plotted as a function of the bias voltage of the LED for the HG channel. This
is shown in figure 5.

Figure 11: Left: Amplitude and charge of the HG pulse as a function of the LED
bias voltage [15]. Right: 4-degrees polynomial fit555The polynomial fit
parameters are: $p_{4}=(1.95\pm 0.21)\cdot 10^{5}$ V-4, $p_{3}=(-1.08\pm
0.12)\cdot 10^{7}$ V-3, $p_{2}=(2.24\pm 0.26)\cdot 10^{7}$ V-2,
$p_{1}=(-2.07\pm 0.24)\cdot 10^{7}$ V-1, and $p_{0}=(7.15\pm 0.86)\cdot
10^{6}$. of the non-saturated charge as a function of the LED bias voltage and
extrapolation of the fit combined with the geometrical saturation formula (1).

The relationship between the bias voltage of the LED and the real light level
is determined by calibrating the LED. The gain for both amplitude and charge
is determined by fitting the multi-photoelectron spectrum, which allows
conversion of amplitude and charge into a number of photoelectrons corrected
for optical crosstalk. The gains obtained from the spectra are 112.08
mV$\cdot$ns for the charge and 3.4771 mV for the amplitude. Combining these
gains with the curves shown in Figure 5 leads to the number of photoelectrons
versus the LED bias voltage for both charge and amplitude. The curve obtained
for the charge is shown in Figure 5-left, where the number of photoelectrons
has been converted to a number of photons by accounting for the 36% PDE (at
$\lambda$LED=525 nm) and for the 8% optical crosstalk. The non-saturated part
of this curve, which ranges up to 1.55 V, is fitted using a 4-degree
polynomial function (which reproduces the LED behavior) in order to extract a
relation between the number of photons emitted by the LED and its bias
voltage.

To account for the geometrical saturation due to the limited number of micro-
cells in the sensor ($N_{cells}$=3584), the 4-degree polynomial fit has to be
convoluted with the following relation between the number of fired cells
$N_{f}$ and that of photons $N_{\gamma}$:

$N_{f}=N_{cells}\quantity[1-\exp(-\frac{N_{\gamma}\cdot
PDE\cdot(1+\mu)}{N_{cells}})]$ (1)

where $N_{cells}$ is the number of micro-cells, PDE is the photodetection
efficiency at the LED wavelength, and $\mu$ is the crosstalk probability.

The resulting relation between the number of fired micro-cells and the bias
voltage of the LED is plotted in Figure 5. This allows converting the x-scale
of the initial plots (Figure 5) into a number of fired micro-cells. This is
shown in Figure 12 for both amplitude and charge. We can see a good linear
behavior until about 300 fired micro-cells for the amplitude. The charge has
also a good behavior in saturation which extends the dynamic range up to 3000
photoelectrons.

Figure 12: Amplitude and charge of the HG pulse as a function
of the light level accounting for the geometrical saturation [15]

The same linearity measurements have been performed for the LCT5 3$\times$3
mm2 sensor with the intermediate PZC configuration (R=4, C=9) and same the
behavior is found for both cases.

##  4 Comparison between simulations and measurements using LCT5 sensors
coupled to the MUSIC ASIC

Some experimental measurements with the SiPMs and the MUSIC ASIC have been
performed and compared to simulations, as shown in Figure 13. In this case,
the MUSIC ASIC is configured in the individual analog readout mode, i.e., only
a single channel is evaluated. Although, the summation response when reading a
single channel behaves similarly. Moreover, the PZC circuit has been tuned
with different configurations or bypassed (indicated with PZC off). The
resistance R and capacitance C controls an RC filter to compensate for the
decay time of the SiPM. As can be seen from Figure 13, the larger the
resistance, the narrower becomes the SiPM response at the expense of larger
attenuation of the signal and the consequent loss in SNR.

It can be seen that simulations and measurements are quite similar,
especially, for the cases when the PZC is applied; the divergence in some
cases might be due to inaccuracies in the simulation with respect to the
experimental test bench. For instance, differences may arise from larger
parasitic components in the manufactured chip compared with the expected
value, packaging process, PCB design or cables.

Figure 13: LCT5 3$\times$3 mm2 pulse shapes comparison between measurements
and simulations for different PZC configurations (adapted from [15]).

##  5 Response prediction of a hexagonal sensor with LCT5 technology

In [5], a SiPM with a hexagonal-like shape was proposed using HPK technology,
referred to as S10943-2832(X) (LCT2). We intend to use the same sensor shape
for the upgraded camera, hence it is relevant to simulate the response of the
readout electronics to a LCT5 sensor of such shape.

In this case study, four LCT5 6$\times$6 mm2 (S13360-6050CS) have been used to
predict the performance of a SiPM with an hexagonal shape using the sum output
of the MUSIC ASIC. For the sake of completeness, the same was performed with
the LVR3 6$\times$6 mm2 (S14520-6050VS) technology. The number of micro-cells
in the S10943-2832(X) has been used to simulate the LCT5 and LVR3 hexagonal
sensors. Moreover, a standard square S13360-6050CS has been also used to
evaluate the performance of an existing sensor. Simulations have been
performed by injecting the signal coming from four identical 6$\times$6 mm2
SiPMs into the MUSIC circuit when PZC (R=4, C=15) is active. Table 1 shows the
SNR that can be achieved using these four types of SiPMs in terms of peak
amplitude and after integrating the first 10 ns of the pulse. As it can be
seen from Table 1, all hexagonal sensors have an equivalent capacitance
smaller than the standard square SiPM and thus electronics noise is lower for
these sensors. As said in the introduction, a large capacitance in general
leads to a small SNR. It should be noted that the SNR for the charge
measurement does not always decrease with capacitance, since it also depends
on how the charge is spread over the signal tail, i.e., how fast the SiPM is
discharged. Lastly, as expected, the SNR is smaller when considering the peak
amplitude instead of the charge.

Parameter | S10943-2832X | S13360-XX50 | S14520-XX50 | S13360-6050CS
---|---|---|---|---
| Hexagonal | Hexagonal | Hexagonal |
Number of Pixels | 9210 | 9210 | 9210 | 14400
Capacitance [pF] | 851 | 944 | 1325 | 1477
SNR (Peak amplitude) | 6.6 | 5.3 | 4.8 | 3.7
SNR (Charge) | 10.8 | 6.0 | 4.6 | 5.7
Table 1: Simulated SNR using the peak amplitude and the charge over 10 ns
integration windows for different proposed sensors. Over-voltage is set to 8
V.

Experimental measurements have been carried out connecting the MUSIC to four
square standard S13360-6050CS SiPMs (LCT5 6$\times$6 mm2). The ASIC is
configured in summation mode with PZC (R=4, C=15). This represents a signal
with a FWHM of the order of 10 ns. Figure 14 shows the analog response of the
SiPM when illuminating the sensors at the single photon regime. At 6 V over-
voltage, the different photons can be clearly distinguished in the peak
amplitude histogram. As expected from simulations, when increasing the over-
voltage to 8 V the peaks can be better identified at the expense of larger
correlated and non-correlated SiPM noise. Figure 15 shows the charge histogram
when integrating the output pulse over 10 ns. The single photo-electron
response can be already identified at 4 V over-voltage, and at 8 V the peaks
are better separated, as expected. To quantify this statement, the SNR can be
measured. It is defined as the ratio of the gain (measured as the distance in
charge between the first peak (pedestal) and the second peak, i.e., the 1st
photo-electron) and the standard deviation of the first peak (electronic
noise). Using the charge, we find two values of the SNR of 3.2 and 5.6 for 4 V
and 8 V over-voltage, respectively. The SNR for 8 V agrees with the
simulation, as shown in Table 1.

Figure 14: Response of the MUSIC connected to four S13360-6050CS SiPMs in
summation mode and with PZC (R=4, C=15). The histogram is based on the
highlighted peak of the summation output. (left) 6 V of over-voltage and
(right) 8 V of over-voltage.

Figure 15: Charge histogram of the MUSIC connected to four S13360-6050CS SiPMs
in summation mode and with PZC (R=4, C=15). (Left) 4 V of over-voltage; and
(Right) 8 V of over-voltage.

##  6 Conclusions and discussion

This work presents a method to foresee the performance of a large area sensor
combined to a preamplifying electronics based on the simple model for SiPM,
the Corsi model. Particularly for our case, where we developed a dedicated
hexagonal sensor which requires a dedicated photo-mask by a company, this
becomes very useful as one can exactly plan the required performance and then
implement the sensor according to such design.

In order to anticipate the response of the pre-amplifying stage, the model of
the sensor and the integrated circuit used for its readout have been validated
on measurements. The excellent agreement between measured and simulated pulses
has allowed anticipating with confidence the response of the MUSIC ASIC to
hexagonal sensors produced in the LCT5 and LVR3 technology. Based on this
work, predictions can be made for any custom sensor. It was shown that the
Corsi model parameters can be carefully extrapolated from an existing smaller
sensor. The properties of the available standard and smaller S13360-6050CS
SiPMs of 6$\times$6 mm2 have been used to experimentally prove that the 1 cm2
hexagonal sensors obtained with LCT5 and LVR3 technologies have photo-electron
identification capabilities. The MUSIC ASIC is perfectly capable of treating
large sensors built out of smaller units thanks to its sum output. For the
tested sensor, it produces an output pulse with approximately 10 ns FWHM and a
power consumption of about 100 mW.

If this performance is enough for small telescopes which work in the highest
energy range, further improvements must be achieved to reach pulse widths as
low as 3 ns FWHM. This will help improving performance in the lowest energy
range of the large-sized telescopes of CTA. A future ASIC, based on the
8-channel FastIC ASIC developed in a 65 nm technology node [17, 18] is planned
for the future upgrade of the LSTs cameras using SiPMs. This new ASIC should
have 16 input channels, the capability to perform an active summation in
groups of four channels, the option to add the signals from these 16 channels
into a single output, a mechanism to capture the level of continuous night sky
background (NSB) noise and the possibility to drive cable loads. The power
consumption of the ASIC should be around 10 mW/ch. As a second version of the
chip, an ADC is planned to be added in order to provide a digital output
directly from the ASIC.

## Acknowledgements

We acknowledge the support of the University of Geneva, the Swiss National
Foundation FLARE Program, grant 20FL21_201539 and of the State Secretariat for
Education Research and Innovation on CTAO. N.D.A. acknowledges the support of
the Swiss National Science Foundation. We acknowledge support by the Spanish
Ministerio de Ciencia, Innovación y Universidades 491 (MICINN), grant
FPA2017-82729-C6-2-R (AEI/FEDER, UE) and grant PID2019-104114RB-C33 / AEI /
10.13039/501100011033. We also acknowledge financial support from the State
Agency for Research of the Spanish Ministry of Science and Innovation through
the “Unit of Excellence María de Maeztu 2020-2023” award to the Institute of
Cosmos Sciences (CEX2019-000918-M).

## References

  * [1] Robert Klanner  “Characterisation of SiPMs”  In _Nuclear Instruments and Methods in Physics Research A_ 926, 2019, pp. 36–56  DOI: 10.1016/j.nima.2018.11.083
  * [2] Giovanni Ambrosi and V. Vagelli  “Applications of silicon photomultipliers in ground-based and spaceborne high-energy astrophysics”  In _The European Physical Journal Plus_ 137, 2022  DOI: 10.1140/epjp/s13360-021-02159-4
  * [3] M. Heller et al.  “An innovative silicon photomultiplier digitizing camera for gamma-ray astronomy”  In _The European Physical Journal C_ 77, 2017, pp. 47  DOI: https://doi.org/10.1140/epjc/s10052-017-4609-z
  * [4] Stefan Gundacker and Arjan Heering  “The silicon photomultiplier: fundamentals and applications of a modern solid-state photon detector”  In _Physics in Medicine and Biology_ 65.17, 2020, pp. 17TR01  DOI: 10.1088/1361-6560/ab7b2d
  * [5] A. Nagai et al.  “Characterization of a large area silicon photomultiplier”  In _Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment_ 948, 2019, pp. 162796  DOI: https://doi.org/10.1016/j.nima.2019.162796
  * [6] T. Montaruli and A. Nagai  “Properties of large SiPM at room temperature”  In _Nucl. Instrum. Meth. A_ 985, 2021, pp. 164609  DOI: 10.1016/j.nima.2020.164609
  * [7] Vittorio Boccone et al.  “Characterization of New Hexagonal Large Area MPPCs”  In _IEEE Transactions on Nuclear Science_ 61.3, 2014, pp. 1474–1482  DOI: 10.1109/TNS.2014.2321339
  * [8] F. Aharonian  “Construction and on-site performance of the LHAASO WFCTA camera”  In _Eur. Phys. J. C_ 81.7, 2021, pp. 657  DOI: 10.1140/epjc/s10052-021-09414-z
  * [9] “Cherenkov Telescope Array Observatory web page”, https://www.cta-observatory.org
  * [10] J.A. Aguilar et al.  “The front-end electronics and slow control of large area SiPM for the SST-1M camera developed for the CTA experiment”  In _Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment_ 830, 2016, pp. 219–232  DOI: https://doi.org/10.1016/j.nima.2016.05.086
  * [11] A. Nagai et al.  “SENSE: A comparison of photon detection efficiency and optical crosstalk of various SiPM devices” New Developments In Photodetection 2017  In _Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment_ 912, 2018, pp. 182–185  DOI: https://doi.org/10.1016/j.nima.2017.11.018
  * [12] M. Heller  “Development of an advanced SiPM camera for the Large Size Telescope of the Cherenkov Telescope Array”  In _37th International Cosmic Ray Conference_ , 2021  arXiv:2108.10112 [astro-ph.IM]
  * [13] Sergio Gómez et al.  “MUSIC: An 8 channel readout ASIC for SiPM arrays.”  In _Proceeding SPIE 9899, Optical Sensing and Detection IV_ , 2016  DOI: https://doi.org/10.1117/12.2231095
  * [14] Sergio Gómez et al.  “Multiple Use SiPM Integrated Circuit (MUSIC) for Large Area and High Performance Sensors”  In _Electronics_ 10.8, 2021  DOI: 10.3390/electronics10080961
  * [15] Nicolas De Angelis  “Studies of readout electronics and optical elements for a gamma-ray telescope. Master thesis. University of Geneva”, 2019  URL: http://cds.cern.ch/record/2683289
  * [16] David Sánchez et al.  “Multimodal simulation of large area silicon photomultipliers for time resolution optimization”  In _Nuclear Instruments and Methods in Physics Research Section A: Accelerators, Spectrometers, Detectors and Associated Equipment_ 1001, 2021, pp. 165247  DOI: https://doi.org/10.1016/j.nima.2021.165247
  * [17] S. Gómez et al.  “FastIC: a fast integrated circuit for the readout of high performance detectors”  In _Journal of Instrumentation_ 17.05 IOP Publishing, 2022, pp. C05027  DOI: 10.1088/1748-0221/17/05/c05027
  * [18] S. Gómez et al.  “FastIC: A Highly Configurable ASIC for Fast Timing Applications”  In _2021 IEEE Nuclear Science Symposium and Medical Imaging Conference (NSS/MIC)_ , 2021, pp. 1–4  DOI: 10.1109/NSS/MIC44867.2021.9875546