value, which includes star forming regions that are not directly involved in
the Ly$\alpha$ radiative transfer. Kimm et al. (2019) perform radiation-
hydrodynamic simulations of Ly$\alpha$ escape from star-forming regions and
find that $f_{\mathrm{esc}}$ increases with stellar age up to $\sim 10$ Myr.
Given that our total apertures contain multiple stellar populations with
different ages and Ly$\alpha$ output, that we do not see this effect in our
sample is expected. The overall lack of correlations with the SFR quantities
could be explained by high SFR being linked to other properties that makes
escape harder. For example, at low redshift the highest (extinction corrected)
SFR is often found in the most dusty galaxies, and while similar systems are
hard to find at high redshift they do exist (e.g., Casey et al., 2014), and
many more are waiting to be detected by JWST. These dusty galaxies are also
more likely to have higher global neutral gas masses. Hence, the increased
mechanical energy injected into the ISM from high SFR may be combined with
higher dust and gas content which will limit the escape of Ly$\alpha$ photons.
Even though mechanical feedback from star formation will increase the escape
probability on scales of H ii regions, this is not seen in global measurements
for the sample.

Both stellar mass and nebular metallicity are quantities that are linked to
each other and to other galaxy properties that may affect Ly$\alpha$ escape,
most notably H i mass (for the stellar mass), and dust extinction (for
metallicity). In the top middle and right panels of Figure 13 we show scatter
plots of $f_{\mathrm{esc}}$ versus these two properties. For stellar mass we
find a weak anti-correlation ($\\{\tau,p_{0}\\}=\\{-0.23,0.029\\}$) with
substantial scatter. The anti-correlation is mainly driven by the lack of
massive galaxies with high $f_{\mathrm{esc}}$. From this data set it seems
clear that Ly$\alpha$ escape is low in galaxies with the highest stellar mass,
perhaps reflecting that these galaxies have larger neutral gas reservoirs
(Pardy et al., 2014). The higher neutral gas content will in turn cause longer
scattering path lengths for the Ly$\alpha$ radiative transfer, and could also
make it harder for the star formation feedback to create clear channels for
escape. For low stellar masses the scatter in $f_{\mathrm{esc}}$ is large,
this could be explained by variations in the neutral gas fraction, e.g., some
of these galaxies may have small stellar masses compared to the neutral gas
mass. In these galaxies $f_{\mathrm{esc}}$ will be low because of high H i
column densities, while low mass galaxies that have used up a greater portion
of their neutral gas will have higher escape fractions. Lower mass galaxies
also have weaker gravitational potential which will make it easier for
feedback to create low column density channels for Ly$\alpha$ escape (Laursen
et al., 2009). This effect would make $f_{\mathrm{esc}}$ higher in low-mass
galaxies on average, but would also show considerable scatter due to
geometrical effects. A more detailed analysis of H i column densities and
kinematics for this sample, and how they relate to Ly$\alpha$ escape will be
given in a forthcoming paper (Le Reste et al., in prep.). The comparison with
nebular metallicity does not show anything conclusive, but we note that there
is a lack of data points in the low metallicity range.

The effect of ionisation state of the main star-forming regions in the
galaxies (as measured by the SDSS O32 ratio) on global $f_{\mathrm{esc}}$ is
shown in the bottom row left panel. We do not find any correlation of O32 with
$f_{\mathrm{esc}}$, but note that there are very few galaxies with low O32
(indicating that the ISM is not strongly ionized) that show strong Ly$\alpha$
escape. A highly ionized ISM in the regions where Ly$\alpha$ is produced thus
seems to be a necessary condition for high ($\gtrsim 0.1$) $f_{\mathrm{esc}}$.
This result is consistent with earlier studies that also show a large scatter
of $f_{\mathrm{esc}}$ at low O32 ratios (e.g., Henry et al., 2015; Jaskot et
al., 2017; Yang et al., 2017b; Izotov et al., 2020), but we do not find any
strong correlation. However, it should be noted that the sample selection is
quite different for some of these studies in that the galaxies are all
Ly$\alpha$ emitters and are selected to have strong emission. Our sample
contains a number of galaxies with moderate to strong Ly$\alpha$ absorption
which will increase the scatter substantially which is also seen in the study
by Izotov et al. (2020).

In the bottom row middle and right panels we show scatter plots of
$f_{\mathrm{esc}}$ versus nebular and stellar extinction. Both of these are
strongly anticorrelated with the escape fraction,
$\\{\tau,p_{0}\\}=\\{-0.37,2.9\cdot 10^{-4}\\}$ for $E(B-V)_{n}$, and
$\\{\tau,p_{0}\\}=\\{-0.31,3.00\cdot 10^{-3}\\}$ for $E(B-V)_{s}$. The
correlation with nebular extinction was previously found by Runnholm et al.
(2020) in a rank correlation study, where $f_{\mathrm{esc}}$ was found to have
a strong dependence on $E(B-V)_{n}$. We will come back to the details of
nebular extinction and escape fraction in Section 6.1. Given that nebular and
stellar extinction are strongly correlated it is not surprising to find that
the escape fraction is also anti-correlated to the stellar extinction.

Figure 14: Scatter plots of $f_{\mathrm{esc}}$ versus petrosian radii in FUV
(left panel) and Ly$\alpha$ (right panel). Note that the galaxy id labels
refer to the left panel, therefore some of the galaxies with limits in the
right panel are not shown as such in the legend. Please refer to Table 5 for
individual galaxy limits. Figure 15: Left panel: $f_{\mathrm{esc}}$ versus the
$\xi_{\mathrm{Ly}\alpha}$ parameter. Right panel: $\xi_{\mathrm{Ly}\alpha}$
versus nebular extinction. Galaxies with lower limits on $r_{\mathrm{p,FUV}}$
(ELARS15 and ELARS22) are excluded from these plots.

###  5.3 The Ly$\alpha$ spatial distribution and escape fractions

As previously mentioned, a more thorough analysis and discussion of how the
spatial distribution of Ly$\alpha$ relates to global observables for the
sample in this study is given in Rasekh et al. (2022). In this paper we
provide two simple observational diagnostics for the Ly$\alpha$ spatial
extent: $r_{\mathrm{p,Ly}\alpha}$, and $\xi_{\mathrm{Ly}\alpha}$ (defined
below). We also investigate the effect of FUV compactness (as measured by
$r_{p}(FUV)$) on the escape fraction.

In Figure 14 we show the dependence of $f_{\mathrm{esc}}$ on the compactness
(as measured by $r_{\mathrm{p}}$) of FUV and Ly$\alpha$ emission. There are no
correlations found for any of these properties, but we note that high escape
fractions are much less common in galaxies with large $r_{\mathrm{p,FUV}}$.
Galaxies with very extended FUV emission thus seem to be less likely to have
high $f_{\mathrm{esc}}$.

We define the Ly$\alpha$ overextension parameter, $\xi_{\mathrm{Ly}\alpha}$,
as the ratio between the Ly$\alpha$ and FUV continuum petrosian radii:

$\xi_{\mathrm{Ly}\alpha}=\frac{r_{\mathrm{p,Ly}\alpha}}{r_{\mathrm{p,FUV}}}.$
(4)

This ratio is high for galaxies which have extensive Ly$\alpha$ halos caused
by radiative transfer of the Ly$\alpha$ photons far away from where they were
produced. The parameter is also high for galaxies with substantial continuum
absorption in the central regions. Lower $\xi_{Ly\alpha}$ values are found in
galaxies in which the Ly$\alpha$ escape happens along the line of sight or in
which resonant scattering is less significant. In Figure 15 we show how
$\xi_{\mathrm{Ly}\alpha}$ relates to $f_{\mathrm{esc}}$ and $E(B-V)_{n}$ for
our sample. We do not find any correlation between $f_{\mathrm{esc}}$ and
$\xi_{\mathrm{Ly}\alpha}$. Naively one might expect higher Ly$\alpha$ escape
for galaxies where the photons escape directly with little or no radiative
transfer (giving low $\xi_{\mathrm{Ly}\alpha}$), but there are likely a number
of competing effects which cause this prediction to fail. For example, in
galaxies which have high escape fractions due to low dust content photons tend
to travel longer paths before being destroyed, which naturally gives rise to
larger Ly$\alpha$ halos. Hayes et al. (2014) found a anti-correlation between
$\xi_{\mathrm{Ly}\alpha}$ and $E(B-V)_{n}$ for the original LARS galaxies,
however, this result is not confirmed with the updated calibration and when
including the 31 additional galaxies in this sample.

##  6 Discussion

Figure 16: Global escape fraction of Ly$\alpha$ versus nebular extinction. The
gray lines show where the points would fall for different ISM extinction
curves. The curve labeled CCM is the galactic extinction curve from Cardelli
et al. (1989) and SMC is from Prevot et al. (1984). Coloured lines show curves
derived from fitting a line to observation from this and other Ly$\alpha$
studies. Note that that Hayes et al. (2011) and Matthee et al. (2016) use
stellar extinction, for this comparison we multiply their $E(B-V)$ values by a
factor 2.27 (Calzetti et al., 2000). The dotted gray lines (labeled CCM
(clumps)) show clumpy dust models with varying numbers of clumps along the
line of sight.

###  6.1 Ly$\alpha$ escape and dust extinction

As shown in Section 5.2 we find a strong anti-correlation between
$f_{\mathrm{esc}}$ and $E(B-V)_{n}$. In Figure 16 we show this relation
together with a number of extinction curves from the literature. When
comparing to dust screen models (CCM and SMC) it is clear that Ly$\alpha$
escape cannot be described by these simple models, both in terms of the
average slope of the curve and in terms of scatter. Similar to Hayes et al.
(2011) we make the assumption that $f_{\mathrm{esc}}$ at $E(B-V)_{n}\approx 0$
can be $<1$ (even very small amounts of dust will affect Ly$\alpha$ photons
due to the scattering), and fit a function of the form:

$f_{\mathrm{esc}}=C\cdot 10^{-0.4k_{\mathrm{Ly}\alpha}E(B-V)_{n}}.$ (5)

We then fit a linear function to $\log$ $f_{\mathrm{esc}}$ (using a censored
bayesian linear regression method, LINMIX
333https://github.com/jmeyers314/linmix and Kelly (2007)) which yields the
parameter values $C=0.075\pm 0.019$ and $k_{\mathrm{Ly}\alpha}=3.30\pm 0.73$
where the slope and intercept uncertainties are derived from the posterior
probabilities. The best fit slope is significantly flatter than the CCM and
SMC curves.

There are multiple causes for the scatter in the relation that can make the
Ly$\alpha$ escape fraction move both up and down (see also Hayes et al.,
2011).

  * •

Resonant scattering of the Ly$\alpha$ photons in the neutral hydrogen gas
cause the effective path lengths (or optical depth) to increase. This means
that more photons are absorbed by dust and $f_{\mathrm{esc}}$ consequently
goes down (e.g., Dijkstra, 2017). For this data set the effect is clearly seen
at low nebular extinction. Even with minimal dust content the effective
optical depth at the Ly$\alpha$ wavelength becomes large and the escape
fraction is strongly decreased compared to the prediction of the dust screen
model.

  * •

In a clumpy ISM (where most of the dust is confined within the dense clumps of
H i gas), the Ly$\alpha$ photons can scatter on the surface of the clumps and
thus avoid entering the regions where the dust content is highest. In this
situation Ly$\alpha$ photons may suffer less from extinction than both the
stellar continuum and the optical nebular line emission (an effect found in
the theoretical work of Neufeld, 1991), causing the escape fraction to be
higher than the dust screen models would imply. However, if this effect was
present for a majority of the galaxies we would also find a positive
correlation of $W_{\mathrm{Ly\alpha}}$ and $E(B-V)_{n}$ (Scarlata et al.,
2009), which we do not find (in fact, they are anti-correlated to the same
degree that $f_{\mathrm{esc}}$ and $E(B-V)_{n}$ is). More recent radiative
transfer simulations have also found that this type of scattering does not
have a large effect on the output Ly$\alpha$ emission (e.g., Laursen et al.,
2013; Duval et al., 2016). It is not excluded that $f_{\mathrm{esc}}$ could be
boosted by this effect in some of the sources, but it does not seem to be in
play for the sample as a whole.

  * •

Pure dust geometrical effects can also cause $f_{\mathrm{esc}}$ to be higher
than for the uniform dust screen models. In Figure 16 we show the effect of a
clumpy dust screen (Natta & Panagia, 1984; Calzetti et al., 1994; Scarlata et
al., 2009) on the Ly$\alpha$ escape fraction. As the dust content inside the
clumps increases the $f_{\mathrm{esc}}$-$E(B-V)_{n}$ relation will move along
these locii (N is the average number of clumps along the line of sight). The
extinction will then become more and more gray — as the individual clump
extinction becomes high the attenuation will tend toward a constant
$e^{-\mathrm{N}}$. Many of the galaxies in our sample (in particular the ones
with high extinction) show escape fractions that are consistent with clumpy
dust geometry.

The overall correlation of extinction and $f_{\mathrm{esc}}$ was found already
in the original LARS publications (Hayes et al., 2014), and also by other
theoretical and observational studies at low and high redshift (Verhamme et
al., 2008; Atek et al., 2009; Scarlata et al., 2009; Kornei et al., 2010;
Hayes et al., 2011; Matthee et al., 2016). In order to make the comparison to
the studies using stellar extinction we use the conversion
$E(B-V)_{n}=2.27E(B-V)_{s}$ (Calzetti et al., 2000), and convert the stellar
extinctions given in the other publications to nebular. This ratio of stellar
to nebular extinction is explained by nebular emission originating from
regions within galaxies with higher density ISM, but was derived for a small
sample of nearby galaxies. The number is thus not necessarily applicable to
the high redshift samples compared with here, but a similar ratio ($2.077\pm
0.088$) was found by Reddy et al. (2020) for sample of $z\sim 2$ galaxies.

Similar to this study these authors find that if Ly$\alpha$ escape were to be
explained by extinction alone, the nebular attenuation curve needs to be
considerably more shallow than the standard dust screen models. In this paper
we include galaxies with substantial extinction, and find a shallower slope
than all but one of the earlier studies. The slope is remarkably similar to
the one found by Matthee et al. (2016), who uses a sample of 488
H$\alpha$-selected galaxies at $z\sim 2$. The sample in their study includes
galaxies with $E(B-V)_{s}\lesssim 0.5$ which, when applying the stellar to
nebular extinction conversion, gives a similar range of extinctions as in this
data set.

The shallow effective attenuation curve for Ly$\alpha$ thus seems to be caused
by a combination of radiative scatter and dust geometry effects. At low
extinction, resonant scattering causes the optical depth for Ly$\alpha$
photons to become longer and $f_{\mathrm{esc}}$ consequently is lowered. The
effect of a clumpy dust geometry on the effective attenuation is small at low
dust content, and the average $f_{\mathrm{esc}}$ thus becomes lower than a
simple dust screen would predict ($0.075$ at $E(B-V)_{n}=0$ for the fit
above). At higher dust content resonant scattering seems to become less
important, and the effects of a clumpy dust geometry results in higher
$f_{\mathrm{esc}}$ than predicted. One interpretation of this result is that
the Ly$\alpha$ radiative transfer in sufficiently dusty galaxies starts to be
less dominated by resonant scattering effects. Instead it is the dust content
and clump distribution that has the largest effect on the escape fraction.
Remarkably, none of the galaxies with $E(B-V)_{n}>0.6$ sits below the dust
screen models in the figure, indicating that optical depth effects from
resonant scattering are less prominent in this regime. The caveat to this
conclusion is that resonant scattering in a clumpy H i +dust medium can also
lead to boosted Ly$\alpha$ escape (Neufeld, 1991). However, as discussed
above, this effect is likely not enough to cause the observed escape
fractions.

While a general prediction of the $f_{\mathrm{esc}}$ evolution with redshift
based on the overall changes to extinction over cosmic time (Hayes et al.,
2011) is outside the scope of this paper we note that the significantly lower
slope found in this study means that $f_{\mathrm{esc}}$ would evolve less with
redshift compared to the earlier study.

###  6.2 Ly$\alpha$ luminosity as a diagnostic of star formation rate

Figure 17: Comparison of the $f_{\mathrm{esc}}$-corrected Ly$\alpha$ SFR
diagnostic with SFR derived from H$\alpha$. The color scale of the circles
indicate the $E(B-V)_{n}$ of the galaxies. Galaxies with a zero or negative
Ly$\alpha$ luminosity are excluded from this plot. The dotted line shows the
1-to-1 relation. The upper plot shows the ratio
SFR${}_{\mathrm{Ly}\alpha}/$SFRHα. The light green band shows the standard
deviation of the ratios for the full sample, while the darker green band shows
the standard deviation for a sub-sample with $E(B-V)_{n}<0.3$.

The use of Ly$\alpha$ luminosity to measure star formation rates in galaxies
is problematic because of the uncertainties in escape fraction. In contrast to
other nebular hydrogen lines (e.g., H$\alpha$ or H$\beta$) a simple dust
correction is not enough to get measurements of similar precision.
Nevertheless earlier studies showed correlations of Ly$\alpha$ luminosity with
FUV (e.g., Wisotzki et al., 2016) or H$\alpha$ luminosity (e.g., Hayes et al.,
2010), albeit with large scatter. Sobral & Matthee (2019) used a sample of 218
Ly$\alpha$ +H$\alpha$ emitters at $z\lesssim 0.3$ and $z\sim$2.2–2.6 to derive
a SFR calibration valid for Ly$\alpha$ luminosity. They find a fairly tight
linear relation betwen $f_{\mathrm{esc}}$ and $W_{\mathrm{Ly\alpha}}$ and use
this to correct the observed Ly$\alpha$ luminosities which can then be used to
define a Ly$\alpha$-based SFR calibration. This calibration does not require
any knowledge of the extinction or escape fraction, but can be used with
$W_{\mathrm{Ly\alpha}}$ and luminosity alone.

In Figure 11 we show a linear fit of $f_{\mathrm{esc}}$ as a function of
$W_{\mathrm{Ly\alpha}}$ (again using LINMIX referenced above):

$f_{\mathrm{esc}}=\left(0.0040\pm 0.00050\cdot W_{\mathrm{Ly\alpha}}\right)\pm
0.013,$ (6)

where the error to the right is the uncertainty of the fitted intercept (the
best-fit intercept value itself is negligible compared to the uncertainty).
While the fit has a quite small statistical uncertainty ($\sim$5–10% for the
range of $W_{\mathrm{Ly\alpha}}$ in this sample), there are also systematic
effects that come into play (see Section 5.1 for a discussion of these), in
particular the effect of extinction. The found relation is significantly
shallower than the one found by Sobral & Matthee (2019), even when comparing
to their low redshift fit. As discussed by these authors, the relation of
$f_{\mathrm{esc}}$ versus $W_{\mathrm{Ly\alpha}}$ is affected strongly by the
dust content of the galaxies, and the difference in slope we find could be due
to the inclusion of more dusty galaxies in our sample (or other selection
effects).

To find an expression for SFR as a function of Ly$\alpha$ luminosity we
rewrite Equation 3 and use the H$\alpha$ SFR calibration of Kennicutt & Evans
(2012). We then replace $f_{\mathrm{esc}}$ with the relation found above and
arrive at:

$\mathrm{SFR}_{\mathrm{Ly}\alpha}=\frac{C_{\mathrm{H}\alpha}^{-1}}{8.7f_{\mathrm{esc}}}L_{\mathrm{Ly}\alpha}=\frac{5.37\cdot
10^{-42}}{0.035\cdot W_{\mathrm{Ly\alpha}}}L_{\mathrm{Ly}\alpha},$ (7)

where $C_{\mathrm{H}\alpha}$ is the SFR calibration constant from Kennicutt &
Evans (2012). Note that this calibration is valid for extinction corrected
H$\alpha$, and SFRLyα thus also provides a dust corrected estimate of the star
formation rate. The uncertainty of the $f_{\mathrm{esc}}$
($W_{\mathrm{Ly\alpha}}$) relation translates into an uncertainty on the SFR
of $\sim$5–10%. In Figure 17 we show how the Ly$\alpha$-based SFR calibration
compares to the H$\alpha$ SFR for our sample. Overall the Ly$\alpha$-based
rate follows the H$\alpha$ calibration well, but the scatter is large. The
average relative residual for the full sample is 5%, but with a scatter of 59%
(0.20 dex). The standard error of mean ($\sigma/\sqrt{N}$) for the calibration
is 9%. The high scatter seems to be driven by the systemic uncertainties of
the $W_{\mathrm{Ly\alpha}}$ to $f_{\mathrm{esc}}$ conversion, and in
particular the high extinction galaxies that have severely underestimated SFR.
If we instead consider only galaxies with low extinction ($E(B-V)_{n}<0.3$)
the average relative residual is $+30$%, the scatter is reduced to 46% (0.16
dex), and the standard error of mean stays the same. The positive offset found
in this case may reflect systematic effects on the $f_{\mathrm{esc}}$ relation
resulting from variations in the ionizing radiation production efficiency.
Galaxies with less bursty star formation (older average stellar ages) will
tend to have over-estimated escape fractions from the fit and will thus get an
overestimated SFR from the diagnostic. The $f_{\mathrm{esc}}$ to
$W_{\mathrm{Ly\alpha}}$ relation stays the same (within the uncertainties)
when excluding high extinction galaxies.

Another systematic uncertainty in the calibration is the assumption on the
value of $C_{\mathrm{H}\alpha}$, which itself depends on the age of the
stellar population generating the ionizing photons (and the ionized ISM
temperature). If the age is significantly lower than the assumption made by
Kennicutt & Evans (2012), with massive star formation already in equilibrium,
the SFR could be overestimated by a factor $\sim 2$ (for ionizing clusters
with age $\leq 3-4$ Myr). Note that if the
$f_{\mathrm{esc}}$–$W_{\mathrm{Ly\alpha}}$ relation is inverted and
$W_{\mathrm{Ly\alpha}}$ at 100% escape calculated, the intrinsic equivalent
widths would fall around 200 – 400 Å. These values are typical of very young
bursts of massive star formation, which are the ones probably driving the
observed correlation. In any case, we want to stress that this does not affect
the comparison shown in Figure 17, given that the SFR indicators there use the
same constant. For simplicity and re-usability we choose to use the standard
calibration in this work.

The scatter found from comparing to SFRHα is comparabe to the one reported by
Sobral & Matthee (2019), who, comparing to dust corrected UV-based SFR find a
scatter of 0.2–0.3 dex. It is also similar to the typical scatter found when
comparing extinction corrected H$\alpha$ SFR to FIR SFR, $\sim 0.3$ dex (e.g.,
Domínguez Sánchez et al., 2012). By applying a cut on the extinction we find a
scatter which is somewhat lower. Variations in star formation histories of the
galaxies (giving rise to variations in the ionizing radiation production
efficiency) are likely responsible for part of the remaining scatter. The
Ly$\alpha$ based SFR calibration presented here can thus provide estimates for
SFR in high redshift LAEs which do not have any optical line measurements.
These estimates will be better if the galaxies have low extinction and similar
star formation history as our sample. Without access to the optical lines it
is hard to constrain the extinction, but a rough selection could be done using
rest-frame UV stellar continuum slopes or SED fitting.

###  6.3 The effect of aperture size on measurements of $f_{\mathrm{esc}}$ and
Ly$\alpha$ luminosity

In Section 4.3 we show and describe how we use stacked growth curves to
estimate the aperture losses of Ly$\alpha$ in the sample and how this
translates to some other observational cases. The quite severe losses
estimated for our sample do not change the results found for
$f_{\mathrm{esc}}$, but the SFR calibration comparison discussed in the
previous section could be affected. However, given the large scatter for this
relation the overall effect is too small to make a difference to the
comparison.

Ly$\alpha$ emitters at $z\sim 0.3$ is an interesting case because it is
possible to observe Lyman continuum (LyC) emission with HST/COS at this
redshift. A number of studies targeting galaxies at this redshift have
detected LyC emission, and also observe the Ly$\alpha$ spectrum with COS
(e.g., Izotov et al., 2018; Flury et al., 2022). Escape of ionizing photons
have been found to correlate with $f_{\mathrm{esc}}$, but if there is a
significant amount of Ly$\alpha$ flux missing from these measurement these
conclusions could be misleading. In this study we find that $27\pm 23$% of the
Ly$\alpha$ flux is outside of the COS aperture. An average offset of this
magnitude for the individual galaxies in these studies will not change the
found correlations, but if there are other dependencies on the amount of
Ly$\alpha$ that is emitted from the outskirts of galaxies there may be
systematics of the relations that are poorly constrained. In addition, Rasekh
et el. (in prep.) finds that the fraction of Ly$\alpha$ emitted far outside
the star-forming regions for Green Peas seems to be considerably smaller for
Lyman continuum leaking galaxies (see also Verhamme et al., 2015). The
aperture losses found for the sample in this paper may thus not be applicable
to these galaxies and are probably lower than the number reported here.

At high redshift ($z\sim 3$) the aperture losses are consistent with zero for
a typical slit width of 1″, and measurements using this setup should provide
good estimates of the total Ly$\alpha$ flux. However, at extreme redshifts and
when using narrow slits (e.g., using the MOS mode with JWST/NIRSPEC and a
single slitlet width) the slit loss is very high. In addition, at high
redshifts, much more of the circumgalactic medium and ISM inside galaxies is
expected to be neutral. Radiative transfer simulations show that this can lead
to more radiative transfer of Ly$\alpha$ photons and thus larger spatial sizes
of Ly$\alpha$ halos (Laursen et al., 2019). Our estimated aperture losses
should thus be treated as lower limits in this case.

When using the estimated average aperture losses to correct $f_{\mathrm{esc}}$
is should be noted that the values presented here require that all of the
H$\alpha$ emission is measured (i.e. the H$\alpha$ aperture covers all of the
star forming regions in the galaxy). In a matched aperture COS (FUV) and SDSS
(optical) setup observing nearby galaxies also H$\alpha$ will be affected by
aperture losses and the simple correction suggested here suffer from
additional systematic uncertainty. This could also be true for the NIRSPEC
setup pictured above where the very narrow slit width will not cover the full
star-forming extent of the galaxies. In the intermediate redshift range all of
the H$\alpha$ emission will be easily measured and our loss values should
hold.

##  7 Summary

In this paper we present Ly$\alpha$ imaging of a sample of 45 star-forming
galaxies at low redshift. We describe the sample selection and data reduction,
including specialized tools to match the FUV PSF to the optical data. We also
provide a detailed description of the pixel SED-fitting code we use to obtain
secure Ly$\alpha$, H$\alpha$, and stellar population property maps, as well as
error maps for these quantities, for the galaxies.

We then perform circular photometry on the resulting maps to obtain global
measurements of fluxes (Ly$\alpha$, H$\beta$, H$\alpha$, and FUV stellar
continuum) and stellar properties (mass, luminosity-weighted stellar
population age, and stellar dust extinction). The measurements are also used
to derive maps and global measurements of Ly$\alpha$ equivalent width,
Ly$\alpha$ escape fraction, star formation rate, and nebular dust extinction.
Different size estimates for the Ly$\alpha$, H$\alpha$, and FUV emission are
obtained using Petrosian radii and growth curve analysis. The main results of
the study are summarized below:

  * •

We detect Ly$\alpha$ emission on a global scale from 39 out of 45 galaxies.
The equivalent widths and escape fractions range from $\sim 1$ Å to $\sim 80$
Å, and 0.1% to 30%, respectively. The total FUV luminosities of the galaxies
in our sample are comparable (but slightly lower) to those in Ly$\alpha$
emitter samples at high redshift, but the Ly$\alpha$ luminosities are
considerably lower.

  * •

We find that Ly$\alpha$ is more extended than the UV continuum (as measured by
Petrosian radius) in 35 out of the 39 galaxies with global Ly$\alpha$
emission. The distribution of Ly$\alpha$ emission in the galaxies tends to be
more clumpy and irregular in the sources with low $f_{\mathrm{esc}}$, but can
be overall asymmetric also in galaxies with high $f_{\mathrm{esc}}$. Comparing
the Ly$\alpha$ maps to nebular extinction maps (H$\alpha/$H$\beta$) we find
that Ly$\alpha$ escape tends to happen along directions in the galaxies where
there is less extinction.

  * •

We investigate how the escape of Ly$\alpha$ depends on a number of galaxy
properties, both stellar (mass, age, dust, and SFR quantities) and ISM
(metallicity, O32, dust). We find that the escape fraction is strongly anti-
correlated with extinction, both for stars and nebular gas. We also find a
weak anti-correlation of $f_{\mathrm{esc}}$ with stellar mass and discuss
various explanations for this. Neither SFR, specific SFR, or SFR surface
density is correlated to $f_{\mathrm{esc}}$ in our study. The ionisation state
of the ISM (as probed by the O32 ratio) in the galaxies seems not to be
related to $f_{\mathrm{esc}}$ in a direct way, but we note that most of the
galaxies with O32$<1$ have $f_{\mathrm{esc}}$ $<10$%.

  * •

The anti-correlation of Ly$\alpha$ escape fraction with nebular extinction is
strong, but shows substantial scatter. Fitting a log-linear relation to the
measurements we find an intercept of $f_{\mathrm{esc}}(E(B-V)_{n}=0)=0.075\pm
0.019$ and a slope of $3.30\pm 0.73$. We compare the observations with the
expected relation for a number of attenuation curves and find that simple dust
screen models are too steep to explain the observations. The shallow slope we
find is comparable to those found by other Ly$\alpha$ studies, in particular
when the samples include also weak and dusty Ly$\alpha$ emitters. We conclude
that the shallow slope found — with the observed $f_{\mathrm{esc}}$ being
lower than expected at low $E(B-V)_{n}$, and higher than expected at high
$E(B-V)_{n}$ — can be explained by a combination of radiative transfer effects
and a clumpy dust geometry.

  * •

The equivalent width of Ly$\alpha$ is strongly correlated to
$f_{\mathrm{esc}}$ as found by multiple previous authors. We fit a linear
function to the observed data points and find a a relatively tight relation:
$f_{\mathrm{esc}}=\left(0.0040\pm 0.00050\cdot W_{\mathrm{Ly\alpha}}\right)\pm
0.013$. However, the scatter around the relation is large and we discuss two
different scenarios (high extinction and a non-standard ionisation production
efficiency) which could cause individual sources to deviate from the relation.
We then use this relation to provide an approximate dust corrected SFR
calibration for Ly$\alpha$ luminosity which uses only $W_{\mathrm{Ly\alpha}}$
to correct the luminosity for the non-zero escape fraction. When comparing the
Ly$\alpha$ SFR estimate to dust corrected H$\alpha$ estimates for our galaxies
we find a good agreement, but with large scatter (0.20 dex for the full
sample, 0.16 dex for a sub-sample of galaxies with $E(B-V)_{n}<0.3$).

  * •

We characterize the aperture losses of Ly$\alpha$ in limited aperture
spectroscopic observations by using stacked growth curves of our sample. We
find that the losses can be quite severe for small slit widths or aperture
sizes at certain redshift ranges, and that this needs to be taken into account
when discussing results or planning new surveys.

The fully reduced data and maps are available to the community along with our
measurements as a High Level Science Product at MAST
(http://dx.doi.org/10.17909/pbe1-m743 (catalog 10.17909/pbe1-m743)). The
Ly$\alpha$ images for our reference sample, and the supporting observations
can thus be used to interpret and compare to future Ly$\alpha$ observations in
both nearby and very distant galaxies.

This work has been supported by the Swedish National Space Agency (SNSA). M.M.
acknowledges the support of the Swedish Research Council, Vetenskapsrådet
(grant 2019-00502). JMMH is funded by Spanish MCIN/AEI/10.13039/501100011033
grant PID2019-107061GB-C61. The Cosmic Dawn Center (DAWN) is funded by the
Danish National Research Foundation under grant No. 140. This research is
based on observations made with the NASA/ESA Hubble Space Telescope obtained
from the Space Telescope Science Institute, which is operated by the
Association of Universities for Research in Astronomy, Inc., under NASA
contract NAS 5–26555. These observations are associated with programs 12310,
13483, 13027, 14923, and 13656. The research has also made use of the
NASA/IPAC Infrared Science Archive, which is funded by the National
Aeronautics and Space Administration and operated by the California Institute
of Technology.

##  Appendix A HST filters used for each galaxy

Table 7: HST imaging for individual galaxies

ID | Ly$\alpha$ (SBC) | H$\alpha$ (UVIS/WFC) | H$\beta$ (UVIS/WFC) | FUV cont. (SBC) | UV/optical cont. (UVIS/WFC)
---|---|---|---|---|---
(1) | (2) | (3) | (4) | (5) | (6)
LARS01 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS02 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS03 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, WFC/F435W, WFC/F814W
LARS04 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS05 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS06 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS07 | F125LP | WFC/FR656N | WFC/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS08 | F125LP | WFC/FR656N | WFC/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS09 | F125LP | WFC/FR716N | WFC/FR505N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS10 | F125LP | WFC/FR716N | WFC/FR505N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS11 | F125LP | WFC/FR716N | WFC/FR551N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS12 | F125LP | WFC/FR716N | WFC/FR551N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
LARS13 | F140LP | WFC/FR782N | WFC/FR551N | F150LP | UVIS/F390W, UVIS/F475W, UVIS/F850LP
LARS14 | F140LP | WFC/FR782N | WFC/FR551N | F150LP | UVIS/F390W, UVIS/F475W, WFC/F850LP
ELARS01 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, WFC/F435W, WFC/F814W
ELARS02 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS03 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS04 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS05 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS06 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS07 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS08 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS09 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS10 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS11 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS12 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS13 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS14 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS15 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS16 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS17 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS18 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS19 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS20 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS21 | F125LP | UVIS/F673N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS22 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS23 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS24 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS25 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS26 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS27 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
ELARS28 | F125LP | UVIS/F680N | UVIS/FQ508N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
Tol1214 | F125LP | UVIS/F680N | UVIS/F502N | F140LP, F150LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
Tol1247 | F125LP | UVIS/FQ674N | WFC/FR505N | F140LP | UVIS/F336W, UVIS/F438W, UVIS/F775W
J1156 | F140LP | WFC/FR782N | WFC/FR601N | F150LP, F165LP | UVIS/F390W, UVIS/F475W, UVIS/F850LP

Note. — The filter columns shows the instruments and filters used.

##  Appendix B Point spread function matching

###  B.1 Building PSF models

Figure 18: PSF models for a selection of filters.

Accurate PSF models are available for most of the filters in HST instruments,
but some of these only cover the inner parts of the PSF. In order to correctly
estimate low surface brightness emission of Ly$\alpha$ using SBC we need match
the PSF over a larger area than what is usually done. For the study reported
in this paper we thus make custom models that reach a radius of at least 1
″(25 pixels for the resampled pixel scale in use). In Figure 18 we show models
for a selection of filters used.

Models for the optical (ACS/WFC and WFC3/UVIS) filters were generated by
adding a synthetic point source to individual reduced data frames, drizzling
these using the same procedure as for the science frames, and then cutting out
a thumbnail image of the PSF which was then properly background subtracted and
normalized. The synthetic point source used for this comes from the TinyTim
software, which is able to produce models out to the needed radius, and each
individual frame has the correct defocus parameter and position on chip (we
use the center of the galaxy) set for the model. We thus obtain one PSF model
per filter and galaxy which has undergone the same treatment (drizzling and
distortion correction) as the real data.

Figure 19: Best fit PSF model for the SBC F125LP filter. Observed data points
come from a stack of stars in NGC 6681. The average residual for the observed
data points is given in the legend.

For ACS/SBC the models from TinyTim are inaccurate in the wings and are
therefore unusable for the matching required in this study. Instead we use HST
calibration data targeting the stellar cluster NGC 6681. We drizzle this data
to the same pixel scale as the one used in our data set (0.04 ″per pixel),
find sufficiently isolated stars (roughly 50 stars) in the frame and fit a 2D
model to them, making sure to subtract background. This is done iteratively
and nearby stars that may affect the fit are removed from the image before the
next pass is done. The model that provides the best fit to the SBC PSF
consists of three components: (i) a Moffat model with seven free parameters,
including the center coordinates and ellipticity, (ii) a Gaussian model with
three free parameters, and (iii) a radially symmetric 4th order polynomial
function. The Moffat function dominates in the center and the polynomial
function in the outer parts. The Gaussian function is needed to correctly
characterize the observed bump in the SBC PSF wings at a radius of $\sim 0.5$
″. An example of a SBC PSF model fit is shown in Figure 19.

A number of different functions have been tried to find a suitable model for
the SBC filters, and this is by far the best option we have found. Another
option would be to use a thumbnail of stacked stars directly, but the
resulting PSF model is then too noisy to be used in the PSF matching code. It
should be noted that the PSF models for SBC are not dependent on defocus, nor
location on chip, they are instead an average model per filter that we use for
all galaxies. The uncertainty resulting from not taking these second order
effects into account is negligible compared to the inaccuracy found for the
SBC TinyTim models. Note that, somewhat unintuitively, the bluest SBC filter
(F125LP for most of the galaxies) has the widest PSF in terms of power
contained in the wings.

From the set of PSF models for a given galaxy (one per filter) we construct a
maximum width 2D PSF model by finding the maximum value per pixel in the set
of normalized models. This is then the PSF model that all of the filters are
later on matched to. We have found that this step is necessary to limit the
appearance of convolution artifacts in the matching kernel. The maximum width
PSF model is generally dominated by the reddest filter in the center regions,
and by the bluest SBC filter in the wings.

###  B.2 The PSF matching technique

Figure 20: Relative residual of PSF models compared to a maximum width PSF
(filter setup and individual models from the ELARS01 data set). This figure
shows residuals without performing any matching. Figure 21: Relative residual
of PSF models after matching them to the maximum width PSF using the technique
described in this paper.

Initial attempts to match the PSFs from FUV to red optical using a single
Gaussian convolution kernel showed very poor results, in particular in the
wings. The reason for this is the non-standard shape of the SBC PSF, which is
very narrow in the center but has strong extended wings. Motivated by this we
have developed a matching tool with a more general kernel base, which has the
freedom to correctly find the optimum matching solution also for complicated
PSFs. We base our algorithm on the work of Becker et al. (2012), who use a
convolution kernel built from a linear combination of multiple delta functions
(one per pixel in the kernel). A kernel of this type is shown to be very
powerful in finding an optimal matching solution, but is very sensitive to
noise (both in the data and numerical), due to the high number of degrees of
freedom in finding the matching kernel. The authors address this issue by
adding regularisation to the kernel solving. They show that moderate
regularisation makes it possible to find robust matching kernels without
significant loss of accuracy.

Following Becker et al. (2012) the expression for the kernel we use is:

$K_{i,j}(u,v)=\delta(u-i)\delta(v-j),$ (B1)

which represents one delta function per kernel pixel. The final convolution
kernel is then modelled as a linear combination of these basis functions,
$K(u,v)=\Sigma_{i,j}a_{i,j}K_{i,j}(u,v)$. The convolution of the input PSF
model ($\mathbf{R}$) with the kernel basis funtions can then be written as
$\mathbf{C}=\mathbf{R}\bigotimes K_{i,j}$. We then want to find the optimum
matching kernel (i.e., solve for $a_{i,j}$) that best transforms the input
model into the maximum width filter PSF ($\mathbf{S}$). This is done by
solving the following matrix equation:

$\mathbf{b}=\mathbf{M}_{\lambda}\,a=(\mathbf{M}+\lambda_{r}\mathbf{H})\,a,$
(B2)

where $\mathbf{b}=\mathbf{C}^{T}\mathbf{\Sigma}^{-1}\mathbf{S}$,
$\mathbf{M}=\mathbf{C}^{T}\mathbf{\Sigma}^{-1}\mathbf{C}$ and
$\lambda_{r}\mathbf{H}$ is the regularisation term. For details of the matrix
equation and definitions refer to Becker et al. (2012). The best value to use
for $\lambda_{r}$ depends on the input models and we run tests to find a good
number that limits the amount of over-fitting artifacts in the matched PSF
models while still keeping the overall matching residuals low. For the models
used in this project we have found that $\lambda_{r}\sim 10^{-5}$ is usually
good, but for some of the galaxies a higher regularisation strength is needed
(up to $10^{-4}$). The algorithm is implemented as a parallelised python code
which can find the optimum matching kernel for PSF models using a given kernel
size. The number of operations performed by the code scale with kernel size to
the power of three, which means that this size should be kept as low as
possible. For our data set we have found that a kernel size of 40 pixels
(1.6″, corresponding to a radius of $\sim$ 0.8″) is adequate to get accurate
matching results while keeping the computing time reasonable.

In Figure 20 we show example relative residual maps
($(\mathbf{R}-\mathbf{S})/\mathbf{S}$) of comparing individual filter models
(using the same filters as the PSF models displayed in Figure 18) to the
maximum width model without performing any PSF matching. As expected the
residuals are very strong in the wings of the optical filter PSFs, and quite
substantial also in the core. The SBC filters look slightly better, but the
core regions are problematic also there. In Figure 21 we show the result of
applying the optimum PSF matching with a kernel size of 40 pixels. This
changes the residual maps drastically, decreasing the average residual by a
factor $\sim$10–20 (depending on the individual filter).

##  Appendix C Surface brightness, $W_{\mathrm{Ly\alpha}}$, and
$f_{\mathrm{esc}}$ profiles

This section contains tables of surface brightness, $W_{\mathrm{Ly\alpha}}$,
and $f_{\mathrm{esc}}$ profiles for the galaxies. Note that the we have
elected to use the same number of significant digits for all the entries, even
if the errors are sometimes large. The dynamic range of the measurements is
very large, and the errors depend both on the measurements and the fitting
uncertainties. This means that the errors are often very small far away from
the bright galaxy centres, and the number of significant digits have been
chosen to be able to show these errors.

Table 8: Surface brightness, $W_{\mathrm{Ly\alpha}}$, and $f_{\mathrm{esc}}$
profiles

Galaxy ID | $r$ | SBFUV | SBHα | SBLyα | $W_{\mathrm{Ly\alpha}}$ | $f_{\mathrm{esc}}$
---|---|---|---|---|---|---
| kpc |  $10^{39}$erg/s/Å/kpc2 |  $10^{40}$erg/s/kpc2 |  $10^{40}$erg/s/kpc2 | Å |
(1) | (2) | (3) | (4) | (5) | (6) | (7)
LARS01 | 0.10 | 70.99200${}_{0.29892}^{0.21275}$ | 38.51040${}_{0.19535}^{0.32427}$ | 158.01700 ${}_{6.68443}^{8.49713}$ | 22.25840 ${}_{0.98618}^{1.10507}$ | 0.42986 ${}_{0.01997}^{0.01997}$
LARS01 | 0.33 | 13.05120${}_{0.04493}^{0.03799}$ | 28.19240${}_{0.06222}^{0.06018}$ | 36.39420 ${}_{1.09275}^{0.82643}$ | 25.86760 ${}_{0.47232}^{0.54082}$ | 0.16459 ${}_{0.00309}^{0.00309}$
LARS01 | 0.56 | 3.83610 ${}_{0.01470}^{0.01055}$ | 10.98370${}_{0.02848}^{0.02638}$ | 12.20720 ${}_{0.29510}^{0.24303}$ | 27.53570 ${}_{0.42431}^{0.45634}$ | 0.13933 ${}_{0.00251}^{0.00251}$
LARS01 | 0.78 | 1.99101 ${}_{0.00719}^{0.00827}$ | 6.78822 ${}_{0.02206}^{0.01633}$ | 8.20741 ${}_{0.18433}^{0.18546}$ | 30.00710 ${}_{0.41020}^{0.37692}$ | 0.13022 ${}_{0.00187}^{0.00187}$
LARS01 | 1.01 | 1.31142 ${}_{0.00568}^{0.00459}$ | 5.01932 ${}_{0.01414}^{0.01185}$ | 5.55087 ${}_{0.11409}^{0.13061}$ | 31.70100 ${}_{0.34369}^{0.38338}$ | 0.12243 ${}_{0.00139}^{0.00139}$
LARS01 | 1.24 | 0.70184 ${}_{0.00446}^{0.00375}$ | 3.16691 ${}_{0.01451}^{0.00720}$ | 3.01222 ${}_{0.07796}^{0.07880}$ | 32.64860 ${}_{0.31905}^{0.37942}$ | 0.11510 ${}_{0.00124}^{0.00124}$
LARS01 | 1.47 | 0.48885 ${}_{0.00377}^{0.00262}$ | 1.29346 ${}_{0.00808}^{0.00739}$ | 2.02582 ${}_{0.07115}^{0.06760}$ | 33.22940 ${}_{0.35021}^{0.31417}$ | 0.11564 ${}_{0.00136}^{0.00136}$
LARS01 | 1.69 | 0.38911 ${}_{0.00263}^{0.00254}$ | 0.76669 ${}_{0.00601}^{0.00729}$ | 1.47740 ${}_{0.05702}^{0.06151}$ | 33.50380 ${}_{0.36378}^{0.35916}$ | 0.11658 ${}_{0.00131}^{0.00131}$
LARS01 | 1.92 | 0.24877 ${}_{0.00191}^{0.00193}$ | 0.56295 ${}_{0.00470}^{0.00503}$ | 1.08880 ${}_{0.03840}^{0.03709}$ | 33.92070 ${}_{0.35802}^{0.37290}$ | 0.11715 ${}_{0.00125}^{0.00125}$
LARS01 | 2.15 | 0.16905 ${}_{0.00192}^{0.00162}$ | 0.42472 ${}_{0.00481}^{0.00466}$ | 0.95145 ${}_{0.03103}^{0.04196}$ | 34.59480 ${}_{0.31824}^{0.34963}$ | 0.11843 ${}_{0.00125}^{0.00125}$
LARS01 | 2.38 | 0.13984 ${}_{0.00167}^{0.00142}$ | 0.56876 ${}_{0.00374}^{0.00442}$ | 0.85177 ${}_{0.03095}^{0.04317}$ | 35.30470 ${}_{0.38651}^{0.29758}$ | 0.11819 ${}_{0.00136}^{0.00136}$
LARS01 | 2.60 | 0.10606 ${}_{0.00127}^{0.00140}$ | 0.43856 ${}_{0.00408}^{0.00355}$ | 0.65158 ${}_{0.03461}^{0.02706}$ | 35.88010 ${}_{0.36342}^{0.29408}$ | 0.11771 ${}_{0.00130}^{0.00130}$
LARS01 | 2.83 | 0.10035 ${}_{0.00142}^{0.00143}$ | 0.26576 ${}_{0.00252}^{0.00410}$ | 0.57028 ${}_{0.02338}^{0.03272}$ | 36.34590 ${}_{0.35660}^{0.34410}$ | 0.11874 ${}_{0.00129}^{0.00129}$
LARS01 | 3.06 | 0.10276 ${}_{0.00162}^{0.00120}$ | 0.28880 ${}_{0.00203}^{0.00460}$ | 0.46350 ${}_{0.02990}^{0.03008}$ | 36.55690 ${}_{0.33116}^{0.36894}$ | 0.11906 ${}_{0.00129}^{0.00129}$
LARS01 | 3.29 | 0.06105 ${}_{0.00107}^{0.00114}$ | 0.16322 ${}_{0.00314}^{0.00358}$ | 0.39626 ${}_{0.01937}^{0.02541}$ | 36.98720 ${}_{0.33296}^{0.33533}$ | 0.12001 ${}_{0.00126}^{0.00126}$
LARS01 | 3.51 | 0.04505 ${}_{0.00081}^{0.00081}$ | 0.11332 ${}_{0.00277}^{0.00277}$ | 0.36863 ${}_{0.02067}^{0.02556}$ | 37.51910 ${}_{0.33064}^{0.35437}$ | 0.12177 ${}_{0.00135}^{0.00135}$
LARS01 | 3.74 | 0.03694 ${}_{0.00103}^{0.00091}$ | 0.09549 ${}_{0.00297}^{0.00212}$ | 0.34906 ${}_{0.02205}^{0.02507}$ | 38.10440 ${}_{0.42559}^{0.35960}$ | 0.12314 ${}_{0.00137}^{0.00137}$
LARS01 | 3.97 | 0.02005 ${}_{0.00082}^{0.00094}$ | 0.05743 ${}_{0.00261}^{0.00252}$ | 0.28614 ${}_{0.01302}^{0.02067}$ | 38.72050 ${}_{0.41239}^{0.36458}$ | 0.12436 ${}_{0.00132}^{0.00132}$
LARS01 | 4.20 | 0.01615 ${}_{0.00079}^{0.00077}$ | 0.04557 ${}_{0.00232}^{0.00257}$ | 0.24280 ${}_{0.01352}^{0.01844}$ | 39.27850 ${}_{0.40120}^{0.37909}$ | 0.12549 ${}_{0.00137}^{0.00137}$
LARS01 | 4.42 | 0.01513 ${}_{0.00073}^{0.00082}$ | 0.03664 ${}_{0.00231}^{0.00274}$ | 0.22110 ${}_{0.01102}^{0.01709}$ | 39.80430 ${}_{0.39713}^{0.40556}$ | 0.12631 ${}_{0.00137}^{0.00137}$
LARS01 | 4.65 | 0.02175 ${}_{0.00071}^{0.00060}$ | 0.03915 ${}_{0.00200}^{0.00217}$ | 0.18782 ${}_{0.01542}^{0.01392}$ | 40.14830 ${}_{0.40845}^{0.47300}$ | 0.12664 ${}_{0.00144}^{0.00144}$
LARS01 | 4.88 | 0.02118 ${}_{0.00061}^{0.00066}$ | 0.03577 ${}_{0.00249}^{0.00185}$ | 0.16818 ${}_{0.01297}^{0.01260}$ | 40.44280 ${}_{0.40367}^{0.47217}$ | 0.12682 ${}_{0.00144}^{0.00144}$
LARS01 | 5.11 | 0.01575 ${}_{0.00068}^{0.00069}$ | 0.03014 ${}_{0.00237}^{0.00234}$ | 0.14806 ${}_{0.01511}^{0.01234}$ | 40.75410 ${}_{0.39275}^{0.44259}$ | 0.12685 ${}_{0.00137}^{0.00137}$
LARS01 | 5.33 | 0.01151 ${}_{0.00055}^{0.00068}$ | 0.01777 ${}_{0.00199}^{0.00249}$ | 0.11454 ${}_{0.01322}^{0.01025}$ | 41.01380 ${}_{0.37782}^{0.41764}$ | 0.12789 ${}_{0.00127}^{0.00127}$
LARS01 | 5.56 | 0.01064 ${}_{0.00053}^{0.00062}$ | 0.01521 ${}_{0.00190}^{0.00209}$ | 0.09811 ${}_{0.01009}^{0.01347}$ | 41.23110 ${}_{0.37352}^{0.42272}$ | 0.12881 ${}_{0.00127}^{0.00127}$
LARS01 | 5.79 | 0.00773 ${}_{0.00052}^{0.00043}$ | 0.01026 ${}_{0.00229}^{0.00176}$ | 0.09727 ${}_{0.01083}^{0.01223}$ | 41.50220 ${}_{0.38337}^{0.45148}$ | 0.12986 ${}_{0.00135}^{0.00135}$
LARS01 | 6.02 | 0.00674 ${}_{0.00046}^{0.00049}$ | 0.00840 ${}_{0.00166}^{0.00181}$ | 0.08448 ${}_{0.01706}^{0.00951}$ | 41.74500 ${}_{0.40574}^{0.42879}$ | 0.13082 ${}_{0.00140}^{0.00140}$
LARS01 | 6.24 | 0.00600 ${}_{0.00048}^{0.00042}$ | 0.00742 ${}_{0.00157}^{0.00189}$ | 0.06429 ${}_{0.00996}^{0.01176}$ | 41.91980 ${}_{0.41548}^{0.41685}$ | 0.13155 ${}_{0.00141}^{0.00141}$
LARS01 | 6.47 | 0.00514 ${}_{0.00047}^{0.00043}$ | 0.00801 ${}_{0.00184}^{0.00200}$ | 0.05201 ${}_{0.01377}^{0.01151}$ | 42.06010 ${}_{0.42472}^{0.42333}$ | 0.13212 ${}_{0.00148}^{0.00148}$
LARS01 | 6.70 | 0.00485 ${}_{0.00037}^{0.00046}$ | 0.00783 ${}_{0.00182}^{0.00208}$ | 0.04873 ${}_{0.01241}^{0.01181}$ | 42.19510 ${}_{0.45770}^{0.39695}$ | 0.13266 ${}_{0.00156}^{0.00156}$
LARS02 | 0.10 | 17.01260${}_{0.12730}^{0.09033}$ | 25.70370${}_{0.17739}^{0.22320}$ | $<2.83704$ | $<1.63848$ | $<0.00958$
LARS02 | 0.38 | 1.83723 ${}_{0.01057}^{0.01042}$ | 5.91272 ${}_{0.03036}^{0.02139}$ | 7.51897 ${}_{0.17952}^{0.28318}$ | 22.34760 ${}_{1.15808}^{1.01958}$ | 0.09356 ${}_{0.00451}^{0.00451}$
LARS02 | 0.65 | 0.96809 ${}_{0.00461}^{0.00567}$ | 1.97144 ${}_{0.01155}^{0.01040}$ | 2.94720 ${}_{0.11656}^{0.11996}$ | 25.59460 ${}_{1.09964}^{0.82924}$ | 0.10887 ${}_{0.00431}^{0.00431}$
LARS02 | 0.93 | 0.52246 ${}_{0.00492}^{0.00399}$ | 0.86781 ${}_{0.00516}^{0.00474}$ | 2.23548 ${}_{0.08832}^{0.09785}$ | 29.88720 ${}_{0.88144}^{0.83957}$ | 0.13106 ${}_{0.00382}^{0.00382}$
LARS02 | 1.21 | 0.41289 ${}_{0.00350}^{0.00235}$ | 0.81836 ${}_{0.00491}^{0.00759}$ | 2.07652 ${}_{0.07493}^{0.05451}$ | 34.17820 ${}_{0.97234}^{0.70548}$ | 0.15393 ${}_{0.00411}^{0.00411}$

Note. — Table 8 is published in its entirety in machine-readable format. A
portion (profiles for LARS01 and part of LARS02) is shown here for guidance
regarding its form and content.

## References

  * Adams et al. (2011) Adams, J. J., Blanc, G. A., Hill, G. J., et al. 2011, ApJS, 192, 5, doi: 10.1088/0067-0049/192/1/5
  * Akritas et al. (1995) Akritas, M. G., Murphy, S. A., & LaValley, M. P. 1995, Journal of the American Statistical Association, 90, 170.  http://www.jstor.org/stable/2291140
  * Aller (1984) Aller, L. H. 1984, Physics of thermal gaseous nebulae, doi: 10.1007/978-94-010-9639-3
  * Anders & Fritze-v. Alvensleben (2003) Anders, P., & Fritze-v. Alvensleben, U. 2003, A&A, 401, 1063, doi: 10.1051/0004-6361:20030151
  * Astropy Collaboration et al. (2013) Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33, doi: 10.1051/0004-6361/201322068
  * Atek et al. (2009) Atek, H., Kunth, D., Schaerer, D., et al. 2009, A&A, 506, L1, doi: 10.1051/0004-6361/200912787
  * Atek et al. (2014) —. 2014, A&A, 561, A89, doi: 10.1051/0004-6361/201321519
  * Baldwin et al. (1981) Baldwin, J. A., Phillips, M. M., & Terlevich, R. 1981, PASP, 93, 5, doi: 10.1086/130766
  * Becker et al. (2012) Becker, A. C., Homrighausen, D., Connolly, A. J., et al. 2012, MNRAS, 425, 1341, doi: 10.1111/j.1365-2966.2012.21542.x
  * Blanton et al. (2001) Blanton, M. R., Dalcanton, J., Eisenstein, D., et al. 2001, AJ, 121, 2358, doi: 10.1086/320405
  * Bradley et al. (2022) Bradley, L., Sipőcz, B., Robitaille, T., et al. 2022, astropy/photutils: 1.5.0, 1.5.0, Zenodo, Zenodo, doi: 10.5281/zenodo.6825092
  * Bridge et al. (2018) Bridge, J. S., Hayes, M., Melinder, J., et al. 2018, ApJ, 852, 9, doi: 10.3847/1538-4357/aa9932
  * Brinchmann et al. (2004) Brinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151, doi: 10.1111/j.1365-2966.2004.07881.x
  * Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682, doi: 10.1086/308692
  * Calzetti et al. (1994) Calzetti, D., Kinney, A. L., & Storchi-Bergmann, T. 1994, ApJ, 429, 582, doi: 10.1086/174346
  * Cardelli et al. (1989) Cardelli, J. A., Clayton, G. C., & Mathis, J. S. 1989, ApJ, 345, 245, doi: 10.1086/167900
  * Casey et al. (2014) Casey, C. M., Narayanan, D., & Cooray, A. 2014, Phys. Rep., 541, 45, doi: 10.1016/j.physrep.2014.02.009
  * Cowie et al. (2011) Cowie, L. L., Barger, A. J., & Hu, E. M. 2011, ApJ, 738, 136, doi: 10.1088/0004-637X/738/2/136
  * Cowie & Hu (1998) Cowie, L. L., & Hu, E. M. 1998, AJ, 115, 1319, doi: 10.1086/300309
  * Davidsen et al. (1977) Davidsen, A. F., Hartig, G. F., & Fastie, W. G. 1977, Nature, 269, 203, doi: 10.1038/269203a0
  * Dijkstra (2014) Dijkstra, M. 2014, PASA, 31, e040, doi: 10.1017/pasa.2014.33
  * Dijkstra (2017) —. 2017, arXiv e-prints, arXiv:1704.03416.  https://arxiv.org/abs/1704.03416
  * Domínguez Sánchez et al. (2012) Domínguez Sánchez, H., Mignoli, M., Pozzi, F., et al. 2012, MNRAS, 426, 330, doi: 10.1111/j.1365-2966.2012.21710.x
  * Duval et al. (2016) Duval, F., Östlin, G., Hayes, M., et al. 2016, A&A, 587, A77, doi: 10.1051/0004-6361/201526876
  * Erb et al. (2014) Erb, D. K., Steidel, C. C., Trainor, R. F., et al. 2014, ApJ, 795, 33, doi: 10.1088/0004-637X/795/1/33
  * Flury et al. (2022) Flury, S. R., Jaskot, A. E., Ferguson, H. C., et al. 2022, ApJS, 260, 1, doi: 10.3847/1538-4365/ac5331
  * Fricke et al. (2001) Fricke, K. J., Izotov, Y. I., Papaderos, P., Guseva, N. G., & Thuan, T. X. 2001, AJ, 121, 169, doi: 10.1086/318016
  * Fruchter & Hook (2002) Fruchter, A. S., & Hook, R. N. 2002, PASP, 114, 144, doi: 10.1086/338393
  * Giavalisco et al. (1996) Giavalisco, M., Koratkar, A., & Calzetti, D. 1996, ApJ, 466, 831, doi: 10.1086/177557
  * Gonzaga et al. (2012) Gonzaga, S., Hack, W., Fruchter, A., & Mack, J. 2012, The DrizzlePac Handbook
  * Gronke & Dijkstra (2016) Gronke, M., & Dijkstra, M. 2016, ApJ, 826, 14, doi: 10.3847/0004-637X/826/1/14
  * Guaita et al. (2015) Guaita, L., Melinder, J., Hayes, M., et al. 2015, A&A, 576, A51, doi: 10.1051/0004-6361/201425053
  * Guseva et al. (2011) Guseva, N. G., Izotov, Y. I., Stasińska, G., et al. 2011, A&A, 529, A149, doi: 10.1051/0004-6361/201016291
  * Hayes (2015) Hayes, M. 2015, PASA, 32, e027, doi: 10.1017/pasa.2015.25
  * Hayes et al. (2016) Hayes, M., Melinder, J., Östlin, G., et al. 2016, ApJ, 828, 49, doi: 10.3847/0004-637X/828/1/49
  * Hayes et al. (2007) Hayes, M., Östlin, G., Atek, H., et al. 2007, MNRAS, 382, 1465, doi: 10.1111/j.1365-2966.2007.12482.x
  * Hayes et al. (2009) Hayes, M., Östlin, G., Mas-Hesse, J. M., & Kunth, D. 2009, AJ, 138, 911, doi: 10.1088/0004-6256/138/3/911
  * Hayes et al. (2005) Hayes, M., Östlin, G., Mas-Hesse, J. M., et al. 2005, A&A, 438, 71, doi: 10.1051/0004-6361:20052702
  * Hayes et al. (2011) Hayes, M., Schaerer, D., Östlin, G., et al. 2011, ApJ, 730, 8, doi: 10.1088/0004-637X/730/1/8
  * Hayes et al. (2010) Hayes, M., Östlin, G., Schaerer, D., et al. 2010, Nature, 464, 562, doi: 10.1038/nature08881
  * Hayes et al. (2013) —. 2013, ApJ, 765, L27, doi: 10.1088/2041-8205/765/2/L27
  * Hayes et al. (2014) Hayes, M., Östlin, G., Duval, F., et al. 2014, ApJ, 782, 6, doi: 10.1088/0004-637X/782/1/6
  * Heckman et al. (2011) Heckman, T. M., Borthakur, S., Overzier, R., et al. 2011, ApJ, 730, 5, doi: 10.1088/0004-637X/730/1/5
  * Henry et al. (2015) Henry, A., Scarlata, C., Martin, C. L., & Erb, D. 2015, ApJ, 809, 19, doi: 10.1088/0004-637X/809/1/19
  * Herenz et al. (2016) Herenz, E. C., Gruyters, P., Orlitova, I., et al. 2016, A&A, 587, A78, doi: 10.1051/0004-6361/201527373
  * Herenz et al. (2017) Herenz, E. C., Urrutia, T., Wisotzki, L., et al. 2017, A&A, 606, A12, doi: 10.1051/0004-6361/201731055
  * Hold-Geoffroy et al. (2014) Hold-Geoffroy, Y., Gagnon, O., & Parizeau, M. 2014, in Proceedings of the 2014 Annual Conference on Extreme Science and Engineering Discovery Environment, XSEDE ’14 (New York, NY, USA: Association for Computing Machinery), doi: 10.1145/2616498.2616565
  * Hunter (2007) Hunter, J. D. 2007, Computing In Science & Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Izotov et al. (2020) Izotov, Y. I., Schaerer, D., Worseck, G., et al. 2020, MNRAS, 491, 468, doi: 10.1093/mnras/stz3041
  * Izotov et al. (2018) Izotov, Y. I., Worseck, G., Schaerer, D., et al. 2018, MNRAS, 478, 4851, doi: 10.1093/mnras/sty1378
  * Jaskot et al. (2017) Jaskot, A. E., Oey, M. S., Scarlata, C., & Dowd, T. 2017, ApJ, 851, L9, doi: 10.3847/2041-8213/aa9d83
  * Jensen et al. (2013) Jensen, H., Laursen, P., Mellema, G., et al. 2013, MNRAS, 428, 1366, doi: 10.1093/mnras/sts116
  * Kelly (2007) Kelly, B. C. 2007, ApJ, 665, 1489, doi: 10.1086/519947
  * Kennicutt & Evans (2012) Kennicutt, R. C., & Evans, N. J. 2012, ARA&A, 50, 531, doi: 10.1146/annurev-astro-081811-125610
  * Kewley & Dopita (2002) Kewley, L. J., & Dopita, M. A. 2002, ApJS, 142, 35, doi: 10.1086/341326
  * Kimm et al. (2019) Kimm, T., Blaizot, J., Garel, T., et al. 2019, MNRAS, 486, 2215, doi: 10.1093/mnras/stz989
  * Kornei et al. (2010) Kornei, K. A., Shapley, A. E., Erb, D. K., et al. 2010, ApJ, 711, 693, doi: 10.1088/0004-637X/711/2/693
  * Krist et al. (2011) Krist, J. E., Hook, R. N., & Stoehr, F. 2011, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 8127, Optical Modeling and Performance Predictions V, ed. M. A. Kahan, 81270J, doi: 10.1117/12.892762
  * Kroupa et al. (1993) Kroupa, P., Tout, C. A., & Gilmore, G. 1993, MNRAS, 262, 545, doi: 10.1093/mnras/262.3.545
  * Kunth et al. (2003) Kunth, D., Leitherer, C., Mas-Hesse, J. M., Östlin, G., & Petrosian, A. 2003, ApJ, 597, 263, doi: 10.1086/378396
  * Kunth et al. (1998) Kunth, D., Mas-Hesse, J. M., Terlevich, E., et al. 1998, A&A, 334, 11.  https://arxiv.org/abs/astro-ph/9802253
  * Laursen et al. (2013) Laursen, P., Duval, F., & Östlin, G. 2013, ApJ, 766, 124, doi: 10.1088/0004-637X/766/2/124
  * Laursen et al. (2009) Laursen, P., Sommer-Larsen, J., & Andersen, A. C. 2009, ApJ, 704, 1640, doi: 10.1088/0004-637X/704/2/1640
  * Laursen et al. (2019) Laursen, P., Sommer-Larsen, J., Milvang-Jensen, B., Fynbo, J. P. U., & Razoumov, A. O. 2019, A&A, 627, A84, doi: 10.1051/0004-6361/201833645
  * Le Reste et al. (2022) Le Reste, A., Hayes, M., Cannon, J. M., et al. 2022, ApJ, 934, 69, doi: 10.3847/1538-4357/ac77ed
  * Leclercq et al. (2017) Leclercq, F., Bacon, R., Wisotzki, L., et al. 2017, A&A, 608, A8, doi: 10.1051/0004-6361/201731480
  * Leitherer et al. (1999) Leitherer, C., Schaerer, D., Goldader, J. D., et al. 1999, ApJS, 123, 3, doi: 10.1086/313233
  * Malhotra & Rhoads (2004) Malhotra, S., & Rhoads, J. E. 2004, ApJ, 617, L5, doi: 10.1086/427182
  * Matthee et al. (2016) Matthee, J., Sobral, D., Oteo, I., et al. 2016, MNRAS, 458, 449, doi: 10.1093/mnras/stw322
  * Meier & Terlevich (1981) Meier, D. L., & Terlevich, R. 1981, ApJ, 246, L109, doi: 10.1086/183565
  * Momose et al. (2014) Momose, R., Ouchi, M., Nakajima, K., et al. 2014, MNRAS, 442, 110, doi: 10.1093/mnras/stu825
  * Montero-Dorta et al. (2016) Montero-Dorta, A. D., Bolton, A. S., Brownstein, J. R., et al. 2016, MNRAS, 461, 1131, doi: 10.1093/mnras/stw1352
  * Nakajima et al. (2018) Nakajima, K., Fletcher, T., Ellis, R. S., Robertson, B. E., & Iwata, I. 2018, MNRAS, 477, 2098, doi: 10.1093/mnras/sty750
  * Natta & Panagia (1984) Natta, A., & Panagia, N. 1984, ApJ, 287, 228, doi: 10.1086/162681
  * Neufeld (1991) Neufeld, D. A. 1991, ApJ, 370, L85, doi: 10.1086/185983
  * Östlin et al. (2009) Östlin, G., Hayes, M., Kunth, D., et al. 2009, AJ, 138, 923, doi: 10.1088/0004-6256/138/3/923
  * Östlin et al. (2014) Östlin, G., Hayes, M., Duval, F., et al. 2014, ApJ, 797, 11, doi: 10.1088/0004-637X/797/1/11
  * Ouchi et al. (2018) Ouchi, M., Harikane, Y., Shibuya, T., et al. 2018, PASJ, 70, S13, doi: 10.1093/pasj/psx074
  * Pardy et al. (2014) Pardy, S. A., Cannon, J. M., Östlin, G., et al. 2014, ApJ, 794, 101, doi: 10.1088/0004-637X/794/2/101
  * Partridge & Peebles (1967) Partridge, R. B., & Peebles, P. J. E. 1967, ApJ, 147, 868, doi: 10.1086/149079
  * Petrosian (1976) Petrosian, V. 1976, ApJ, 210, L53, doi: 10.1086/182301
  * Prevot et al. (1984) Prevot, M. L., Lequeux, J., Maurice, E., Prevot, L., & Rocca-Volmerange, B. 1984, A&A, 132, 389
  * Puschnig et al. (2017) Puschnig, J., Hayes, M., Östlin, G., et al. 2017, MNRAS, 469, 3252, doi: 10.1093/mnras/stx951
  * Rasekh et al. (2022) Rasekh, A., Melinder, J., Östlin, G., et al. 2022, A&A, 662, A64, doi: 10.1051/0004-6361/202140734
  * Reddy et al. (2020) Reddy, N. A., Shapley, A. E., Kriek, M., et al. 2020, ApJ, 902, 123, doi: 10.3847/1538-4357/abb674
  * Rivera-Thorsen et al. (2015) Rivera-Thorsen, T. E., Hayes, M., Östlin, G., et al. 2015, ApJ, 805, 14, doi: 10.1088/0004-637X/805/1/14
  * Rodighiero et al. (2011) Rodighiero, G., Daddi, E., Baronchelli, I., et al. 2011, ApJ, 739, L40, doi: 10.1088/2041-8205/739/2/L40
  * Runnholm et al. (2020) Runnholm, A., Hayes, M., Melinder, J., et al. 2020, ApJ, 892, 48, doi: 10.3847/1538-4357/ab7a91
  * Salim et al. (2007) Salim, S., Rich, R. M., Charlot, S., et al. 2007, ApJS, 173, 267, doi: 10.1086/519218
  * Scarlata et al. (2009) Scarlata, C., Colbert, J., Teplitz, H. I., et al. 2009, ApJ, 704, L98, doi: 10.1088/0004-637X/704/2/L98
  * Schlafly & Finkbeiner (2011) Schlafly, E. F., & Finkbeiner, D. P. 2011, ApJ, 737, 103, doi: 10.1088/0004-637X/737/2/103
  * Smith et al. (2015) Smith, A., Safranek-Shrader, C., Bromm, V., & Milosavljević, M. 2015, MNRAS, 449, 4336, doi: 10.1093/mnras/stv565
  * Sobral & Matthee (2019) Sobral, D., & Matthee, J. 2019, A&A, 623, A157, doi: 10.1051/0004-6361/201833075
  * Sobral et al. (2018) Sobral, D., Matthee, J., Darvish, B., et al. 2018, MNRAS, 477, 2817, doi: 10.1093/mnras/sty782
  * Steidel et al. (2011) Steidel, C. C., Bogosavljević, M., Shapley, A. E., et al. 2011, ApJ, 736, 160, doi: 10.1088/0004-637X/736/2/160
  * Terlevich et al. (1993) Terlevich, E., Diaz, A. I., Terlevich, R., & Vargas, M. L. G. 1993, MNRAS, 260, 3, doi: 10.1093/mnras/260.1.3
  * Terlevich et al. (1991) Terlevich, R., Melnick, J., Masegosa, J., Moles, M., & Copetti, M. V. F. 1991, A&AS, 91, 285
  * Trainor et al. (2019) Trainor, R. F., Strom, A. L., Steidel, C. C., et al. 2019, ApJ, 887, 85, doi: 10.3847/1538-4357/ab4993
  * van der Walt et al. (2014) van der Walt, S., Schönberger, J. L., Nunez-Iglesias, J., et al. 2014, PeerJ, 2, e453, doi: 10.7717/peerj.453
  * van Dokkum (2001) van Dokkum, P. G. 2001, PASP, 113, 1420, doi: 10.1086/32389410.48550/arXiv.astro-ph/0108003
  * Verhamme et al. (2012) Verhamme, A., Dubois, Y., Blaizot, J., et al. 2012, A&A, 546, A111, doi: 10.1051/0004-6361/201218783
  * Verhamme et al. (2015) Verhamme, A., Orlitová, I., Schaerer, D., & Hayes, M. 2015, A&A, 578, A7, doi: 10.1051/0004-6361/201423978
  * Verhamme et al. (2008) Verhamme, A., Schaerer, D., Atek, H., & Tapken, C. 2008, A&A, 491, 89, doi: 10.1051/0004-6361:200809648
  * Wisotzki et al. (2016) Wisotzki, L., Bacon, R., Blaizot, J., et al. 2016, A&A, 587, A98, doi: 10.1051/0004-6361/201527384
  * Yang et al. (2017a) Yang, H., Malhotra, S., Rhoads, J. E., et al. 2017a, ApJ, 838, 4, doi: 10.3847/1538-4357/aa6337
  * Yang et al. (2017b) Yang, H., Malhotra, S., Gronke, M., et al. 2017b, ApJ, 844, 171, doi: 10.3847/1538-4357/aa7d4d
  * Yin et al. (2007) Yin, S. Y., Liang, Y. C., Hammer, F., et al. 2007, A&A, 462, 535, doi: 10.1051/0004-6361:20065798


# Trajectory planning with a dynamic obstacle clustering strategy using Mixed-
Integer Linear Programming*

Vinícius Antonio Battagello1, Nei Yoshihiro Soma1 and Rubens Junqueira
Magalhães Afonso2,3 *This work was supported by Conselho Nacional de
Desenvolvimento Científico e Tecnológico CNPq (grant 141007/2016-8), São Paulo
Research Foundation - FAPESP (grant 2016/01860-1), CNPq PQ, and CAPES (grant
88881.145490/2017-01) and the German Ministry of Education and Research
through the Alexander von Humboldt Foundation 1Vinícius A. Battagello is with
the Divisão de Engenharia de Computação, Instituto Tecnológico de Aeronáutica,
São José dos Campos, SP 12228-900, Brazil<EMAIL_ADDRESS>Y. Soma is
with the Divisão de Engenharia de Computação, Instituto Tecnológico de
Aeronáutica, São José dos Campos, SP 12228-900, Brazil<EMAIL_ADDRESS>J.
M. Afonso is with the Institute of Flight System Dynamics, Department of
Aerospace and Geodesy, Technical University of Munich (TUM), Garching bei
München, Bavaria, Germany<EMAIL_ADDRESS>J. M. Afonso is on leave
from the Divisão de Engenharia Eletrônica, Instituto Tecnológico de
Aeronáutica, São José dos Campos, SP 12228-900, Brazil<EMAIL_ADDRESS>
###### Abstract

In this paper we propose a technique that assigns obstacles to clusters used
for collision avoidance via Mixed-Integer Programming. This strategy enables a
reduction in the number of binary variables used for collision avoidance, thus
entailing a decrease in computational cost, which has been a hindrance to the
application of Model Predictive Control approaches with Mixed-Integer
Programming formulations in real-time. Moreover, the assignment of obstacles
to clusters and the sizes of the clusters are decided within the same
optimization problem that performs the trajectory planning, thus yielding
optimal cluster choices. Simulation results are presented to illustrate an
application of the proposal.

##  I INTRODUCTION

Motion planning in the presence of obstacles is of paramount importance for
autonomous agents. Generally, motion planning may be split into three
different layers, namely (i) path planning, (ii) trajectory planning, and
(iii) low-level control.

Path planning methods aim at finding routes in the euclidean space without
time parametrization, whereas trajectory planning deals with the search for
routes that obey dynamics constraints, usually assigning a time tag to each
position in space.

Some proposals combine both path and trajectory planning in a single layer
with advantages regarding the optimal cost of the ensuing trajectory. Current
successful approaches involve using Model Predictive Control (MPC)
formulations based on the solution of a Mixed-Integer Linear Program (MILP) to
cope with the nonconvexity brought about by the presence of obstacles. In [1]
the authors propose an approximate model of aircraft dynamics using linear
constraints to the trajectory planning of airplanes with the application of a
MILP approach. The formulation ensures collision avoidance for each aircraft
and guarantees the desired hard constraints fulfillment. This approach is
extended in [2], which applies a constraint tightening strategy to obtain a
robust solution that guarantees finite-time arrival into an arbitrary target
set and fulfillment of output/input constraints, in spite of unknown
disturbance. The central idea is to hold a “border” for feedback action as
time goes by.

However, MILP is known to be NP-hard [3], thus such approaches may demand a
prohibitive computational effort for complex environments. The computation
time demanded in trajectory planning grows strongly with the number of binary
variables used for obstacle avoidance, therefore this has been the aim of
several recent contributions both using pre-processing and alternative MILP
encodings.

In [4] a three-stage algorithm is used to perform pre-processing: (i) a
collision-free path is computed via graph search, (ii) then convex polytopes
that cover the entire path are generated by using a triangulation of the free
space, and finally (iii) a MILP formulation using the polytopes as a “tunnel”
from the starting position to the goal is solved to determine the optimal
trajectory. Since it relies on a graph search to prune the free space in the
first step, optimality of the planned trajectory cannot be guaranteed, as the
cost is not taken into account in this initial phase.

Another example of pre-processing is proposed in [5], where ellipsoidal
regions of convex obstacle-free space are calculated to provide a collision-
free path with fewer integer variables.

More recently, several techniques using both pre-processing of the environment
as well as taking advantage of the current available solution of the MPC
problem were proposed in [6]. Regarding pre-processing, the proposal involved
clustering the obstacles according to their distances to the agent and among
themselves, replacing several obstacles and the binary variables associated to
them with fewer clusters and their corresponding binary variables. Moreover,
as the MPC is run in a receding horizon fashion, the current available
solution was explored to determine which obstacles can be replaced by convex
constraints as they are not circled by the agent.

On the other hand, in [7] an alternative encoding of the exclusion regions in
terms of the binary variables associated to them is proposed, leading to a
decrease in the number of binary variables. In fact, a linear dependency of
the number of binary variables with the number of sides of the obstacles was
reduced to a logarithmic growth. The number of additional inequalities to
allow that was shown to be one per obstacle in [8]. For a more thorough
discussion of the MILP encoding schemes used in control the reader is referred
to [9].

In the present paper, we propose a novel scheme to perform the clustering of
the obstacles within the same optimization problem that is solved to yield the
optimal trajectory. Therefore, a reduction in the online computational time is
provided, whereas the clustering is directly coupled with the cost
optimization, these two characteristics being the main advantages of our
proposal.

The remainder of this paper is divided as follows: the problem is formally
presented in Sec. II along with a detailed description of how the obstacle
avoidance is performed in a MILP framework; the proposed clustering strategy
and the constraints to avoid collisions with the clusters are developed in
Sec. III; the simulation scenario is presented in Sec. IV; Sec. V contains the
simulation results and their discussion; concluding remarks and suggestions
for future work are given in Sec. VI.

##  II PROBLEM STATEMENT

Here we study the problem of maneuvering an agent in a plane with a coordinate
system $(r_{x},r_{y})$ with state $\textbf{x}\in\mathbb{R}^{4}$ and position
$\textbf{r}=[r_{x}\ r_{y}]^{T}\in\mathbb{R}^{2}$ into a target region in a
two-dimensional environment, while collisions against obstacles must be
avoided. Extending the proposal to a three-dimensional scenario is
straightforward and is not treated in the present paper for simplicity of
presentation.

The dynamics of the agent in discrete-time is given as:

$\textbf{x}[k+1]=\textbf{A}\textbf{x}[k]+\textbf{B}\textbf{u}[k].$ (1)

MPC is used to perform the maneuver of the system in (1). The task of the
control problem is to choose optimally the predictions
$\hat{\textbf{x}}[k+j|k]$ with the corresponding $\hat{\textbf{u}}[k+j|k]$ in
each timestep $j$. The maneuvering is performed in the literature [1, 10, 2]
by solving the following optimization problem:

$J\left[\textbf{x}[k]\right]=\min_{\hat{\textbf{u}},\hat{\textbf{x}},N[k]}\sum_{j=0}^{N[k]}\left(1+\gamma\mathinner{\\!\left\lVert\hat{\textbf{u}}[k+j|k]\right\rVert}_{1}\right)$
(2)

subject to

$\displaystyle\hat{\textbf{x}}[k|k]=\textbf{x}[k],$ (3a)
$\displaystyle\hat{\textbf{x}}[k+j|k]\in\mathcal{X},\ j=1,2,\ldots,N[k]+1$
(3b) $\displaystyle\hat{\textbf{u}}[k+j|k]\in\mathcal{U},\ j=0,1,\ldots,N[k]$
(3c) $\displaystyle\hat{\textbf{r}}[k+j|k]=\textbf{C}\hat{\textbf{x}}[k+j|k],\
j=1,2,\ldots,N[k]+1$ (3d)
$\displaystyle\hat{\textbf{r}}[k+j|k]\notin\mathcal{O}_{i},\
j=1,2,\ldots,N[k]+1,\ i=1,2,\ldots,N_{o}$ (3e)
$\displaystyle\hat{\textbf{x}}[k+N[k]+1|k]\in\mathcal{Q},$ (3f)
$\displaystyle\hat{\textbf{x}}[k+j+1|k]\ =\
\textbf{A}\hat{\textbf{x}}[k+j|k]+\textbf{B}\hat{\textbf{u}}[k+j|k],$ (3g)
$\displaystyle j=0,1,\ldots,N[k],$ (3h)

where $N[k]$ is a variable horizon, $\mathcal{X}\subset\mathbb{R}^{4}$ is a
polytope representing the pointwise-in-time state constraints,
$\mathcal{U}\subset\mathbb{R}^{2}$ is a polytope representing the pointwise-
in-time control constraints, $\mathcal{Q}\subset\mathbb{R}^{4}$ is a polytope
representing the terminal constraints, and $\textbf{C}\in\mathbb{R}^{2\times
4}$ is a matrix that extracts position information from the state vector. The
sets $\mathcal{O}_{i}\subset\mathbb{R}^{2}$, $i=1,2,\ldots,N_{o}$, represent
the $N_{o}$ polytopic obstacles that must be avoided.

For details on how to rewrite this optimization problem as a MILP with a
finite horizon $N_{s}$ the reader is referred to [1, 10, 2]. In the present
paper we present only the obstacle avoidance constraints in detail as they are
the object of our contribution.

###  II-A Obstacle avoidance

The presence of obstacles may render the problem of optimizing a trajectory
not convex. In this section we describe the commonly adopted remodeling used
in the MPC-MILP literature with binary optimization variables and the “big-M”
approach [11].

For the sake of simplicity, the obstacles are assumed to be rectangles with
sides parallel to the coordinate axes. Generalizing to polytopes of any number
of sides is straightforward and our simplification does not represent a
limitation of the proposal. Bearing that in mind, the obstacles
$\mathcal{O}_{i},\,i=1,2,\dots,\,N_{o}$ can be represented by:

$\mathcal{O}_{i}=\\{\textbf{z}\in\mathbb{R}^{2}|[o_{i}^{x,min}\
o_{i}^{y,min}]^{T}\leq\textbf{z}\leq[o_{i}^{x,max}\ o_{i}^{y,max}]^{T}]\\},$
(4)

where $o_{i}^{x,min}$, $o_{i}^{y,min}$, $o_{i}^{x,max}$, and $o_{i}^{y,max}$
represent the coordinates of the left, lower, right, and upper sides of the
i-th obstacle, respectively.

Collisions can be avoided by imposing that the position of the agent is in one
of the outer half-planes defined by the support lines of the sides of the
obstacle. This is done by using the binary variable $b_{j,i,f}^{\mathcal{O}}$,
which marks at time step $k+j$ if the agent is outside side $f$ of obstacle
$\mathcal{O}_{i}$. Such constraints may be written using the “big-M” method as

$\displaystyle\hat{r}_{x}[k+j|k]$ $\displaystyle\leq\,\ \
o_{i}^{x,min}+b_{j,i,1}^{\mathcal{O}}M_{\mathcal{O}},$ (5a)
$\displaystyle\hat{r}_{y}[k+j|k]$ $\displaystyle\leq\,\ \
o_{i}^{y,min}+b_{j,i,2}^{\mathcal{O}}M_{\mathcal{O}},$ (5b)
$\displaystyle-\hat{r}_{x}[k+j|k]$ $\displaystyle\leq-
o_{i}^{x,max}+b_{j,i,3}^{\mathcal{O}}M_{\mathcal{O}},$ (5c)
$\displaystyle-\hat{r}_{y}[k+j|k]$ $\displaystyle\leq-
o_{i}^{y,max}+b_{j,i,4}^{\mathcal{O}}M_{\mathcal{O}},$ (5d)
$\displaystyle\sum_{f=1}^{4}b_{j,i,f}^{\mathcal{O}}$ $\displaystyle\leq 3,\
j=1,2,\ldots,N_{s},\ i=1,2,\ldots,N_{o},$ (5e)

where $M_{\mathcal{O}}\in\mathbb{R}^{+}$ is chosen large enough to render the
inequalities (5a), (5b), (5c) and (5d) inactive for all achievable values of
$\hat{r}_{x}[k+j|k]$ and $\hat{r}_{y}[k+j|k]$, according to the “big-M”
method.

If we use this scheme to model environments with increasingly more obstacles,
as each obstacle partitions the search space into four disjoint regions, the
computation time in environments with tens of obstacles or more becomes
prohibitive. Therefore, in this paper we propose to group the obstacles in
clusters and perform the collision avoidance with the clusters instead of the
obstacles, _i.e._ using binary variables $b_{j,\ell,f}^{\mathcal{C}}$ for
cluster avoidance rather than $b_{j,i,f}^{\mathcal{O}}$ for obstacle
avoidance.

##  III PROPOSED CLUSTER-OBSTACLE ASSIGNMENT AND COLLISION AVOIDANCE

In this section we will present our proposal to enable avoiding clusters of
obstacles, with the difference regarding the literature that the definition of
the clusters is also an object of optimization within the MILP problem.

Let the $N_{c}$ clusters be represented as
$\mathcal{C}_{\ell}\subset\mathbb{R}^{2}$, $\ell=1,2,\ldots,N_{c}$. As in
section II-A, we consider rectangular clusters with sides parallel to the
coordinate axes for simplicity. Again, this simplification serves the purpose
of ease of presentation, but represents no limitation of the proposal, as it
can seamlessly be adapted for general polytopic sets. The first requirement is
that each obstacle $\mathcal{O}_{i}$ must be a subset of at least one of the
clusters $\mathcal{C}_{\ell}$, so that avoiding collisions with the clusters
entails avoiding collisions with the obstacles. This is formally stated as the
“or” constraints:

$\mathcal{O}_{i}\subset\mathcal{C}_{1}\lor\mathcal{O}_{i}\subset\mathcal{C}_{2}\lor\ldots\lor\mathcal{O}_{i}\subset\mathcal{C}_{N_{c}},\
i=1,2,\ldots,N_{o}.$ (6)

It is easy to see that, if (6) holds, then avoiding collisions with the
obstacles as required by (3e) can be done by avoiding collisions with the
clusters, since

$\hat{\textbf{r}}\notin\mathcal{C}_{\ell},\
\ell=1,2,\ldots,N_{c}\implies\hat{\textbf{r}}\notin\mathcal{O}_{i},\
i=1,2,\ldots,N_{o}.$ (7)

Since the constraints in (6) are of the logical type “or”, again binary
variables and the “big-M” are one tool to encode them within the MILP
formulation [11]. Then, using binary variables $b^{R}_{\ell,i}$ we can
implement (6) through the following inequalities:

$\displaystyle c^{x,min}_{\ell}$ $\displaystyle\leq\,\ \
o_{i}^{x,min}+(1-b^{R}_{\ell,i})M_{c},$ (8a) $\displaystyle c^{y,min}_{\ell}$
$\displaystyle\leq\,\ \ o_{i}^{y,min}+(1-b^{R}_{\ell,i})M_{c},$ (8b)
$\displaystyle-c^{x,max}_{\ell}$ $\displaystyle\leq-
o_{i}^{x,max}+(1-b^{R}_{\ell,i})M_{c},$ (8c) $\displaystyle-c^{y,max}_{\ell}$
$\displaystyle\leq-o_{i}^{y,max}+(1-b^{R}_{\ell,i})M_{c},$ (8d)
$\displaystyle\sum_{\ell=1}^{N_{c}}b^{R}_{\ell,i}$ $\displaystyle=1,\
i=1,2,\ldots,N_{o},$ (8e)

where $c^{x,min}_{\ell}$, $c^{y,min}_{\ell}$, $c^{x,max}_{\ell}$, and
$c^{y,max}$ are decision variables that represent the coordinates of the left,
lower, right, and upper sides of cluster $\mathcal{C}_{\ell}$,
$\ell=1,2,\ldots,N_{c}$, respectively. Constraint (8e) ensures that for each
obstacle $\mathcal{O}_{i}$, at least one of the binary variables
$b^{R}_{\ell,i}$ is set to one, thus assigning it to the corresponding cluster
$\mathcal{C}_{\ell}$.

Fig. 1 depicts examples with the same obstacles $\mathcal{O}_{i}$,
$i=1,2,\ldots,5$ and two clusters $\mathcal{C}_{1}$ and $\mathcal{C}_{2}$. The
corresponding values of the cluster-obstacle assignment binary variables
$b^{R}_{\ell,i}$ are shown in Tab. I for the cases “a” and “b” depicted in
Fig. 1a and Fig. 1b, respectively. The difference in the clustering choices
seen in Fig. 1 can be explained by the changes in the value of the binary
variables $b_{1,3}^{R}$ and $b_{2,3}^{R}$ between cases “a” and “b” in Tab. I.
In case “a”, $b_{2,3}^{R}=1$, effectively imposing that
$\mathcal{O}_{3}\subset\mathcal{C}_{2}$. On the other hand, in case “b”,
$b_{1,3}^{R}=1$, thus entailing $\mathcal{O}_{3}\subset\mathcal{C}_{1}$.

Figure 1: Examples of cluster-obstacle assignment for $N_{c}=2$ and $N_{o}=5$. | Binary variables
---|---
Case | $b^{R}_{1,1}$ | $b^{R}_{1,2}$ | $b^{R}_{1,3}$ | $b^{R}_{1,4}$ | $b^{R}_{1,5}$ | $b^{R}_{2,1}$ | $b^{R}_{2,2}$ | $b^{R}_{2,3}$ | $b^{R}_{2,4}$ | $b^{R}_{2,5}$
a | 1 | 1 | 0 | 0 | 0 | 0 | 0 | 1 | 1 | 1
b | 1 | 1 | 1 | 0 | 0 | 0 | 0 | 0 | 1 | 1
TABLE I: Cluster-obstacle assignment binary variables $b^{R}_{\ell,i}$ for the
example in Fig. 1.

Cluster avoidance can then be ensured by imposing that
$\textbf{r}\notin\mathcal{C}_{\ell},\ell=1,2,\ldots,N_{c}$, which can be
written in MILP form using cluster avoidance binary variables
$b_{j,\ell,f}^{\mathcal{C}}$ and the “big-M” method as

$\displaystyle\hat{r}_{x}[k+j|k]$ $\displaystyle\leq\,\ \
c^{x,min}_{\ell}+b_{j,\ell,1}^{\mathcal{C}}M_{\mathcal{C}}$ (9a)
$\displaystyle\hat{r}_{y}[k+j|k]$ $\displaystyle\leq\,\ \
c^{y,min}_{\ell}+b_{j,\ell,2}^{\mathcal{C}}M_{\mathcal{C}}$ (9b)
$\displaystyle-\hat{r}_{x}[k+j|k]$
$\displaystyle\leq-c^{x,max}_{\ell}+b_{j,\ell,3}^{\mathcal{C}}M_{\mathcal{C}}$
(9c) $\displaystyle-\hat{r}_{y}[k+j|k]$
$\displaystyle\leq-c^{y,max}_{\ell}+b_{j,\ell,4}^{\mathcal{C}}M_{\mathcal{C}}$
(9d) $\displaystyle\sum_{f=1}^{4}b_{j,\ell,f}^{\mathcal{C}}$
$\displaystyle\leq 3,\ j=1,2,\ldots,N_{s},\ \ell=1,2,\ldots,N_{c},$ (9e)

where $M_{\mathcal{C}}\in\mathbb{R}^{+}$ is chosen large enough to render the
inequalities (9a), (9b), (9c) and (9d) inactive for all achievable values of
$\hat{r}_{x}[k+j|k]$ and $\hat{r}_{y}[k+j|k]$, according to the “big-M”
method.

Whereas usually the constraints in (5) are employed in the literature to
implement (3e), requiring $4N_{s}N_{o}$ binary variables
$b^{\mathcal{O}}_{j,i,f}$, in our proposal (9) is used to impose
$\hat{\textbf{r}}[k+j|k]\notin\mathcal{C}_{\ell}$, requiring $4N_{s}N_{c}$
binary variables $b^{\mathcal{C}}_{j,\ell,f}$. Moreover, in view of (7), one
must append constraints (8) to the problem in our proposal to ensure that (6)
holds, which in turn demands $N_{o}N_{c}$ additional binary variables
$b^{R}_{\ell,i}$. Therefore, the total number of of binary variables to avoid
collisions with our proposal amounts to
$4N_{s}N_{c}+N_{o}N_{c}=(4N_{s}+N_{o})N_{c}$. Notice that in the state-of-the-
art the product between the number of obstacles $N_{o}$ and the maximal
horizon $N_{s}$ dictates the growth in the number of binary variables. In
contrast, with our proposal, this product does not appear and both $N_{s}$ and
$N_{o}$ are multiplied by the number of clusters $N_{c}$. Since the number of
obstacles and the maximal horizon are linked to the scenario itself, they
cannot be controlled by the user to achieve a lower computational burden. In
opposition, with our approach the user can decide the number $N_{c}$ and
effectively limit the growth of the number of binary variables.

Clearly, the lowest $N_{c}$, the farthest from the optimal cost the solution
yielded with our proposal, as the possibility of passing between many
obstacles is hindered by clustering them. Therefore, one expects a compromise
between the computational workload and the optimality of the planned
trajectory. However, it is important to emphasize that although our proposal
entails this compromise, one key advantage is that the clusters are decided
within the same optimization problem that minimizes the cost, as opposed to
techniques that try to eliminate the obstacles or to produce a collision-free
tunnel before calculating the optimal trajectory, as mentioned in the
literature review. Furthermore, it is interesting to remark that the
controller runs in a receding horizon fashion, therefore the optimization
problem is solved at each timestep with feedback from the current measured
state. As a consequence, the clustering can be changed throughout the
maneuver.

##  IV SIMULATION SCENARIO

We consider in this work the model of an agent as that of a particle moving in
a plane with axes $r_{x}$ and $r_{y}$ orthogonal to each other, as in [1, 10,
2]. The position vector is defined as $\textbf{r}=[r_{x}\ r_{y}]^{T}$, the
inputs are the accelerations $a_{x}$ and $a_{y}$, which yield velocities
$v_{x}$ and $v_{y}$ aligned respectively with the $r_{x}$ and $r_{y}$ axes.
The continuous-time model in state-space is given by
$\dot{\textbf{x}}=\textbf{A}_{\textbf{c}}\textbf{x}+\textbf{B}_{\textbf{c}}\textbf{u}$,
with the state x and control u vectors given by $\textbf{x}^{T}=[r_{x}\ v_{x}\
r_{y}\ v_{y}]$ and $\mathbf{u}^{T}=[a_{x}\ a_{y}]$.

Since the MPC controller is implemented in discrete-time, we must discretize
the plant model to use it as an internal controller, which is done by Zero-
Order Hold discretization [12]. This generates a model in the form
${\textbf{x}[k+1]=\textbf{A}\textbf{x}[k]+\textbf{B}\textbf{u}[k]}$, and by
choosing a sample period of $0.8$ time units, matrices A and B are given in
(10). Here, units of length and time are assumed to be $m$ and $s$.
$\mathbf{A}=\begin{bmatrix}1&0.8&0&0\\\ 0&1&0&0\\\ 0&0&1&0.8\\\
0&0&0&1\end{bmatrix},$ $\mathbf{B}=\begin{bmatrix}0.32&0\\\ 0.8&0\\\ 0&0.32\\\
0&0.8\end{bmatrix}.$ (10)

The controller parameters are

  * •

the weight of the sum of the absolute values of the accelerations $\gamma=1$,

  * •

the maximal horizon $N_{s}=18$,

  * •

the constraints over velocities and accelerations at every time step: $-10\
m/s\leq v_{x},v_{y}\leq 10\ m/s$ and $-3\ m/s^{2}\leq a_{x},a_{y}\leq 3\
m/s^{2}$,

  * •

and the terminal set $\mathcal{Q}=\left\\{\textbf{x}\in\mathbb{R}^{4}|[14\
-v_{\epsilon}\ 0\ -v_{\epsilon}]^{T}\leq\textbf{x}\leq[15\ v_{\epsilon}\ 1\
v_{\epsilon}]^{T}\right\\}$, for $v_{\epsilon}=5\times 10^{-3}\ m/s$.

The initial state is $\textbf{x}[0]=[0\ 0\ 10\ 0]^{T}$and the scenario has
$N_{c}=45$ obstacles, which are depicted as brown rectangles in the plots.

The CPLEX toolbox from IBM ILOG was used for solving the MILP problem in
Matlab environment.

##  V RESULTS AND DISCUSSION

Three different simulations are run for comparison: (i) state-of-the-art
collision avoidance employing the constraints in Sec. II-A, herein deemed
“unclustered”, collision avoidance using the proposal in Sec. III with (ii)
$N_{c}=2$ and (iii) $N_{c}=3$ clusters, herein deemed “clustering”.

The results of the execution of the unclustered collision avoidance strategy
can be found in the trajectory of Fig. 2, with the evolution of the cost
function $J^{\star}[k]$ and its prediction $\hat{J}[k+1|k]$ depicted in Fig.
3.

Figure 2: Trajectory of the Unclustered avoidance strategy. Figure 3:
Evolution of the cost function and its prediction for the Unclustered
avoidance strategy.

Fig. 4 contains the evolution of the clustering strategy as the agent follows
towards the target considering $N_{c}=2$ clusters. Initially, in Fig. 4(a) the
clustering strategy groups $39$ obstacles into a large cluster and the
remaining $6$ obstacles into a small cluster, which leaves a small corridor
for the agent to proceed towards the target. The agent follows the initially
planned path along this corridor until $k=6$, then the obstacles are assigned
differently among clusters, freeing previously clustered space that becomes
available for maneuvering. The agent then moves downward, optimizing a shorter
path to the target, as can be seen in Fig. 4(b), and after overcoming the
obstacles a few steps later is able to reach the target set, clustering the
obstacles into a single cluster.

(a) $k=6$

(b) $k=7$

(c) $k=13$

Figure 4: Trajectory and clusters at three different time instants with the
clustering strategy with $N_{c}=2$.

The predicted $\hat{J}$ and optimal $J^{\star}$ costs at each timestep can be
seen in Fig. 5. It can be noticed that
$J^{\star}(k+1)=\hat{J}(k+1|k)=J^{\star}(k)-\|\textbf{u}(k|k)\|_{1}-1,1\leq
k\leq 6\ \text{and}\ 7<k\leq 13$, which coincides with instants where the
cluster-obstacle assignment did not change. Otherwise, if a previously
clustered space is cleared for maneuvering due to a change in the cluster-
obstacle assignment, room for shortcuts not initially considered becomes
available, which enables the agent to obtain a shorter path to the target, as
can be seen by comparing the clusters in Fig. 4(a) and Fig. 4(b). In turn,
this yields a lower cost than was predicted, as can be confirmed by
$J^{\star}(7)<\hat{J}(7|6)$ in Fig. 5. For comparison, notice that in Fig. 3
the actual optimal cost and the predicted always coincided. Indeed, the final
trajectory executed by the agent in closed-loop seen in Fig. 2 is the one
predicted at the first time instant, as neither changes in the environment nor
disturbances were considered.

Figure 5: Evolution of the solution cost and its prediction for the clustering
strategy with $N_{c}=2$ clusters.

As the number of clusters available is limited, the clustering algorithm
adjusts the cluster configuration according to the position of the agent to
allow a minimum-cost trajectory to the target, which causes the decrease in
the value of the cost function regarding what was predicted in the preceding
step. This illustrates the advantage of our proposal regarding methods that
rely on pre-processing the environment, such as the alternatives found in the
literature [4, 5, 6]. The embedding of the clustering strategy within the
optimization problem enables to decide the clusters considering the
minimization of the cost.

The results for the clustering strategy considering $N_{c}=3$ are shown in
Fig. 6. One can see that the existence of an additional cluster allows to
obtain in Fig. 6(a) a configuration that already enables a closer
circumvention of the rightmost cluster during the last eight steps to reach
the target. The previous large cluster of Fig. 4(a) is split into two
different clusters, one comprising the obstacles in the left hand-side, and
another consisting of the remaining obstacles, spreading in the bottom right-
hand side of the environment. The final trajectory and the last cluster
configuration are shown Fig. 6(b).

(a) $k=4$

(b) $k=13$

Figure 6: MILP Clustering strategy with $N_{c}=3$.

In Fig. 7 there is a representation of the cost function evolution along the
trajectory. The clustering procedure obtains clusters for the obstacles such
that there is no cost reduction due to changes in the cluster configuration
along the agent movement, as opposed to what happened previously for
$N_{c}=2$.

Figure 7: Evolution of the solution cost and its prediction for the MILP
Clustering strategy for $N_{c}=3$ clusters.

###  V-A Computation time

In accordance to what was pointed in Sec. III, the MILP Clustering strategy
has two main consequences, the first one being its impact upon the optimized
trajectory, as discussed in the last section. The second consequence is the
reduction of the computation times. Table II depicts the sum of the
computation times of every iteration until the vehicle reached the target set
and the cost in closed-loop of the final trajectories for the three cases
compared in the previous section. In the comparisons that follow, the
computation time and the cost of the unclustered strategy will be used as
benchmark, corresponding to $100\%$. The computation time with $N_{c}=2$ is
reduced to $5.5\%$ of that obtained without clustering. On the other hand, the
cost increases to $106\%$. Moreover, a compromise can be seen by the
clustering with $N_{c}=3$, whose computation time is $67\%$ of the unclustered
strategy, but the cost is only $102\%$. We conjecture that the growth in
computation time observed between the scenarios with $N_{c}=2$ and $N_{c}=3$
clusters is caused by symmetry in the optimization problem. This symmetry is
mainly due to different possibilities of assignment of obstacles to clusters
that yield the same cost.

Strategy | Computation time ($s$) | Cost in closed-loop
---|---|---
Unclustered | 337.4 | 22.4
Clustering, $N_{c}=2$ | 18.4 | 24.8
Clustering, $N_{c}=3$ | 226.2 | 22.9
TABLE II: Simulation times in a computer with $2.8$ GHz Intel i7 processor, 8
GB DDR3 RAM, and cost in closed-loop.

##  VI CONCLUSIONS

A novel strategy was proposed in this paper that clusters the obstacles within
the optimization problem solved to determine the trajectory. As a result, the
choice of the clusters considers directly the optimality of the resulting
trajectory. Moreover, this strategy was shown to greatly reduce the
computation time of the online solution of the optimization problem, while
yielding moderately higher cost. A compromise between the optimality and the
computation cost is achieved by the choice of the number of clusters to be
used.

Future work encompasses (i) applying symmetry breaking techniques to reduce
the solution search time, particularly when the number of clusters is
increased, and (ii) extending the approach to deal with cooperating multiple
agents.

## References

  * [1] A. Richards and J. How, “Aircraft trajectory planning with collision avoidance using mixed integer linear programming,” in _Proceedings of the 2002 American Control Conference (IEEE Cat. No.CH37301)_. IEEE, 2002.
  * [2] A. Richards and J. P. How, “Robust variable horizon model predictive control for vehicle maneuvering,” _International Journal of Robust and Nonlinear Control_ , vol. 16, no. 7, pp. 333–351, 2006.
  * [3] M. Garey, _Computers and intractability : a guide to the theory of NP-completeness_. San Francisco: W.H. Freeman, 1979.
  * [4] M. Vitus, V. Pradeep, G. Hoffmann, S. Waslander, and C. Tomlin, “Tunnel-MILP: Path planning with sequential convex polytopes,” in _AIAA Guidance, Navigation and Control Conference and Exhibit_. American Institute of Aeronautics and Astronautics, aug 2008.
  * [5] R. Deits and R. Tedrake, “Efficient mixed-integer planning for UAVs in cluttered environments,” in _2015 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, may 2015.
  * [6] V. A. Battagello, N. Y. Soma, and R. J. M. Afonso, “Computational load reduction of the agent guidance problem using mixed integer programming,” _PLOS ONE_ , vol. 15, no. 6, p. e0233441, jun 2020.
  * [7] I. Prodan, F. Stoican, S. Olaru, and S.-I. Niculescu, “Enhancements on the hyperplanes arrangements in mixed-integer programming techniques,” _Journal of Optimization Theory and Applications_ , vol. 154, no. 2, pp. 549–572, mar 2012.
  * [8] R. J. M. Afonso and R. K. H. Galvão, “Comments on “enhancements on the hyperplanes arrangements in mixed-integer programming techniques”,” _Journal of Optimization Theory and Applications_ , vol. 162, no. 3, pp. 996–1003, nov 2013\.
  * [9] I. Prodan, F. Stoican, S. Olaru, and S.-I. Niculescu, _Mixed-Integer Representations in Control Design_. Springer-Verlag GmbH, 2015.
  * [10] J. Bellingham, A. Richards, and J. How, “Receding horizon control of autonomous aerial vehicles,” in _Proceedings of the 2002 American Control Conference (IEEE Cat. No.CH37301)_. IEEE, 2002.
  * [11] A. Bemporad and M. Morari, “Control of systems integrating logic, dynamics, and constraints,” _Automatica_ , vol. 35, no. 3, pp. 407–427, mar 1999.
  * [12] G. Franklin, _Digital control of dynamic systems_. Menlo Park, Calif: Addison-Wesley, 1998.