††thanks: These authors contributed equally to this work.††thanks: These
authors contributed equally to this work.††thanks: These authors contributed
equally to this work.††thanks: These authors contributed equally to this
work.††thanks: These authors contributed equally to this work.

# A hybrid integrated dual-microcomb source

Nikita Yu. Dmitriev  Russian Quantum Center, Moscow, 143026, Russia  Moscow
Institute of Physics and Technology (MIPT), Dolgoprudny, Moscow Region,
141701, Russia  Sergey N. Koptyaev  Samsung R&D Institute Russia, SAIT-Russia
Laboratory, Moscow, 127018, Russia  Andrey S. Voloshin  Institute of Physics,
Swiss Federal Institute of Technology Lausanne (EPFL), CH-1015 Lausanne,
Switzerland  Nikita M. Kondratiev  Russian Quantum Center, Moscow, 143026,
Russia  Valery E. Lobanov  Russian Quantum Center, Moscow, 143026, Russia
Kirill N. Min’kov  Russian Quantum Center, Moscow, 143026, Russia  Maxim V.
Ryabko  Samsung R&D Institute Russia, SAIT-Russia Laboratory, Moscow, 127018,
Russia  Stanislav V. Polonsky  Samsung R&D Institute Russia, SAIT-Russia
Laboratory, Moscow, 127018, Russia  Igor A. Bilenko  Russian Quantum Center,
Moscow, 143026, Russia  Faculty of Physics, M.V. Lomonosov Moscow State
University, 119991 Moscow, Russia

## Abstract

Dual-comb interferometry is based on self-heterodyning two optical frequency
combs, with corresponding mapping of the optical spectrum into the radio-
frequency domain. The dual-comb enables diverse applications, including
metrology, fast high-precision spectroscopy with high signal-to-noise ratio,
distance ranging, and coherent optical communications. However, current dual-
frequency-comb systems are designed for research applications and typically
rely on scientific equipment and bulky mode-locked lasers. Here we demonstrate
for the first time a fully integrated power-efficient dual-microcomb source
that is electrically driven and allows turnkey operation. Our implementation
uses commercially available components, including distributed-feedback and
Fabry–Perot laser diodes, and silicon nitride photonic circuits with
microresonators fabricated in commercial multi-project wafer runs. Our devices
are therefore unique in terms of size, weight, power consumption, and cost.
Laser-diode self-injection locking relaxes the requirements on microresonator
spectral purity and Q-factor, so that we can generate soliton microcombs
resilient to thermal frequency drift and with pump-to-comb sideband efficiency
of up to 40% at mW power levels. We demonstrate down-conversion of the optical
spectrum from 1400 nm to 1700 nm into the radio-frequency domain, which is
valuable for fast wide-band Fourier spectroscopy, which was previously not
available with chip-scale devices. Our findings pave the way for further
integration of miniature microcomb-based sensors and devices for high-volume
applications, thus opening up the prospect of innovative products that
redefine the market of industrial and consumer mobile and wearable devices and
sensors.

##  I Introduction

Over the past few decades, optical frequency combs have become a versatile
tool for addressing scientific and technical challenges Udem  _et al._ (2002);
Fortier and Baumann (2019); Diddams  _et al._ (2020). One of the most
promising and widely employed applications is the use of a double (or, dual)
optical combs for the efficient transfer of signals from the optical domain
into the radio frequency (RF) range, thus greatly simplifying data acquisition
and subsequent processing. For instance, dual-comb spectroscopy is a
remarkable form of Fourier spectroscopy, enabling ultra-fast measurement of
broadband optical absorption spectra that provide fingerprints of specific
materials or their quantity in the sample using a single photodetector,
without the need for moving components, and only a few basic optical
components Suh  _et al._ (2016); Coddington  _et al._ (2016); Ideguchi (2017);
Picqué and Hänsch (2019); Baumann  _et al._ (2019); Koptyaev  _et al._ (2019).
The basic idea of the dual-comb technique is to combine two coherent optical
frequency combs with shifted pump lines ($f_{1}$,$f_{2}$) and slightly
different line spacing in the frequency domain ($f_{rep1}$,$f_{rep2}$) (Fig.
1(a)). Thereby, the optical spectrum of the combs transmitted through the
matter is down-converted at the photodetector into the RF band for
measurement. The resulting signal is a RF frequency comb with
$\delta=|f_{rep1}-f_{rep2}|$ line spacing, a central line located at
$\Delta=|f_{1}-f_{2}|$, and line amplitudes uniquely defined by the
corresponding lines of the optical combs. These characertistics make dual-comb
techniques highly attractive for myriad practical applications, such as ultra-
broadband near-IR spectroscopy Okubo  _et al._ (2014, 2015), near-field
microscopy for sub-wavelength spatial resolution Brehm  _et al._ (2006); von
Ribbeck  _et al._ (2008), precision metrology of molecular-line center
frequencies Zolot  _et al._ (2013), greenhouse-gas monitoring Baumann  _et
al._ (2011); Rieker  _et al._ (2014); Zhu  _et al._ (2015), combustion
diagnosis Schroeder  _et al._ (2017), and distance ranging (LIDAR) Coddington
_et al._ (2009); Suh and Vahala (2018); Trocha  _et al._ (2018); Nürnberg  _et
al._ (2021); Riemensberger  _et al._ (2020).

There exist a number of rapidly developing approaches for implementing dual-
comb techniques. Various dual-comb systems are based on conventional fibre
mode-locked lasers (MLLs) Coddington  _et al._ (2010); Hoghooghi  _et al._
(2021); Kayes  _et al._ (2018); Link  _et al._ (2016). The use of mode-locked
integrated external-cavity surface emitting lasers (MIXSELs) potentially
allows dual-comb signal generation with a single cavity, exploiting different
polarisation states Link  _et al._ (2015); Nürnberg  _et al._ (2019). Hybrid
THz dual-comb spectrometers based on quantum cascade lasers (QCLs) Consolino
_et al._ (2020) and coherently averaged dual-comb spectrometers Komagata  _et
al._ (2021) demonstrated fast and high-accuracy spectroscopy measurements in
the mid-wavelength infrared (MWIR) and long-wavelength infrared (LWIR) ranges.
These results have been obtained in laboratory settings, which means that many
of these dual-comb systems are only partially integrated and rather bulky and
complicated, requiring a range of auxiliary equipment and technical expertise.
As a consequence, they are not suitable for industrial and consumer
applications, despite outstanding performance in laboratory environments.
Therefore, integration and device miniaturisation are pressing issues that
need to be addressed before any of the approaches can be utilized in
industrial-grade device. Fully integrated dual-comb systems hold the promise
to unlock major applications, including airborne and spaceborne sensors,
distance ranging with unprecedented speed and resolution, and compact
spectroscopy sensors. In turn, such devices could become a core technology for
consumer and wearable applications, including non-invasive spectroscopic
sensors.

The most promising platform for a fully integrated dual-comb source is
silicon-based integrated photonics. In recent years, this platform has
experienced major advances and reached considerable maturity Pfeiffer  _et
al._ (2016); Liu  _et al._ (2018, 2021). Today, the performance of low-loss
silicon-based photonic systems has become comparable with that of free-space
optic systems. In addition, high-level compatibility with CMOS fabrication
processes Hochberg and Baehr-Jones (2010); Hochberg  _et al._ (2013); Weimann
_et al._ (2017) and also with the III–V semiconductor platform Xiang  _et al._
(2021a, b); Cuyvers  _et al._ (2021) have been demonstrated. Recent progress
in silicon photonics, combined with the self-injection locking (SIL) effect
Dahmani  _et al._ (1987); Laurent  _et al._ (1989); Oraevsky  _et al._ (2001);
Liang  _et al._ (2010, 2015a, 2015b); Kondratiev  _et al._ (2017); Galiev  _et
al._ (2020); Voloshin  _et al._ (2021), enabled optical microcomb generation
using semiconductor laser diodes (LDs) instead of bulky narrow-linewidth
lasers, thereby greatly simplifying the process of microcomb generation and
paving the way for the development of fully integrated chip-scale single
microcomb sources based on high-Q microresonators (MRs) Xiang  _et al._
(2021a, b); Shen  _et al._ (2020); Raja  _et al._ (2019); Stern  _et al._
(2018). Proof-of-concept experiments based on passive high-Q MRs achieved
broadband optical spectrum down-conversion to the RF range, generating ultra-
wide frequency combs with widely variable (from GHz to THz) line spacing in
the near-infrared (NIR), SWIR, and visible-wavelength ranges in bulk and on-
chip structures Kippenberg  _et al._ (2011); Lin  _et al._ (2017); Pasquazi
_et al._ (2018); Gaeta  _et al._ (2019); Kovach  _et al._ (2020); Lin  _et
al._ (2020); Bao  _et al._ (2019); Dutt  _et al._ (2018); Wang  _et al._
(2018); Yu  _et al._ (2018); Pavlov  _et al._ (2017). Also, the recently
demonstrated scanning dual-comb spectroscopy (SDCS) technique allows
additionally to increase the resolution of spectroscopic systems based on
high-Q MRs Lin  _et al._ (2020).

Figure 1: Principle of the hybrid integrated dual-microcomb source a. Sketch
of the experimental setup illustrating integrated dual-comb signal generation
based on the two separate soliton microcomb sources (marked with green and
blue rectangles): LD – semiconductor laser diode; SiN – photonic chip with
high-Q silicon nitride microresonator; TC and CC – temperature and current
controllers; FPD – fast photo detector; OSA – optical spectrum analyzer; ESA –
electrical spectrum analyzer.b. Photograph of the portable turnkey dual-comb
source comprising two standalone matched integrated soliton microcomb sources
c. Evolution of the experimentally observed dual-comb signal based on two
matched soliton microcombs during the search for the operating point.
Demonstration of the dual-comb signal transition from a noisy (Area I) to a
soliton state (Area II) by fine-tuning the LD current. d. Various concepts of
dual-comb spectrometers: two separate microcomb-generating photonic chips
pumped with two LDs; single photonic chip generating two microcombs pumped
with two LDs; and single photonic chip generating two microcombs pumped with a
single LD. FPGA is a field-programmable gate array.

Here, we report a feasibility study of dual-comb integration and introduce the
first hybrid integrated dual-microcomb source for the SWIR range based on
commercially available low-cost components (Fig. 1(b)). With the assembled
prototype we have successfully down-converted a 300-nm wide optical spectrum
to a 600-MHz wide RF signal. Our findings establish that electrically driven
soliton microcombs comprising integrated SiN high-Q MRs combined with SIL
semiconductor LDs (Fig. 1(a)) are a promising technology platform for highly
integrated energy-efficient dual-comb sources covering wide spectral ranges.
Specifically, we demonstrate that SIL provides up to 40% pump-to-comb sideband
conversion efficiency ($\eta_{p2c}$) for bright solitons. In the light of such
a high conversion value we consider in detail the $\eta_{p2c}$ efficiency of
SIL-enabled generation of bright dissipative Kerr solitons in photonic chip-
based microresonators, and its dependence on key parameters. In this way we
found that SIL relaxes the requirements on MR properties such as Q-factor,
spectral purity, the number of mode crossings, the width of a so-called
‘soliton step’, making it possible to generate microcombs with the majority of
commercial photonic chip-based MRs featuring moderate Q-factors. SIL greatly
facilitates tuning to the soliton regime, due to the compensation of the
thermal effects inevitable in systems where a free-running laser is used as a
pump source. Consequently, SIL enables soliton microcomb generation in cases
where it is otherwise not possible with optically isolated external-cavity
diode lasers (ECDLs).

The here-presented prototype of the integrated dual-microcomb source, based on
laser diode self-injection locked to a microresonator, allows to combine all
the benefits of Fourier-transform infrared broadband spectroscopy in a chip-
scale spectroscopic sensor and looks promising as a platform for future mobile
and wearable devices. The demonstrated versatile approach to dual-microcomb
source integration provides various design options (Fig1(d)), offering
interesting perspectives in terms of device miniaturization and performance,
in particular with a view to broadband infrared dual-comb sensors for high-
volume applications. With further integration, there is a clear route to
satisfying the so-called SWaP-C (Size, Weight, Power and Cost) requirements
that are of central importance for industrial, airborne, space, and consumer
applications.

##  II Results

###  II.1 Hybrid integrated platform for optical dual-microcomb source

Our hybrid integrated dual-microcomb source comprises specially matched (see
Methods and Supplementary Note 1) microcomb sources. Each consists of a
thermally stabilized SiN photonic chip with a high-Q MR, a butt-coupled
semiconductor LD, and an output lensed fiber (Fig. 1 (a,b)). This versatile
approach enables fast prototyping by testing different photonic-chip designs
and various types of LDs.

The MRs based on CMOS-compatible SiN photonic chips that we used in our
experiments were fabricated in commercial multi-project wafer (MPW) runs. We
use two sets of chips with MRs of two diameters, corresponding to $\sim$150
GHz and to $\sim$1 THz free spectral range (FSR) with integrated microheaters
enabling grid matching of the eigenfrequencies of different MRs by tuning
their FSR (spacing between fundamental modes in the frequency domain). Also,
each chip has edge waveguide couplers, ensuring insertion losses as low as
-1.1 dB on both sides of a chip for coupling light in and out.

Fabry–Perot (FP) and distributed feedback (DFB) laser diodes have been used
for experiments and were compared in terms of performance for microcomb
generation (see Methods and Supplementary Note 2). FP diodes have a single
spatial mode, 35-GHz longitudinal mode spacing, 1535-nm central wavelength,
and $\sim$ 200 mW optical power at 500 mA of injection current. DFB diodes
have 1545 nm wavelength and optical power of $\sim$ 100 mW at 400 mA.

During the experiment we simultaneously monitored the optical spectra of the
microcombs and the resulting RF dual-comb signal using an optical spectrum
analyzer (OSA) and an electrical spectrum analyzer (ESA). Owing to the SIL
effect, the assembled microcomb sources offer turnkey operation. This kind of
turnkey operation was described in Shen  _et al._ (2020) and also demonstrated
in Jin  _et al._ (2021). The spectrogram presented in Fig. 1(c) shows dual-
comb signal evolution while selecting the turnkey operating point by slow
manual tuning of the LD injection current. Away from the operating point, the
optical spectrum consists of one coherent soliton microcomb and one chaotic
non-coherent microcomb (MI state). The heterodyne beatnote signal of this
state is noisy (regime I in Fig. 1(c)). When the LD-current value reaches the
operating point, the system locks to the state with both soliton microcombs,
featuring high mutual coherence and low-noise RF beatnotes of the optical
components (regime II in Fig. 1(c)). The SIL mechanism compensates thermal
effects and the microcomb sources lock to the comb states without additional
manipulations, which are inevitable when pumping with a free-running laser. At
the operating point, our dual-microcomb source quickly transits to a coherent
state. However, it should be noted that due to the independent thermal
stabilization of the two photonic chips, a relative thermal drift of the RF
beat frequencies arises for observation periods exceeding several minutes.

Figure 2: Soliton microcomb generation using self-injection-locked laser-diode
pumping in the case of an inaccessible conventional soliton step a. Photograph
of the silicon nitride photonic chip used in the experiment b. Nonlinear-
resonance shape for the different pump powers of the external optically
isolated laser, demonstrating inaccessibility of the conventional soliton
step. The measured nonlinear threshold of this microresonator is approximately
13 mW. c. Photograph of the assembled prototype of the microcomb source, based
on the same photonic chip pumped with a butt-coupled Fabry–Perot laser diode.
Inset: Inside view the prototype. d. Nonlinear-resonance shape for 35-mW pump
power for the same resonance shown in (b) for the case of the self-injection
locking effect (laser diode pumping). The area where the soliton exists is
highlighted in green colour. e. Output of the microcomb source (blue line) and
theoretically predicted envelope of the single soliton state microcomb
spectrum (green line) based on the measured microresonator parameters (see
Supplementary Note 4). Left inset: RF spectrum of the output signal. Right
inset: Beatnote signal of the generated microcomb line and a Toptica CTL-1550
tunable laser recorded with a resolution bandwidth (RBW) of 10kHz (blue
points), together with a Voigt-profile fit (red line) with 0.5 kHz and 68 kHz
Lorentzian and Gaussian widths, respectively.

###  II.2 Highly efficient soliton microcombs for dual-microcomb source

Current technology used in the fabrication of high-Q silicon nitride MRs,
especially for commercially available runs, does not guarantee Q-factors
higher than one million and a spectral purity of MRs sufficient high for
sustainable soliton microcomb generation using external pumping with single-
frequency narrow-linewidth lasers. The procedure of soliton comb excitation by
means of an external pump using ECDL requires accessible soliton steps and is
complicated by the need for additional equipment to achieve a high tuning rate
and to overcome thermal instabilities Herr  _et al._ (2014a). We note that in
our experiments only high-FSR MRs (FSR=1 THz) provide easily accessible
soliton steps and support soliton generation using ECDL.

Figure 3: High-power and broadband microcombs in 150 GHz and 1 THz microchips.
a,c. Measured and fitted anomalous dispersion landscape in a $Si_{3}N_{4}$
microresonator with FSR = 143.6 and 999.8 GHz, and estimated second-order
dispersion coefficient $D_{2}/2\pi$ $\approx$ 1.38 and $\approx$ 14.3 MHz,
respectively. b,d. Optical spectra of generated high-power microcomb with
$\sim$150 GHz (b) and $\sim$1 THz (d) repetition rate (blue lines) and values
of key parameters: total comb power ($P_{comb}$), laser diode power
($P_{pump}$), the power of the microcomb central line (pump line) ($P_{pump\
line}$), and pump-to-comb sideband conversion efficiency ($\eta_{p2c}$).
Theoretically predicted envelope for corresponding microcombs based on
measured microresonator parameters (see Supplementary Note 4) for single-
soliton state (green line) and multi-soliton state (purple line). e,f.
Experimentally measured resonances in the linear regime (blue line) and loaded
Lorentzian profile fits (red dashed line) for $\sim$150 GHz (e) and $\sim$1
THz (f) microresonators. The estimated intrinsic loss ($\kappa_{0}$), coupling
rate ($\kappa_{c}$) and backward-wave coupling rate ($\gamma$) correspond to a
coupling coefficient $\eta$ = 0.42 and 0.59, and Q-factor $\sim 1.8\times
10^{6}$ and $\sim 1\times 10^{6}$ for 150 GHz and 1 THz microresonators,
respectively.

We have not been able to observe soliton generation with 150 GHz MRs using the
external amplified ECDL providing more than 150 mW of in-chip pump power. This
power level is more than ten times higher than the parametric instability
threshold, corresponding to the normalized pump amplitude
$f=\sqrt{P_{pump}/P_{th}}>3$, where $P_{pump}$ is the optical power of the
amplified ECDL reduced by losses for the butt-coupling of the lensed fiber
with the chip (power in the bus waveguide) and $P_{th}$ is the nonlinearity
threshold power Herr  _et al._ (2014a). Apparently, the influence of thermal
processes, high-order dispersion (Fig. 3(a)) and avoided mode crossing points
shortens a soliton step and makes it inaccessible (see Fig. 2(b)), Herr  _et
al._ (2014b); Kim  _et al._ (2021); Li  _et al._ (2017); Joshi  _et al._
(2016); Kuse  _et al._ (2020).

However, the fully integrated SIL scheme (Fig. 2(c)) provides outstanding
turnkey operation without any additional equipment. The same MR pumped by the
self-injection locked LD allows to clearly observe a soliton step in the
locked state, even for $f\sim 1.6$ (Fig. 2(d,e)). Indeed, most of the thermal
effects are suppressed as the laser frequency is locked to the MR and the
laser–microresonator detuning $\zeta_{\rm eff}$ is fixed. If the MR frequency
fluctuates due to the thermal effects, the generation frequency also changes,
keeping the comb-generation regime stable Kondratiev  _et al._ (2021). In the
SIL regime the laser–microresonator detuning becomes fixed to the value
$\zeta_{\rm eff}^{0}\approx-3(f^{2}/2)^{1/3}+(f^{2}/2)^{-1/3}$ (for small
backscattering and $f>1$), which lies inside the soliton-generation region
$\zeta_{\rm eff}\in[-3(f^{2}/4)^{1/3}+(f^{2}/4)^{-1/3}/4;-\pi^{2}f^{2}/8]$
(see Supplementary Note 3).

In addition, the detuning control is much more robust in the SIL regime as the
effective detuning does not change while the LD is tuned within the locking
range. More precisely, the speed of the laser-frequency tuning is effectively
reduced by the factor of the stabilization coefficient $K_{0}$ (see
Supplementary Note 3).

Figure 4:  Dual-microcomb signals Top row: Optical spectra of the generated
microcombs that were combined for dual-microcomb signal generation. Inset:
Zoom-in on the central line area, with the distance between them ($\Delta$)
indicated. Bottom row: The dual-microcomb signal is an RF comb with central
line at $\Delta$ and repetition rate $\delta=|f_{rep1}-f_{rep2}|$, where
$FSR_{i}$ are the repetition rates of the combined optical combs. The top
full-span spectrum (from 0 to 15 GHz) demonstrates the absence of low-
frequency noise, showing a strong dual-microcomb signal. The bottom spectrum
is a zoom-in on the central (Frequency - $\Delta$) dual-microcomb signal area.
The red line is the experimentally observed signal and the grey line the
predicted signal based on the optical spectra above. Inset: estimation of the
RF-comb linewidth from a Voigt-profile fit. a,b,c. Dual-comb signal obtained
by combining two optical frequency combs with $f_{rep1}\sim f_{rep2}\sim$1 THz
and $\delta=|f_{rep1}-f_{rep2}|\approx$20 MHz. The distance between the pump
lines of the optical microcombs is $\Delta=$7.93 GHz. d,e,f. Dual-microcomb
signal obtained by combining two optical frequency combs with repetition rates
$f_{rep1}\sim$1 THz, $f_{rep2}\sim$2 THz and
$\delta=|2f_{rep1}-f_{rep2}|\approx$95 MHz. The distance between the pump
lines of the optical microcombs is $\Delta=$1.70 GHz. g,h,i. Dual-microcomb
signal obtained by combining two optical frequency combs with $f_{rep1}\sim
f_{rep2}\sim$150 GHz and $\delta=|f_{rep1}-f_{rep2}|\approx$18 MHz. The
distance between the pump lines of the optical microcombs is $\Delta=$14.06
GHz.

The assembled prototype works as a turnkey device, and thereby greatly
simplifies the process of comb generation and significantly improves its
stability. Once calibrated to define the operating point, this device can
generate microcombs immediately after being turned on.

Microcomb spectra and the microresonator parameters are shown in Fig. 3.
Notably, there were no optical elements filtering the pump line. Despite the
moderate Q-factor ($\sim 10^{6}$), the generated combs are broadband. The
spectrum width of the microcombs with 1-THz line spacing reaches 500 nm, with
20 lines of power > -20 dBm. The width of the 150-GHz microcomb spectrum
exceeds 200 nm, with 30 lines of power > -20 dBm. The generated microcombs
also feature high optical power per line and therefore a high signal-to-noise
ratio, and demonstrates high pump-to-comb sideband power conversion efficiency
($\eta_{\rm p2c}$). The latter can be expressed as $\eta_{\rm
p2c}=(P_{comb}-P_{pump\ line})/P_{pump}$, where $P_{comb}$ is the total
generated microcomb power in the output fiber, $P_{pump\ line}$ is the power
of the microcomb central line (pump line) in the output fibre, and $P_{pump}$
is the optical power of the free-running laser diode for the same injection
current value reduced by the coupling losses (power in the bus waveguide). For
the 150-GHz microcomb, $P_{comb}\approx$ 20 mW, $P_{pump\ line}\approx$ 11 mW,
and $P_{pump}\approx$ 35 mW; for the 1-THz microcomb, $P_{comb}\approx$ 20 mW,
$P_{pump\ line}\approx$ 3 mW, and $P_{pump}\approx$ 40 mW. Evaluated
$\eta_{\rm p2c}$ efficiency values for the 150-GHz and 1-THz FSR microcombs
are 25$\%$ and 40$\%$, respectively. Before this work, comparable efficiency
has been demonstrated for the generation of dark microcombs only Xue  _et al._
(2017); Kim  _et al._ (2019). The high efficiencies obtained are in a good
agreement with our estimates (see Supplementary Note 4) and are compatible
with Jang  _et al._ (2021).

###  II.3 Spectral characteristics of dual-microcomb source

Combined microcomb optical spectra and the resulting RF dual-microcomb signals
are shown in Figs 4 (a,d,g) and in 4(b,c,e,f,h,i), respectively. Insets in the
bottom row give information about the linewidth of the lines of the generated
dual-microcomb signal, estimated with Voigt-profile fits ($\mathrm{W_{L}}$ is
the Lorentzian linewidth and $\mathrm{W_{G}}$ the Gaussian linewidth). With
1-THz FSR microresonators we successfully down-converted $\sim$300-nm-wide
optical spectra to 600-MHz-wide spectra in the RF range (Fig. 4 (a-f)); with
150-GHz FSR microresonators we achieved down-conversion from $\sim$100 nm to
800 MHz (Fig. 4(g-i)).

By applying voltage to the microheaters we can adjust the microresonator
temperature and thereby control the microcomb line spacing. One optical
microcomb can be shifted relatively to the other to control the dual-microcomb
signal, changing its central line position $\Delta$ and repetition rate
$\delta$. This capability is illustrated in Fig. 4(a–f). These two dual-
microcomb signals were observed for the same pair of MRs, and using
microheaters we changed the central-frequency difference $\Delta$ from 7.93
GHz to 1.70 GHz.

Verification of the measurement data is conducted by comparing the
experimental data with the theoretically predicted beatnote signal of the two
combined optical microcombs. Knowing the MR parameters, including dispersion
profile, and with the measured optical frequency combs, we can calculate the
expected dual-comb spectrum profile in the RF domain. These simulated dual-
microcomb spectra are shown as grey lines in Fig. 4(c,f,i).

##  III Discussion

We found that SIL relaxes the requirements on the microresonator Q-factor and
its spectral purity, and makes widely available microresonators commercially
fabricated in MPW runs suitable for being used in on-chip dual-microcomb
sources. Note that soliton-comb excitation with such microresonators pumped by
a tunable optically isolated ECDL (not SIL) was unsuccessful, because of the
required soliton steps could not be accessed (Fig. 2(b)). We connect this
failure with the influence of the thermal processes, high-order dispersion
effects, and avoided mode crossing points, which makes the soliton step
shorter and therefore inaccessible. When exploing the SIL effect, the
frequencies of the laser diode and the microresonator are connected and
fluctuate in a correlated manner, resulting in microcomb generation with
higher tolerance to thermal drift (Fig. 2(d)).

The demonstrated pump-to-comb sideband conversion efficiency of up to 40% can
be explained by the matching between Q-factor and pump power, taking into
account the following considerations. First, the comb power does not depend
directly on the pump power (see Supplementary Note 4). Although both the
maximal and the locked detunings do depend on the pump power, these
dependencies are weaker, and the pump-to-soliton (and pump-to-comb sideband)
power conversion efficiency decreases with it. This creates the illusion that
using a low-power laser and reducing the threshold we can make a highly
energy-efficient device. However, the second point is that the comb power does
depend on the threshold power of the parametric instability. Therefore, using
the strategy described above, we end up with a negligible output signal. This
brings us to the rather counter-intuitive conclusion that for best
performance, the threshold power should be increased (meaning lower Q-factors
or less nonlinearity, for example) and matched with the used laser pump power.
This immediately brings up a trade-off problem, as the lower Q-factor means
less stabilization and wider beatnote or even the the self-injection locking
regime cannor be reached at all. Another way to increase the pump-to-comb
sideband power conversion efficiency is increasing the second-order dispersion
coefficient. This presents, however, a trade-off problem for the comb width,
which should be solved separately for the desired application.

We also can see an additional mechanism to increase the comb power (and
subsequently the power efficiency). In our system, multi-soliton states can
emerge; however, the number of solitons has a maximum for a given dispersion
value Karpov  _et al._ (2019) and the pump-to-comb sideband conversion
efficiency ($\eta_{p2c}$) also saturates with it (see Supplementary Note 4).
The multi-soliton option was realized for the 1-THz microresonator (see Fig.
3d), where we found good correspondence of the measured spectrum with the
theoretical prediction for 3-soliton states (the soliton positions were
optimized to fit the data) for the experimentally estimated parameters ($f$,
$D_{2}$,$\kappa$) and SIL detuning value ($\zeta_{\rm eff}^{0}$). We note that
while the form of the spectrum is highly dependent on the intersoliton
distances, the total comb power (and the conversion efficiency) does not. For
the 150-GHz comb we saw a slightly different picture. By increasing the number
of solitons we were able to match the total comb output power, but the comb
envelope at the sides show much smoother behaviour than the multi-soliton
state can provide (see green curve in Fig. 3b). At the same time, the comb
power is too high for a single-soliton state (see purple curve in Fig. 3d).
Such comb enhancement can be attributed to the comb-line amplification inside
the active medium of the laser or so-far unexplored effects of the multi-
frequency locking, while the non-smooth envelope near the pump – to the
dispersion distortions from the mode-crossings Herr  _et al._ (2014b).

A comparison of the semiconductor FP and DFB laser diodes highlights the
benefits of DFB in terms of predictable wavelength of locking and more
convenient matching of two combs, while the FP is much more powerful and
cheaper, and hence more promising with a view to practical applications.

The early integrated dual-microcomb source prototype presented here is still
affected by the relative thermal drift of the microresonators on the two
separate photonic chips. Various design options (Fig. 1(d)) could help to
overcome the drift, paving the way to more compact devices. A first
improvement could be the combination of the two microresonators on the same
photonic chip positioned on a common temperature-stabilized substrate with two
pumping laser diodes. This option would simplify comb matching owing to small
fabrication errors for the two microresonators at close distance, thus
enhancing the stability of the dual-comb beatnotes. The next improvement for
higher stability and smaller size is using the same microresonator for the
generation of two soliton combs propagating in the same or opposite directions
along the ring. Two laser diodes self-injection locked into the same
microresonator provide the highest mutual coherence and generated microcombs
would lead to dual-comb beatnotes with the lowest phase noise. To reduce the
number of components that need to be aligned during the integration, keeping
small size and high stability, the design option with just one laser diode and
one photonic chip with two microresonator could be considered. A laser diode
locked to one of the microresonators would provide the first microcomb and the
second resonator tuned relative to the frequency of the locked laser would
provide the second one. The discussed design options based on the microcomb
generation using self-injection locking should improve the performance of the
integrated dual-microcomb source, and will be explored in future research.

Taking into account recent advances of III–V heterogeneous integration with
silicon nitride photonic waveguides Xiang  _et al._ (2021a), deeper
integration of the proposed dual-microcomb source can lead to reaching tiny
chip-scale sizes. Also, the demonstrated possibility to use FP laser diodes
might allow using just a gain section without an additional laser cavity in
future generations of the chip-based dual-comb sources.

In conclusion, obtained results demonstrate that SIL Kerr microcombs based on
silicon photonics can successfully compete with other on-chip optical-comb
sources and outperform them owing to the unique combination of power
efficiency with mWs comb power, a wider spectrum, and low phase noise.

Methods

Silicon nitride chips characterisation: The SiN photonic chip-based
microresonators used in our experiments were fabricated by Ligentec SA,
Switzerland. The pumped microresonator resonance is measured using an external
tunable laser Toptica CTL1550, and fitted taking backscattering into account
Gorodetsky  _et al._ (2000) to obtain the intrinsic loss
$\kappa_{\mathrm{0}}$, the coupling rate $\kappa_{c}$, and the backward-wave
coupling rate $\gamma$ (mode splitting). Based on these data, we evaluated the
full resonance linewidth $\kappa=\kappa_{\mathrm{0}}+\kappa_{c}$, the pump
coupling efficiency $\eta={\kappa_{c}}/{\kappa}$, and the normalized
backscattering coefficient $\beta=\gamma/\kappa$.

For nonlinearity-threshold estimation the experimental setup with the external
laser Toptica CTL1550 and a booster NKT Photonics Koheras Boostik HP E15 was
used. Gradually increasing the pump power, we simultaneously monitor the
optical spectrum and the resonance shape using an OSA (Yokogawa AQ6370D) and
an oscilloscope (Keysight DSO-X 3024A). First, we reach the thermal
nonlinearity threshold ($P_{\mathrm{th}}^{\mathrm{Thermal}}$), where the
resonance shape becomes triangular. By further increasing the pump power, we
reach the parametric-instability threshold ($P_{\mathrm{th}}$), accompanied by
sideband generation.

Dispersion characteristics of microresonators were measured using the original
experimental setup based on the tunable laser Toptica CTL1550 and a calibrated
fibre Mach–-Zehnder interferometer (MZI) with a FSR of 102 MHz.

The results obtained for the 150-GHz and 1-THz FSR microresonators,
respectively, were as follows: $\kappa_{\mathrm{0}}/2\pi=$ 106 and 188 MHz,
$\kappa_{c}/2\pi=$ 75 and 269 MHz, $\gamma/2\pi=$ 34 and 331 MHz, $\eta=$ 0.42
and 0.59, Q-factor = $1.8\times 10^{6}$ and $1\times 10^{6}$,
$P_{\mathrm{th}}=$ 14 and 11 mW, $P_{\mathrm{th}}^{\mathrm{Thermal}}=$3.7 and
2.5 mW, $D_{1}/2\pi=$ 143.6 and 999.8 GHz, and $D_{2}/2\pi=$ 1.38 and 14.3
MHz.

Microcomb source: After chip characterisation we have successfully generated
microcombs with all available chips and compared the reproducibility of their
parameters (see Supplementary Note 1). As a result, we have matched photonic
chips for RF dual-microcomb signal generation and defined the operating points
for all laser diodes and photonic chip pairs, which were used for the assembly
process. The operating point is defined by the following parameters: laser-
diode temperature, injection current, and microheater voltage.

In addition, we have compared two types of laser diodes, DFB and FP (Seminex
Corp., USA), in terms of their suitability for being used for microcomb
generation (see Supplementary Note 2).

We have also estimated the linewidth of the microcomb components using the
heterodyne technique. The beatnote signal between the external laser and the
microcomb line is presented in the inset of Fig. 3(e, f). The Lorentzian width
of the beatnote signal (estimated using a Voigt-profile fit) is 0.5 kHz. The
Gaussian width is about 68 kHz, which can be improved by further reduction of
the technical noises.

Funding Information: This work was supported by the Russian Science Foundation
(grant 20-12-00344).

Acknowledgments: V.E.L. acknowledges personal support from the Foundation for
the Advancement of Theoretical Physics and Mathematics “BASIS.” A.S.V. is
supported by the EU H2020 research and innovation programme under the Marie
Sklodowska-Curie grant agreement No 101033663 (RaMSoM). N.Yu.D. was partially
supported by the Samsung Research Center in Moscow. The authors thank H.-S.
Lee and Y.-G. Roh from the Samsung Advanced Institute of Technologies for help
in establishing the project and its further support. Authors express special
thanks to Alexander Gorodnitskiy, graduate student of MIPT and Sofya
Agafonova, student of MIPT for assistance at an early stage of experiments.

Author contribution: A.S.V., S.N.K., N.Yu.D. conducted the experiment. N.M.K.,
V.E.L. developed a theoretical model and performed numerical simulations. All
authors analyzed the data and prepared the manuscript. I.A.B., M.V.R. and
S.V.P. initiated the collaboration and supervised the project.

Data Availability Statement: Data underlying the results presented in this
paper are not publicly available at this time but may be obtained from the
authors upon reasonable request.

Competing interests: S.N.K., A.S.V., V.E.L., M.V.R., S.V.P. and I.A.B. are
listed as co-authors in joint Samsung Electronics and Russian Quantum Center
US patent no. US10224688B2, which is related to the technology reported in
this article.

The authors declare no competing interests.

## References

  * Udem  _et al._ (2002) T. Udem, R. Holzwarth, and T. W. Hänsch, Nature 416, 233 (2002).
  * Fortier and Baumann (2019) T. Fortier and E. Baumann, Communications Physics 2, 153 (2019).
  * Diddams  _et al._ (2020) S. A. Diddams, K. Vahala, and T. Udem, Science 369, eaay3676 (2020).
  * Suh  _et al._ (2016) M.-G. Suh, Q.-F. Yang, K. Y. Yang, X. Yi, and K. J. Vahala, Science 354, 600 (2016).
  * Coddington  _et al._ (2016) I. Coddington, N. Newbury, and W. Swann, Optica 3, 414 (2016).
  * Ideguchi (2017) T. Ideguchi, Opt. Photon. News 28, 32 (2017).
  * Picqué and Hänsch (2019) N. Picqué and T. W. Hänsch, Nature Photonics 13, 146 (2019).
  * Baumann  _et al._ (2019) E. Baumann, E. V. Hoenig, E. F. Perez, G. M. Colacion, F. R. Giorgetta, K. C. Cossel, G. Ycas, D. R. Carlson, D. D. Hickstein, K. Srinivasan, _et al._ , Optics express 27, 11869 (2019).
  * Koptyaev  _et al._ (2019) S. N. Koptyaev, G. V. Lihachev, N. G. Pavlov, A. A. Shchekin, I. A. Bilenko, M. V. Riabko, M. L. Gorodetsky, S. V. Polonsky, A. S. Voloshin, A. D. Lantsov, _et al._ , “Optical dual-comb source apparatuses including optical microresonator,”  (2019), uS Patent 10,224,688.
  * Okubo  _et al._ (2014) S. Okubo, K. Iwakuni, H. Inaba, K. Hosaka, A. Onae, H. Sasada, and F.-L. Hong, in  _CLEO: 2014_ (Optical Society of America, 2014) p. STh1N.6.
  * Okubo  _et al._ (2015) S. Okubo, K. Iwakuni, H. Inaba, K. Hosaka, A. Onae, H. Sasada, and F.-L. Hong, Applied Physics Express 8, 082402 (2015).
  * Brehm  _et al._ (2006) M. Brehm, A. Schliesser, and F. Keilmann, Opt. Express 14, 11222 (2006).
  * von Ribbeck  _et al._ (2008) H.-G. von Ribbeck, M. Brehm, D. van der Weide, S. Winnerl, O. Drachenko, M. Helm, and F. Keilmann, Opt. Express 16, 3430 (2008).
  * Zolot  _et al._ (2013) A. Zolot, F. Giorgetta, E. Baumann, W. Swann, I. Coddington, and N. Newbury, Journal of Quantitative Spectroscopy and Radiative Transfer 118, 26 (2013).
  * Baumann  _et al._ (2011) E. Baumann, F. R. Giorgetta, W. C. Swann, A. M. Zolot, I. Coddington, and N. R. Newbury, Phys. Rev. A 84, 062513 (2011).
  * Rieker  _et al._ (2014) G. B. Rieker, F. R. Giorgetta, W. C. Swann, J. Kofler, A. M. Zolot, L. C. Sinclair, E. Baumann, C. Cromer, G. Petron, C. Sweeney, P. P. Tans, I. Coddington, and N. R. Newbury, Optica 1, 290 (2014).
  * Zhu  _et al._ (2015) F. Zhu, A. Bicer, R. Askar, J. Bounds, A. A. Kolomenskii, V. Kelessides, M. Amani, and H. A. Schuessler,  12, 095701 (2015).
  * Schroeder  _et al._ (2017) P. Schroeder, R. Wright, S. Coburn, B. Sodergren, K. Cossel, S. Droste, G. Truong, E. Baumann, F. Giorgetta, I. Coddington, N. Newbury, and G. Rieker, Proceedings of the Combustion Institute 36, 4565 (2017).
  * Coddington  _et al._ (2009) I. Coddington, W. C. Swann, L. Nenadovic, and N. R. Newbury, Nature Photonics 3, 351 (2009).
  * Suh and Vahala (2018) M.-G. Suh and K. J. Vahala, Science 359, 884 (2018).
  * Trocha  _et al._ (2018) P. Trocha, M. Karpov, D. Ganin, M. H. P. Pfeiffer, A. Kordts, S. Wolf, J. Krockenberger, P. Marin-Palomo, C. Weimann, S. Randel, W. Freude, T. J. Kippenberg, and C. Koos, Science 359, 887 (2018).
  * Nürnberg  _et al._ (2021) J. Nürnberg, B. Willenberg, C. R. Phillips, and U. Keller, Opt. Express 29, 24910 (2021).
  * Riemensberger  _et al._ (2020) J. Riemensberger, A. Lukashchuk, M. Karpov, W. Weng, E. Lucas, J. Liu, and T. J. Kippenberg, Nature 581, 164 (2020).
  * Coddington  _et al._ (2010) I. Coddington, W. Swann, and N. Newbury, Physical Review A 82, 043817 (2010).
  * Hoghooghi  _et al._ (2021) N. Hoghooghi, R. K. Cole, and G. B. Rieker, Applied Physics B 127, 1 (2021).
  * Kayes  _et al._ (2018) M. I. Kayes, N. Abdukerim, A. Rekik, and M. Rochette, Opt. Lett. 43, 5809 (2018).
  * Link  _et al._ (2016) S. M. Link, A. Klenner, and U. Keller, Opt. Express 24, 1889 (2016).
  * Link  _et al._ (2015) S. Link, C. Zaugg, A. Klenner, M. Mangold, M. Golling, B. Tilma, and U. Keller, in  _Vertical External Cavity Surface Emitting Lasers (VECSELs) V_, Vol. 9349 (International Society for Optics and Photonics, 2015) p. 93490Q.
  * Nürnberg  _et al._ (2019) J. Nürnberg, C. G. Alfieri, D. Waldburger, L. Krüger, M. Golling, and U. Keller, in  _2019 Conference on Lasers and Electro-Optics Europe & European Quantum Electronics Conference (CLEO/Europe-EQEC)_ (IEEE, 2019) pp. 1–1.
  * Consolino  _et al._ (2020) L. Consolino, M. Nafa, M. De Regis, F. Cappelli, K. Garrasi, F. P. Mezzapesa, L. Li, A. G. Davies, E. H. Linfield, M. S. Vitiello, S. Bartalini, and P. De Natale, Communications Physics 3, 69 (2020).
  * Komagata  _et al._ (2021) K. Komagata, A. Shehzad, G. Terrasanta, P. Brochard, R. Matthey, M. Gianella, P. Jouy, F. Kapsalidis, M. Shahmohammadi, M. Beck, _et al._ , Optics Express 29, 19126 (2021).
  * Pfeiffer  _et al._ (2016) M. H. Pfeiffer, A. Kordts, V. Brasch, M. Zervas, M. Geiselmann, J. D. Jost, and T. J. Kippenberg, Optica 3, 20 (2016).
  * Liu  _et al._ (2018) J. Liu, A. S. Raja, M. Karpov, B. Ghadiani, M. H. Pfeiffer, B. Du, N. J. Engelsen, H. Guo, M. Zervas, and T. J. Kippenberg, Optica 5, 1347 (2018).
  * Liu  _et al._ (2021) J. Liu, G. Huang, R. N. Wang, J. He, A. S. Raja, T. Liu, N. J. Engelsen, and T. J. Kippenberg, Nature communications 12, 1 (2021).
  * Hochberg and Baehr-Jones (2010) M. Hochberg and T. Baehr-Jones, Nature photonics 4, 492 (2010).
  * Hochberg  _et al._ (2013) M. Hochberg, N. C. Harris, R. Ding, Y. Zhang, A. Novack, Z. Xuan, and T. Baehr-Jones, IEEE Solid-State Circuits Magazine 5, 48 (2013).
  * Weimann  _et al._ (2017) C. Weimann, M. Lauermann, F. Hoeller, W. Freude, and C. Koos, Optics express 25, 30091 (2017).
  * Xiang  _et al._ (2021a) C. Xiang, J. Liu, J. Guo, L. Chang, R. N. Wang, W. Weng, J. Peters, W. Xie, Z. Zhang, J. Riemensberger, _et al._ , arXiv preprint arXiv:2103.02725 (2021a).
  * Xiang  _et al._ (2021b) C. Xiang, J. Guo, W. Jin, L. Wu, J. Peters, W. Xie, L. Chang, B. Shen, H. Wang, Q.-F. Yang, _et al._ , Nature Communications 12, 1 (2021b).
  * Cuyvers  _et al._ (2021) S. Cuyvers, B. Haq, C. Op de Beeck, S. Poelman, A. Hermans, Z. Wang, A. Gocalinska, E. Pelucchi, B. Corbett, G. Roelkens, _et al._ , Laser & Photonics Reviews , 2000485 (2021).
  * Dahmani  _et al._ (1987) B. Dahmani, L. Hollberg, and R. Drullinger, Opt. Lett. 12, 876 (1987).
  * Laurent  _et al._ (1989) P. Laurent, A. Clairon, and C. Breant, IEEE Journ. Quant. El. 25, 1131 (1989).
  * Oraevsky  _et al._ (2001) A. N. Oraevsky, A. V. Yarovitsky, and V. L. Velichansky,  31, 897 (2001).
  * Liang  _et al._ (2010) W. Liang, V. S. Ilchenko, A. A. Savchenkov, A. B. Matsko, D. Seidel, and L. Maleki, Opt. Lett. 35, 2822 (2010).
  * Liang  _et al._ (2015a) W. Liang, D. Eliyahu, V. S. Ilchenko, A. A. Savchenkov, A. B. Matsko, D. Seidel, and L. Maleki, Nature Communications 6, 7957 (2015a).
  * Liang  _et al._ (2015b) W. Liang, V. S. Ilchenko, D. Eliyahu, A. A. Savchenkov, A. B. Matsko, D. Seidel, and L. Maleki, Nature Communications 6, 7371 (2015b).
  * Kondratiev  _et al._ (2017) N. M. Kondratiev, V. E. Lobanov, A. V. Cherenkov, A. S. Voloshin, N. G. Pavlov, S. Koptyaev, and M. L. Gorodetsky, Opt. Express 25, 28167 (2017).
  * Galiev  _et al._ (2020) R. R. Galiev, N. M. Kondratiev, V. E. Lobanov, A. B. Matsko, and I. A. Bilenko, Phys. Rev. Applied 14, 014036 (2020).
  * Voloshin  _et al._ (2021) A. S. Voloshin, N. M. Kondratiev, G. V. Lihachev, J. Liu, V. E. Lobanov, N. Y. Dmitriev, W. Weng, T. J. Kippenberg, and I. A. Bilenko, Nature Communications 12, 235 (2021).
  * Shen  _et al._ (2020) B. Shen, L. Chang, J. Liu, H. Wang, Q.-F. Yang, C. Xiang, R. N. Wang, J. He, T. Liu, W. Xie, J. Guo, D. Kinghorn, L. Wu, Q.-X. Ji, T. J. Kippenberg, K. Vahala, and J. E. Bowers, Nature 582, 365 (2020).
  * Raja  _et al._ (2019) A. S. Raja, A. S. Voloshin, H. Guo, S. E. Agafonova, J. Liu, A. S. Gorodnitskiy, M. Karpov, N. G. Pavlov, E. Lucas, R. R. Galiev, A. E. Shitikov, J. D. Jost, M. L. Gorodetsky, and T. J. Kippenberg, Nature Communications 10, 680 (2019).
  * Stern  _et al._ (2018) B. Stern, X. Ji, Y. Okawachi, A. L. Gaeta, and M. Lipson, Nature 562, 401 (2018).
  * Kippenberg  _et al._ (2011) T. J. Kippenberg, R. Holzwarth, and S. A. Diddams, Science 332, 555 (2011).
  * Lin  _et al._ (2017) G. Lin, A. Coillet, and Y. K. Chembo, Adv. Opt. Photon. 9, 828 (2017).
  * Pasquazi  _et al._ (2018) A. Pasquazi, M. Peccianti, L. Razzari, D. J. Moss, S. Coen, M. Erkintalo, Y. K. Chembo, T. Hansson, S. Wabnitz, P. Del’Haye, X. Xue, A. M. Weiner, and R. Morandotti, Physics Reports 729, 1 (2018).
  * Gaeta  _et al._ (2019) A. L. Gaeta, M. Lipson, and T. J. Kippenberg, Nature Photonics 13, 158 (2019).
  * Kovach  _et al._ (2020) A. Kovach, D. Chen, J. He, H. Choi, A. H. Dogan, M. Ghasemkhani, H. Taheri, and A. M. Armani, Adv. Opt. Photon. 12, 135 (2020).
  * Lin  _et al._ (2020) T. Lin, A. Dutt, C. Joshi, X. Ji, C. T. Phare, Y. Okawachi, A. L. Gaeta, and M. Lipson, arXiv preprint arXiv:2001.00869 (2020).
  * Bao  _et al._ (2019) C. Bao, M.-G. Suh, and K. Vahala, Optica 6, 1110 (2019).
  * Dutt  _et al._ (2018) A. Dutt, C. Joshi, X. Ji, J. Cardenas, Y. Okawachi, K. Luke, A. L. Gaeta, and M. Lipson, Science Advances 4, e1701858 (2018).
  * Wang  _et al._ (2018) W. Wang, W. Zhang, Z. Lu, S. T. Chu, B. E. Little, Q. Yang, L. Wang, and W. Zhao, Photon. Res. 6, 363 (2018).
  * Yu  _et al._ (2018) M. Yu, Y. Okawachi, A. G. Griffith, N. Picqué, M. Lipson, and A. L. Gaeta, Nature communications 9, 1 (2018).
  * Pavlov  _et al._ (2017) N. Pavlov, G. Lihachev, S. Koptyaev, E. Lucas, M. Karpov, N. Kondratiev, I. Bilenko, T. Kippenberg, and M. Gorodetsky, Optics letters 42, 514 (2017).
  * Jin  _et al._ (2021) W. Jin, Q.-F. Yang, L. Chang, B. Shen, H. Wang, M. A. Leal, L. Wu, M. Gao, A. Feshali, M. Paniccia, K. J. Vahala, and J. E. Bowers, Nature Photonics 15, 346 (2021).
  * Herr  _et al._ (2014a) T. Herr, V. Brasch, J. D. Jost, C. Y. Wang, N. M. Kondratiev, M. L. Gorodetsky, and T. J. Kippenberg, Nat. Photon. 8, 145 (2014a).
  * Herr  _et al._ (2014b) T. Herr, V. Brasch, J. D. Jost, I. Mirgorodskiy, G. Lihachev, M. L. Gorodetsky, and T. J. Kippenberg, Phys. Rev. Lett. 113, 123901 (2014b).
  * Kim  _et al._ (2021) C. Kim, K. Yvind, and M. Pu, Opt. Lett. 46, 3508 (2021).
  * Li  _et al._ (2017) Q. Li, T. C. Briles, D. A. Westly, T. E. Drake, J. R. Stone, B. R. Ilic, S. A. Diddams, S. B. Papp, and K. Srinivasan,  4, 193 (2017).
  * Joshi  _et al._ (2016) C. Joshi, J. K. Jang, K. Luke, X. Ji, S. A. Miller, A. Klenner, Y. Okawachi, M. Lipson, and A. L. Gaeta, Opt. Lett. 41, 2565 (2016).
  * Kuse  _et al._ (2020) N. Kuse, T. Tetsumoto, G. Navickaite, M. Geiselmann, and M. E. Fermann, Opt. Lett. 45, 927 (2020).
  * Kondratiev  _et al._ (2021) N. M. Kondratiev, R. R. Galiev, and V. E. Lobanov, in  _Nonlinear Optics and Applications XII_, Vol. 11770, edited by M. Bertolotti, A. V. Zayats, and A. M. Zheltikov, International Society for Optics and Photonics (SPIE, 2021) pp. 53 – 58\.
  * Xue  _et al._ (2017) X. Xue, P.-H. Wang, Y. Xuan, M. Qi, and A. M. Weiner, Laser & Photonics Reviews 11, 1600276 (2017).
  * Kim  _et al._ (2019) B. Y. Kim, Y. Okawachi, J. K. Jang, M. Yu, X. Ji, Y. Zhao, C. Joshi, M. Lipson, and A. L. Gaeta, Opt. Lett. 44, 4475 (2019).
  * Jang  _et al._ (2021) J. K. Jang, Y. Okawachi, Y. Zhao, X. Ji, C. Joshi, M. Lipson, and A. L. Gaeta, Opt. Lett. 46, 3657 (2021).
  * Karpov  _et al._ (2019) M. Karpov, M. H. P. Pfeiffer, H. Guo, W. Weng, J. Liu, and T. J. Kippenberg, Nature Physics 15, 1071 (2019).
  * Gorodetsky  _et al._ (2000) M. L. Gorodetsky, A. D. Pryamikov, and V. S. Ilchenko, J. Opt. Soc. Am. B 17, 1051 (2000).


# On the energy of the protons producing
the very high-energy astrophysical neutrinos

Esteban Roulet(a) and Francesco Vissani(b,c)
(a)Centro Atómico Bariloche, Comisión Nacional de Energía Atómica
Consejo Nacional de Investigaciones Científicas y Técnicas (CONICET)
Av. Bustillo 9500, R8402AGP, Bariloche, Argentina
(b) INFN, Laboratori Nazionali del Gran Sasso, 67100 Assergi, L’Aquila, Italy
(c) Gran Sasso Science Institute (GSSI), Viale F. Crispi 7, 67100 L’Aquila,
Italy

###### Abstract

We study the distribution of the energy of the protons that can produce an
astrophysical neutrino with a given observed energy, in the TeV–PeV range,
both through $pp$ or $p\gamma$ interactions. Due to the increasing
multiplicity of the pion production at high center of mass energies, the
resulting average proton energies can be much larger than the often used
approximate value $E_{p}\simeq 20E_{\nu}$. Also the threshold of the $p\gamma$
process can lead to a pronounced increase in the values of $\langle
E_{p}\rangle/E_{\nu}$ for decreasing neutrino energies. The results depend
sensitively on the assumed proton spectrum, since steeper spectra give less
weight to the lower energy neutrino tail resulting from the decays of the
abundant low-energy pions. In the $p\gamma$ scenarios they also depend
sensitively on the spectrum of the target photons. The results are in
particular relevant to relate possible characteristics of the neutrino
spectrum to those of the corresponding cosmic rays that generated them. We
also discuss the associated production of gamma rays at the sources.

##  1 Introduction

Very high energy ($>10$ TeV) neutrinos of non-atmospheric origin have been
observed in the last decade with the IceCube detector [1, 2, 3, 4, 5, 6]. They
are believed to be produced in the astrophysical sources where the ultra-high
energy cosmic rays (CRs) are accelerated, through the interactions of these
with the gas or the radiation present in the acceleration sites or their
surroundings. These interactions can produce unstable particles, mostly pions,
which through their weak decays give rise to neutrinos, with the decays of the
neutral pions giving rise to comparable amounts of gamma rays (for recent
reviews see [7, 8]). Also the interactions of the CRs with the background
radiation that they traverse as they propagate through the Universe, such as
the cosmic microwave background (CMB) and the extragalactic background light
(EBL), can give rise to high-energy cosmogenic neutrinos [9]. Cosmogenic
neutrinos are however likely subdominant with respect to the astrophysical
ones produced at the sources, except eventually at energies well beyond 10 PeV
[10]. In spite of the major progress being done, the sources of the observed
neutrinos are however still largely uncertain, with only one or two potential
identifications having been done up to now through the coincidence between the
arrival directions of well determined muon neutrino induced tracks and the
location of the blazar TXS-0506+056 [11] or the Seyfert active galaxy NGC 1068
[3]. The lack of correlations with source candidates of the remaining several
tens of astrophysical neutrinos is also a relevant piece of the puzzle.

Besides the searches for angular correlations, also the detailed study of the
energy distribution of the neutrinos should give important information about
the mechanism producing the neutrinos, and hence on the nature of the CR
sources. However, analyses using different neutrino datasets (muon tracks [6],
cascade events [4] or high energy starting events [5]) don’t lead to a
definite picture about the astrophysical neutrino spectral shape. In
particular, the slope $\gamma$ of a power-law fit $\phi_{\nu}\propto
E_{\nu}^{-\gamma}$ which is obtained from the different studies actually
covering the range of values between 2 and 3, so that more data would be
required in order to clarify this issue.

The two main classes of production mechanisms for the astrophysical neutrinos
are the $pp$ and the $p\gamma$ interactions. The $pp$ mechanism takes place
when a CR proton, or more in general a nucleon in a CR nucleus, interacts with
the gas in the acceleration region (such as when a jet punches through the
envelope of a star) or in the surrounding medium (as for instance in cluster
of galaxies or starburst galaxies), while the $p\gamma$ one when it interacts
with dense radiation fields, such as the synchrotron radiation from co-
accelerated electrons in the jets of gamma ray bursts (GRB) or active galactic
nuclei (AGN) or the strong UV emission from the supermassive black-hole
accretion disks in AGNs, etc.

Due to the approximate scaling invariance of the hadronic interactions, which
implies for instance that the pion production yield of the $pp$ interactions
mostly depends on the fraction of the proton energy carried by the pion,
$x_{\pi}\equiv E_{\pi}/E_{p}$, rather than on the two energies separately,
together with the small (logarithmic) energy dependence of the $pp$ cross
section and the low threshold for pion production, one expects that in the
$pp$ mechanism the spectrum of the neutrinos from the pion decays will have a
spectral shape following closely that of the parent CR protons. In particular,
a CR power-law differential spectrum $\phi_{\rm CR}(E_{p})\sim
E_{p}^{-\alpha}$ should approximately lead to $\phi_{\nu}(E)\sim E^{-\gamma}$,
with $\gamma\simeq\alpha$. On the other hand, for the $p\gamma$ mechanism the
proton energy threshold for the pion production to be possible is high,
requiring that $E_{p}>m_{p}m_{\pi}(1+m_{\pi}/2m_{p})/(2\varepsilon)\simeq
70\,{\rm TeV}/(\varepsilon/{\rm keV})$, with $\varepsilon$ being the target
photon energy. In addition, at least near the threshold where the $\Delta$
resonant production dominates, the cross section has a significant energy
dependence. As a result, the resulting neutrino spectrum depends not only on
the spectral shape of the CRs, but also on the background photon spectrum.

In order to relate the energy of the neutrinos with that of the parent
protons, one usually considers that the typical pion energy is $E_{\pi}\simeq
E_{p}/5$ and that the four particles resulting from the decays
$\pi^{+}\to\mu^{+}\nu_{\mu}$ and $\mu^{+}\to e^{+}\nu_{e}\bar{\nu}_{\mu}$ (and
similarly for $\pi^{-}$) share similar amounts of energy, so that
$E_{\nu}\simeq E_{\pi}/4\simeq E_{p}/20$. For the neutral pions one has that
$\pi^{0}\to\gamma\gamma$, so that one similarly finds the naive estimate
$E_{\gamma}\simeq E_{p}/10$. These estimates are however very simplistic,
since the pions have actually a distribution of energies and similarly the
energies of the neutrinos from muon decays are continuously distributed, and
hence the neutrinos resulting from the interaction of a proton with a given
energy are expected to have a quite broad distribution of energies. Moreover,
the inelasticity of the interaction and the multiplicity of the produced
particles depend on energy, with the multipion production becoming relevant at
high energies and giving rise to large numbers of lower energy pions. The
impact of the neutrinos carrying different fractions of the proton energy in
the final neutrino spectrum will also depend on the shape of the CR spectrum,
since steeper (i.e. softer) spectra will make less relevant the contribution
from the low energy tail of the neutrinos from pion decays, enhancing the
relative contribution of those with large values of $x_{\nu}=E_{\nu}/E_{p}$.
On the other hand, flatter (i.e. harder) CR spectra will lead to neutrino
fluxes that depend more sensitively on the low energy tails and on the
multipion production processes. The relevance of multipion production channels
for the astrophysical neutrino production was emphasized e.g. in [12] for the
$pp$ scenarios, and in [13, 14] for the $p\gamma$ scenarios.

In this work we will analyze the distribution of the energies of the protons
that could have produced neutrinos with a given observed energy (in the
TeV–PeV range). This is in particular useful to understand what are the CR
energies being probed with the observed astrophysical neutrinos. We will also
discuss the resulting shape of the spectra of the neutrinos for the different
scenarios. We will base the study on the analytical fits to the neutrino
yields that parametrize the results of numerical simulations obtained in [12],
based on the SIBYLL 2.1 model [15] to account for the $pp$ hadronic
interactions, and in [16], relying on the SOPHIA code [17] to account for the
$p\gamma$ interactions. Further studies of the photon spectra in $pp$
interactions were performed in [18], and the neutrino and photon spectra
resulting in some post-LHC models were considered in [19], obtaining results
which are essentially consistent with those we will adopt from [12, 16].

##  2 The $pp$ scenarios

Let us consider the interactions of a flux of protons with energy $E_{p}$ with
a target of hydrogen gas with density $n_{\rm H}$ (heavier target nuclei can
approximately be accounted for by associating to them an effective H density,
while for the generalization to the case of a flux of heavy nuclei see e.g.
[18, 20]). The neutrinos produced can be computed as

$\frac{{\rm d}N_{\nu_{i}}}{{\rm d}E_{\nu}}=n_{\rm
H}D\int_{E_{\nu}}^{\infty}\frac{{\rm
d}E_{p}}{E_{p}}\sigma_{pp}F_{\nu_{i}}(x_{\nu},E_{p})\frac{{\rm d}N_{p}}{{\rm
d}E_{p}},$ (1)

where $\sigma_{pp}$ is the inelastic $pp$ cross section, with an approximate
expression for its energy dependence being
$\sigma_{pp}\simeq(34.9+1.985L+0.18L^{2})$ mb, with $L={\rm ln}(E_{p}/{\rm
TeV})$ [18]. The quantity $D$ is the distance traversed by the protons through
the target gas, so that $N_{\rm H}\equiv n_{\rm H}D$ is the associated column
density, and we have assumed that the source medium is thin, i.e. that
$\sigma_{pp}N_{\rm H}\ll 1$. Given that $\sigma_{pp}\sim 10^{-25}$ cm2, the
source will be thin as long as $N_{\rm H}\ll 10^{25}$ cm-2. The function
$F_{\nu_{i}}$, depending on $x_{\nu}\equiv E_{\nu}/E_{p}$ and on $E_{p}$,
describes the neutrino yield in the $pp$ interaction, and is in principle
different for the different neutrino and antineutrino flavors. At low energies
there is a noticeable excess of positive pions being produced in $pp$
interactions, leading in particular to a net excess of $\nu_{e}$ over
$\bar{\nu}_{e}$. However, when considering very high energies it is a good
approximation to assume that the $\pi^{+}$ and $\pi^{-}$ multiplicities are
similar. Hence, for the energies that we will explore here the fluxes of
neutrinos and antineutrinos are expected to be similar for the $pp$ scenarios,
so that we will consider for each flavor the sum of neutrinos and
antineutrinos that were obtained in [12], and just refer to them as neutrinos.
Note also that neutrino telescopes, such as IceCube, cannot distinguish
neutrinos from antineutrinos, with the rare exception of the $\bar{\nu}_{e}$
induced $W$-resonance events at 6.3 PeV neutrino energies resulting from the
interaction with atomic electrons ($\bar{\nu}_{e}e\to W\to anything$). The
index $i$ in Eq. (1) identifies the type of neutrino produced, with
$\nu_{\mu}^{(1)}$ being those from the pion decay, while $\nu_{\mu}^{(2)}$ and
$\nu_{e}$ those from the subsequent muon decay. The spectrum of $\nu_{e}$
turns out to be, to an accuracy better than 5%, equal to that of
$\nu_{\mu}^{(2)}$, and hence the two fluxes are just assumed to be the same.
To compute these functions one has to account for the fact that the muon
produced in the pion decay is polarized [21], since this affects the energy
distribution of its decay products. These functions have been parameterized in
[12] using simulations based on the SIBYLL 2.1 code to account for the details
of the $pp$ hadronic interactions (finding comparable results also for
simulations based on the QGSJet hadronic interaction model). There is also a
small contribution to the neutrino fluxes from charged kaon decays (which is
below the 10% level, see e.g. [22]), but following [12] we do not include it
since this should not affect significantly our results.

Analogous expressions also hold for the flux of secondary $\gamma$ rays, which
arises mostly from $\pi^{0}$ decays, with a smaller contribution from the
isosinglet $\eta$ decays. The different functions $F_{i}$ are plotted in Fig.
1 for the illustrative case with $E_{p}=10$ PeV.

Figure 1: Yields of neutrinos and photons, $x^{2}F_{i}(x,E_{p})$ vs. $x$, for
an incident proton energy $E_{p}=10$ PeV. The spectrum of $\nu_{e}$ is very
similar to that of the $\nu_{\mu}^{(2)}$ produced in the $\mu$ decays.

The probability that a neutrino of type $i$ (or a photon) and energy $E_{i}$
originated from protons in a given logarithmic energy interval is then

$\frac{{\rm d}P_{i}}{{\rm d\,ln}E_{p}}(E_{i},E_{p})={\cal
N}\sigma_{pp}F_{i}(x_{i},E_{p})\frac{{\rm d}N_{p}}{{\rm d}E_{p}},$ (2)

where the normalization is obtained from ${\cal
N}^{-1}=\int_{E_{\nu}}^{\infty}\sigma_{pp}F_{i}[{\rm d}N_{p}/{\rm
d}E_{p}]\,{\rm d}E_{p}/E_{p}$. The resulting distributions are shown in Fig.
2, adopting a power-law CR spectrum with spectral index $\alpha=2$ and no
cutoff. One can see that the distributions peak at different values for the
different types of secondaries, they are quite broad and have significant
high-energy tails.

Figure 2: Probability distribution for the parent proton energies that produce
neutrinos or photons with different energies in $pp$ interactions, adopting
$\alpha=2$.

A useful quantity is the median energy of the protons, $\bar{E}_{p}$, leading
to a given neutrino energy, which is defined such that the integral of the
probability up to $\bar{E}_{p}$ gives 0.5. In the left panel of Fig. 3 we plot
the value of $\bar{E}_{p}/E_{i}$ for different representative values of the
neutrino energies and for different values of the spectral index $\alpha$ of
the assumed power-law proton distribution. We see that the ratio
$\bar{E}_{p}/E_{i}$ depends on the kind of neutrino considered and on the
value of $\alpha$. Its value increases with increasing neutrino energy, due to
the larger contribution from multipion production channels at higher energies.
The values obtained differ from the usual assumption that $E_{p}/E_{\nu}\simeq
20$, specially when the proton spectrum is hard since in this case the
contribution from neutrinos with small values of $x_{\nu}$ becomes more
relevant. In particular, we see that for $\alpha\simeq 2$ a neutrino of 1 PeV
will originate from protons with $\bar{E}_{p}\simeq 50$ PeV, while for
$\alpha\simeq 2.5$ it would originate from protons with $\bar{E}_{p}\simeq 20$
PeV. On the other hand, due to the fact that the distribution of parent proton
energies extends up to very large values, see Fig. 2, the average proton
energy $\langle E_{p}\rangle$ leading to a given value of $E_{\nu}$ turns out
to be much higher than the median energy $\bar{E}_{p}$, and can be few hundred
times larger than the neutrino energy if the spectrum is hard (and assuming no
cutoff), as is shown in the right panel of Fig. 3.

If the sources were at cosmological distances, the observed neutrino energies
would appear redshifted by a factor $(1+z)^{-1}$ with respect to the energy
produced at the source111If one were to compute the cumulative flux from all
such sources, it would also be necessary to integrate over the contribution
from sources at all redshifts accounting for the cosmological evolution of
their density., but for simplicity we ignore here those effects, which anyway
do not affect the ratios $\bar{E}_{p}/E_{i}$ significantly, just shifting the
energy in the horizontal axis in Fig. 3 by that factor.

Figure 3: Left panel: Ratio $\bar{E}_{p}/E_{i}$ involving the median proton
energy $\bar{E}_{p}$ giving rise to a neutrino (or gamma) energy $E_{i}$, for
different CR spectral indices $\alpha$, as a function of $E_{i}$. Right panel:
similar plot but considering the average proton energy.

The energy dependence of $\bar{E}_{p}/E_{i}$, as well as the increase of the
cross section and pion multiplicity with energy, also implies that the
resulting neutrino spectral shape will slightly differ from that of the
original proton spectrum. While the first effect would tend to make it
steeper, the other two would tend to make it harder. Detailed computation
shows that the neutrino spectrum tends indeed to be slightly harder than the
proton spectrum, due to the increase in the average number of neutrinos
produced in $pp$ interactions with increasing energies (see e.g. [23, 18,
19]). Another potential implication of these results is that they affect the
relation between the cutoff energy of the CR spectrum and the associated
cutoff in the neutrino spectrum. As an example of these effects, we show in
Fig. 4 the total differential neutrino flux summed over all flavors (and
multiplied by $E_{\nu}^{2}$), as well as the individual contributions from
$\nu_{\mu}^{(1)}$ and $\nu_{\mu}^{(2)}$, resulting from an $E^{-2}$ CR proton
spectrum with a sharp cutoff at $E_{c}=100$ PeV (left panel) and one with a
smooth cutoff of the form d$N_{p}/{\rm d}E\propto
E^{-2}/\cosh((E/E_{c})^{\delta})$, with $\delta=1$ (right panel).222We use the
hyperbolic cosine suppression rather than an exponential one because with it
the spectrum below $E_{c}$ is closer to the power-law $E^{-\alpha}$ and the
determination of the hardening becomes clearer. One can see that the neutrino
spectrum is harder than the original proton spectrum and it has a very smooth
suppression starting at about 1 PeV. A good fit to the total neutrino spectrum
can be obtained with the expression d$N_{\nu}/{\rm d}E\propto
E^{-\gamma}/\cosh((E/E_{\star})^{\delta^{\prime}})$, as shown in the plot. We
obtain $\gamma\simeq 1.91$, $E_{\star}\simeq E_{c}/107$ and
$\delta^{\prime}\simeq 0.45$ for the sharp cutoff case (left panel), while
$\gamma\simeq 1.90$, $E_{\star}\simeq E_{c}/107$ and $\delta^{\prime}\simeq
0.39$ for the softer cutoff case with $\delta=1$ (right panel). Note that due
to the wide spread of the energies of the neutrinos produced in $pp$
interactions, one indeed expects that the falloff in the neutrino spectrum
should be significantly less pronounced than that of the original CR spectrum.
It is worth noting that these results cannot be obtained if one computes the
neutrino production in the so-called delta function approximation, i.e.
setting the pion and neutrino energies in a $pp$ interaction to their average
values (see e.g. [24]).

Figure 4: Total $\nu$ spectrum (from all flavors), as well as the individual
contributions from $\nu_{\mu}^{(1)}$ and $\nu_{\mu}^{(2)}$, from $pp$
scenarios with a power-law proton spectrum with $\alpha=2$. The left panel
corresponds to a sharp cutoff at $E_{c}=100$ PeV while the right panel is for
a smoother cutoff scaling as $\cosh^{-1}(E/E_{c})$, whose shapes are shown
with black dashes, arbitrarily normalized to unity at 1 TeV. Also shown is the
fitting function to the total neutrino spectrum discussed in the text.

Let us mention that the spectra of the different neutrinos at the source could
be affected by the presence of strong (larger than 1 kG) magnetic fields in
the source interaction region, since this could lead to significant
synchrotron losses of the high-energy muons before they decay, corresponding
to the so-called muon damped neutrino scenarios. In this hypothetical case,
above a certain critical energy the neutrino fluxes at the source will be
dominated by the $\nu_{\mu}^{(1)}$ component, with the $\nu_{\mu}^{(2)}$ and
$\nu_{e}$ neutrinos from the $\mu$ decay being strongly suppressed. This
critical energy can be estimated equating the muon decay time, $t_{\rm
dec}^{\mu}=(E_{\mu}/m_{\mu})\tau_{\mu}\simeq 20(E_{\mu}/{\rm PeV})$ s, with
the muon synchrotron loss time, $t_{\rm
sync}^{\mu}=(9/4)m_{\mu}^{4}/(e^{4}B^{2}E_{\mu})\simeq 3\times 10^{7}\,{\rm
s}/[(B/{\rm G})^{2}(E_{\mu}/{\rm PeV})]$, and is $E_{\mu}\simeq{\rm
PeV}/(B/{\rm kG})$ (the associated neutrino energy would be a factor of about
3 smaller). At energies about a factor 20 larger, the pion themselves will
lose their energy by synchrotron emission before they can decay, and hence a
strong suppression in the total neutrino flux will result above the
corresponding neutrino energies, with the spectrum becoming steeper by two
extra powers of the energy, corresponding to the fraction of the mesons that
decay before loosing significant amounts of energy, until eventually the
neutrinos from kaon decays become the dominant ones [25, 26, 27]. The non-
observation of a strong suppression in the astrophysical neutrino fluxes at or
below PeV energies hence suggests that in the $pp$ scenarios the magnetic
fields in the interaction region are not much larger than a kG. One should
also keep in mind that the superposition of sources with different magnetic
fields, or even the emission of neutrinos from different regions of a source
having different magnetic field values, could soften the suppression effect.

The neutrino fluxes at the Earth are obtained after accounting for the
averaged neutrino flavor oscillations, and would then have a different flavor
ratio below and above the threshold associated to the muon damping. Below the
muon damping threshold, the original source flavor ratio
$(\nu_{e}:\nu_{\mu}:\nu_{\tau})\simeq(1:2:0)$ will be converted by
oscillations to an approximately equal flavor ratio $(1:1:1)$, while above the
muon damping threshold (but below the pion damping regime), the source flavor
ratio $(0:1:0)$ will be converted approximately to the $(0.2:0.4:0.4)$ flavor
ratio. Note also that when the muons are damped, the ratio
$\bar{E}_{p}/E_{\nu}$ will be that associated to $\nu_{\mu}^{(1)}$ rather than
the value averaged over the different neutrino species. To describe in detail
the neutrino spectra above the onset of the muon damping effects, one would
need to consider the yield of pions from the proton interactions, obtain the
yield of muons from pion decays and then account for their energy losses
before considering the neutrino production in their decays (see e.g. [28,
29]).

For completeness, we also showed in Fig. 3 the values of
$\bar{E}_{p}/E_{\gamma}$ and $\langle{E}_{p}\rangle/E_{\gamma}$ obtained for
different gamma ray energies, which are found to differ from the usually
assumed value of 10 especially for hard proton spectra. One has to keep in
mind that the high-energy photons can cascade down to lower energies either by
interactions with gas or with the photons at the source (synchrotron or
inverse Compton photons) or with the background of optical/infrared and
microwave photons that they traverse in their trip up to us. Hence, only for
nearby sources, such as Galactic ones, that are thin to escaping photons the
very-high energy gamma rays may reach the Earth unattenuated, while the weakly
interacting neutrinos will always arrive unscathed.

##  3 The $p\gamma$ scenarios

In the $p\gamma$ scenarios the CR protons can interact with the dense
radiation fields present in the sources, such as for instance the synchrotron
photons from the co-accelerated electrons in the AGN or GRB jets or the
thermal optical/UV photons of the blue bump produced in the inner torus near
the AGN cores. We will compute the production of secondaries in the frame in
which the photon background is isotropic, denoting the quantities in this
frame with a prime. In the case of production in relativistic jets, such as in
GRBs or blazars, this will be the shock rest frame (SRF), and hence the energy
in the comoving frame will typically appear boosted with respect to the energy
$E^{\prime}$ in the SRF by the Lorentz factor $\Gamma$ of the jet (and then
eventually redshifted in the observer’s frame, so that $E\simeq\Gamma
E^{\prime}/(1+z)$).

The neutrino (and $\gamma$) fluxes from the decays of mesons produced in
$p\gamma$ interactions can be computed as [16]

$\frac{{\rm d}N_{i}}{{\rm d}E^{\prime}_{i}}(E^{\prime}_{i})=\int\frac{{\rm
d}E^{\prime}_{p}}{E^{\prime}_{p}}\,{\rm
d}\varepsilon^{\prime}\,f_{p}(E^{\prime}_{p})f_{\rm
ph}(\varepsilon^{\prime})\Phi_{i}(\eta,x_{i}),$ (3)

where $\varepsilon^{\prime}$ is the energy of the isotropic target photons,
with $f_{\rm ph}(\varepsilon^{\prime}){\rm d}\varepsilon^{\prime}$ and
$f_{p}(E^{\prime}_{p}){\rm d}E^{\prime}_{p}$ being the number densities of
background photons and CR protons in the energy intervals
d$\varepsilon^{\prime}$ and d$E^{\prime}_{p}$ respectively. The variable
$\eta\equiv 4\varepsilon^{\prime}E^{\prime}_{p}/m_{p}^{2}$ characterizes the
center of mass total energy squared of the interaction and
$x_{i}=E^{\prime}_{i}/E^{\prime}_{p}$ is the fraction of the proton energy
carried by the secondary, with $i$ describing the (anti) neutrino type or the
photons. The threshold condition to produce a pion is that
$(p_{\gamma}+p_{p})^{2}>(m_{p}+m_{\pi})^{2}$, and one then has that

$\eta>\frac{2m_{\pi}}{m_{p}}+\frac{m_{\pi}^{2}}{m_{p}^{2}}\equiv\eta_{0}\simeq
0.313,$ (4)

where we considered for simplicity $m_{\pi^{+}}\simeq m_{\pi^{0}}\simeq 137$
MeV. Hence, for a given photon energy $\varepsilon^{\prime}$ this requires
that $E^{\prime}_{p}>70\,{\rm TeV}/(\varepsilon^{\prime}/{\rm keV})$. The
different functions $\Phi_{i}(\eta,x_{i})$ where obtained in [16] using SOPHIA
based [17] simulations and were conveniently parameterized there, with the
parameters of the fits being tabulated for different values of
$\eta/\eta_{0}$.333We used a linear interpolation between the tabulated
values, and for $\eta/\eta_{0}>100$ we used for it the values of
$\Phi_{i}(100\eta_{0},x)$ in order to be able to approximately extend the
computation up to the highest energies that we considered. These functions are
plotted in Fig. 5 for the different neutrino and antineutrino flavors as well
as for the photons (although additional photons may result for instance from
cascading of secondary electrons).

Figure 5: Distribution functions of the secondaries resulting from $p\gamma$
interactions for different values of $\rho\equiv\eta/\eta_{0}$. The left panel
is for muon neutrinos (blue) and antineutrinos (red) and the right panel for
electron neutrinos (blue), electron antineutrinos (red) and photons (black).

Near the threshold of the $p\gamma$ interaction, the pion production is
dominated by the $\Delta(1232)$ resonance (leading to a ratio
$\pi^{0}/\pi^{+}\simeq 2$ because of isospin, or equivalently because
$\pi^{0}$ can result from the production of a $u\bar{u}$ or a $d\bar{d}$ pair
while $\pi^{+}$ only from $d\bar{d}$ pair production), with also some
contribution to $\pi^{+}$ production from direct $t$-channel $\pi^{+}$
exchange (since a charged particle is required to couple to the photon), and
then other resonances (such as $\Delta^{\prime}$, $\rho$, etc.) contribute
above the $\Delta$ resonance. For $\eta>2.14\eta_{0}$ the multipion production
becomes possible and it actually becomes the dominant channel for
$\sqrt{s}\gg{\rm GeV}$, leading to comparable amounts of $\pi^{+}$, $\pi^{0}$
and $\pi^{-}$. Other mesons get also produced, but in smaller amounts, with
the neutrino production from $K^{\pm}$ decays being only potentially relevant
at very high energies in scenarios in which the magnetic fields at the source
are so strong that the charged pions get damped before they can decay. Note
that for beamed emission from relativistic jets, the condition for the muons
to get damped before decaying translates into $E_{\mu}\simeq(\Gamma/100)\,{\rm
PeV}/(B^{\prime}/100\,{\rm kG})$, with $B^{\prime}$ being the magnetic field
in the SRF (where $\Gamma=O(10)$ for AGN jets, while it is
$O(10^{2}$–$10^{3}$) for GRB jets). Before the onset of multipion production,
no $\pi^{-}$ get produced in $p\gamma$ interactions, and hence fluxes of
${\bar{\nu}}_{e}$ appear at the source only above the threshold for two pion
production, i.e. for $\eta>2.14\eta_{0}$ (we neglect the $\bar{\nu}_{e}$ from
the decays of the produced neutrons, since their energies are very small). One
should keep in mind however that a flux of $\bar{\nu}_{e}$ will anyhow result
at the Earth from the oscillations of the $\bar{\nu}_{\mu}$ produced at the
source.

We will consider in the following two specific $p\gamma$ scenarios to study
their characteristic features. One is a typical GRB model, in which the target
photon spectrum consists of a broken power-law with a break at
$\varepsilon^{\prime}_{\rm b}=1$ keV (corresponding to the observed GRB
spectral breaks at photon energies of few hundred keV), with $f_{\rm
ph}=N_{\gamma}(\varepsilon^{\prime}/\varepsilon^{\prime}_{\rm b})^{-\beta}$,
where $\beta=1$ for $\varepsilon^{\prime}<\varepsilon^{\prime}_{\rm b}$ and
$\beta=2$ for $\varepsilon^{\prime}>\varepsilon^{\prime}_{\rm b}$, and
considering photon energies in the range $2\,{\rm
eV}<\varepsilon^{\prime}<300$ keV. The other is a thermal photon spectrum with
$kT=10$ eV, typical of the AGN’s blue bump.444One has to keep in mind however
that in this model the photon source would be external to the jet where CRs
get accelerated, and hence in reality would not be exactly isotropic in the
interaction region.

Figure 6: Probability distribution for the parent proton energies that produce
neutrinos or photons with different energies in GRB scenarios, adopting
$\alpha=2$.

Considering first the GRB scenario, we show in Fig. 6 the probability
distribution of the proton energies giving rise to different neutrino energies
for the different flavors. The left panel is for $\nu_{\mu}$ and
$\bar{\nu}_{\mu}$, while the right panel is for $\nu_{e}$, $\bar{\nu}_{e}$ and
$\gamma$. In the left-panel of Fig. 7 we show the median proton energies
$\bar{E}_{p}$, normalized to the corresponding neutrino (or photon) energies
$E_{i}$, as a function of the neutrino (or $\gamma$) energies. The main
distinguishing feature observed is the rise in the ratio $\bar{E}_{p}/E_{i}$
below $E^{\prime}_{i}\simeq 100$ TeV, due to the fact that in this case the
threshold condition for $\Delta$ production requires median proton energies in
excess of $20E_{i}$ (see also the shift in the distributions in Fig. 6 for
decreasing $E_{i}$). The results for the electron flavor are not very
different from those shown for $\bar{\nu}_{\mu}$. The right panel is similar
but in terms of the average energy of the protons $\langle
E^{\prime}_{p}\rangle$ that give rise to a given secondary energy
$E^{\prime}_{i}$, in which case the effect is much more pronounced.

Figure 7: Left panel: ratio $\bar{E}^{\prime}_{p}/E^{\prime}_{i}$ as a
function of $E^{\prime}_{i}$ in $p\gamma$ GRB scenarios, adopting $\alpha=2$.
Right panel: the same in terms of the average proton energy $\langle
E_{p}^{\prime}\rangle$.

Figure 8: Total all flavor $\nu$ spectrum in GRB scenarios (solid line),
adopting a power-law proton spectrum with $\alpha=2$ having a sharp cutoff at
100 PeV (left panel) or a smooth one with $\delta=1$ (right panel). Also shown
are the different flavor contributions (at the source) arising from the
interactions with the photons with
$\varepsilon^{\prime}<\varepsilon^{\prime}_{\rm b}$ and with
$\varepsilon^{\prime}>\varepsilon^{\prime}_{\rm b}$. Black dashed lines
illustrate the shape of the proton spectrum adopted (with arbitrary
normalization).

We then show in Fig. 8 the resulting neutrino spectrum in the SRF, for an
assumed power-law proton spectrum with $\alpha=2$ and a sharp cutoff at
$E^{\prime}_{p}=100$ PeV (left panel) or a smooth one with the hyperbolic
cosine suppression with $\delta=1$ (right panel). In order to better
understand the origin of the neutrinos with different energies, we also show
the separate contributions to the flux arising from the interactions with
photons with $\varepsilon^{\prime}<\varepsilon^{\prime}_{\rm b}$ and
$\varepsilon^{\prime}>\varepsilon^{\prime}_{\rm b}$, for the different
flavors555Note that SOPHIA does not account for the polarization of the
decaying muons, but the associated effects are small in comparison with the
uncertainties on the assumed photon spectra, the CR spectral cutoff shape,
etc.. Most of the high-energy neutrinos in the flatter region above 10 TeV
turn out to originate from interactions with the photons with the harder
spectrum present below $\varepsilon^{\prime}_{\rm b}$, and actually the bulk
of the contribution arises from photons with energies within an order of
magnitude of $\varepsilon^{\prime}_{\rm b}$. The spectral cutoff in the
neutrino spectrum is found to be much softer than the cutoff adopted for the
CR spectrum. We also note that often the non-thermal radiation in AGN jets is
modelled as a photon background scaling as $\varepsilon^{\prime-1}$ which
extends up to about 1 keV, and hence in this case the resulting neutrino
spectrum would be similar to the GRB component computed for
$\varepsilon^{\prime}<\varepsilon^{\prime}_{\rm b}$ that is shown in the
figure.

Considering now the case of a thermal radiation background with $kT=10$ eV,
characteristic of the AGN blue bump, we show in Fig. 9 the distribution of
probabilities for the different neutrino flavors and for the gamma rays. In
Fig. 10 we show the corresponding median and average proton energies giving
rise to secondaries with energies $E_{i}$. The main difference with the GRB
example is that, given that the photon energies are on average much smaller
and that they have a narrower distribution, as a consequence of the effect of
the threshold for $\Delta$ production the required proton energies to give
rise to a given neutrino energy can be much larger than in the previous
example. In particular, the pronounced increase in $\bar{E}_{p}/E_{i}$ for
values of $E_{i}$ smaller than few hundred TeV can be explained because the
pions are produced mostly at the $\Delta$ resonance and the main contribution
to the neutrinos comes from the low energy tail of the distribution of decay
products. One can see that in this case the median energy of the protons
remains at a value of about 5 PeV for neutrinos with energies below 200 TeV
(for the $\bar{\nu}_{e}$ the median proton energy is a factor about 2.5 times
higher, since their production requires the multipion production to be
allowed). Since the low energy tail of the neutrino distribution is quite
suppressed, one then expects a strong suppression to appear in the neutrino
spectrum at low energies. The resulting spectrum is shown in Fig. 11, where
indeed one can see that it steeply increases with energy and only flattens at
energies above few hundred TeV. The spectral cutoff in the neutrino spectrum
is also found to be much smoother than the cutoff adopted for the CR spectrum.

Figure 9: Probability distribution for the parent proton energies that produce
neutrinos or photons with different energies $E_{i}$ in scenarios with a
thermal photon background with $kT=10$ eV, adopting $\alpha=2$.

Figure 10: Left panel: ratio $\bar{E}_{p}/E_{i}$ as a function of $E_{i}$ in
$p\gamma$ scenarios with a thermal photon background with $kT=10$ eV, adopting
$\alpha=2$. The results for $\nu_{e}$ are comparable to those of
$\bar{\nu}_{\mu}$. Right panel: the same in terms of the average proton energy
$\langle E_{p}\rangle$.

Figure 11: Total source $\nu$ spectrum as well as the different flavor
contributions in $p\gamma$ scenarios with a thermal photon background with
$kT=10$ eV, adopting a power-law proton spectrum with $\alpha=2$ having a
sharp cutoff at 100 PeV (left panel) or a smooth one with $\delta=1$ (right
panel). Black dashed lines illustrate the shape of the proton spectrum adopted
(with arbitrary normalization).

##  4 Discussion

We have studied the distribution of the energies of the CR protons giving rise
to astrophysical neutrinos with a given energy, both for the $pp$ and the
$p\gamma$ scenarios. In the $pp$ scenarios, given the low energy threshold for
pion production and the slow (logarithmic) rise in the inelastic cross section
and charged particle multiplicities, for power-law CR spectra the neutrino
spectral shape is quite similar, though slightly harder. On the other hand,
for increasing energies the median energy of the protons giving rise to
neutrinos steadily increases beyond the usual estimate of $E_{p}\simeq
20E_{\nu}$, and depends on the production channel (pion or muon decay). This
increase depends on the CR spectral slope, and for instance it is about
$\bar{E}_{p}\simeq 50E_{\nu}$ for $\alpha=2$, becoming smaller for softer CR
spectra. This has an impact also on the relationship between the CR spectral
cutoff and the resulting neutrino cutoff, which becomes much smoother than the
proton one once the broad distribution of the energies of the produced
neutrinos is taken into account. In this sense, if a sharp suppression in the
astrophysical neutrino spectrum were to be observed, this could be an
indication in favor of the presence of strong magnetic fields in the
interaction region that damp the decaying charged mesons as a consequence of
their synchrotron losses. However, for the effect to be sharp the magnetic
field values in the neutrino producing regions would be required to be similar
in most of the sources. On the other hand, the observation of an astrophysical
neutrino spectrum with an approximately power-law shape could be an indication
of an underlying $pp$ scenario with a large proton cutoff energy. A
distinguishing feature of the $pp$ scenarios is that the spectrum of both the
generated neutrinos and the photons keep their approximate power law shape
down to lower energies, and through the cascading of the photons by
interactions with the EBL radiation this could overshoot the observed diffuse
gamma ray fluxes below TeV energies if the CR spectrum is soft ($\alpha\geq
2.2$) [30, 31]. This tension may be overcome for hidden neutrino sources, i.e.
those that are opaque to the photons.

In the $p\gamma$ scenarios, the threshold for the $\Delta$ production has a
strong impact on the neutrino production, leading to a steep suppression of
the neutrino fluxes below an energy $E_{\nu}\sim 70\,{\rm
TeV}/(\varepsilon/{\rm keV})$ for target photons with energy $\varepsilon$.
This also leads to an increase in the values of $\bar{E}_{p}/E_{\nu}$ for
lower energies, and in particular in the case considered of the relatively
narrow thermal distribution with $kT=10$ eV we found that due to the effect of
the $\Delta$-resonance threshold the value of $\bar{E}_{p}$ remains about 5
PeV for $E_{\nu}<200$ TeV (being about 12 PeV for the case of
$\bar{\nu}_{e}$). For target photons with a broken power-law energy
distribution, as the one that we considered with a break at
$\varepsilon^{\prime}_{\rm b}=1$ keV characteristic of GRB scenarios, the high
energy neutrinos with $E^{\prime}_{\nu}>10$ TeV get produced predominantly
through interactions with the photons below the break, which have a harder
spectrum, but due to the much broader photon spectral distribution the changes
in the values of $\bar{E}_{p}/E_{\nu}$ are less pronounced than in the case of
a thermal photon distribution.

To conclude, the research field of high-energy neutrinos is today in a very
exciting state with great promise for the future, particularly in terms of
what it will teach us about the sources of cosmic rays. We have quantitatively
discussed the connection between neutrinos and cosmic rays, examining the
correspondence between the energies of the primaries and of the secondaries
for different source scenarios. Knowing which kind of sources produce the
neutrinos will be fundamental for a proper interpretation of the results and
the identification of specific features in the neutrino spectrum will be of
crucial importance for this.

## Acknowledgments

The work of E.R. was supported by CONICET (PIP 2015-0369) and ANPCyT (PICT
2016-0660). The work of F.V. was partially supported by the grant 2017W4HA7S
“NAT-NET: Neutrino and Astroparticle Theory Network” under the program PRIN
2017 funded by the Italian Ministero dell’Istruzione, dell’Università e della
Ricerca (MIUR).

## References

  * [1] M.G. Aartsen et al. (IceCube Collaboration), Phys. Rev. Lett. 113 (2014) 101101
  * [2] M.G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 833 (2016) 3
  * [3] M.G. Aartsen et al. (IceCube Collaboration), Phys. Rev. Lett. 124 (2020) 051103
  * [4] M.G. Aartsen et al. (IceCube Collaboration), Phys. Rev. Lett. 125 (2020) 121104
  * [5] R. Abbasi et al. (IceCube Collaboration), (2020) arXiv:2011.03545
  * [6] J. Stettner (for the IceCube Collaboration), PoS-ICRC2019-1017
  * [7] P. Mészáros, Ann. Rev. Nucl. Part. Sci. 67 (2017) 45
  * [8] M. Ahlers and F. Halzen, Prog. Part. Nucl. Phys. 102 (2018) 73
  * [9] V. Berezinsky and G. Zatsepin, Phys. Lett. 28B (1969) 423
  * [10] E. Roulet, G. Sigl, A. van Vliet and S. Mollerach, JCAP 01 (2013) 028
  * [11] M.G. Aartsen et al. (IceCube Collaboration), Science 361 (2018) 147
  * [12] S.R. Kelner, F.A. Aharonian and V.V. Bugayov, Phys. Rev. D 74 (2006) 034018
  * [13] A. Müecke, J.P. Rachen, R. Engel, R.J. Protheroe and T. Stanev, Nucl. Phys. B. Proc. Suppl. 80 (2000) 0810
  * [14] K. Murase and S. Nagataki, Phys. Rev. D 73 (2006) 063002
  * [15] R. Engel, F. Riehn, A. Fedynitch, T.K. Gaisser and T. Stanev, EPJ Web of Conferences 145 (2017) 08001
  * [16] S.R. Kelner and F.A. Aharonian, Phys. Rev. D 78 (2008) 034013
  * [17] A. Mücke, R. Engel, J.P. Rachen, R.J. Protheroe and T. Stanev, Comput. Phys. Commun. 124 (2000) 290
  * [18] E. Kafexhiu, F. Aharonian, A. Taylor and G.S. Vila, Phys. Rev. D 90 (2014) 123014
  * [19] D. Supanitsky, Phys. Rev. D 93 (2016) 043001
  * [20] J.C. Joshi, W. Winter and N. Gupta, MNRAS 439 (2014) 3414
  * [21] P. Lipari, Astropart. Phys. 1 (1993) 195
  * [22] F.L. Villante and F. Vissani, Phys. Rev. D 78 (2008) 103007
  * [23] T. Kamae, N. Karlsson, T. Mizuno, T. Abe and T. Koi, Astrophys. J 647 (2006) 692
  * [24] J. Becker Tjus, B. Eichmann, F. Halzen, A. Kheirandish and S.M. Saba, Phys. Rev. D 89 (2014) 123005
  * [25] J.P. Rachen and P. Mészáros, Phys. Rev. D 58 (1998) 123005
  * [26] E. Waxman and J. Bahcall, Astrophys. J. 541 (2000) 707
  * [27] T. Kashti and E. Waxman, Phys. Rev. Lett. 95 (2005) 181101
  * [28] S. Hümmer, M. Rüger, F. Spanier and W. Winter, Astrophys. J. 721 (2010) 630
  * [29] W. Winter, Adv. High Energy Phys. 2012 (2012) 586413
  * [30] K. Murase, M. Ahlers and B.C. Lacki, Phys. Rev. D 88 (2013) 121301
  * [31] K. Bechtol, M. Ahlers, M. Di Mauro, M. Ajello and J. Vandenbroucke, Astrophys. J. 836 (2017) 47