# Improved Measurements of the IXPE Crab Polarization

Josephine Wong Kavli Institute for Particle Astrophysics and Cosmology,
Stanford University, 452 Lomita Mall, Stanford, CA 94305, USA  Department of
Physics, Stanford University, 382 Via Pueblo Mall, Stanford CA 94305
<EMAIL_ADDRESS>Roger W. Romani Kavli Institute for Particle Astrophysics
and Cosmology, Stanford University, 452 Lomita Mall, Stanford, CA 94305, USA
Department of Physics, Stanford University, 382 Via Pueblo Mall, Stanford CA
94305  Jack T. Dinsmore Kavli Institute for Particle Astrophysics and
Cosmology, Stanford University, 452 Lomita Mall, Stanford, CA 94305, USA
Department of Physics, Stanford University, 382 Via Pueblo Mall, Stanford CA
94305

###### Abstract

X-ray polarization from the Imaging X-ray Polarimetry Explorer (IXPE) provides
an important new probe of the geometry of the pulsar emission zone and of
particle acceleration in the surrounding pulsar wind nebula (PWN). However,
with IXPE’s modest $\sim\,20-30^{\prime\prime}$ spatial resolution, separation
of the pulsar signal from the nebula is a challenge. Conventional analysis
defines an “off” phase window as pure nebular emission and subtracts its
polarization to isolate the phase-varying pulsar (“on-off fitting”). We
present a more sensitive scheme that uses external measurements of the nebula
structure and pulsar light curve to isolate their contributions to the phase-
and spatially-varying polarization via least-squares regression (“simultaneous
fitting”). Tests with simulation data show $\sim$30% improvement in pulse
phase polarization uncertainties, decreased background systematics, and
substantially improved nebular polarization maps. Applying “simultaneous
fitting” to early IXPE Crab data extracts additional phase bins with
significant polarization. These bins show interesting departures from the
well-known optical polarization sweeps, although additional exposure will be
needed for precise model confrontation.

X-ray polarimetry, pulsar wind nebula, data reduction

##  1 Introduction

The Imaging X-ray Polarimetry Explorer (IXPE) (Weisskopf et al., 2022) is
NASA’s first satellite dedicated to measuring X-ray polarization. Since it was
launched in Dec 2021, it has observed several pulsar wind nebula (PWN),
including Crab, Vela, MSH 15-5(2), and B0540. In these sources, a central
pulsar is embedded in a bright halo of relativistic charged particles
accelerated by the pulsar’s spin-down. The pulsed emission is generally
thought to originate from the magnetically-controlled flow within and just
beyond the light cylinder. Far from the light cylinder, particles and fields
form a randomized wind that powers the nebula (the PWN). The nature of
particle acceleration in both regions is still not fully understood although
many PSR/PWN models have been proposed. For PWN reviews, see Gaensler & Slane
(2006) and Amato (2020). IXPE measurements of synchrotron X-ray polarization
allows us to probe the underlying magnetic structure of these zones and test
these models.

Due to IXPE’s modest HPD (Half Power Diameter)$\sim 20-30^{\prime\prime}$
spatial resolution, isolating the pulsar signal from the surrounding nebular
emission is a challenge. The conventional method is to assume an “off” phase
of pure nebular emission and subtract it to isolate the varying pulsar
contribution (Bucciantini et al., 2023b). In this paper, we introduce a more
sensitive method using Chandra measurements of the nebula map and pulsar light
curve, with high spatial and time resolution, to weight their respective
contributions to the phase- and spatially-varying polarization, extracting the
polarization properties via least-squares regression. In Section 2, we
describe this method and demonstrate, using simulated Crab observations, that
it produces an improved nebula polarization map and more accurate pulsar
polarization curve. In Section 3, we apply this technique to the IXPE Feb/Mar
2022 observation of the Crab and comment on new features visible in the data.
In Section 4, we discuss notable differences between our X-ray measurements
and previous optical polarization measurements and further extensions of our
work. We conclude with brief remarks in Section 5.

##  2 Simulation Analysis

We used IXPE’s simulation and analysis software (IXPEsim V14.3.1 and
IXPEobssim V28.1.0, using V11 response functions) (Baldini et al., 2022) to
generate a mock 100ks observation of the Crab PWN. In the simulation, the
pulsar light curve and phase varying spectrum are derived from Chandra (CXO)
HRC/LETG observations (Weisskopf et al., 2004) and the nebula image comes from
a contemporaneous CXO ACIS image (OBSID 23539, PI Slane) obtained to support
the IXPE study.

ACIS CCD limitations produce two artifacts in the image. First, photons
collected during readout produce trails along the charge-transfer direction: a
bright narrow narrow streak from the pulsar and a diffuse rectangular patch
from the nebula. We excised the streak, replacing it with a random sample of
events, half each from two streak-sized regions above and below. The streak is
close to the nebula symmetry axis, so this gives a fairly continuous flux and
spectrum in the corrected image. The diffuse flux from the nebula was
corrected by defining a nebula boundary, marking the zones extended to either
side along the readout axis and replacing events in these zones with events
sampled from a readout-free background region. Second, CCD pile-up eliminates
counts in the pulsar core and leaving an annulus of pulsar counts in the PSF
wings. We excised a circular region of radius $\sim$3.2′′ and replaced it with
photons drawn from a $\sim$18′′-wide surrounding hemisphere (excluding flux at
the base of the jet).

The IXPE detector reads out pixels directly without charge transfer (and with
small $\sim 1.1$ ms dead time/event). Thus our corrected image better matches,
in morphology and spectrum, the nebular flux that IXPE records. To this
corrected PWN image, we add the phase-variable pulsar point source, with
photons spatially distributed according to the IXPE PSF.

Figure 1: On-off and simultaneous fits for the pulsar polarization of a
simulated 100ks Crab observation. The black line represents the polarization
model. The light curve is displayed in the background (gray) with the off-
pulse phase region shaded in purple. By definition, $Q_{\rm pulsar}=U_{\rm
pulsar}=0$ in these phase bins for on-off analysis.

To complete the simulation, we assembled a plausible model of the X-ray
polarization. For the pulsar, we used OPTIMA optical polarization measurements
(Słowikowska et al., 2009) but compressed the sweep so that the position
relative to the optical peaks was mapped to the X-ray peaks. For the nebula,
we assumed that magnetic field lines follow elongated features (e.g. toroidal
near the central wisps, arced with the limb at the outskirts, and parallel to
the jet) to construct a polarization map.

In the on-off method, we divided the data into 29 equal-spaced phase bins,
with $0.897-1.103$ (6 phase bins, see Fig. 1) as the off-pulse phase. For the
pulsar polarization, we restricted our analysis to a circular 20′′ radius
aperture centered on the pulsar to minimize nebula contamination and a single
energy bin between 2-8 keV (IXPE optimal energy range). To determine the
nebula polarization, we used 5 energy bins and a 150′′ $\times$ 150′′ area,
divided into a $15\times 15$, 10′′ pixel grid, mapping the flux from the off-
pulse phases. This sub-HPD pixel grid gives a more detailed map of the
polarization morphology of the nebula, although the polarization value
assigned to each pixel will be slightly influenced by adjacent pixels. This
may be mitigated by PSF-based deconvolution of the final maps.

In the simultaneous method, we used the latter binning scheme and solved for
the desired $q$ and $u$ of the pulsar phase-resolved and nebula spatially-
resolved spectra using the observed $Q$ and $U$ fluxes as:

$\begin{split}Q_{ijklm}&=I_{{\rm psr},ijklm}\times q_{{\rm psr},ij}+I_{{\rm
neb},ijklm}\times q_{{\rm neb},jkl}\\\ U_{ijklm}&=I_{{\rm psr},ijklm}\times
u_{{\rm psr},ij}+I_{{\rm neb},ijklm}\times u_{{\rm neb},jkl}\end{split}$ (1)

where the indices $i$, $j$, $k$, and $l$ represent the phase, energy, and
spatial position of the bin and $m=1-3$ refers to the three IXPE telescopes.
Assuming equally-spaced phase bins,

$\begin{split}I_{{\rm psr},ijklm}&=\mathscr{I}_{{\rm
psr},ijm}\times\mathrm{PSF}_{jklm}\\\ I_{{\rm neb},ijklm}&=\mathscr{I}_{{\rm
neb},jklm}\ /\ i_{\rm max}.\end{split}$ (2)

where $\mathscr{I}_{{\rm psr},ijm}$ and $\mathscr{I}_{{\rm neb},jklm}$ are the
expected counts determined from CXO measurements of the phase-dependent pulsar
spectrum and the energy-resolved image of the nebula, passed through
IXPEobssim to account for the instrument response. The three IXPE telescopes
have significantly different PSFs (and slightly different effective areas),
hence the $m$ dependence. In practice, we use a long 1Ms simulation to predict
the counts for a shorter observation to reduce statistical errors.

We want to find the parameters that minimize the Gaussian error, where the
variance of $Q_{ijklm}$ and $U_{ijklm}$ are given by Kislat et al. 2015:

$\begin{split}\mathrm{var}(Q)&=N(\frac{2}{\mu^{2}}-\bar{q}^{2})\\\
\mathrm{var}(U)&=N(\frac{2}{\mu^{2}}-\bar{u}^{2})\\\
\mathrm{cov}(QU)&=-N\bar{q}\bar{u}\end{split}$ (3)

We have now constructed an over-determined least-squares problem for which a
best-fit solution must exist. The fitting was performed using the scipy
function, lsq_linear, which allows for the specification of parameter bounds.
We found that it was sometimes helpful to introduce physical limits
$\\{-1,1\\}$ on the Stokes parameters to obtain a good fit to the model.
Because bounded-value least squares is an iterative optimization algorithm,
when values reach the physical $q,u\in\\{-1,1\\}$ bounds, the error bars could
not be obtained analytically. To handle such cases, we recovered error bars
with bootstrap analysis. We found that 500 bootstrap iterations was enough for
the uncertainties to converge to within 1.5%, and confirmed that, when bounded
value least-squares converged without hitting the bounds, standard error
propagation gives accurate fit errors.

Figure 2: On-off and simultaneous fits for nebula polarization from a
simulated 100ks Crab observation. Left column is the model maps of Stokes
$q=Q/I$ and $u=U/I$. Their color bar is common to all panels. Beside each
reconstructed map (middle and right columns) is a histogram of the
uncertainties of the polarization values. Note that simultaneous fitting
produces maps that more cleanly resemble the model with smaller uncertainties.

Figures 1 and 2 show the reconstructed polarization of the two fitting methods
as well as the input model. To compare the two methods, we used three summary
statistics elaborated below: the median error bar size, the GoF, and the
number of measurements. For the pulsar measurements, the median error bars
decrease by $\times 1.20$, averaged between q and u; there are also
$(29/23=1.26)\times$ more measurements in the simultaneous case. We can
further characterize the systematic errors in the recovery by a ‘Goodness of
Fit’

$\mathrm{GoF}_{q}=\\{\sum_{n}[(q-q_{\rm mod})/\sigma_{q}]^{2}/n\\}^{1/2}$

and similarly for $u$, with an average $\times 1.08$ improvement in recovery
of the original model, even relative to the improved errors. Thus, in total,
simultaneous fitting can be considered to provide a $1.20\times 1.26\times
1.08\approx 1.63\times$ improvement in recovery of the pulsar polarization.

For the nebula, since simultaneous fitting uses nearly all the data away from
the pulsar, especially off the brightest portion of the peaks, the effective
nebula exposure is larger by $1/\Delta\phi_{\rm off}\approx 4.9\times$ than in
the simple “off” portion of the pulse phase used to map the nebula in the on-
off method. Further, the method takes account of the small expected pulsar
flux in off phases to provide a cleaner measurement of the true nebulae
structure very close to the pulsar. The polarization maps in Figure 2 show
substantially improved source recovery and decreased typical uncertainty in
the recovered pixels. Quantitatively, the median error bars decrease by $\sim
2\times$ for both q and u. Although the GoF for simultaneous fitted nebula is
slightly larger than the on-off nebula, since the uncertainty is halved, this
actually means that simultaneous fitting is able to recover polarization
values closer to truth.

##  3 Application to First IXPE Crab Observation

In February and March 2022, IXPE made a $\sim$91 ks observation of the Crab,
its first PWN source. At the time, the mirrors were misaligned from the
pointing axis by $\sim$3 arcmin, and we do not have energy-dependent response
functions calibrated for this offset. This (and the incursion of Poisson
statistics in low count bins, see below) led us to analyze in a single 2–8 keV
IXPE energy bin.

Moreover, it has been discovered that errors in the present track
reconstruction of IXPE photon conversion points are correlated with the
initial direction of the photoelectron track, and hence, with the inferred
event polarization. This leads to “polarization leakage,” which induces
polarized fringes about sharply localized X-ray sources (point source and
compact nebulae). As described by Bucciantini et al. (2023a), these fringes
average away for point sources, and hence do not affect aperture polarization
measurements, but do affect the edges of extended sources, such as the Crab
Nebula. The paper also gives a prescription for correcting these fringes,
assuming a smooth Gaussian blurring of the point sources.

We implement here an improved version of this correction, using more detailed
IXPE PSFs for each mirror assembly, derived from ground calibration data and
on-sky observations of point sources. The mirror PSF are lightly smoothed to
suppress numerical noise. We then compute maps of the Hessian terms $H_{xx}$,
$H_{yy}$ and $H_{xy}$. On-sky images have residual blur beyond the mirror
PSFs, produced by incomplete aspect correction and imperfect estimation of the
photon conversion points. We treat these as simple Gaussian blurs with
$\sigma_{\rm G}=2.1^{\prime\prime},\,1.4^{\prime\prime},\,1.2^{\prime\prime}$
for detector units (DUs) 1, 2, 3 respectively. In addition, the correlation
between the conversion point and the EVPA induces a ‘leakage’ blur, which is
energy dependent and grows at large photon energies. We model the effect as
$\sigma_{\rm L}=(10+3\Delta E)^{1/2}$ arcsec, with $\Delta E=E_{\rm keV}-4$keV
for photon energy above $4$ keV and $\Delta E=0$ below. This is common to all
detectors. The effective PSFs for unpolarized sources are

$\begin{split}P_{I}^{*}(\vec{\mathbf{x}})&=(P_{M}\star
G(\sigma_{G}))(\vec{\mathbf{x}})\\\
&\qquad+\frac{\sigma_{L}^{2}}{4}(H_{xx}(\vec{\mathbf{x}})+H_{yy}(\vec{\mathbf{x}}))\\\
P_{Q}^{*}(\vec{\mathbf{x}})&=\frac{\sigma_{L}^{2}}{4}(H_{xx}(\vec{\mathbf{x}})-H_{yy}(\vec{\mathbf{x}}))\\\
P_{U}^{*}(\vec{\mathbf{x}})&=\frac{\sigma_{L}^{2}}{2}(H_{xy}(\vec{\mathbf{x}}))\end{split}$
(4)

with $P_{M}$ the appropriate mirror PSF and $\star\ G$ the symmetric Gaussian
convolution. For a polarized source, $P_{I,Q,U}$ are mixed by the blurring
effect as

$\begin{split}P_{I}(\vec{\mathbf{x}})&=P_{I}^{*}(\vec{\mathbf{x}})+{1\over
2}[q_{\rm src}P_{Q}^{*}(\vec{\mathbf{x}})+u_{\rm
src}P_{U}^{*}(\vec{\mathbf{x}})]\\\ P_{Q}(\vec{\mathbf{x}})&=q_{\rm
src}P_{I}^{*}(\vec{\mathbf{x}})+P_{Q}^{*}(\vec{\mathbf{x}})\\\
P_{U}(\vec{\mathbf{x}})&=u_{\rm
src}P_{I}^{*}(\vec{\mathbf{x}})+P_{U}^{*}(\vec{\mathbf{x}}),\end{split}$ (5)

where the $q_{\rm src}$ and $u_{\rm src}$ are the estimates of the true source
polarization.

In our implementation, we iterate, starting by fitting to raw IXPE $I$, $Q$,
$U$ data, to arrive at leakage-corrected polarization fits in $\sim$5 steps.
In the method, we form leakage-corrected $q$ and $u$ maps for each phase and
energy bin and detector. These are fed to the simultaneous fitting minimizer,
which separates the nebula and pulsar signals. In our case, during fitting,
the energy bins are collapsed. The results are final derived nebula $q$ and
$u$ maps and pulsar phased $q$ and $u$ values, corrected for the spatially-
dependent polarization leakage. In simulation, our prescription has shown to
improve recovery of the original $q_{src}$ and $u_{src}$. For details, see
Dinsmore et al. (in prep 2023). Our analysis below includes this correction,
which can make substantial modifications to the nebula polarization map near
the pulsar and at the outer edges. Near-pulsar corrections, in turn, feed back
to changes in the phase-resolved polarization.

Using the same procedure as for simulated data, we performed on-off and
simultaneous fits to the Crab observation. For on-off, we used the same phase
bins, area, and off-pulse range as the initial IXPE discovery paper on the
Crab (Bucciantini et al. 2023b) and replicated those results.

For simultaneous fitting, we used a 150′′ $\times$ 150′′ area, binned into a
15 $\times$ 15 pixel grid. Using the phase bins of Bucciantini et al. (2023b),
we also obtain significant polarization detection in the main peak (P1) phase
bin $[0.12,0.14]$ with a higher significance of PD = $15.1\pm 2.1\%$. The
second peak (P2) phase bin $[0.515,0.545]$ was also found to be significantly
polarized with PD = $8.8\pm 2.8\%$; these errors decrease from the on-off
values.

Figure 3: IXPE Crab pulsar Stokes Q/I and U/I, after leakage correction, from
on-off and simultaneous fitting. Statistically significant detections
($>3\sigma$) are circled red and blue, respectively. With simultaneous
fitting, due to the extra signal we extract from the data, we can afford to
use finer (and more) phase bins, which helps us capture the polarization
sweep. Both methods have $>5\sigma$ detection at the first peak. The X-ray
light curve is displayed in gray.

Seeing that we can measure polarization with higher significance, we decided
to use smaller bins around the peak phases, where we can expect fast sweeps in
polarization angle, to better measure the polarization curve. Here we employ
16 phase bins, as compared to the 11(+off) phase bins used in the on-off
analysis. With these new bins and corrections, we obtain refined Stokes q and
u phase points (Figures 3, 4 and Table 1), which may be compared with the
standard on-off results.

Polarization angle (PA) and degree (PD) are useful for comparison with models,
but these do not have simple Gaussian errors. We can however plot the
confidence regions for the phase bins around the two peaks (Figure 4), where
three phase bins have above $3\sigma$ significance, and four more are above
$2\sigma$. In P1, we see a counterclockwise (CCW) sweep from $\sim
100^{\circ}-150^{\circ}$. P2 hints at a CCW sweep as well, though since it
only has one significant bin, more data is required to see this.

Figure 4: IXPE Crab pulsar PD-PA contour plots showing
$1\sigma,2\sigma,3\sigma$ (red, green, blue) contours for select simultaneous-
fitted bins around the two peaks. Phases with greater than $5\sigma$ detection
are the brightest in color and fade for the $>2\sigma$ bins. Numbers
correspond with the phase bins numbered in Figure 3. For P1, a CCW sweep can
be inferred; P2 requires more bins to show a polarization pattern.

Figure 5 shows the reconstructed nebula PD map, cut at 4.7$\sigma$
significance, with green bars showing PA and magnetic field direction. As seen
in Bucciantini et al. 2023b, the PD is reduced at the sides of the torus where
the PA sweeps rapidly through 180∘. With leakage correction, these feature are
enhanced, appearing as PD holes on opposite sides of the nebula. In general,
the PD is largest towards the northern and southern edges of the nebula; in
part, this may be because such regions have toroidal PA at a similar angle
across the PSF, minimizing polarization beam dilution. A highly polarized PD
$>$ 50% region west of the jet is also present. Generally, the magnetic field
lines follow the filamentary structure partly visible in the background image.

Figure 5: IXPE Crab nebula polarization map, after leakage correction. Green lines indicate magnetic field direction (perpendicular to the electric field PA) and the tile color indicates PD. Contour intensity levels (red) derived from the contemporaneous Chandra Crab observation (OBSID 23539), both displayed in the background, are included for reference; the pulsar position is marked with a red cross. On-Off | Simultaneous
---|---
| Phase | q | u | PD (%) | PA (∘) |  | Phase | q | u | PD (%) | PA (∘)
P1 | $0.120-0.140$ | $-0.132\pm 0.025$ | $-0.079\pm 0.025$ | $15.4\pm 2.5$ | $105\pm 18$ | 2 | $0.085-0.124$ | $-0.087\pm 0.029$ | $-0.022\pm 0.029$ | $9\pm 3$ | $97\pm 9$
3 | $0.124-0.138$ | $-0.132\pm 0.027$ | $-0.065\pm 0.027$ | $15\pm 3$ | $103\pm 5$
4∗ | $0.138-0.142$ | $-0.048\pm 0.053$ | $-0.100\pm 0.053$ | $11\pm 5$ | $122\pm 14$
5∗ | $0.142-0.179$ | $0.028\pm 0.032$ | $-0.059\pm 0.032$ | $6\pm 3$ | $145\pm 14$
P2 | $0.515-0.545$ | $-0.013\pm 0.030$ | $-0.06\pm 0.030$ | $6.3\pm 3.0$ | $129\pm 14$ | 12∗ | $0.495-0.530$ | $-0.066\pm 0.031$ | $-0.044\pm 0.031$ | $8\pm 3$ | $107\pm 11$
13 | $0.530-0.545$ | $0.036\pm 0.039$ | $-0.158\pm 0.039$ | $16\pm 4$ | $141\pm 7$
14∗ | $0.545-0.595$ | $-0.006\pm 0.040$ | $-0.090\pm 0.040$ | $9\pm 4$ | $133\pm 13$
Table 1:  Numerical comparison of published on-off measurements (Bucciantini
et al., 2023b) and simultaneous fit measurement of pulsar peak polarizations.
∗Note: These low significance bins have substantial PD-PA covariance; see
Figure 4 for 2D error contours.

##  4 Discussion

In this paper, we show that simultaneous fitting has several merits. As tested
with simulated data, it provides improved recovery of the pulsar polarization
sweep with smaller statistical uncertainties on both pulsar and nebula
measurements. This allows us to use smaller phase bins to more finely resolve
the pulsar polarization. Our method essentially uses externally measured, high
precision data on the time and spatial intensity structure to provide a
weighting that improves extraction of the polarization signal. Perhaps the
greatest limitation in the present implementation is our assumption of
Gaussian statistics. While this allows an efficient linear algebra solution
for the many pulse phase and nebula image polarization values, it does break
down when the counts in a given bin are too low. We have noticed that such
Poisson effects may be important in the low count outskirts of the PWN,
although in the computations presented here $<2.5\%$ of the $16\times 15\times
15\times 3=10,800$ data bins have $<10$ counts. However, this does limit our
ability to extend our decomposition to be fully energy-dependent, as the
higher energy IXPE bins often have low count rate. With additional Crab
exposure it will be straight-forward to extend the analysis to a modest number
of spectral bins.

Figure 6: IXPE X-ray polarization measurements, after leakage correction,
compared with optical values. Optical $q=Q/I$ and $u=U/I$ from off-pulse
subtracted OPTIMA observations are shown as black dots (with the optical light
curve in green for phase reference). The average OPTIMA $q$ and $u$ fluxes
will be zero in the dashed line-delimited off-pulse region. Red error bars
show the IXPE measurement uncertainties. Phase bins whose PD exceeds $3\sigma$
significance are circled in blue, those exceeding $2\sigma$ in dark green.
Note that while optical and X-ray $q$ values at the second peak are
consistent, there are substantial differences for the other peak
polarizations.

The polarization of the optical Crab pulsed emission has been measured
repeatedly in the 50 years since its discovery. The measurements of
Słowikowska et al. (2009) with the OPTIMA fiber photometer are of particularly
high quality. The central fiber used for these measurements had a diameter of
$2.3^{\prime\prime}$ on the sky and thus included both the Crab pulsar and
nebular emission such as the “inner knot”. Accordingly these authors define a
“background” phase of $\phi=0.78-0.84$ (dashed lines in Figure 6) and subtract
the Q and U fluxes from this phase to get the pulsed Crab emission. This turns
the optical Crab curves into “on-off” measurements. In Figure 6, we compare
this pulsed optical polarization signal (measurements kindly supplied by G.
Kanbach) with our new IXPE measurements. Alignment of the optical and X-ray
curves to the radio phase convention was checked via their light curves, with
the X-ray peak at $\phi=0.99$ and the optical peak at $\phi=0.994$ (Enoto et
al., 2021).

There have been several attempts to compare the OPTIMA polarization data with
theoretical pulsar models (e.g. Dyks et al., 2004), none particularly
satisfactory. Clearly some ingredients are missing in our basic understanding
the pulse emission. In particular, absorption effects or contributions from
multiple emission regions may introduce complications beyond the essentially
geometric models that have been applied to date. Such effects should differ
between the X-ray and optical bands. So it is encouraging that we do see
statistically significant differences between the optical OPTIMA and IXPE
polarization curves. Most notably, the P1 $q$ sweep starts from larger
negative values and is delayed in phase. Negative $u$ values also persist
later in P1. Similarly $u$ is much more negative in the core of P2. We do not
find any significant measurements in the pulse minimum bin which is
substantially nebula dominated, with a maximum of $\sim 20\%$ pulsar counts in
the central pixel (and steeply tapering in adjacent pixels). Additional
signal-to-noise (S/N) can allow better separation of the pulse signal near
this phase. Near the peaks the polarization values have significant
sensitivity to our choice of bin boundaries, so we suspect that unresolved
rapid sweeps still suppress the polarization. With additional exposure, IXPE
can measure several more phase bins with good precision. Comparison of the
optical and X-ray signals can then give new insights into the geometry of the
pulse emission zones. By September 2023, IXPE will conclude a follow-up 300 ks
observation of the Crab, which, combined with the current dataset and event
quality weighting, should more than double the S/N presented here.

##  5 Conclusion

Analysis of the first IXPE Crab observation using simultaneous fitting has
improved our measurements of the Crab polarization. Compared to the original
paper (Bucciantini et al. 2023b), we recover more bins of significant pulsar
polarization and are able to use a finer phase resolution to see departures
from the well-measured optical polarization. Nebula features, such as the PD
holes at the edge of the torus, are better recovered. With the additional S/N
from planned Crab follow-up exposure, we can substantially extend these gains.
The method is of course general and can be applied to any phase varying source
embedded in extended emission. We anticipate that application to other IXPE
sources, such as MSH 15-5(2) and possibly B0540 and Vela, can provide improved
measurements as well.

This work was supported by NASA under grant NNM17AA26C. The Imaging X-ray
Polarimetry Explorer (IXPE) is a joint US and Italian mission. The US
contribution is supported by the National Aeronautics and Space Administration
(NASA) and led and managed by its Marshall Space Flight Center (MSFC), with
industry partner Ball Aerospace (contract NNM15AA18C). The Italian
contribution is supported by the Italian Space Agency (Agenzia Spaziale
Italiana, ASI) through contract ASI-OHBI-2017-12-I.0, agreements ASI-
INAF-2017-12-H0 and ASI-INFN-2017.13-H0, and its Space Science Data Center
(SSDC) with agreements ASI-INAF-2022-14-HH.0 and ASI-INFN 2021-43-HH.0, and by
the Istituto Nazionale di Astrofisica (INAF) and the Istituto Nazionale di
Fisica Nucleare (INFN) in Italy. This research used data products provided by
the IXPE Team (MSFC, SSDC, INAF, and INFN) and distributed with additional
software tools by the High-Energy Astrophysics Science Archive Research Center
(HEASARC), at NASA Goddard Space Flight Center (GSFC).

## References

  * Amato (2020) Amato, E. 2020, arXiv e-prints, arXiv:2001.04442
  * Baldini et al. (2022) Baldini, L., Bucciantini, N., Lalla, N. D., et al. 2022, SoftwareX, 19, 101194
  * Bucciantini et al. (2023a) Bucciantini, N., Di Lalla, N., Romani, R. W. R., et al. 2023a, arXiv e-prints, arXiv:2302.00346
  * Bucciantini et al. (2023b) Bucciantini, N., Ferrazzoli, R., Bachetti, M., et al. 2023b, Nature Astronomy
  * Dinsmore et al. (in prep 2023) Dinsmore, J. T., Peirson, A. L., & Romani, R. W. in prep 2023
  * Dyks et al. (2004) Dyks, J., Harding, A. K., & Rudak, B. 2004, ApJ, 606, 1125
  * Enoto et al. (2021) Enoto, T., Terasawa, T., Kisaka, S., et al. 2021, Science, 372, 187
  * Gaensler & Slane (2006) Gaensler, B. M., & Slane, P. O. 2006, ARA&A, 44, 17
  * Kislat et al. (2015) Kislat, F., Clark, B., Beilicke, M., & Krawczynski, H. 2015, Astroparticle Physics, 68, 45
  * Słowikowska et al. (2009) Słowikowska, A., Kanbach, G., Kramer, M., & Stefanescu, A. 2009, MNRAS, 397, 103
  * Weisskopf et al. (2004) Weisskopf, M. C., O’Dell, S. L., Paerels, F., et al. 2004, ApJ, 601, 1050
  * Weisskopf et al. (2022) Weisskopf, M. C., Soffitta, P., Baldini, L., et al. 2022, Journal of Astronomical Telescopes, Instruments, and Systems, 8, 026002


# Calibration of the ICARUS cryogenic photo-detection system at FNAL

M. Bonesini<EMAIL_ADDRESS>R. Benocci  R. Bertoni  A. Chatterjee
M. Diwan  A. Menegolli  G. Raselli  M. Rossella  A. Scarpelli  N. Suarez  on
behalf of the Icarus Collaboration  University and Sezione INFN Milano
Bicocca, Milano, Italy  University of Pittsburg, PA, USA  Brookhaven National
Laboratory, NY, USA  University of Pavia, Pavia, Italy  Sezione INFN Pavia,
Pavia, Italy

###### Abstract

The calibration of the ICARUS photo-detection system is based on a low power
laser diode at 405 nm. Laser pulses arrive to one optical switch and then are
sent to 36 UHV flanges, by 20 meters long optical patches. Light is then
delivered to the ten PMTs connected to a single flange, by 7m long injection
optical patches. Extensive tests of the used components and care in the design
of the optical system have guaranteed to each PMT a sizeable calibration
signal with minimal distortion, with respect to the original one. Gain
equalization of PMTs has reached a 1 $\%$ resolution. In this procedure data
from background photons were also used. The distribution of the PMTs’ signal
arrival time has a distribution with a resolution less than 1 ns, thus
allowing a good determination of the absolute event timing. The status of the
laser calibration system with its possible upgrades will be reported.

###### keywords:

Laser diodes , PMTs , calibration

††journal: Journal of LaTeX Templates

##  1 Introduction

ICARUS T600 [1] is presently used as far detector of the Short Baseline
Neutrino (SBN) program at Fermilab (USA) to search for a possible LSND-like
sterile neutrino signal at $\Delta m^{2}\sim$ O(eV2 ) [2]. It is made of two
identical cryostats, filled with about 760 tons of ultra-pure liquid Argon.
Each cryostat houses two TPCs with 1.5 m maximum drift path, sharing a common
central cathode. Charged particles interacting in liquid argon produce both
scintillation light at 128 nm and ionization electrons. ICARUS is placed at
shallow depth on the Booster Neutrino beam (BNB). To reduce the cosmic ray
background, in addition to a full coverage cosmic ray tagger (CRT), a system
based on 360 large area Hamamatsu R5912-MOD photomultipliers (PMTs) [3],
directly immersed in liquid Argon, is used to detect scintillation light. To
exploit the BNB bunch structure, for futher rejection of out-of-bunch cosmic
events, an overall time resolution $\sim 1$ ns is needed. This and the trigger
system require an accurate calibration in gain and time of each PMT.

##  2 The laser calibration system

The PMTs equalization may be performed by using fast laser pulses, as done in
previous experiments such as Borexino at LNGS [4] and HARP at CERN PS [5]. The
light pulse from a PLP10 Hamamatsu laser diode is sent to each PMT via a
distribution system including an Arden Photonics Mode Scrambler (MS), a DD-100
OZ/Optics programmable attenuator(ATT), a 10 m armed fiber patch cable,
followed by one 1x46 Agiltron optical switch (where 10 channels are spares).
From this 20 m fiber patch cords go to VACOM UHV optical feedthroughs (36) on
CF40 flanges, mounted on CF40-CF20 nibbles, as shown in figure 1 . Inside the
T600 tank, 1x10 Lightel fused fiber splitters attached to each nibble deliver
the input laser signal to the window of each PMT. While an optical switch
delivers the laser signal to a single output channel, a 1xN optical switch
distributes the input signal among N output channels. The light delivery
system must have a minimal spread in channel-to-channel total delay ($\Delta
t$) and delivered signal power in front of each PMT.

Figure 1: Left panel: schematic layout of the laser system realized for
calibration of the 360 PMTs of the Icarus T600 detector. Right panel: image of
one UHV flange (a) with the 1x10 splitter (b) and the FC/FC patch panel (c).

In addition laser pulses have to be delivered to the PMT’s photocathode with a
minimal attenuation and without a deterioration of the original timing
characteristics.

##  3 Optical components characterization

To keep costs low, commercial optical components (fibers, optical switches,
fused optical splitters, …) were used. As these components are easily
available only at Telecom wavelengths ($\sim 850\div 1300$ nm), we had to
characterize their behaviour at 400 nm with a dedicated test bench at INFN
Milano Bicocca. To minimize alignment problems multimode fibers (MM) instead
of single mode ones (SM) were used. At this point, the remaining main task was
to found a suitable optical switch and fused 1x10 optical splitters with
channel-to-channel minimal spread in insertion losses (IL) and delay. The
choice was to use an Agiltron 1x46 optical switch and Lightel 1x10 fused
optical splitters. In the optical switch, the FWHM of the laser pulses is
increased at most by $3\%$ and the insertion loss (IL) is less than 0.5 dB,
with a cross talk $\sim 50$ dB. Channel-to-channel output uniformity is better
than $10\%$ and the RMS of internal time delays is $\sim$ 10 ps. Instead, a
typical splitter increases laser pulse FWHM of less than $4\%$ and has an IL
less than 3 dB. For each splitter, the signal uniformity in the 10 output
channels is less than $5\%$ aside one leading channel, by construction, that
has to be attenuated. Average internal delays of the 36 splitters are within
80 ps (aside two, belonging to a different batch).

UHV optical feedthroughs convey the laser pulse inside the Icarus cryogenic
tank. The adopted solution from VACOM, with a MM 50 $\mu$m core fiber, had a
measured transmission around 80 %, introduced an additional delay $\sim 100$
ps and a negligible additional time dispersion (FWHM). To adapt these flanges
to the existing T600 chimneys CF40-CF200 nibbles are used, see figure 1 for
details. As the splitters have to work, inside the T600 chimneys, at a
slightly lower temperature than the T600 roof and the 7m injection patches at
cryogenic temperatures, inside the LAr bath, tests were done to assess their
dependence from temperature, using a Lauda thermal machine (precision $0.1$
0C) in the range $0\div 50\ ^{0}$C and with $LN_{2}$ bath inside a dewar, see
reference [6].

##  4 System performances

At the PMTs’ front face, a total delay of the laser signal ($\Delta t\sim 250$
ns) with a channel-to-channel spread estimated at less than 200 ps was
measured both in situ and in laboratory. Up to the UHV flanges a $4.59\pm
0.16$ dB optical signal attenuation was measured. The attenuation of the 7 m
injection patches was measured as $0.61\pm 0.16$ dB for the full 410 sample
(best 360 used).

Figure 2: Left panel: time difference $\Delta T$ distribution. Right panel:
PMTs gain distribution after equalization.

For calibration purposes, charge spectra from PMTs were recorded integrating
laser signals over a 100 ns window. For each PMT data were taken at different
voltages to allow the determination of the gain vs voltage curve and thus
provide the voltage needed for a nominal gain set at $5\times 10^{6}$. A fine
tuning with background single photoelectrons is then applied to reduce the
relative gain spread to less than $1\%$, as shown in figure 2. In the same
figure the distribution of the time difference $\Delta T$ between the trigger
and the PMT signal is shown. Timing of PMTs signals is given by
$T_{PMT}=T_{TRIG}+\Delta T$, where the trigger pulse is synchronized with the
beam RF. Laser pulses permit the monitoring of possible variation of $\Delta
T$ along the data taking and between different channels. With the obtained
resolution of 1 ns, it is possible to determine the absolute timing of
collected events, for trigger and cosmic rays rejection purposes.

As the PMTs’ response to calibration pulses is sizeable, as measured in situ
using the optical attenuator, the replacement of the optical switch (10 PMTs
calibrated in a single run) with a custom 1x36 optical splitter (360 PMTs
calibrated in the same run) is foreseen in the next future. This will allow a
much faster laser calibration procedure.

##  5 Conclusions

The light detection system of ICARUS T600 is fully operational allowing a
regular data taking. Calibration of the PMTs’ gain has a spread less then
$1\%$ while PMTs timing has a resolution better than 1 ns. Some upgrades are
foreseen to improve the automation and reduce the time needed for the laser
calibration procedures.

## Acknowledgements

This work was supported by EU Horizon 2020 grant Agreements n. 74303, 822185,
858189 and 101003460.

## References

  * [1] M. Bonesini, A. Menegolli, PoS (EPS-HEP2019) (2019).
  * [2] R. Acciarri, et al., arxiv.org/abs/arXiv:1503.01520.
  * [3] B. Ali-Mohammadzadeh, et al., JINST 15 (10) (2020) T10007.
  * [4] B. Caccianiga, et al., Nucl. Instr. Meth A496 (2003) 353.
  * [5] M. Bonesini, et al., IEEE Trans. Nucl. Phys. 50 (4) (2003).
  * [6] M. Bonesini, et al., JINST 15 (05) (2020) C05042.


Consider a centrality measure $\mathsf{C}$ for digraphs. The following
statements are equivalent:

  1. (1)

$\mathsf{C}$ is a monotonic subgraph motif centrality measure for digraphs.

  2. (2)

For every $n>0$, $|\text{\rm Val}^{n}(\mathsf{C})|\leq n+1-|\text{\rm
BVal}^{n}(\mathsf{C})|$.

###### Theorem G.8.

Consider a centrality measure $\mathsf{C}$ for digraphs. The following
statements are equivalent:

  1. (1)

$\mathsf{C}$ is a monotonic subgraph motif centrality measure for digraphs
relative to the induced ranking.

  2. (2)

For every digraph $G$ and $n>0$, $|\text{\rm Val}_{G}^{n}(\mathsf{C})|\leq
n+1-|\text{\rm BVal}_{G}^{n}(\mathsf{C})|$.

As for the non-monotonic case, the proofs of the above results are identical
with those for Theorems 6.3 and 6.5, respectively. This is because the proofs
only exploit the sizes of the sets $\text{\rm Val}_{G}^{n}(\mathsf{C})$,
$\text{\rm Val}^{n}(\mathsf{C})$, $\text{\rm BVal}_{G}^{n}(\mathsf{C})$ and
$\text{\rm BVal}^{n}(\mathsf{C})$ rather than the structure of their members.

We can finally show a result, analogous to Theorem 5.11, for directed graphs,
which essentially provides a unifying framework for all centrality measures
over digraphs as long as we focus on connected graphs, and we are only
interested in the induced ranking. The notation $\mathsf{ConC}$ for a
centrality measures $\mathsf{C}$ for digraphs is defined in the same way as
for undirected graphs.

###### Theorem G.9.

For a centrality measure $\mathsf{C}$ for digraphs it holds that
$\mathsf{ConC}$ is a monotonic subgraph motif measure for digraphs relative to
the induced ranking.

###### Proof.

The proof is along the lines of the proof for Theorem 5.11. Given a connected
digraph $G$, we build the equivalence relation $\equiv_{G}$ over $V(G)$ in the
same way as in the proof of Theorem 5.11. Due to the fact that $G$ is a
connected digraph, for every vertex $v\in V(G)$, $|\mathsf{A}(v,G)|\geq|V(G)|$
since between $v$ and any other node $u\in V(G)$ there exists a connected
subgraph containing $v$ and $u$. Such a subgraph can be seen as an undirected
path in the multi-graph $\text{UN}(G)$. Now, we can define a subgraph family
$\mathsf{F}$ such that, for every $(v,G)$ with $[v]_{\equiv_{G}}=C_{i}$,
$|\mathsf{F}(v,G)|=i-1$. Finally, with the filtering function $f$ such that
$f(i)=i+1$, it is easy to verify that $\text{\rm Rank}(\mathsf{C})=\text{\rm
Rank}(\mathsf{C}\langle\mathsf{F},f\rangle)$.

### Classification

We proceed to perform a classification as the one performed in Section 7 for
undirected graphs. It turned out that the classification for digraphs is
identical to the one for undirected graphs (see Tables 1 and 2) apart from one
case already discussed at the end of Section 7, that is, the negative result
$\times[\mathit{trees}]$ for $\mathsf{Eigenvector}$ in both Tables 1 and 2.
Although we can show that for directed graphs, $\mathsf{Eigenvector}$ is not a
(monotonic) subgraph motif centrality measure, it remains open whether these
negative results hold even for directed trees. Despite the fact that the
classification for digraphs is the same as the one for undirected graphs, the
underlying proofs are different since now we need to work with directed
graphs, and the versions of the centrality measures in question that operate
in digraphs.

#### Stress.

We define the set of shortest paths from a node $v$ to a node $u$ in a digraph
$G$ as $S_{G}(v,u)=\\{\pi_{vu}\subseteq G\mid|\pi_{vu}|=d_{G}(v,u)\\}$. For a
node $w\in V(G)$, the set of shortest paths from $u$ to $v$ going through $w$
is the set $S_{G}^{w}(v,u)=\\{\pi_{vu}\in S_{G}(v,u)\mid w\in\pi_{vu}\\}$. The
stress centrality of a vertex $v\in V(G)$ is

$\mathsf{Stress}(v,G)\ =\
\displaystyle\sum_{(u,w)\in(V(G)/\\{v\\})^{2}}|S_{G}^{v}(u,w)|.$

###### Proposition G.10.

$\mathsf{Stress}$ is a monotonic subgraph motif measure for directed graphs.

###### Proof.

In the same way as we did for undirected graph, we define the family of
shortest paths going through $v$ in $G$ as

$\mathsf{S}(v,G)\ =\ \bigcup_{u,w\in
V(G)\setminus\\{v\\}}\left\\{G_{\pi}\mid\pi\in S_{G}^{v}(u,w)\right\\},$

$G_{\pi}$ is the subgraph induced by the path $\pi$ in $G$. We also consider
the filtering function $id$, i.e., the identity function over $\mathbb{N}$,
which is clearly monotonic. It is easy to see that, for very pair $(v,G)$,

$\mathsf{Stress}(v,G)\ =\mathsf{C}\langle\mathsf{S},f\rangle(v,G),$

and the claim follows.

#### All-Subgraphs.

The all-subgraphs centrality measure can be defined for digraphs in the
obvious way. It is easy to see that:

###### Proposition G.11.

$\mathsf{All\text{\rm-}Subgraphs}$ is a monotonic subgraph motif measure for
directed graphs.

###### Proof.

The $\mathsf{All\text{\rm-}Subgraphs}$ measure for digraphs simply counts
weakly connected subgraphs, and the filtering function is $\log_{2}$.
Therefore,
$\mathsf{All\text{\rm-}Subgraphs}(v,G)=\log_{2}(|\mathsf{A}(v,G)|)$, for every
pair $(v,G)$, and the claim follows.

#### Degree.

Given a digraph $G$ and a vertex $v\in V(G)$, the degree centrality of $v$ in
$G$ is

$\mathsf{Degree}(v,G)\ =\ |\text{OUT}(v,G)|+|\text{IN}(v,G)|.$

We can show the following result:

###### Proposition G.12.

$\mathsf{Degree}$ is a monotonic subgraph motif measure for directed graphs.

###### Proof.

Consider the subgraph family $\mathsf{D}$ such that, for every
$(v,G)\in\mathbf{VG}$,

$\mathsf{D}(v,G)\ =\ \\{(\\{v,u\\},\\{(v,u)\\})\mid
u\in\text{OUT}(v,G)\\}\cup\\{(\\{v,u\\},\\{(u,v)\\})\mid
u\in\text{IN}(v,G)\\}$

and the filtering function $id$ that is the identity over $\mathbf{N}$, i.e.,
$id(x)=x$, which is monotonic. It is clear that
$\mathsf{Degree}(v,G)=id(|\mathsf{D}(v,G)|)$, and the claim follows.

#### Cross-Clique.

The cross-clique centrality measure (Everett and Borgatti, 1998) was
introduced in the context of analyzing virus and information propagation. It
is defined as the amount of cliques containing a certain vertex. Given a
digraph $G$ and a vertex $v\in V(G)$, we define the set
$\mathit{Cl}(v,G)=\\{S\in\mathsf{A}(v,G)\mid S\text{ is a clique}\\}$. The
cross-clique centrality of $v$ in $G$ is

$\mathsf{Cross\text{\rm-}Clique}(v,G)\ =\ |\mathit{Cl}(v,G)|.$

It is straightforward to see that the following holds:

###### Proposition G.13.

$\mathsf{Cross\text{\rm-}Clique}$ is a subgraph motif measure for directed
graphs.

###### Proof.

We consider the subgraph family $\mathsf{K}(v,G)=\\{S\in\mathsf{A}(v,G)\mid
S\text{ is a clique}\\}$ and the identity filtering function $id$. Clearly,
$\mathsf{Cross\text{\rm-}Clique}(v,G)=id(|\mathsf{K}(v,G)|)$, for every pair
$(v,G)$, and the claim follows.

#### Closeness.

For a digraph $G$ and vertices $u,v\in V(G)$, the distance from $u$ to $v$ in
$G$ is $d_{G}(u,v)=\displaystyle\min\\{|\pi_{uv}|\mid\pi_{uv}\text{ is a path
from }u\text{ to v}\\}$. Then, we define the sum of distances for v in $G$ as
$\text{SD}(v,G)=\displaystyle\sum_{u\in K_{v}(G)/\\{v\\}}d_{G}(v,u)$. Finally,
we define closeness centrality as

$\mathsf{Closeness}(v,G)\ =\ \displaystyle\frac{1}{\text{SD}(v,G)}.$

By definition, the $\mathsf{Closeness}$ centrality for a pair $(v,G)$ such
that $|K_{v}(G)|=1$ is $0$.

###### Proposition G.14.

$\mathsf{Closeness}$ is not a subgraph motif measure (relative to the induced
ranking) for directed graphs.

###### Proof.

As we did for Proposition 4.3, we will show that $\mathsf{Closeness}$ does not
satisfy the bounded value property nor the non-uniform bounded value property.
For the former, we consider a set $S\subsetneq\mathbf{VG}$ such that, with
$A_{n}=\\{\mathsf{Closeness}(v,G)\mid(v,G)\in S\text{ and
}|\mathsf{A}(v,G)|\leq n\\}$, it holds that $|A_{n}|>n+1$. This is enough to
prove that $|\text{\rm Val}^{n}(\mathsf{Closeness})|>n+1$ since
$A_{n}\subseteq\text{\rm Val}^{n}$. Consider the set
$\Sigma_{n}=\\{\text{SD}(v,G)\mid(v,G)\in S\text{ and }|\mathsf{A}(v,G)|\leq
n\\}$. Observe that $|\text{\rm Val}^{n}(\mathsf{Closeness})|\leq n+1$ iff
$|\Sigma_{n}|\leq n+1$ because
$\mathsf{Closeness}(v,G)=\displaystyle\frac{1}{\text{SD}(v,G)}$ for every
$(v,G)\in\mathbf{VG}$. Therefore, it suffices to show that $|\Sigma_{n}|>n+1$.
The rest of the proof proceed in two main steps:

  1. (1)

For every $j>0$, it is easy to see that for the extreme of the line graph
$\text{SD}(1,L_{j})=\displaystyle\frac{j\cdot j-1}{2}$ while
$|\mathsf{A}(1,L_{j})|=j$. Therefore, it is clear that for $i,j\in[n]$ such
that $i\not=j$, $\text{SD}(1,L_{i})\not=\text{SD}(1,L_{j})$. Therefore,
$|\Sigma_{n}|>n$ if we add the isolated vertex and the pairs $\\{(1,L_{j})\mid
j\in\\{1,...,n\\}\\}$ to our set $S$.

  2. (2)

It remains to show that there exists $n>0$ for which $|\Sigma_{n}|>n+1$.
Therefore, we just need to find one pair $(v^{\prime},G^{\prime})$ such that
$|\mathsf{A}(v^{\prime},G^{\prime})|\leq n$ and
$\textsc{SD}(v^{\prime},G^{\prime})\not\in\\{\textsc{SD}(1,L_{i})\mid
i\in\\{1,...,n\\}\\}$. This is the case for $n=5$. Observe that for the center
of the star graph with $3$ nodes $(1,S_{3})$ we have that
$|\mathsf{A}(1,S_{3})|=4$, whereas $\textsc{SD}(2,S_{3})=2$. Consider now the
digraph $H=\\{\\{1,2,3,4\\},\\{(1,2),(2,3),(2,4)\\}\\}$. For the leaf node $1$
in $H$ we have $|\mathsf{A}(1,H)|=5$, whereas $\textsc{SD}(1,H)=5$. At the
same time, there is no $i\in\\{1,...,5\\}$ such that
$\textsc{SD}(1,L_{i})\in\\{2,5\\}$. This means for
$S=\\{(1,L_{1}),(1,L_{2}),(1,L_{3}),(1,L_{4}),(2,S_{3}),(1,H)\\}$, $|\text{\rm
Val}^{4}(\mathsf{Closeness})|=|\Sigma_{4}|\geq|\\{\textsc{SD}(v,G)\mid(v,G)\in
S\\}|>5$.

To show that $\mathsf{Closeness}$ does not satisfy the non-uniform bounded
value property we exploit the fact that it is local and closed under
isomorphism. Gather all the pairs of $S$ inside one disconnected graph, where
each connected component is isomorphic to one of the graphs used for the
previous proof. In such a graph, $|\text{\rm Val}_{G}^{5}|=6$, which leads to
a contradiction due to Proposition G.6.

#### Harmonic.

The harmonic centrality for directed graphs is defined as follows:

$\mathsf{Harmonic}(v,G)\ =\ \begin{cases}0&\text{ if }|K_{v}(G)|=1\\\ &\\\
\displaystyle\sum_{u\in
V(G)\setminus\\{v\\}}\frac{1}{d_{G}(v,u)}&\text{otherwise}.\end{cases}$

Notice that whenever $u\not\in K_{v}(G)$, $d_{G}(v,u)=\infty$, and thus, we
assume that $\frac{1}{d_{G}(v,u)}=0$.

###### Proposition G.15.

The following hold:

  1. (1)

$\mathsf{Harmonic}$ is not a subgraph motif measure for directed graphs.

  2. (2)

$\mathsf{Harmonic}$ is not a subgraph motif measure relative to the induced
ranking for directed graphs.

###### Proof.

Item (1). We follow the same strategy as in the proof of Proposition F.3. We
first observe that
$\mathsf{Harmonic}(1,L_{i})=\sum_{j=1}^{i-1}\displaystyle\frac{1}{j}$.
Therefore, for $i\geq 1$, we have that
$\mathsf{Harmonic}(1,L_{i+1})=\mathsf{Harmonic}(1,L_{i})+\displaystyle\frac{1}{i}$.
This means $\mathsf{Harmonic}(1,L_{i})\not=\mathsf{Harmonic}(1,L_{j})$ for
every $i,j$ with $i\neq j$. Hence, $|\text{\rm
Val}^{n}(\mathsf{Harmonic})|\geq n-1$. Since $\mathsf{Harmonic}(1,L_{1})=0$,
we get that $|\text{\rm Val}^{n}(\mathsf{Harmonic})|\geq n$. Now, let $n=6$.
Consider the following pairs: $(2,H_{1})$ and $(2,H_{2})$ where
$H_{1}=\\{\\{1,2,3\\},\\{(2,1),(2,3)\\}\\}$ and
$H_{2}=\\{\\{1,2,3,4\\},\\{(2,1),(2,3),(3,4)\\}\\}$. We first observe that,
for $n>0$,

$\mathsf{Harmonic}(2,H_{1})\ =\ 2\ >\
\sum_{i=1}^{n-1}\displaystyle\frac{1}{i}\ =\ \mathsf{Harmonic}(1,L_{n}).$

Moreover, we have that

$\mathsf{Harmonic}(2,H_{2})\ =\ 1+1+\frac{1}{2}\ =\ 2.5.$

Therefore, $|\text{\rm Val}^{6}(\mathsf{Harmonic})|\geq 7$, which means
$\mathsf{Harmonic}$ violates the bounded value property for $n=6$. Thus, by
Theorem G.4, we conclude that $\mathsf{Harmonic}$ is not a subgraph motif
measure for directed graphs. This hold even for trees since every graph used
above is a tree.

Item (2). By Proposition G.6, it suffices to show that $\mathsf{Harmonic}$
does not enjoy the non-uniform bounded value property. To this end, we follow
the same idea as in the proof of Proposition 5.10. Roughly speaking, we gather
all the graphs employed above to prove that $|\text{\rm
Val}^{6}(\mathsf{Harmonic})|\geq 7$ in one disconnected graph $G$. Then, by
exploiting the fact that $\mathsf{Harmonic}$ is closed under isomorphism (see
Definition A.1) and local (see Definition A.2), we can show that $|\text{\rm
Val}_{G}^{6}(\mathsf{Harmonic})|\geq 7$.

#### PageRank.

Given a graph $G$ with $V(G)=\\{v_{1},\ldots,v_{|V(G)|}\\}$, we define the
adjacency matrix of $G$ as the matrix $A_{G}$ such that ${A_{G}}_{ij}=1$ if
$(v_{i},v_{j})\in E(G)$. We then define P as the column-stochastic matrix such
that:

$\textbf{P}_{ij}\ =\ \displaystyle\frac{{A_{G}}_{ij}}{|N_{G}(i)|}.$

In other words, P contains the probabilities of jumping from node $i$ to node
$j$ during a random walk over $G$. Let $\mathbf{v}$ be a stochastic vector
$(\mathbf{1}^{T}\mathbf{v}=1)$, and let $0<\alpha<1$ be a teleportation
parameter. Then, the ($\mathbf{v},\alpha$)-PageRank vector is defined as the
vector $x_{\mathbf{v},\alpha}$ that solves the equation
$(I-\alpha\textbf{P})x=(1-\alpha)\mathbf{v}$. Now, given a graph $G$ with
$V(G)=\\{v_{1},\ldots,v_{|V(G)|}\\}$, the PageRank centrality of $v_{i}$ in
$G$ is defined as

$\mathsf{PageRank}(v_{i},G)=x_{0.8,\mathbf{u}}[i],$

where $\mathbf{u}=\mathbf{1}\cdot\frac{1}{|V(G)|}$. Several variations of
$\mathsf{PageRank}$ have been proposed in the literature. In particular, if we
do not force $x_{0.8,\mathbf{u}}$ to be a probability distribution, we can ask
that the sum of centrality values inside a connected component is equal to the
number of nodes in that component. Therefore, this give rise to the local
version of $\mathsf{PageRank}$, denoted $\mathsf{LPR}$. We call it local since
the values in a connected component are always the same regardless of the
graph.

###### Proposition G.16.

The following hold:

  1. (1)

$\mathsf{PageRank}$ is not a subgraph motif measure for directed graphs.

  2. (2)

$\mathsf{PageRank}$ is not a subgraph motif measure relative to induced
ranking for directed graphs.

###### Proof.

For both items, it suffices to consider the local version of
$\mathsf{PageRank}$ since $\text{\rm Rank}(\mathsf{PageRank})=\text{\rm
Rank}({\mathsf{LPR}})$.

Item (1). It can verified that $\text{\rm
Val}_{L_{1}}^{5}(\mathsf{LPR})=\\{1\\}$, $\text{\rm
Val}_{L_{2}}^{5}(\mathsf{LPR})=\\{1,1.80\\}$, $\text{\rm
Val}_{L_{3}}^{5}(\mathsf{LPR})=\\{1,1.80,2.44\\}$, $\text{\rm
Val}_{L_{4}}^{5}(\mathsf{LPR})=\\{1,1.80,2.44,2.95\\}$, $\text{\rm
Val}_{L_{5}}^{5}(\mathsf{LPR})=\\{1,1.80,2.44,2.95,3.36\\}$ and $\text{\rm
Val}_{S_{3}}^{5}(\mathsf{LPR})=\\{1,2.60\\}$. Moreover, for the digraph
$H=(\\{1,2,3,4\\},\\{(2,1),(2,3),(2,4)\\})$, we have that $\text{\rm
Val}_{H}^{5}(\mathsf{LPR})=\\{1.0,1.27\\}$ Therefore,

$\\{1,1.80,2.44,2.95,3.36,2.60,1.27\\}\ \subseteq\ \text{\rm
Val}^{5}(\mathsf{LPR}).$

Thus, $|\text{\rm Val}^{5}(\mathsf{LPR})|\geq 7$, which implies that
$\mathsf{LPR}$ does not enjoy the bounded value property. Hence, by Theorem
G.4, we conclude that $\mathsf{LPR}$ is not a subgraph motif measure. This
holds even if we focus on trees since every graph considered is a tree.

Item (2). By Proposition G.6, it suffices to show that $\mathsf{LPR}$ does not
enjoy the non-uniform bounded value property. To this end, we follow the same
idea as in the proof of Proposition 5.10. Roughly speaking, we gather all the
graphs employed above to prove that $|\text{\rm Val}^{5}(\mathsf{LPR})|>7$ in
one disconnected graph $G$. Then, by exploiting the fact that $\mathsf{LPR}$
is closed under isomorphism (see Definition G.1) and local (see Definition
G.2), we can show that $|\text{\rm Val}_{G}^{5}(\mathsf{LPR})|\geq 7$, and the
claim follows.

#### EigenVector.

Given a graph $G$ with $V(G)=\\{u_{1},...,u_{|V(G)|}\\}$, we define the
adjacency matrix of $G$ as the matrix $A_{G}$ such that ${A_{G}}_{ij}=1$ if
$(u_{i},u_{j})\in E(G)$. Let $\lambda_{\max}(G)$ be the greatest eigenvalue of
$A_{G}$. We define the unique $L^{2}$ normalized eigenvector associated with
$\lambda_{\max}(G)$ as $\mathbf{v}_{\max}$. The eigenvector centrality
(Bonacich, 1987) of $u_{i}$ in $G$ is the $i$-th entrance of
$\mathbf{v}_{\max}$:

$\mathsf{Eigenvector}(u_{i},G)\ =\mathbf{v}_{\max}^{i}.$

###### Proposition G.17.

$\mathsf{Eigenvector}$ is not a subgraph motif measure for directed graphs.

###### Proof.

We first prove the following: given a disconnected graph $G$ with
$n=|\text{\rm Comp}(G)|$, for each component $S$ of the set
$M=\\{S\in\text{\rm Comp}(G)\mid\lambda_{\max}(G[V(S)])=\lambda_{\max}(G)\\}$,
we have that $|\\{\mathsf{Eigenvector}(v,G)\mid v\in S\\}|=|\text{\rm
Val}_{G[V(S)]}(\mathsf{Eigenvector})|.$ In other words, each component with
the same maximum eigenvalue as $G$ will have the same amount of values as the
component seen as a connected graph. On the other hand, we know that
$\mathsf{Eigenvector}$ must be a unitary vector. Therefore, if $G$ has
multiple connected components with the maximum eigenvalue $\lambda_{\max}$, in
order to give values to every vertex in those components, the vector is forced
to give different values than those given to the component seen as a connected
graph.

Given a graph $G$, we can order its vertices $V(G)=\\{v_{1},\ldots,v_{n}\\}$
in such way that, if $\text{\rm Comp}(G)=\\{S_{1},\ldots,S_{k}\\}$, then, for
every $m\in\\{1,\ldots,k\\}$, there exists
$I_{m}=\\{i_{m},i_{m}+1,i_{m}+2,\ldots,i_{m}+|S_{m}|-1\\}\subseteq\\{1,\ldots,n\\}$
such that $\\{v_{s}\\}_{s\in I_{m}}=S_{m}$. Therefore, the rows and columns of
the adjacency matrix of $G$, $A_{G}$, are ordered such that there are square
sections in the diagonal in $A_{G}$ where the $i$-th square section is the
adjacency matrix of $G[V(S_{i})]$, $A_{G[V(S_{i})]}$. Out of these square
sections the matrix $A_{G}$ just have values $0$ since vertices in different
components have no connections between them. Now, for each $S_{m}\in M$, we
know that if $b_{m}^{\max}$ is the eigenvector associated with
$\lambda_{\max}$ in $G[V(S_{m})]$, then the vector
$v_{m}^{\star}=(0,0,\ldots,b_{m}^{\max}(1),\ldots,b_{m}^{\max}(|S_{m}|),0,\ldots,0)$
is a an eigenvector with eigenvalue $\lambda_{\max}$ for $A_{G}$ such that
$v_{m}^{\star}(i_{m}+j)=b_{m}^{\max}(j-1)$ for $j\in\\{0,\ldots,|S_{m}|-1\\}$.
This implies that any linear combination of the vectors $v_{m}^{\star}$ will
by an eigenvector for $A_{G}$ with eigenvalue $\lambda_{\max}$. Finally, we
have that the unique vector
$v^{\star}=c_{1}v_{1}^{\star}+c_{2}v_{2}+\ldots+v_{|M|}^{\star}$ such that
$\|v^{\star}\|=1$ is the eigenvector for $\mathsf{Eigenvector}$ centrality in
$G$. In other words, for each node $v_{i}\in V(G)$, then
$\mathsf{Eigenvector}(v,G)=v^{\star}(i)$. As we showed, the vector $v^{\star}$
assign values for every connected component $S_{m}\in M$ which are
proportional to $b_{m}^{\max}$. Therefore, we can conclude that for every $n$,
$|\text{\rm Val}_{G}^{n}(\mathsf{Eigenvector})|\geq|\bigcup_{S\in M}\text{\rm
Val}_{G[V(S)]}^{n}(\mathsf{Eigenvector})|$. Now, we will construct a sequence
of graphs $G_{i}$ in order to prove that $|\text{\rm
Val}^{4}(\mathsf{Eigenvector})|>5$. For $i>1$, define the digraph
$G_{i}=(\\{0,1,\ldots,i\\},\\{(i,i+1),(i+1,i)\mid i\mod 2=0\\})$. For $i\geq
2$, we have that $G_{i}$ contains $\lceil i/2\rceil$ connected components.
Each one is isomorphic to cycle with two nodes $C_{2}$ whenever $i\mod 2=0$.
Since the cycle with two nodes has maximum eigenvalue $\lambda_{\max}=1$ the
eigenvector $v^{\star}=(1/\sqrt{i},1/\sqrt{i},\ldots,1/\sqrt{i})$ is the
eigenvector centrality vector for $G_{i}$. This fact implies that $|\text{\rm
Val}^{4}(\mathsf{Eigenvector})|\geq|\bigcup_{j=1}^{i}\text{\rm
Val}_{G_{i}}^{n}(\mathsf{Eigenvector})|=i$ for every $i>0$. In other words,
the set $\text{\rm Val}^{4}(\mathsf{Eigenvector})$ has infinitely many values
which contradicts the bounded value property. Thus, by Theorem G.4, we
conclude that $\mathsf{Eigenvector}$ is not a subgraph motif measure for
directed graphs.

By Theorem G.9, $\mathsf{Eigenvector}$ is a subgraph motif measure relative to
the induced ranking over connected graphs. We proceed to strengthen this
result by showing that it holds over a broader class of graphs than the class
of connected graphs. Let

$\mathbf{K}_{1}\ =\ \\{G\in\mathbf{G}\mid\text{for every
}G_{1},G_{2}\in\text{\rm
Comp}(G),\lambda_{\max}(G[V(G_{1})])=\lambda_{\max}(G[V(G_{2})])=\lambda_{\max}(G)\text{
implies }G_{1}\cong G_{2}\\}.$

It is easy to verify that every connected graph belongs to $\mathbf{K}_{1}$,
and that $\mathbf{K}_{1}$ may contain disconnected graphs. We show that:

###### Proposition G.18.

$\mathsf{Eigenvector}$ is a subgraph motif measure relative to the induced
ranking when restricted to $\mathbf{K_{1}}$.

###### Proof.

From the proof of Proposition G.17, we know that for every component of
$M(G)=\\{S\in\text{\rm
Comp}(G)\mid\lambda_{\max}(G[V(S)])=\lambda_{\max}(G)\\}$, the eigenvector
centrality assigns the same amount of values to the connected component as it
seen as a connected graph. On the other hand, if $S\not\in M(G)$, then, for
every vertex $v\in V(S)$, we have that $\mathsf{Eigenvector}(v,G)=0$. Note
that, whenever two graphs are isomorphic they have the same values in every
eigenvector. Values are just permuted depending on the isomorphism you map
each graph to the other. Therefore, we have that, for every graph
$G\in\mathbf{K}_{1}$, there exists $S\in M(G)$ such that $\text{\rm
Val}_{G}^{n}(\mathsf{Eigenvector})=\text{\rm
Val}_{G[V(S)]}^{n}(\mathsf{Eigenvector})$. In other words, when
$\mathsf{Eigenvector}$ is restricted to $\mathbf{K}_{1}$ it behaves like
$\mathsf{Con}\mathsf{Eigenvector}$. Hence, we can use the same subgraph family
and filtering function employed in the proof of Theorem G.9 with the addition
of value $0$ for vertices outside $S$ in $M(G)$.

#### Betweenness.

Betweenness centrality for $v$ in $G$ is defined as

$\mathsf{Betweenness}(v,G)\ =\
\displaystyle\sum_{(u,w)\in(V(G)/\\{v\\})^{2}}\displaystyle\frac{|S_{G}^{v}(u,w)|}{|S_{G}(u,w)|}.$

If for some pair $(u,w)$ the set $|S_{G}(u,w)|=0$, then
$\displaystyle\frac{|S_{G}^{v}(u,w)|}{|S_{G}(u,w)|}=0$ by definition.

###### Proposition G.19.

$\mathsf{Betweenness}$ is a subgraph motif measure for directed graphs when
restricted to trees.

###### Proof.

Proof is analogous to the one of Proposition F.7. For any tree $T$, and any
pair of nodes $u\not=v$ in $V(T)$, $|S_{T}(u,v)|\leq 1$. Therefore, we can use
the family of subgraphs defined in Proposition G.10 and the identity filtering
function over $\mathbb{N}$ to prove that $\mathsf{Betweenness}$ is a subgraph
motif measure when restricted to trees.

###### Proposition G.20.

The following hold:

  1. (1)

$\mathsf{Betweenness}$ is not a monotonic subgraph motif measure for directed
graphs, even if we focus on connected graphs.

  2. (2)

$\mathsf{Betweenness}$ is not a monotonic subgraph motif measure relative to
the induced ranking for directed graphs.

###### Proof.

Item (1). We are going to show that $\mathsf{Betweenness}$ does not enjoy the
monotonic bounded value property, and thus, by Theorem G.7, it is not a
monotonic subgraph motif measure. For $n\geq 1$, consider the connected graph
$H_{n}=(V_{n},E_{n})$, where

$V_{n}\ =\ \\{1,2,...,n+2\\}\qquad\text{and}\qquad E_{n}=\\{(1,i)\mid
i\in\\{3,4,\ldots,n+2\\}\\}\cup\\{(2,i)\mid i\in\\{3,4,.\ldots,n+2\\}\\}.$

It can be verified that $\mathsf{Betweenness}(3,H_{n})=\frac{2}{n}$. However,
for $n>1$, $|\mathsf{A}(3,H_{n})|>4$. Therefore, we have that

$|\text{\rm Val}^{4}(\mathsf{Betweenness})|+|\text{\rm
BVal}^{4}(\mathsf{Betweenness})|\ >\ 5,$

which implies that $\mathsf{Betweenness}$ does not enjoy the monotonic bounded
value property, as needed.

Item(2). By Theorem G.8, it suffices to show that $\mathsf{Betweenness}$ does
not enjoy the non-uniform bounded value property. Since $\mathsf{Betweenness}$
is closed under isomorphism and local, we will use the same strategy as we did
in analogous proofs. In particular, we define the graph $G$ as the minimal
graph with the following property: for every graph
$G^{\prime}\in\\{H_{1},H_{2},\ldots,H_{5},H_{6}\\}$, there exists a component
$G^{\prime\prime}\in\text{\rm Comp}(G)$ such that $G^{\prime\prime}$ is
isomorphic to $G^{\prime}$. We can then show that

$|\text{\rm Val}_{G}^{4}(\mathsf{Betweenness})|+|\text{\rm
BVal}_{G}^{4}(\mathsf{Betweenness})|\ >\ 5,$

which implies that $\mathsf{Betweenness}$ does not enjoy the non-uniform
monotonic bounded value property, and the claim follows.