# BUFFALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1

Jose M Diego, Sung Kei Li, Ashish K Meena, Anna Niemiec, Ana Acebron, Mathilde Jauzac, Mitchell F Struble, Alfred Amruth, Tom J Broadhurst, Catherine Cerny, et al.

<figure><img src="image_1.png" /></figure>

## To cite this version:

Jose M Diego, Sung Kei Li, Ashish K Meena, Anna Niemiec, Ana Acebron, et al. `BUF-FALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1`. Astronomy and Astrophysics - A&A, 2024, 681, pp.A124. 10.1051/0004-6361/202346761. `hal-04092779`

HAL Id: `hal-04092779`

https://hal.science/hal-04092779

Submitted on 5 May 2024

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# BUFFALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1

Jose M Diego, Sung Kei Li, Ashish K Meena, Anna Niemiec, Ana Acebron, Mathilde Jauzac, Mitchell F Struble, Alfred Amruth, Tom J Broadhurst, Catherine Cerny, et al.

► To cite this version:

Jose M Diego, Sung Kei Li, Ashish K Meena, Anna Niemiec, Ana Acebron, et al. `BUF-FALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1`. Astronomy & Astrophysics - A&A, 2024, 681, pp.A124. 10.1051/0004-6361/202346761. `hal-04092779`

HAL Id: hal-04092779

https://hal.science/hal-04092779v1

Submitted on 5 May 2024

**HAL** is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire **HAL**, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# BUFFALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1

Jose M. Diego<fnref n="1" />, ∗ Sung Kei Li<fnref n="2" />, † Ashish K. Meena<fnref n="3" />, ‡ Anna Niemiec<fnref n="4" />, Ana Acebron<fnref n="5" />, <fnref n="6" />, ‖ Mathilde Jauzac<fnref n="7" />, <fnref n="8" />, <fnref n="9" />, <fnref n="10" />, Mitchell F. Struble<fnref n="11" />, Alfred Amruth<fnref n="2" />, Tom J. Broadhurst<fnref n="13" />, <fnref n="14" />, <fnref n="15" />, Catherine Cerny<fnref n="7" />, <fnref n="8" />, Harald Ebeling<fnref n="16" />,

Alexei V. Filippenko<fnref n="17" />, Eric Jullo<fnref n="18" />, Patrick Kelly<fnref n="19" />, Anton M. Koekemoer<fnref n="20" />, ‡ David Lagattuta<fnref n="7" />, <fnref n="8" />, Jeremy Lim<fnref n="2" />, ‡ Marceau Limousin<fnref n="18" />, Guillaume Mahler<fnref n="7" />, <fnref n="8" />, ‡ Nency Patel<fnref n="8" />, Juan Remolina<fnref n="21" />, Johan Richard<fnref n="22" />, Keren Sharon<fnref n="21" />, Charles Steinhardt<fnref n="23" />, Keiichi Umetsu<fnref n="24" />, ‡ Liliya Williams<fnref n="25" />, ‡ Adi Zitrin<fnref n="3" />, ‡ Jose María Palencia<fnref n="1" />, Liang Dai<fnref n="26" />, ‡ Lingyuan Ji<fnref n="26" />, and Massimo Pascale<fnref n="17" />, <fnref n="12" />

Received 27 April 2023 / Accepted 1 November 2023

## ABSTRACT

In this work, we present a constraint on the abundance of supergiant (SG) stars at redshift $z \approx 1$, based on recent observations of a strongly lensed arc at this redshift. First we derived a free-form model of MACS J0416.1-2403 using data from the Beyond Ultra-deep Frontier Fields and Legacy Observations (BUFFALO) program. The new lens model is based on 72 multiply lensed galaxies that produce 214 multiple images, making it the largest sample of spectroscopically confirmed lensed galaxies on this cluster. The larger coverage in BUFFALO allowed us to measure the shear up to the outskirts of the cluster, and extend the range of lensing constraints up to ~1 Mpc from the central region, providing a mass estimate up to this radius. As an application, we make predictions for the number of high-redshift multiply lensed galaxies detected in future observations with the James Webb Space Telescope (JWST). Then we focus on a previously known lensed galaxy at $z = 1.0054$, nicknamed Spock, which contains four previously reported transients. We interpret these transients as microcaustic crossings of SG stars and explain how we computed the probability of such events. Based on simplifications regarding the stellar evolution, we find that microlensing (by stars in the intracluster medium) of SG stars at $z = 1.0054$ can fully explain these events. The inferred abundance of SG stars is consistent with either (1) a number density of stars with bolometric luminosities beyond the Humphreys-Davidson (HD) limit ($L_{max} \approx 6 \times 10^5 L_{\odot}$ for red stars), which is below ~400 stars kpc$^{-2}$, or (2) the absence of stars beyond the HD limit but with a SG number density of ~9000 kpc$^{-2}$ for stars with luminosities between $10^5 L_{\odot}$ and $6 \times 10^5 L_{\odot}$. This is equivalent to one SG star per $10 \times 10$ pc$^2$. Finally, we make predictions for future observations with JWST’s NIRcam. We find that in observations made with the F200W filter that reach 29 mag AB, if cool red SG stars exist at $z \approx 1$ beyond the HD limit, they should be easily detected in this arc.

**Key words.** gravitation – gravitational lensing: strong – supergiants

<footnote n="1">¹ Instituto de Física de Cantabria (CSIC-UC), Avda. Los Castros s/n, 39005 Santander, Spain
e-mail: jdiego@ifca.unican.es</footnote>

<footnote n="2">² Department of Physics, The University of Hong Kong, Pokfulam Road, Hong Kong</footnote>

<footnote n="3">³ Physics Department, Ben-Gurion University of the Negev, PO Box 653, Be’er-Sheva 84105, Israel</footnote>

<footnote n="4">⁴ LPNHE, CNRS/IN2P3, Sorbonne Université, Université Paris-Cité, Laboratoire de Physique Nucléaire et de Hautes Énergies, 75005 Paris, France</footnote>

<footnote n="5">⁵ Department of Physics of the University of Milano, via Celoria 16, 20133 Milano, Italy</footnote>

<footnote n="6">⁶ INAF - IASF Milano, via A. Corti 12, 20133 Milano, Italy</footnote>

<footnote n="7">⁷ Center for Extragalactic Astronomy, Durham University, South Rd, Durham DH1 3LE, UK</footnote>

<footnote n="8">⁸ Institute for Computational Cosmology, Durham University, South Road, Durham DH1 3LE, UK</footnote>

<footnote n="9">⁹ Astrophysics Research Center, University of KwaZulu-Natal, Westville Campus, Durban 4041, South Africa</footnote>

<footnote n="10">¹⁰ School of Mathematics, Statistics & Computer Science, University of KwaZulu-Natal, Westville Campus, Durban 4041, South Africa</footnote>

<footnote n="11">¹¹ Department of Physics and Astronomy, University of Pennsylvania, 209 South 33rd Street, Philadelphia, PA 19104, USA</footnote>

<footnote n="12">¹² Department of Physics, University of Hong Kong, Hong Kong, Hong Kong SAR</footnote>

<footnote n="13">¹³ Department of Physics, University of the Basque Country UPV/EHU, 48080 Bilbao, Spain</footnote>

<footnote n="14">¹⁴ DIPC, Basque Country UPV/EHU, 48080 San Sebastian, Spain</footnote>

<footnote n="15">¹⁵ Ikerbasque, Basque Foundation for Science, 48011 Bilbao, Spain</footnote>

<footnote n="16">¹⁶ Institute for Astronomy, University of Hawaii, 640 N Aohoku Pl, Hilo, HI 96720, USA</footnote>

<footnote n="17">¹⁷ Department of Astronomy, University of California, Berkeley, CA 94720-3411, USA</footnote>

<footnote n="18">¹⁸ Aix-Marseille Univ., CNRS, CNES, LAM, 13388 Marseille, France</footnote>

<footnote n="19">¹⁹ Minnesota Institute for Astrophysics, University of Minnesota, 116 Church Street SE, Minneapolis, MN 55455, USA</footnote>

<footnote n="20">²⁰ Space Telescope Science Institute, 3700 San Martin Dr., Baltimore, MD 21218, USA</footnote>

<footnote n="21">²¹ Department of Astronomy, University of Michigan, 1085 S. University Ave, Ann Arbor, MI 48109, USA</footnote>

<footnote n="22">²² Univ. Lyon, ENS de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon UMR5574, 69230 Saint-Genis-Laval, France</footnote>

<footnote n="23">²³ Niels Bohr Institute, University of Copenhagen, Lyngbyvej 2, København Ø 2100, Denmark</footnote>

<footnote n="24">²⁴ Institute of Astronomy and Astrophysics, Academia Sinica (ASIAA), AS/NTU Astronomy-Mathematics Building, No. 1, Sec. 4, Roosevelt Rd., Taipei 10617, Taiwan</footnote>

<footnote n="25">²⁵ School of Physics and Astronomy, University of Minnesota, 116 Church Street, Minneapolis, MN 55455, USA</footnote>

<footnote n="26">²⁶ Department of Physics, University of California, 366 Physics North MC 7300, Berkeley, CA 94720, USA</footnote>
# BUFFALO/Flashlights: Constraints on the abundance of lensed supergiant stars in the Spock galaxy at redshift 1

Jose M. Diego<fnref n="1" />, ∗, Sung Kei Li<fnref n="2" />, †, Ashish K. Meena<fnref n="3" />, ‡, Anna Niemiec<fnref n="4" />, Ana Acebron<fnref n="5" />, <fnref n="6" />, ‖, Mathilde Jauzac<fnref n="7" />, <fnref n="8" />, <fnref n="9" />, <fnref n="10" />, Mitchell F. Struble<fnref n="11" />, Alfred Amruth<fnref n="2" />, Tom J. Broadhurst<fnref n="13" />, <fnref n="14" />, <fnref n="15" />, Catherine Cerny<fnref n="7" />, <fnref n="8" />, Harald Ebeling<fnref n="16" />,

Alexei V. Filippenko<fnref n="17" />, Eric Jullo<fnref n="18" />, Patrick Kelly<fnref n="19" />, Anton M. Koekemoer<fnref n="20" />, †, David Lagattuta<fnref n="7" />, <fnref n="8" />, Jeremy Lim<fnref n="2" />, ‡, Marceau Limousin<fnref n="18" />, Guillaume Mahler<fnref n="7" />, <fnref n="8" />, ‖, Nency Patel<fnref n="8" />, Juan Remolina<fnref n="21" />, Johan Richard<fnref n="22" />, Keren Sharon<fnref n="21" />, Charles Steinhardt<fnref n="23" />, Keiichi Umetsu<fnref n="24" />, †, Liliya Williams<fnref n="25" />, ‡, Adi Zitrin<fnref n="3" />, ‡, Jose María Palencia<fnref n="1" />, Liang Dai<fnref n="26" />, ‡, Lingyuan Ji<fnref n="26" />, and Massimo Pascale<fnref n="17" />, <fnref n="12" />

Received 27 April 2023 / Accepted 1 November 2023

## ABSTRACT

In this work, we present a constraint on the abundance of supergiant (SG) stars at redshift $z \approx 1$, based on recent observations of a strongly lensed arc at this redshift. First we derived a free-form model of MACS J0416.1-2403 using data from the Beyond Ultra-deep Frontier Fields and Legacy Observations (BUFFALO) program. The new lens model is based on 72 multiply lensed galaxies that produce 214 multiple images, making it the largest sample of spectroscopically confirmed lensed galaxies on this cluster. The larger coverage in BUFFALO allowed us to measure the shear up to the outskirts of the cluster, and extend the range of lensing constraints up to ~1 Mpc from the central region, providing a mass estimate up to this radius. As an application, we make predictions for the number of high-redshift multiply lensed galaxies detected in future observations with the James Webb Space Telescope (JWST). Then we focus on a previously known lensed galaxy at $z = 1.0054$, nicknamed Spock, which contains four previously reported transients. We interpret these transients as microcaustic crossings of SG stars and explain how we computed the probability of such events. Based on simplifications regarding the stellar evolution, we find that microlensing (by stars in the intracluster medium) of SG stars at $z = 1.0054$ can fully explain these events. The inferred abundance of SG stars is consistent with either (1) a number density of stars with bolometric luminosities beyond the Humphreys-Davidson (HD) limit ($L_{max} \approx 6 \times 10^5 L_{\odot}$ for red stars), which is below ~400 stars kpc$^{-2}$, or (2) the absence of stars beyond the HD limit but with a SG number density of ~9000 kpc$^{-2}$ for stars with luminosities between $10^5 L_{\odot}$ and $6 \times 10^5 L_{\odot}$. This is equivalent to one SG star per $10 \times 10$ pc$^2$. Finally, we make predictions for future observations with JWST's NIRcam. We find that in observations made with the F200W filter that reach 29 mag AB, if cool red SG stars exist at $z \approx 1$ beyond the HD limit, they should be easily detected in this arc.

**Key words.** gravitation – gravitational lensing: strong – supergiants

<footnote n="1">¹ Instituto de Física de Cantabria (CSIC-UC), Avda. Los Castros s/n, 39005 Santander, Spain
e-mail: jdiego@ifca.unican.es</footnote>

<footnote n="2">² Department of Physics, The University of Hong Kong, Pokfulam Road, Hong Kong</footnote>

<footnote n="3">³ Physics Department, Ben-Gurion University of the Negev, PO Box 653, Be'er-Sheva 84105, Israel</footnote>

<footnote n="4">⁴ LPNHE, CNRS/IN2P3, Sorbonne Université, Université Paris-Cité, Laboratoire de Physique Nucléaire et de Hautes Énergies, 75005 Paris, France</footnote>

<footnote n="5">⁵ Department of Physics of the University of Milano, via Celoria 16, 20133 Milano, Italy</footnote>

<footnote n="6">⁶ INAF - IASF Milano, via A. Corti 12, 20133 Milano, Italy</footnote>

<footnote n="7">⁷ Center for Extragalactic Astronomy, Durham University, South Rd, Durham DH1 3LE, UK</footnote>

<footnote n="8">⁸ Institute for Computational Cosmology, Durham University, South Road, Durham DH1 3LE, UK</footnote>

<footnote n="9">⁹ Astrophysics Research Center, University of KwaZulu-Natal, Westville Campus, Durban 4041, South Africa</footnote>

<footnote n="10">¹⁰ School of Mathematics, Statistics & Computer Science, University of KwaZulu-Natal, Westville Campus, Durban 4041, South Africa</footnote>

<footnote n="11">¹¹ Department of Physics and Astronomy, University of Pennsylvania, 209 South 33rd Street, Philadelphia, PA 19104, USA</footnote>

<footnote n="12">¹² Department of Physics, University of Hong Kong, Hong Kong, Hong Kong SAR</footnote>

<footnote n="13">¹³ Department of Physics, University of the Basque Country UPV/EHU, 48080 Bilbao, Spain</footnote>

<footnote n="14">¹⁴ DIPC, Basque Country UPV/EHU, 48080 San Sebastian, Spain</footnote>

<footnote n="15">¹⁵ Ikerbasque, Basque Foundation for Science, 48011 Bilbao, Spain</footnote>

<footnote n="16">¹⁶ Institute for Astronomy, University of Hawaii, 640 N Aohoku Pl, Hilo, HI 96720, USA</footnote>

<footnote n="17">¹⁷ Department of Astronomy, University of California, Berkeley, CA 94720-3411, USA</footnote>

<footnote n="18">¹⁸ Aix-Marseille Univ., CNRS, CNES, LAM, 13388 Marseille, France</footnote>

<footnote n="19">¹⁹ Minnesota Institute for Astrophysics, University of Minnesota, 116 Church Street SE, Minneapolis, MN 55455, USA</footnote>

<footnote n="20">²⁰ Space Telescope Science Institute, 3700 San Martin Dr., Baltimore, MD 21218, USA</footnote>

<footnote n="21">²¹ Department of Astronomy, University of Michigan, 1085 S. University Ave, Ann Arbor, MI 48109, USA</footnote>

<footnote n="22">²² Univ. Lyon, ENS de Lyon, CNRS, Centre de Recherche Astrophysique de Lyon UMR5574, 69230 Saint-Genis-Laval, France</footnote>

<footnote n="23">²³ Niels Bohr Institute, University of Copenhagen, Lyngbyvej 2, København Ø 2100, Denmark</footnote>

<footnote n="24">²⁴ Institute of Astronomy and Astrophysics, Academia Sinica (ASIAA), AS/NTU Astronomy-Mathematics Building, No. 1, Sec. 4, Roosevelt Rd., Taipei 10617, Taiwan</footnote>

<footnote n="25">²⁵ School of Physics and Astronomy, University of Minnesota, 116 Church Street, Minneapolis, MN 55455, USA</footnote>

<footnote n="26">²⁶ Department of Physics, University of California, 366 Physics North MC 7300, Berkeley, CA 94720, USA</footnote>
# 1. Introduction

The galaxy cluster MACS J0416.1-2403 (MACS0416 hereinafter; Mann & Ebeling 2012) at redshift $z = 0.396$ is one of the six clusters studied with exquisite detail by the Hubble Space Telescope (HST) as part of the Hubble Frontier Field (HFF) program. It was selected based on its prominent lensing signature (Zitrin et al. 2013). During this program, HST observed ~10 arcmin² around the central portion of the cluster in wavelengths ranging from 0.45 µm to 1.6 µm and to a depth of ~28.5 mag in the optical and infrared (IR) bands. The HFF data from MACS0416 have been extensively used by lens modelers (Jauzac et al. 2014, 2015; Johnson et al. 2014; Diego et al. 2015; Kawamata et al. 2016; Hoag et al. 2016; Caminha et al. 2017; Richard et al. 2021; Bergamini et al. 2023) in the search for high-redshift galaxies (McLeod et al. 2015; Oesch et al. 2015; Kawamata et al. 2016; Ishigaki et al. 2018), and for the discovery of some of the first lensed stars beyond $z = 0.9$ (Rodney et al. 2018; Chen et al. 2019; Kaurov et al. 2019).

One of the goals of the HFF program was to find the most distant galaxies, taking advantage of the magnifying boost given by the gravitational lensing effect. The frequency coverage of HST allows one to potentially detect galaxies up to $z \approx 12$ thanks to the dropout technique (Steidel et al. 1999). The census of these high-redshift galaxies has provided useful information about the evolution of the ultraviolet (UV) luminosity function with redshift. However, the number of galaxies at the highest redshifts remains low, resulting in large uncertainties in the luminosity function. In order to increase the statistics, the Beyond Ultra-deep Frontier Fields and Legacy Observations (BUFFALO) program doubles the area around these clusters that is being covered by the IR detectors. In addition, the extended area around the central region of the clusters allows us to complement the strong-lensing constraints with weak-lensing measurements. The addition of weak lensing in the external parts of the cluster is important because it helps to reduce the uncertainty in the magnification, which is vital to properly derive the luminosity functions of distant lensed galaxies observed around this cluster.

MACS0416 is so far the most remarkable cluster in terms of the number of transient events, allegedly due to microlensing. This is mostly owed to the existence of two galaxies at $z \approx 1$ that are crossing the cluster caustic. Earlier HST observations have revealed at least six transient events in these two galaxies (Rodney et al. 2018; Chen et al. 2019; Kaurov et al. 2019; Kelly et al. 2022), with all of them consistent with being microlensing events<fnref n="1" />. The low redshift of these two caustic-crossing galaxies is one of the reasons behind the wealth of microlensing events of distant stars detected in this cluster. More abundant but fainter stars can be promoted by lensing beyond the detection limit with typical microlensing magnification factors of several thousand. To detect similar stars in more distant galaxies at $z \approx 2$ or beyond would require more extreme (and less likely) magnification factors in order to compensate for the increase in luminosity distance. Because of this, MACS0416 can be referred to as the “Christmas tree galaxy cluster” owing to the large number of flickering microlensing events detected behind this cluster. The main focus of this paper is to study these microlensing events in one of the two galaxies, the Spock galaxy at $z = 1.0054$ (Rodney et al. 2018; Kelly et al. 2022) which shows four microlensing events in HST data.

Given the interconnection between the microlensing events and the macromodel magnification, we divide this paper into several parts. First we discuss a new lens macromodel derived using the latest HST data from the BUFFALO program, and then we study the population of luminous stars in the Spock galaxy in more detail. Finally, we explore how microlenses in the cluster help the transients observed in this arc to be interpreted as microlensing events of the most luminous stars in the Spock galaxy.

Two of the transient events were found as part of the FrontierSN program (HST-PID 13386; PI S. Rodney), and the other two as part of the Flashlights program (Kelly et al. 2022). The FrontierSN program (Rodney 2014) is designed to produce light curves of supernovae (SNe), so it has a natural good cadence of ~1 observation per day during a period of ~1 month. The Flashlights program (Kelly et al. 2019) is designed to find flux variations of unresolved objects, and to accomplish this, it makes two deep observations with HST with the same position angle (so subtracted images minimize instrumental effects) using very wide filters (to reach fainter objects). The Flashlights program has proven to be the most successful at finding lensed star candidates with HST, with over two dozen candidates found in the HFF clusters (see Kelly et al. 2022, for a detailed description of the program and main results). The transients in the Spock arc have been interpreted as possible microlensing events of distant stars in the galaxy at $z = 1.0054$, which are temporarily being magnified by stars in the cluster (at $z = 0.396$). Other possible interpretations for the transients are considered by Rodney et al. (2018) – for instance, luminous blue variable (LBV) stars or recurrent novae –, but in this work we adopt the hypothesis that the transient nature of these events is due to changes in the magnification produced by microlensing. This choice is motivated by the fact that we do not know if these stars are intrinsically variable (LBV or nova), although many luminous stars are. However, we do know that microlenses in the lens plane are relatively numerous in this portion of the lens (from the observed intracluster light), so microlensing events are expected (Diego et al. 2018).

As shown in Fig. 1, the observed photometry during the maximum observed luminosity is also consistent with the hypothesis that the two events of Rodney et al. (2018) are very luminous stars undergoing microlensing episodes near the critical curve of the cluster. In particular, we find that one of the events is well described by a red supergiant (SG) star, while the other one is better described by a hotter SG star. Data points shown in this figure correspond to the measured photometry around the position of maximum luminosity. The photometry in the F606W filter (0.6 micron) has been corrected to account for the fact that this data point corresponds to 1 day (a half day in the rest frame) before maximum luminosity. The two events from the Flashlights program were obtained with a wide filter, so we do not have color information that can be used to test the star hypothesis by comparing the photometry with stellar models.

Under the microlensing hypothesis, the stars that one expects to detect are very luminous ($L > 10^4 L_\odot$; see for instance Kelly et al. 2018). Among these, hot blue stars can have luminosities exceeding $10^6 L_\odot$. Red (cooler) stars are typically less luminous. Local measurements show that the most luminous red SG stars, such as UY Scuti or Stephenson 2–18, have luminosities below $6 \times 10^5 L_\odot$. Similar results are found in nearby galaxies (McDonald et al. 2022). This maximum observed luminosity is known as the Humphreys-Davidson (HD) limit (Humphreys 1978).

It is unclear whether brighter red stars can exist beyond this limit, but earlier work shows that they may (Chen et al. 2015;

<footnote n="1">¹ After submission of this paper, Yan et al. (2023) expand this number with 11 additional transients found in new JWST data in these two galaxies.</footnote>
# 1. Introduction

The galaxy cluster MACS J0416.1-2403 (MACS0416 hereinafter; Mann & Ebeling 2012) at redshift $z = 0.396$ is one of the six clusters studied with exquisite detail by the Hubble Space Telescope (HST) as part of the Hubble Frontier Field (HFF) program. It was selected based on its prominent lensing signature (Zitrin et al. 2013). During this program, HST observed ~10 arcmin² around the central portion of the cluster in wavelengths ranging from 0.45 µm to 1.6 µm and to a depth of ~28.5 mag in the optical and infrared (IR) bands. The HFF data from MACS0416 have been extensively used by lens modelers (Jauzac et al. 2014, 2015; Johnson et al. 2014; Diego et al. 2015; Kawamata et al. 2016; Hoag et al. 2016; Caminha et al. 2017; Richard et al. 2021; Bergamini et al. 2023) in the search for high-redshift galaxies (McLeod et al. 2015; Oesch et al. 2015; Kawamata et al. 2016; Ishigaki et al. 2018), and for the discovery of some of the first lensed stars beyond $z = 0.9$ (Rodney et al. 2018; Chen et al. 2019; Kaurov et al. 2019).

One of the goals of the HFF program was to find the most distant galaxies, taking advantage of the magnifying boost given by the gravitational lensing effect. The frequency coverage of HST allows one to potentially detect galaxies up to $z \approx 12$ thanks to the dropout technique (Steidel et al. 1999). The census of these high-redshift galaxies has provided useful information about the evolution of the ultraviolet (UV) luminosity function with redshift. However, the number of galaxies at the highest redshifts remains low, resulting in large uncertainties in the luminosity function. In order to increase the statistics, the Beyond Ultra-deep Frontier Fields and Legacy Observations (BUFFALO) program doubles the area around these clusters that is being covered by the IR detectors. In addition, the extended area around the central region of the clusters allows us to complement the strong-lensing constraints with weak-lensing measurements. The addition of weak lensing in the external parts of the cluster is important because it helps to reduce the uncertainty in the magnification, which is vital to properly derive the luminosity functions of distant lensed galaxies observed around this cluster.

MACS0416 is so far the most remarkable cluster in terms of the number of transient events, allegedly due to microlensing. This is mostly owed to the existence of two galaxies at $z \approx 1$ that are crossing the cluster caustic. Earlier HST observations have revealed at least six transient events in these two galaxies (Rodney et al. 2018; Chen et al. 2019; Kaurov et al. 2019; Kelly et al. 2022), with all of them consistent with being microlensing events<fnref n="1" />. The low redshift of these two caustic-crossing galaxies is one of the reasons behind the wealth of microlensing events of distant stars detected in this cluster. More abundant but fainter stars can be promoted by lensing beyond the detection limit with typical microlensing magnification factors of several thousand. To detect similar stars in more distant galaxies at $z \approx 2$ or beyond would require more extreme (and less likely) magnification factors in order to compensate for the increase in luminosity distance. Because of this, MACS0416 can be referred to as the “Christmas tree galaxy cluster” owing to the large number of flickering microlensing events detected behind this cluster. The main focus of this paper is to study these microlensing events in one of the two galaxies, the Spock galaxy at $z = 1.0054$ (Rodney et al. 2018; Kelly et al. 2022) which shows four microlensing events in HST data.

Given the interconnection between the microlensing events and the macromodel magnification, we divide this paper into several parts. First we discuss a new lens macromodel derived using the latest HST data from the BUFFALO program, and then we study the population of luminous stars in the Spock galaxy in more detail. Finally, we explore how microlenses in the cluster help the transients observed in this arc to be interpreted as microlensing events of the most luminous stars in the Spock galaxy.

Two of the transient events were found as part of the FrontierSN program (HST-PID 13386; PI S. Rodney), and the other two as part of the Flashlights program (Kelly et al. 2022). The FrontierSN program (Rodney 2014) is designed to produce light curves of supernovae (SNe), so it has a natural good cadence of ~1 observation per day during a period of ~1 month. The Flashlights program (Kelly et al. 2019) is designed to find flux variations of unresolved objects, and to accomplish this, it makes two deep observations with HST with the same position angle (so subtracted images minimize instrumental effects) using very wide filters (to reach fainter objects). The Flashlights program has proven to be the most successful at finding lensed star candidates with HST, with over two dozen candidates found in the HFF clusters (see Kelly et al. 2022, for a detailed description of the program and main results). The transients in the Spock arc have been interpreted as possible microlensing events of distant stars in the galaxy at $z = 1.0054$, which are temporarily being magnified by stars in the cluster (at $z = 0.396$). Other possible interpretations for the transients are considered by Rodney et al. (2018) – for instance, luminous blue variable (LBV) stars or recurrent novae –, but in this work we adopt the hypothesis that the transient nature of these events is due to changes in the magnification produced by microlensing. This choice is motivated by the fact that we do not know if these stars are intrinsically variable (LBV or nova), although many luminous stars are. However, we do know that microlenses in the lens plane are relatively numerous in this portion of the lens (from the observed intracluster light), so microlensing events are expected (Diego et al. 2018).

As shown in Fig. 1, the observed photometry during the maximum observed luminosity is also consistent with the hypothesis that the two events of Rodney et al. (2018) are very luminous stars undergoing microlensing episodes near the critical curve of the cluster. In particular, we find that one of the events is well described by a red supergiant (SG) star, while the other one is better described by a hotter SG star. Data points shown in this figure correspond to the measured photometry around the position of maximum luminosity. The photometry in the F606W filter (0.6 micron) has been corrected to account for the fact that this data point corresponds to 1 day (a half day in the rest frame) before maximum luminosity. The two events from the Flashlights program were obtained with a wide filter, so we do not have color information that can be used to test the star hypothesis by comparing the photometry with stellar models.

Under the microlensing hypothesis, the stars that one expects to detect are very luminous ($L > 10^4 L_\odot$; see for instance Kelly et al. 2018). Among these, hot blue stars can have luminosities exceeding $10^6 L_\odot$. Red (cooler) stars are typically less luminous. Local measurements show that the most luminous red SG stars, such as UY Scuti or Stephenson 2–18, have luminosities below $6 \times 10^5 L_\odot$. Similar results are found in nearby galaxies (McDonald et al. 2022). This maximum observed luminosity is known as the Humphreys-Davidson (HD) limit (Humphreys 1978).

It is unclear whether brighter red stars can exist beyond this limit, but earlier work shows that they may (Chen et al. 2015;

<footnote n="1">¹ After submission of this paper, Yan et al. (2023) expand this number with 11 additional transients found in new JWST data in these two galaxies.</footnote>
Gilkis et al. 2021). In this case they would be detectable in arcs such as that of Spock, in regions with a sufficiently high magnification. At the redshift of the Spock galaxy, the distance modulus is 42.63 mag. A very luminous star with an absolute magnitude of -10 would have an apparent magnitude ~33–34 (depending on the star temperature and filter and ignoring specific color corrections). If this star is being magnified by a factor $\mu \approx 400$, its apparent brightness would increase by 6.5 mag, making it detectable in standard HST observations. A lack of detection of such very luminous stars implies that they cannot exist in the volume being amplified by these factors. Since no stars are observed in the Spock arc in several of the observations carried out by HST, one can use this fact to set an upper limit on the luminosity of the stars in that portion of the arc. During the more frequent periods when the star is not observed, it is logical to assume that the star is being magnified by the most likely value of the magnification. Meanwhile, the star is briefly being observed during relatively short periods when the star is crossing one of the microcaustics in the source plane. The main objective of this paper is to quantify the probability that a star is being observed during those short periods. However, in order to reach that point, we first need to start with deriving the mass distribution in the cluster (the macrolens model), followed by the microlens model, then characterizing the morphology and stellar properties of the Spock arc, and finally combining all of the previous steps to compute the probability that a star is detected above some established limiting magnitude and for a given filter.

The paper is organized as follows. In Sect. 2 we briefly discuss the observations carried out during the BUFFALO program as well as the ancillary Multi Unit Spectroscopic Explorer (MUSE) data used to derive the spectroscopic redshifts. The lensing constraints and data used to derive the lens model are discussed in Sect. 3. In Sect. 4, we present the lens model derived using only lensed systems with known redshifts (spectroscopic). Section 5 presents a simple application of the lens model, where we make predictions as to the number of observed lensed galaxies around the cluster at $z = 6$ and $z = 9$. Section 6 discusses the physical parameters derived for the Spock galaxy. In Sect. 7 we briefly describe the microlensing simulations that introduce modifications to the macromodel magnification. We combine all of the elements from the previous sections and explain how we computed the probability of seeing microlensing events in the Spock galaxy (our main result) in Sect. 8. Finally, Sects. 9 and 10 discuss these results and present our conclusions, respectively.

We adopted a standard flat cosmological model with $\Omega_m = 0.3$ and $h = 0.7$. At the redshift of the lens ($z = 0.396$), and for this cosmological model, 1 arcsec corresponds to 5.39 kpc, while at $z_s = 1.0054$, 1 arcsec corresponds to 8.13 kpc. Unless otherwise noted, magnitudes are given in the AB system (Oke & Gunn 1983). Common definitions used throughout the paper are the following. We use the term macromodel when we refer to the global lens model derived in Sect. 4. The main causatic and main critical curves refer to this model. We use the term micromodel for the perturbations introduced by microlenses in the lens plane, and the term microcaustics for the regions of divergent magnification in the source plane associated with these microlenses. Finally, when referring to the effective temperature of a star, we simply denote this temperature by $T$.

## 2. Observations

In this section, we briefly summarize the high-resolution HST imaging and the follow-up spectroscopy of MACS0416 that are used to develop the lens model presented in this work.

<figure><img src="image_3.png" /><figcaption>**Fig. 1.** Photometry of the two events (SE and NW) in the Spock galaxy reported by Rodney et al. (2018) vs. stellar models. For each event, the magnitudes in this plot correspond to the moment of maximum observed luminosity in Rodney et al. (2018) and are obtained from that work. The colored curves are synthetic models from Coelho (2014) for different star temperatures. The SE event is well described by a red SG with $T \approx 3800$ K, while the NW event is better described by a hot-ter star system with temperature $T \approx 18000$ K. Blue dots correspond to the optical bands while red dots are for the IR channels. The open circle is the magnitude measured in F606W 1 day (a half day in the rest frame) before maximum. We corrected this magnitude by 0.5 mag from the relation $dm/dt \approx 1$ mag day$^{-1}$ found in Rodney et al. (2018) where time is expressed in the rest frame of the Spock galaxy (see their Supplementary Fig. 9).</figcaption></figure>

### 2.1. Hubble Space Telescope imaging

The first HST observations of MACS0416 were performed with WFPC2 (F606W/F814W) in 2007/2008 (GO-11103; PI H. Ebeling). These early high-resolution images of MACS0416 established this cluster as a powerful lens, leading to its selection as one of the six targets in the HFF program (Lotz et al. 2017). Thanks to the HFF program, MACS0416 is among the galaxy clusters with the deepest HST observations obtained to date. The HFF observations have been extended through the “Beyond the Ultra-deep Frontier Fields And Legacy Observations” program (BUFFALO; GO-15117; PIs C. Steinhardt & M. Jauzac, Steinhardt et al. 2020).

The BUFFALO HST mosaics include the original HFF and archival data, as described by Steinhardt et al. (2020). The mosaics are produced following the approaches described by Koekemoer et al. (2011). In Fig. 2, we show the BUFFALO color image of MACS0416. The yellow rectangle highlights the area covered by the IR filters in the original HFF program. The red rectangle shows the larger area covered in the IR bands by the BUFFALO program. Our lens model is constructed over the entire $6 \times 6$ arcmin$^2$ region shown in this figure.

### 2.2. Spectroscopy

Redshifts for both member galaxies and lensed galaxies are obtained from observations made with MUSE (Bacon et al. 2010). The MUSE observations of MACS0416 used in this work were originally presented as part of the MUSE Lensing Clusters data release (Richard et al. 2021) of 12 massive, strong-lensing clusters. Data for MACS0416 were obtained through the MUSE Guaranteed Time Observations (GTO) consortium (PID 094.A-0115, PI J. Richard) and two additional General Observation (GO) proposals (PID 094.A-0525, PI F. Bauer; PID
Gilkis et al. 2021). In this case they would be detectable in arcs such as that of Spock, in regions with a sufficiently high magnification. At the redshift of the Spock galaxy, the distance modulus is 42.63 mag. A very luminous star with an absolute magnitude of -10 would have an apparent magnitude ~33–34 (depending on the star temperature and filter and ignoring specific color corrections). If this star is being magnified by a factor $\mu \approx 400$, its apparent brightness would increase by 6.5 mag, making it detectable in standard HST observations. A lack of detection of such very luminous stars implies that they cannot exist in the volume being amplified by these factors. Since no stars are observed in the Spock arc in several of the observations carried out by HST, one can use this fact to set an upper limit on the luminosity of the stars in that portion of the arc. During the more frequent periods when the star is not observed, it is logical to assume that the star is being magnified by the most likely value of the magnification. Meanwhile, the star is briefly being observed during relatively short periods when the star is crossing one of the microcaustics in the source plane. The main objective of this paper is to quantify the probability that a star is being observed during those short periods. However, in order to reach that point, we first need to start with deriving the mass distribution in the cluster (the macrolens model), followed by the microlens model, then characterizing the morphology and stellar properties of the Spock arc, and finally combining all of the previous steps to compute the probability that a star is detected above some established limiting magnitude and for a given filter.

The paper is organized as follows. In Sect. 2 we briefly discuss the observations carried out during the BUFFALO program as well as the ancillary Multi Unit Spectroscopic Explorer (MUSE) data used to derive the spectroscopic redshifts. The lensing constraints and data used to derive the lens model are discussed in Sect. 3. In Sect. 4, we present the lens model derived using only lensed systems with known redshifts (spectroscopic). Section 5 presents a simple application of the lens model, where we make predictions as to the number of observed lensed galaxies around the cluster at $z = 6$ and $z = 9$. Section 6 discusses the physical parameters derived for the Spock galaxy. In Sect. 7 we briefly describe the microlensing simulations that introduce modifications to the macromodel magnification. We combine all of the elements from the previous sections and explain how we computed the probability of seeing microlensing events in the Spock galaxy (our main result) in Sect. 8. Finally, Sects. 9 and 10 discuss these results and present our conclusions, respectively.

We adopted a standard flat cosmological model with $\Omega_m = 0.3$ and $h = 0.7$. At the redshift of the lens ($z = 0.396$), and for this cosmological model, 1 arcsec corresponds to 5.39 kpc, while at $z_s = 1.0054$, 1 arcsec corresponds to 8.13 kpc. Unless otherwise noted, magnitudes are given in the AB system (Oke & Gunn 1983). Common definitions used throughout the paper are the following. We use the term macromodel when we refer to the global lens model derived in Sect. 4. The main causatic and main critical curves refer to this model. We use the term micromodel for the perturbations introduced by microlenses in the lens plane, and the term microcaustics for the regions of divergent magnification in the source plane associated with these microlenses. Finally, when referring to the effective temperature of a star, we simply denote this temperature by $T$.

## 2. Observations

In this section, we briefly summarize the high-resolution HST imaging and the follow-up spectroscopy of MACS0416 that are used to develop the lens model presented in this work.

<figure><img src="image_2.png" /><figcaption>Fig. 1. Photometry of the two events (SE and NW) in the Spock galaxy reported by Rodney et al. (2018) vs. stellar models. For each event, the magnitudes in this plot correspond to the moment of maximum observed luminosity in Rodney et al. (2018) and are obtained from that work. The colored curves are synthetic models from Coelho (2014) for different star temperatures. The SE event is well described by a red SG with $T \approx 3800$ K, while the NW event is better described by a hot-ter star system with temperature $T \approx 18000$ K. Blue dots correspond to the optical bands while red dots are for the IR channels. The open circle is the magnitude measured in F606W 1 day (a half day in the rest frame) before maximum. We corrected this magnitude by 0.5 mag from the relation $dm/dt \approx 1$ mag day$^{-1}$ found in Rodney et al. (2018) where time is expressed in the rest frame of the Spock galaxy (see their Supplementary Fig. 9).</figcaption></figure>

### 2.1. Hubble Space Telescope imaging

The first HST observations of MACS0416 were performed with WFPC2 (F606W/F814W) in 2007/2008 (GO-11103; PI H. Ebeling). These early high-resolution images of MACS0416 established this cluster as a powerful lens, leading to its selection as one of the six targets in the HFF program (Lotz et al. 2017). Thanks to the HFF program, MACS0416 is among the galaxy clusters with the deepest HST observations obtained to date. The HFF observations have been extended through the “Beyond the Ultra-deep Frontier Fields And Legacy Observations” program (BUFFALO; GO-15117; PIs C. Steinhardt & M. Jauzac, Steinhardt et al. 2020).

The BUFFALO HST mosaics include the original HFF and archival data, as described by Steinhardt et al. (2020). The mosaics are produced following the approaches described by Koekemoer et al. (2011). In Fig. 2, we show the BUFFALO color image of MACS0416. The yellow rectangle highlights the area covered by the IR filters in the original HFF program. The red rectangle shows the larger area covered in the IR bands by the BUFFALO program. Our lens model is constructed over the entire $6 \times 6$ arcmin$^2$ region shown in this figure.

### 2.2. Spectroscopy

Redshifts for both member galaxies and lensed galaxies are obtained from observations made with MUSE (Bacon et al. 2010). The MUSE observations of MACS0416 used in this work were originally presented as part of the MUSE Lensing Clusters data release (Richard et al. 2021) of 12 massive, strong-lensing clusters. Data for MACS0416 were obtained through the MUSE Guaranteed Time Observations (GTO) consortium (PID 094.A-0115, PI J. Richard) and two additional General Observation (GO) proposals (PID 094.A-0525, PI F. Bauer; PID
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** BUFFALO data around the cluster. The image shown is 6 arcmin on a side. The yellow rectangle marks the region where previous IR data from the HFF program were taken. The red rectangle shows the expanded area covered by BUFFALO in the IR filters. This image corresponds to the combination of filters F435W, F814W, and F160W. Strong-lensing constraints are present only in the central region (yellow square) while weak-lensing constraints are spread over the entire image.</figcaption></figure>

0100.A-0763, PI E. Vanzella). To maximize the observational depth, all exposures were combined into a master dataset, and for consistency each individual frame was reduced using a common pipeline. Although this procedure largely followed the standard prescription given in the MUSE Data Reduction Pipeline User Manual<fnref n="2" /> based on the method of Weilbacher et al. (2020), some modifications were made owing to the crowded nature of cluster fields. A full description of the method and its modifications is presented in Sect. 2.5 of Richard et al. (2021) and Sect. 2.2 of Lagattuta et al. (2022). After final combination, the master

dataset has a total exposure time of 30 hr, covering an area of ~2.5 arcmin<fnref n="2" />.

Spectroscopic candidates were identified in the field in two ways. The first method (prior targets) selected all objects in the combined multiband HST imaging using SOURCE EXTRACTOR (Bertin & Arnouts 1996), with the size of each detected object slightly broadened (i.e., convolved with the MUSE point-spread function (PSF)) to account for resolution differences. In the second method (muselet sources), strong emission-line features were detected in the MUSE data directly, using the MUSELET routine included in the MUSE Python Data Analysis Framework (MPDAF; Piqueras et al. 2019). This allowed for a more

<footnote n="2">² https://www.eso.org/sci/software/pipelines/muse/</footnote>
<figure><img src="image_3.png" /><figcaption>**Fig. 2.** BUFFALO data around the cluster. The image shown is 6 arcmin on a side. The yellow rectangle marks the region where previous IR data from the HFF program were taken. The red rectangle shows the expanded area covered by BUFFALO in the IR filters. This image corresponds to the combination of filters F435W, F814W, and F160W. Strong-lensing constraints are present only in the central region (yellow square) while weak-lensing constraints are spread over the entire image.</figcaption></figure>

0100.A-0763, PI E. Vanzella). To maximize the observational depth, all exposures were combined into a master dataset, and for consistency each individual frame was reduced using a common pipeline. Although this procedure largely followed the standard prescription given in the MUSE Data Reduction Pipeline User Manual<fnref n="2" /> based on the method of Weilbacher et al. (2020), some modifications were made owing to the crowded nature of cluster fields. A full description of the method and its modifications is presented in Sect. 2.5 of Richard et al. (2021) and Sect. 2.2 of Lagattuta et al. (2022). After final combination, the master

dataset has a total exposure time of 30 hr, covering an area of ~2.5 arcmin<fnref n="2" />.

Spectroscopic candidates were identified in the field in two ways. The first method (prior targets) selected all objects in the combined multiband HST imaging using SOURCE EXTRACTOR (Bertin & Arnouts 1996), with the size of each detected object slightly broadened (i.e., convolved with the MUSE point-spread function (PSF)) to account for resolution differences. In the second method (muselet sources), strong emission-line features were detected in the MUSE data directly, using the MUSELET routine included in the MUSE Python Data Analysis Framework (MPDAF; Piqueras et al. 2019). This allowed for a more

<footnote n="2">² https://www.eso.org/sci/software/pipelines/muse/</footnote>
complete search of multiple-image systems, identifying objects
that were not directly visible in the HST data alone. Regard-
less of the detection method used, a spectrum of each candi-
date was extracted using an optimally weighted (Horne 1986)
algorithm. An initial redshift for each extraction was estimated
using the MARZ software (Hinton et al. 2016), though these were
later refined by human inspection. A finalized, master list of
all redshifts in the MACS0416 MUSE footprint is included in
Richard et al. (2021).

**3. Strong- and weak-lensing constraints**

The strong- and weak-lensing catalogs have been compiled
within the BUFFALO collaboration following the procedure
outlined by Niemiec et al. (2023). Here we only summarize
the main aspects of the process, and refer the reader to that
publication for a detailed step-by-step description. The strong-
lensing constraints were compiled from the existing litera-
ture (Jauzac et al. 2014; Johnson et al. 2014; Diego et al. 2015;
Kawamata et al. 2016; Caminha et al. 2017; Richard et al. 2021)
to produce a self-consistent and homogeneous catalog. All multi-
ple images are revoted within the collaboration and are assigned
a grade (gold, silver, or bronze) according to their reliability,
with gold systems being the most reliable (spectroscopically
confirmed), silver systems being less reliable (they lack spec-
troscopic confirmation but are good candidates based on con-
sistency in color, morphology, and photometric redshift), and
bronze being the least reliable but still potentially correct sys-
tems (as for silver systems, bronze systems lack spectroscopic
confirmation and are usually unresolved, lacking morphologi-
cal information). The position of each image is then visually
inspected by several members of the collaboration, to ensure
that the chosen positions (e.g., the emission knot, brightest pixel,
etc.) are consistent within each multiple-image system. In this
analysis, we only use the constraints labeled as gold to derive
the lens model, which corresponds to 72 multiply lensed galax-
ies, resulting in 214 images when lensed by the cluster. This can
be compared to the recent analysis of Bergamini et al. (2023)
which used 60 individual galaxies that are multiply lensed. Two
of the Bergamini et al. (2023) systems in the northeast (NE) por-
tion of the cluster are excluded from our analysis, not because
they are dubious systems, but simply because our lens mod-
els were derived before the results of Bergamini et al. (2023)
were published. Adding the two Bergamini et al. (2023) systems
would certainly improve our lens model, but not substantially
since the region of the lens plane where these two new systems
are found is already well constrained by other systems nearby.
The locations of our constraints (and the new ones discovered
by Bergamini et al. 2023) are shown in Fig. 3. As seen there,
the new systems of Bergamini et al. (2023; marked as red cir-
cles) are not close to the Spock arc, so not adding them in our
lens model has negligible impact on our main results. However,
future lens models should definitely include these systems, so
we add them in our curated list of constraints in the Appendix,
and for the convenience of future lens modeling efforts with this
cluster.

The weak-lensing (WL) catalog is produced from the
reduced BUFFALO data as follows. The galaxy shapes are
measured from their second- and fourth-order normalized mul-
tipole moments, using the publicly available PYRRG code
(Harvey et al. 2021). All objects (background and foreground
galaxies, stars), are first extracted from the mosaic images using
two runs of SOURCE EXTRACTOR (Bertin & Arnouts 1996) with
different detection threshold and kernel sizes, in order to best

extract both the small and large sources, with the so-called
"hot/cold" method. Galaxies and stars are then separated from
the image artifacts with a maximum surface brightness (μ<sub>max</sub>)
vs. magnitude diagram: the shape of galaxies is what carries
the weak-lensing signal, while the stars are used to model the
PSF and correct the measured galaxy moments. Cuts are then
applied to clean the obtained galaxy-shape catalog, and (i) mask
the edges and sources contaminated by the light of massive clus-
ter members or saturated stars, (ii) remove the galaxies that are
too small or too faint to have a proper shape measurement, and
(iii) remove the sources with ill-measured ellipticities or errors.

The last step is to select galaxies that are located behind the
lensing cluster, as foreground galaxies are not weakly lensed and
would dilute the measured lensing signal. To perform this sepa-
ration, we apply color–color or color–magnitude cuts, depending
on the number of photometric bands available for each source.
Galaxies located within both the BUFFALO ACS and WFC3
fields of view are selected in a (m<sub>814W</sub>−m<sub>160W</sub> vs. m<sub>606W</sub>−m<sub>160W</sub>)
diagram, while galaxies only covered by ACS observations in
a (m<sub>606W</sub> − m<sub>814W</sub>) vs. m<sub>814W</sub> diagram. The cuts are calibrated
using the subsample of galaxies that have spectroscopically mea-
sured redshifts. The final background galaxy catalog contains
3235 sources. The number density of weak-lensing measure-
ments varies across the field of view. We bin these galaxies
in 1 × 1 arcmin<sup>2</sup> bins and obtain a mean surface density of
43.6 galaxies per arcmin<sup>2</sup> with dispersion σ = 18.9 galaxies per
arcmin<sup>2</sup>. The central region of the cluster contains no WL mea-
surements owing to contamination from the cluster itself, but this
region is well constrained by the numerous strong-lensing arcs.
We consider a region of 6 × 6 arcmin<sup>2</sup>, hence the number of WL
constraints is 72 (6 × 6 = 36 bins, each with two constraints γ<sub>1</sub>
and γ<sub>2</sub>).

**4. The derived lens model**

We use the free-form code WSLAP+ (Diego et al. 2005, 2007)
to combine the WL and SL constraints from the previous section
and derive the mass distribution. This code has been described in
detail in earlier work, so we refer the reader to earlier references
for details. A basic description, including all relevant references,
is given in Appendix A. The model is available through the BUF-
FALO website<fnref n="3" />.

The resulting model shows the characteristic bimodal distri-
bution found in earlier work, with two peaks of mass centered in
each of the BCG galaxies. The critical curves for *z* = 1, 2, 4, and
10 are shown in Fig. 3. The critical curves are very elongated and
intersect multiple arcs at their respective redshifts. These inter-
section points offer unique opportunities to zoom in the highly
magnified background galaxies, as discussed in more detail in
Sect. 6.

Regarding the mass profile, we find that the two main clumps
of mass are very similar to each other in terms of density pro-
files. These are shown in Fig. 4. We compute one profile cen-
tered in the NE BCG galaxy and a second profile centered in
the southwest (SW) BCG galaxy. The profiles represent the pro-
jected mass along the line of sight in a cylindrical aperture. The
most remarkable difference between the two profiles is that the
SW BCG is more massive within the central kpc, but beyond a
few kpc both profiles are very similar. We fit a generic Navarro-
Frenk-White profile (gNFW), finding a reasonable fit for a model
with parameters γ = 0.5, α = 2.5, β = 3, and scale radius
r_s = 80 kpc, where γ, α, and β are the slopes of the density in

<footnote n="3">³ https://buffalo.ipac.caltech.edu/</footnote>
complete search of multiple-image systems, identifying objects
that were not directly visible in the HST data alone. Regard-
less of the detection method used, a spectrum of each candi-
date was extracted using an optimally weighted (Horne 1986)
algorithm. An initial redshift for each extraction was estimated
using the MARZ software (Hinton et al. 2016), though these were
later refined by human inspection. A finalized, master list of
all redshifts in the MACS0416 MUSE footprint is included in
Richard et al. (2021).

**3. Strong- and weak-lensing constraints**

The strong- and weak-lensing catalogs have been compiled
within the BUFFALO collaboration following the procedure
outlined by Niemiec et al. (2023). Here we only summarize
the main aspects of the process, and refer the reader to that
publication for a detailed step-by-step description. The strong-
lensing constraints were compiled from the existing litera-
ture (Jauzac et al. 2014; Johnson et al. 2014; Diego et al. 2015;
Kawamata et al. 2016; Caminha et al. 2017; Richard et al. 2021)
to produce a self-consistent and homogeneous catalog. All multi-
ple images are revoted within the collaboration and are assigned
a grade (gold, silver, or bronze) according to their reliability,
with gold systems being the most reliable (spectroscopically
confirmed), silver systems being less reliable (they lack spec-
troscopic confirmation but are good candidates based on con-
sistency in color, morphology, and photometric redshift), and
bronze being the least reliable but still potentially correct sys-
tems (as for silver systems, bronze systems lack spectroscopic
confirmation and are usually unresolved, lacking morphologi-
cal information). The position of each image is then visually
inspected by several members of the collaboration, to ensure
that the chosen positions (e.g., the emission knot, brightest pixel,
etc.) are consistent within each multiple-image system. In this
analysis, we only use the constraints labeled as gold to derive
the lens model, which corresponds to 72 multiply lensed galax-
ies, resulting in 214 images when lensed by the cluster. This can
be compared to the recent analysis of Bergamini et al. (2023)
which used 60 individual galaxies that are multiply lensed. Two
of the Bergamini et al. (2023) systems in the northeast (NE) por-
tion of the cluster are excluded from our analysis, not because
they are dubious systems, but simply because our lens mod-
els were derived before the results of Bergamini et al. (2023)
were published. Adding the two Bergamini et al. (2023) systems
would certainly improve our lens model, but not substantially
since the region of the lens plane where these two new systems
are found is already well constrained by other systems nearby.
The locations of our constraints (and the new ones discovered
by Bergamini et al. 2023) are shown in Fig. 3. As seen there,
the new systems of Bergamini et al. (2023; marked as red cir-
cles) are not close to the Spock arc, so not adding them in our
lens model has negligible impact on our main results. However,
future lens models should definitely include these systems, so
we add them in our curated list of constraints in the Appendix,
and for the convenience of future lens modeling efforts with this
cluster.

The weak-lensing (WL) catalog is produced from the
reduced BUFFALO data as follows. The galaxy shapes are
measured from their second- and fourth-order normalized mul-
tipole moments, using the publicly available PYRRG code
(Harvey et al. 2021). All objects (background and foreground
galaxies, stars), are first extracted from the mosaic images using
two runs of SOURCE EXTRACTOR (Bertin & Arnouts 1996) with
different detection threshold and kernel sizes, in order to best

extract both the small and large sources, with the so-called
"hot/cold" method. Galaxies and stars are then separated from
the image artifacts with a maximum surface brightness (μ<sub>max</sub>)
vs. magnitude diagram: the shape of galaxies is what carries
the weak-lensing signal, while the stars are used to model the
PSF and correct the measured galaxy moments. Cuts are then
applied to clean the obtained galaxy-shape catalog, and (i) mask
the edges and sources contaminated by the light of massive clus-
ter members or saturated stars, (ii) remove the galaxies that are
too small or too faint to have a proper shape measurement, and
(iii) remove the sources with ill-measured ellipticities or errors.

The last step is to select galaxies that are located behind the
lensing cluster, as foreground galaxies are not weakly lensed and
would dilute the measured lensing signal. To perform this sepa-
ration, we apply color–color or color–magnitude cuts, depending
on the number of photometric bands available for each source.
Galaxies located within both the BUFFALO ACS and WFC3
fields of view are selected in a (m<sub>814W</sub>−m<sub>160W</sub> vs. m<sub>606W</sub>−m<sub>160W</sub>)
diagram, while galaxies only covered by ACS observations in
a (m<sub>606W</sub> − m<sub>814W</sub>) vs. m<sub>814W</sub> diagram. The cuts are calibrated
using the subsample of galaxies that have spectroscopically mea-
sured redshifts. The final background galaxy catalog contains
3235 sources. The number density of weak-lensing measure-
ments varies across the field of view. We bin these galaxies
in 1 × 1 arcmin<sup>2</sup> bins and obtain a mean surface density of
43.6 galaxies per arcmin<sup>2</sup> with dispersion σ = 18.9 galaxies per
arcmin<sup>2</sup>. The central region of the cluster contains no WL mea-
surements owing to contamination from the cluster itself, but this
region is well constrained by the numerous strong-lensing arcs.
We consider a region of 6 × 6 arcmin<sup>2</sup>, hence the number of WL
constraints is 72 (6 × 6 = 36 bins, each with two constraints γ<sub>1</sub>
and γ<sub>2</sub>).

**4. The derived lens model**

We use the free-form code WSLAP+ (Diego et al. 2005, 2007)
to combine the WL and SL constraints from the previous section
and derive the mass distribution. This code has been described in
detail in earlier work, so we refer the reader to earlier references
for details. A basic description, including all relevant references,
is given in Appendix A. The model is available through the BUF-
FALO website<fnref n="3" />.

The resulting model shows the characteristic bimodal distri-
bution found in earlier work, with two peaks of mass centered in
each of the BCG galaxies. The critical curves for *z* = 1, 2, 4, and
10 are shown in Fig. 3. The critical curves are very elongated and
intersect multiple arcs at their respective redshifts. These inter-
section points offer unique opportunities to zoom in the highly
magnified background galaxies, as discussed in more detail in
Sect. 6.

Regarding the mass profile, we find that the two main clumps
of mass are very similar to each other in terms of density pro-
files. These are shown in Fig. 4. We compute one profile cen-
tered in the NE BCG galaxy and a second profile centered in
the southwest (SW) BCG galaxy. The profiles represent the pro-
jected mass along the line of sight in a cylindrical aperture. The
most remarkable difference between the two profiles is that the
SW BCG is more massive within the central kpc, but beyond a
few kpc both profiles are very similar. We fit a generic Navarro-
Frenk-White profile (gNFW), finding a reasonable fit for a model
with parameters γ = 0.5, α = 2.5, β = 3, and scale radius
r_s = 80 kpc, where γ, α, and β are the slopes of the density in

<footnote n="3">³ https://buffalo.ipac.caltech.edu/</footnote>
<figure><img src="image_5.png" /><figcaption>Fig. 3. Central region with arc positions. Yellow circles are the 72 systems identified in BUFFALO data and with spectroscopic redshifts. Red circles are systems not included in the previous sample but listed by Bergamini et al. (2023). The critical curves for different redshifts are shown (z = 1, 2, 4, and 10; blue, green, yellow, and red, respectively). This image is a composite made with the seven HST filters (F435W, F606W, F814W, F105W, F125W, F140W, and F160W). The white rectangle marks the position of the Spock arc. The small white square marks the position of the third counterimage of the Spock arc. A zoomed-in version is shown at top right, where the counterimage is contained within the yellow circle. The radius of this circles is 0.25".</figcaption></figure>

the innermost, intermediate, and outermost regions of the cluster, respectively (see Nagai et al. 2007, for a description of this model). The fit is done centering on each halo so it does not account for superposition effects of one halo over the other. This profile fails at reproducing the inner ~20 kpc region, but this is where baryons are expected to cool more efficiently than dark matter and form a baryonic cusp. At larger radii, our lens model benefits from the addition of weak-lensing data that extends the range of constraints up to the boundaries of the field of view shown in Fig. 2.

In terms of integrated mass as a function of radius, we show the cumulative mass profile in the inset of the same figure. The blue and red curves are again for the NE and SW clumps, while the dashed line is for the gNFW profile. In this plot we add a fourth curve (black solid line) that corresponds to a third profile centered in the middle point between the two BCG galaxies (the separation between these galaxies is ~40"). We find that the total mass contained within 1 Mpc for the black continuous line is $4.9 \times 10^{14} M_{\odot}$. Bergamini et al. (2023) derive a parametric model of this cluster and show how the integrated mass for each halo is also very similar to each other. In Fig. 4, we show a direct comparison of the profiles (centered in each halo) for our model and theirs. The models agree well in the intermediate region (30–200 kpc) where strong lensing constraints are available. Our model traces the light in the very central region while the model of Bergamini et al. (2023) is steeper. The opposite trend is found in the outskirts of the cluster. Both regimes are unconstrained by strong-lensing data, so it is natural to see the largest deviations in the innermost and outermost regions of the cluster. Quantitatively, from their Fig. 7, they demonstrate how the integrated cylindrical mass around each halo, measured at a radius of 93 kpc, is ~4.8 $\times 10^{13} M_{\odot}$. For the same

radius, we find masses $5.6 \times 10^{13} M_{\odot}$ and $5.2 \times 10^{13} M_{\odot}$ for the NE and SW halos (respectively), which is 5%–10% higher. At larger radii, Bergamini et al. (2023) find a mass for each halo of $\sim 3.4 \times 10^{14} M_{\odot}$ at R = 400 kpc. At this distance we find $3.6 \times 10^{14} M_{\odot}$ for each, the NE and SW halos.

## 5. Expected number of multiply lensed galaxies at high redshift in JWST images

In this section, we estimate the expected number of lensed galaxies observed at z = 6 and z = 9, based on the magnification prediction from the lens model. For sources at z = 6 the corresponding curve $A(>\mu)$ is almost identical to that shown in Fig. 5, with the amplitude of this curve being just ~14% smaller than the curve for z = 9 (based on the scaling with z of the lensing strength factor, Eq. (2.4) in Seitz & Schneider 1997). Here we consider only the case of multiple images, which takes place when the total magnification is above $\mu \approx 3$. Outside the caustic region the total magnification is below 2.5, so this assumption is well justified.

At a given redshift, the number of lensed galaxies can be expressed in the form

$$ \frac{dN}{d\mu dM_{UV}} = \Phi(z, \tilde{M}_{UV}) * dV(z)N_{ci} \frac{dA}{d\mu}, \qquad (1) $$

where $M_{UV}$ is the absolute UV luminosity of the galaxy and $\tilde{M}_{UV} = M_{UV} - 2.5 \log_{10}(\mu)$ is the magnitude before magnification. $\Phi(z, M_{UV})$ is the luminosity function at redshift z (expressed in the usual units of number per magnitude per Mpc$^3$), and $dA/d\mu$ is the area in the source plane (expressed in arcseconds squared) with total magnification $(\mu - 0.5) < \mu < (\mu + 0.5)$. Since
<figure><img src="image_4.png" /><figcaption>Fig. 3. Central region with arc positions. Yellow circles are the 72 systems identified in BUFFALO data and with spectroscopic redshifts. Red circles are systems not included in the previous sample but listed by Bergamini et al. (2023). The critical curves for different redshifts are shown (z = 1, 2, 4, and 10; blue, green, yellow, and red, respectively). This image is a composite made with the seven HST filters (F435W, F606W, F814W, F105W, F125W, F140W, and F160W). The white rectangle marks the position of the Spock arc. The small white square marks the position of the third counterimage of the Spock arc. A zoomed-in version is shown at top right, where the counterimage is contained within the yellow circle. The radius of this circles is 0.25".</figcaption></figure>

the innermost, intermediate, and outermost regions of the cluster, respectively (see Nagai et al. 2007, for a description of this model). The fit is done centering on each halo so it does not account for superposition effects of one halo over the other. This profile fails at reproducing the inner ~20 kpc region, but this is where baryons are expected to cool more efficiently than dark matter and form a baryonic cusp. At larger radii, our lens model benefits from the addition of weak-lensing data that extends the range of constraints up to the boundaries of the field of view shown in Fig. 2.

In terms of integrated mass as a function of radius, we show the cumulative mass profile in the inset of the same figure. The blue and red curves are again for the NE and SW clumps, while the dashed line is for the gNFW profile. In this plot we add a fourth curve (black solid line) that corresponds to a third profile centered in the middle point between the two BCG galaxies (the separation between these galaxies is ~40"). We find that the total mass contained within 1 Mpc for the black continuous line is $4.9 \times 10^{14} M_{\odot}$. Bergamini et al. (2023) derive a parametric model of this cluster and show how the integrated mass for each halo is also very similar to each other. In Fig. 4, we show a direct comparison of the profiles (centered in each halo) for our model and theirs. The models agree well in the intermediate region (30–200 kpc) where strong lensing constraints are available. Our model traces the light in the very central region while the model of Bergamini et al. (2023) is steeper. The opposite trend is found in the outskirts of the cluster. Both regimes are unconstrained by strong-lensing data, so it is natural to see the largest deviations in the innermost and outermost regions of the cluster. Quantitatively, from their Fig. 7, they demonstrate how the integrated cylindrical mass around each halo, measured at a radius of 93 kpc, is ~4.8 $\times 10^{13} M_{\odot}$. For the same

radius, we find masses $5.6 \times 10^{13} M_{\odot}$ and $5.2 \times 10^{13} M_{\odot}$ for the NE and SW halos (respectively), which is 5%–10% higher. At larger radii, Bergamini et al. (2023) find a mass for each halo of $\sim 3.4 \times 10^{14} M_{\odot}$ at R = 400 kpc. At this distance we find $3.6 \times 10^{14} M_{\odot}$ for each, the NE and SW halos.

## 5. Expected number of multiply lensed galaxies at high redshift in JWST images

In this section, we estimate the expected number of lensed galaxies observed at z = 6 and z = 9, based on the magnification prediction from the lens model. For sources at z = 6 the corresponding curve $A(>\mu)$ is almost identical to that shown in Fig. 5, with the amplitude of this curve being just ~14% smaller than the curve for z = 9 (based on the scaling with z of the lensing strength factor, Eq. (2.4) in Seitz & Schneider 1997). Here we consider only the case of multiple images, which takes place when the total magnification is above $\mu \approx 3$. Outside the caustic region the total magnification is below 2.5, so this assumption is well justified.

At a given redshift, the number of lensed galaxies can be expressed in the form

$$ \frac{dN}{d\mu dM_{UV}} = \Phi(z, \tilde{M}_{UV}) * dV(z)N_{ci} \frac{dA}{d\mu}, \quad (1) $$

where $M_{UV}$ is the absolute UV luminosity of the galaxy and $\tilde{M}_{UV} = M_{UV} - 2.5 \log_{10}(\mu)$ is the magnitude before magnification. $\Phi(z, M_{UV})$ is the luminosity function at redshift z (expressed in the usual units of number per magnitude per Mpc³), and $dA/d\mu$ is the area in the source plane (expressed in arcseconds squared) with total magnification $(\mu - 0.5) < \mu < (\mu + 0.5)$. Since
<figure><img src="image_6.png" /></figure>

**Fig. 4.** Mass profile of the lens model. The blue and red solid curves show the profiles when adopting the NE and SW BCGs as central points, respectively. The mass profiles are given in dimensionless $\kappa$ units assuming a source at $z = 3$. The dashed line is a gNFW profile with parameters $\gamma = 0.5$, $\alpha = 2.5$, $\beta = 3$, and scale radius $r_s = 80$ kpc. The inset shows the integrated mass as a function of radius in units of $10^{14} M_{\odot}$, and for the same three profiles shown in the larger plot. The black dashed line is the integrated profile when the center is chosen as the middle point between the two BCG galaxies (RA = 4h 16m 08.428s, Dec = -24° 04' 21.0"). For comparison, we show as light blue and orange the profiles from Bergamini et al. (2023).

<figure><img src="image_7.png" /></figure>

**Fig. 5.** Area above magnification $\mu$ at $z = 9$. The dashed line is computed in the source plane from the caustic map. The dotted line is the area computed in the image plane from the magnification map corrected by the factor $\mu$. The solid line is similar to the dotted line but divided by a factor of 2 (in order to account for the contribution from the double image at large magnification factors). At large magnification factors, the solid curve falls as the expected $\mu^{-2}$ power law. For comparison, we show as an orange dotted line the same lensing probability for the model of Bergamini et al. (2023) at $z = 9$.

the area $dA/d\mu$ is given in arcsec$^2$, the volume $dV(z)$ corresponds to the volume within 1 arcsec$^2$ and ($z - 0.5$) < $z$ < ($z + 0.5$). The chosen binning in redshift is arbitrary. Smaller values of $\delta z$ provide more accurate results, but for our purposes a redshift binning of $\delta z = 1$ will suffice. Finally, $N_{\text{ci}}$ is the number of detectable counterimages produced by each lensed galaxy. We make the simple approximation that the total magnification $\mu$ of the galaxy is distributed into three counterimages. This is a general statement which is valid for most lensed galaxies. When the galaxy is in the central valley of the cluster caustic region, the

three counterimages have similar magnification. The central valley for MACS0416 has typical magnification $4 < \mu < 12$. Hence, we assume that galaxies falling in this central-valley region form three counterimages, each having magnification $\mu/3$. To differentiate from the valley region, we define the region with magnification $\mu > 12$ as the peak region. For typical clusters, the peak region is just a few kpc away from the caustics. For galaxies in the peak region, there are still three counterimages, but for galaxies near fold caustics (the vast majority), one of the counterimages typically has smaller magnification factors than counterimages produced in the valley region, while the other two counterimages carry most of the magnification. Hence, we do the additional approximation that galaxies located in the valley region form three detectable images ($N_{\text{ci}} = 3$), with moderate magnification factors, while galaxies in the peak region form two detectable images ($N_{\text{ci}} = 2$) with large magnification factors.

For the luminosity function of background galaxies, $\Phi(z, M_{\text{UV}})$, we follow Bowler et al. (2015) for $z = 6$ and Harikane et al. (2023), Pérez-González et al. (2023) for $z = 9$. Using these luminosity functions, we consider three different types of observations: (i) shallow, Where observations can detect galaxies brighter than $M_{\text{UV}} = -17$ mag; (ii) medium, Where galaxies with $M_{\text{UV}} = -16$ mag can be detected; and (iii) deep, reaching down to $M_{\text{UV}} = -15$ mag. These absolute magnitudes correspond to those before correcting for amplification. That is, if a galaxy is inferred to have $M_{\text{UV}} = -15$ mag, but is magnified by a factor $\mu = 10$, the magnification-corrected magnitude would be $M_{\text{UV}} = -12.5$ mag. The distance modulus for a source at $z = 6$ and $z = 9$ is 44.63 and 44.87 mag, respectively. Hence, a galaxy with $M_{\text{UV}} = -15$ mag at $z = 9$ would have an apparent magnitude (in the IR) of 29.87, within reach of deep observations made by JWST (ignoring color corrections and extinction).

The results are summarized in Fig. 6. In the larger plot, we show the number of lensed galaxies as a function of their total magnification – that is, the magnification experienced by the galaxy in the source plane. As mentioned earlier, in the image plane, each lensed galaxy usually produces between two and three detectable counterimages$^4$.

In our case, the transition between two and three detectable counterimages takes place when the total magnification is $\mu = 12$. This transition at $\mu = 12$ can be easily appreciated in Fig. 6, which shows the results obtained for the two redshifts ($z = 6$ and 9) and the three different depths (from bottom to top: $M_{\text{UV}} = -17, -16, \text{and } -15$ mag). As expected, most lensed galaxies have relatively small magnification factors simply because smaller magnification factors correspond to larger areas in the source plane. For the deepest observations ($M_{\text{UV}} = -15$ mag), we expect ~1 strongly lensed galaxy with total magnification $\mu \approx 30$. This is a galaxy near a caustic or even with parts of it crossing a caustic. Large magnification factors are more likely for high-redshift galaxies owing to their smaller intrinsic size. Recent results from JWST reveal very small sizes for faint galaxies at $z \approx 9$ with estimated sizes a few tens of parsecs (Williams et al. 2023). These galaxies appear unresolved even in JWST images. With tangential magnification factors of order 10, these galaxies can be resolved, so we expect that one such galaxy may be detectable in deep observations of MACS0416 with JWST. Surprisingly, from the same figure, this is the same number of galaxies we expect to observe with the same magnification factor, but at $z = 6$. This is due to the fact that the adopted

<footnote n="4">⁴ More counterimages may appear owing to substructure or when the galaxy is also producing radial images, but this situation is not considered here.</footnote>
<figure><img src="image_5.png" /></figure>

**Fig. 4.** Mass profile of the lens model. The blue and red solid curves show the profiles when adopting the NE and SW BCGs as central points, respectively. The mass profiles are given in dimensionless $\kappa$ units assuming a source at $z = 3$. The dashed line is a gNFW profile with parameters $\gamma = 0.5$, $\alpha = 2.5$, $\beta = 3$, and scale radius $r_s = 80$ kpc. The inset shows the integrated mass as a function of radius in units of $10^{14} M_{\odot}$, and for the same three profiles shown in the larger plot. The black dashed line is the integrated profile when the center is chosen as the middle point between the two BCG galaxies (RA = 4h 16m 08.428s, Dec = -24° 04' 21.0"). For comparison, we show as light blue and orange the profiles from Bergamini et al. (2023).

<figure><img src="image_6.png" /></figure>

**Fig. 5.** Area above magnification $\mu$ at $z = 9$. The dashed line is computed in the source plane from the caustic map. The dotted line is the area computed in the image plane from the magnification map corrected by the factor $\mu$. The solid line is similar to the dotted line but divided by a factor of 2 (in order to account for the contribution from the double image at large magnification factors). At large magnification factors, the solid curve falls as the expected $\mu^{-2}$ power law. For comparison, we show as an orange dotted line the same lensing probability for the model of Bergamini et al. (2023) at $z = 9$.

the area $dA/d\mu$ is given in arcsec$^2$, the volume $dV(z)$ corresponds to the volume within 1 arcsec$^2$ and ($z - 0.5$) < $z$ < ($z + 0.5$). The chosen binning in redshift is arbitrary. Smaller values of $\delta z$ provide more accurate results, but for our purposes a redshift binning of $\delta z = 1$ will suffice. Finally, $N_{\text{ci}}$ is the number of detectable counterimages produced by each lensed galaxy. We make the simple approximation that the total magnification $\mu$ of the galaxy is distributed into three counterimages. This is a general statement which is valid for most lensed galaxies. When the galaxy is in the central valley of the cluster caustic region, the

three counterimages have similar magnification. The central valley for MACS0416 has typical magnification $4 < \mu < 12$. Hence, we assume that galaxies falling in this central-valley region form three counterimages, each having magnification $\mu/3$. To differentiate from the valley region, we define the region with magnification $\mu > 12$ as the peak region. For typical clusters, the peak region is just a few kpc away from the caustics. For galaxies in the peak region, there are still three counterimages, but for galaxies near fold caustics (the vast majority), one of the counterimages typically has smaller magnification factors than counterimages produced in the valley region, while the other two counterimages carry most of the magnification. Hence, we do the additional approximation that galaxies located in the valley region form three detectable images ($N_{\text{ci}} = 3$), with moderate magnification factors, while galaxies in the peak region form two detectable images ($N_{\text{ci}} = 2$) with large magnification factors.

For the luminosity function of background galaxies, $\Phi(z, M_{\text{UV}})$, we follow Bowler et al. (2015) for $z = 6$ and Harikane et al. (2023), Pérez-González et al. (2023) for $z = 9$. Using these luminosity functions, we consider three different types of observations: (i) shallow, Where observations can detect galaxies brighter than $M_{\text{UV}} = -17$ mag; (ii) medium, Where galaxies with $M_{\text{UV}} = -16$ mag can be detected; and (iii) deep, reaching down to $M_{\text{UV}} = -15$ mag. These absolute magnitudes correspond to those before correcting for amplification. That is, if a galaxy is inferred to have $M_{\text{UV}} = -15$ mag, but is magnified by a factor $\mu = 10$, the magnification-corrected magnitude would be $M_{\text{UV}} = -12.5$ mag. The distance modulus for a source at $z = 6$ and $z = 9$ is 44.63 and 44.87 mag, respectively. Hence, a galaxy with $M_{\text{UV}} = -15$ mag at $z = 9$ would have an apparent magnitude (in the IR) of 29.87, within reach of deep observations made by JWST (ignoring color corrections and extinction).

The results are summarized in Fig. 6. In the larger plot, we show the number of lensed galaxies as a function of their total magnification – that is, the magnification experienced by the galaxy in the source plane. As mentioned earlier, in the image plane, each lensed galaxy usually produces between two and three detectable counterimages$^4$.

In our case, the transition between two and three detectable counterimages takes place when the total magnification is $\mu = 12$. This transition at $\mu = 12$ can be easily appreciated in Fig. 6, which shows the results obtained for the two redshifts ($z = 6$ and 9) and the three different depths (from bottom to top: $M_{\text{UV}} = -17, -16, \text{ and } -15$ mag). As expected, most lensed galaxies have relatively small magnification factors simply because smaller magnification factors correspond to larger areas in the source plane. For the deepest observations ($M_{\text{UV}} = -15$ mag), we expect ~1 strongly lensed galaxy with total magnification $\mu \approx 30$. This is a galaxy near a caustic or even with parts of it crossing a caustic. Large magnification factors are more likely for high-redshift galaxies owing to their smaller intrinsic size. Recent results from JWST reveal very small sizes for faint galaxies at $z \approx 9$ with estimated sizes a few tens of parsecs (Williams et al. 2023). These galaxies appear unresolved even in JWST images. With tangential magnification factors of order 10, these galaxies can be resolved, so we expect that one such galaxy may be detectable in deep observations of MACS0416 with JWST. Surprisingly, from the same figure, this is the same number of galaxies we expect to observe with the same magnification factor, but at $z = 6$. This is due to the fact that the adopted

<footnote>$^4$ More counterimages may appear owing to substructure or when the galaxy is also producing radial images, but this situation is not considered here.</footnote>
<figure><img src="image_8.png" /></figure>

**Fig. 6.** Expected number of lensed galaxies in MACS0416. The larger plot shows the number of lensed galaxies behind MACS0416 at $z = 6$ (blue curves) and $z = 9$ (red curves), as a function of magnification (larger plot) or observed absolute magnitude (smaller plot). We assumed a standard luminosity function at each redshift. The different lines are for different limiting magnitudes: -17 (solid lines), -16 (dotted lines), or -15 (dashed lines). In the smaller plot the solid lines are the number of lensed galaxies for each observed magnitude, while the dashed lines are the number of galaxies one would have observed in the same magnitude bins but without lensing magnification.

luminosity function is steeper at $z = 9$ than at $z = 6$. This is better appreciated in the inset of Fig. 6, where we plot the number of detected galaxies (now corrected for multiplicity) as a function of the survey depth. Our model predicts that ~100 galaxies should be detectable in deep JWST observations both at $z = 6$ and $z = 9$.

In the inset of Fig. 6, we show as a dashed line the expected number of galaxies in a blank field (that is, with magnification $\mu = 1$) of similar area to the region where counterimages are expected to form in MACS0416. For galaxies at $z = 6$ we predict the classic depletion effect, where one is less likely to find faint galaxies in lensing fields than in blank fields. In this case, lensing increases the number of detections only at the bright end. However, for $z = 9$ galaxies, one is expected to be more successful (by a factor ~3) in the cluster field than in the parallel field. A similar conclusion was reached by Mahler et al. (2019) for galaxies at $z = 10$. These predictions can be tested with ongoing observations of MACS0416 that already cover the cluster region and several parallel fields. If the number density of galaxies found in the parallel fields is greater than in the cluster field, this would be a clear indication that the slope of the luminosity function must be shallower than 2 at the faint end (that is, for galaxies fainter than $M_{UV} = -15$ mag).

## 6. The Spock galaxy

In this section, we focus on one of the lensed galaxies behind MACS0416. Nicknamed Spock by Rodney et al. (2018), it is at relatively low redshift ($z = 1.0054$), but crosses one of the fold caustics of the cluster and hence is being significantly stretched and amplified. Relevant for the main results of this work are its size and geometry, relative position in relation to the main caustic, and stellar population derived from spectral energy distribution (SED) fitting. (These parameters and variables are all discussed in this section.)

The elongated arc from the Spock galaxy is shown in the center of Fig. 7. The lensed image forms a long arc perpendicu-

<figure><img src="image_9.png" /></figure>

**Fig. 7.** Spock arc. The radial arc at $z = 1.0054$ shows the location of the two Spock transients (S1 and S2 in yellow circles) reported by Rodney et al. (2018). Two additional transients are also marked with light-blue crosses (F1 and F2). These were discovered in the same arc as part of the Flashlights program (Kelly et al. 2022). Redshifts of nearby galaxies are indicated in light gray. The white line is the critical curve at $z = 1.0054$ assuming the cluster is at $z = 0.396$. The red curve assumes the cluster is at a slightly larger redshift of 0.4 and consistent with these nearby galaxies. The magenta line is the $z_s = 1.0054$ critical curve predicted by the model of Bergamini et al. (2023).

lar to the main axis and the critical curve of the cluster. A third image (with a much smaller magnification $\mu \approx 3.5$ and shown in the upper-right corner of Fig. 3) is found in the northern part of the cluster (and outside the field of view shown in Fig. 7). The elongated arc is in fact a merging image showing only a small portion of the Spock galaxy. Several lens models are presented by Rodney et al. (2018) showing how this merged image is the superposition of at least two counterimages, but possibly more than two.

The existence of the transients S1 and S2 (labeled in the image and originally reported by Rodney et al. 2018) supports the existence of multiple caustic crossings. This is further supported by the discovery of two new transients by Kelly et al. (2022) as part of the Flashlights program (and marked F1 and F2 in the figure) near the position of the two original transients. Kelly et al. (2022) suggest that S1–F1 form an image pair of the same star, while S2–F2 form a different pair of another star. This would imply that one critical curve must pass through S1–F1 and a second one through S2–F2, making the arc a superposition of at least three counterimages. The critical curve at $z = 1.0054$ from our new lens model is shown in white in Fig. 7. As expected, this critical curve intersects the arc, but it does so only once, and very close to the point between S1 and F1. The points S2 and F2 do not have a critical curve passing through them, but they are <1" away from one. Small changes in the mass of the nearby member galaxies, an invisible millilens (for instance, a dwarf galaxy in the cluster), or the large-scale mass distribution, can modify this portion of the critical curve by a small amount, making a second crossing through S2–F2 possible. We show in red an alternative lens model where the redshift of the cluster is increased slightly from the originally adopted $z = 0.396$ to $z = 0.4$. The increase in redshift in this portion of the lens plane is motivated by the existence of member galaxies nearby at redshifts close to (or even larger than) $z = 0.4$, including the BCG to the north which is also slightly above $z = 0.4$. Interestingly, this alternative model places the critical curve closer to S2–F2, making the hypothesis of a second caustic crossing more plausible. Yet another alternative model is shown as a magenta curve, corresponding to the
<figure><img src="image_7.png" /></figure>

**Fig. 6.** Expected number of lensed galaxies in MACS0416. The larger plot shows the number of lensed galaxies behind MACS0416 at $z = 6$ (blue curves) and $z = 9$ (red curves), as a function of magnification (larger plot) or observed absolute magnitude (smaller plot). We assumed a standard luminosity function at each redshift. The different lines are for different limiting magnitudes: -17 (solid lines), -16 (dotted lines), or -15 (dashed lines). In the smaller plot the solid lines are the number of lensed galaxies for each observed magnitude, while the dashed lines are the number of galaxies one would have observed in the same magnitude bins but without lensing magnification.

luminosity function is steeper at $z = 9$ than at $z = 6$. This is better appreciated in the inset of Fig. 6, where we plot the number of detected galaxies (now corrected for multiplicity) as a function of the survey depth. Our model predicts that ~100 galaxies should be detectable in deep JWST observations both at $z = 6$ and $z = 9$.

In the inset of Fig. 6, we show as a dashed line the expected number of galaxies in a blank field (that is, with magnification $\mu = 1$) of similar area to the region where counterimages are expected to form in MACS0416. For galaxies at $z = 6$ we predict the classic depletion effect, where one is less likely to find faint galaxies in lensing fields than in blank fields. In this case, lensing increases the number of detections only at the bright end. However, for $z = 9$ galaxies, one is expected to be more successful (by a factor ~3) in the cluster field than in the parallel field. A similar conclusion was reached by Mahler et al. (2019) for galaxies at $z = 10$. These predictions can be tested with ongoing observations of MACS0416 that already cover the cluster region and several parallel fields. If the number density of galaxies found in the parallel fields is greater than in the cluster field, this would be a clear indication that the slope of the luminosity function must be shallower than 2 at the faint end (that is, for galaxies fainter than $M_{UV} = -15$ mag).

## 6. The Spock galaxy

In this section, we focus on one of the lensed galaxies behind MACS0416. Nicknamed Spock by Rodney et al. (2018), it is at relatively low redshift ($z = 1.0054$), but crosses one of the fold caustics of the cluster and hence is being significantly stretched and amplified. Relevant for the main results of this work are its size and geometry, relative position in relation to the main caustic, and stellar population derived from spectral energy distribution (SED) fitting. (These parameters and variables are all discussed in this section.)

The elongated arc from the Spock galaxy is shown in the center of Fig. 7. The lensed image forms a long arc perpendicu-

<figure><img src="image_8.png" /></figure>

**Fig. 7.** Spock arc. The radial arc at $z = 1.0054$ shows the location of the two Spock transients (S1 and S2 in yellow circles) reported by Rodney et al. (2018). Two additional transients are also marked with light-blue crosses (F1 and F2). These were discovered in the same arc as part of the Flashlights program (Kelly et al. 2022). Redshifts of nearby galaxies are indicated in light gray. The white line is the critical curve at $z = 1.0054$ assuming the cluster is at $z = 0.396$. The red curve assumes the cluster is at a slightly larger redshift of 0.4 and consistent with these nearby galaxies. The magenta line is the $z_s = 1.0054$ critical curve predicted by the model of Bergamini et al. (2023).

lar to the main axis and the critical curve of the cluster. A third image (with a much smaller magnification $\mu \approx 3.5$ and shown in the upper-right corner of Fig. 3) is found in the northern part of the cluster (and outside the field of view shown in Fig. 7). The elongated arc is in fact a merging image showing only a small portion of the Spock galaxy. Several lens models are presented by Rodney et al. (2018) showing how this merged image is the superposition of at least two counterimages, but possibly more than two.

The existence of the transients S1 and S2 (labeled in the image and originally reported by Rodney et al. 2018) supports the existence of multiple caustic crossings. This is further supported by the discovery of two new transients by Kelly et al. (2022) as part of the Flashlights program (and marked F1 and F2 in the figure) near the position of the two original transients. Kelly et al. (2022) suggest that S1–F1 form an image pair of the same star, while S2–F2 form a different pair of another star. This would imply that one critical curve must pass through S1–F1 and a second one through S2–F2, making the arc a superposition of at least three counterimages. The critical curve at $z = 1.0054$ from our new lens model is shown in white in Fig. 7. As expected, this critical curve intersects the arc, but it does so only once, and very close to the point between S1 and F1. The points S2 and F2 do not have a critical curve passing through them, but they are <1" away from one. Small changes in the mass of the nearby member galaxies, an invisible millilens (for instance, a dwarf galaxy in the cluster), or the large-scale mass distribution, can modify this portion of the critical curve by a small amount, making a second crossing through S2–F2 possible. We show in red an alternative lens model where the redshift of the cluster is increased slightly from the originally adopted $z = 0.396$ to $z = 0.4$. The increase in redshift in this portion of the lens plane is motivated by the existence of member galaxies nearby at redshifts close to (or even larger than) $z = 0.4$, including the BCG to the north which is also slightly above $z = 0.4$. Interestingly, this alternative model places the critical curve closer to S2–F2, making the hypothesis of a second caustic crossing more plausible. Yet another alternative model is shown as a magenta curve, corresponding to the
<figure><img src="image_10.png" /><figcaption>**Fig. 8.** Portion of the caustic region with the magnification shown as gray scale. Left: small region of 4'' × 4'' showing the caustics for a source at $z = 1.0054$. The yellow circle marks the predicted position of the Spock galaxy in relation to the caustics. The green square marks the region that is enlarged in the right panel. Right: zoom-in of a small portion in the caustics obtained after interpolating the deflection field to a resolution of 2 milliarcseconds per pixel and applying the inverse ray shooting method. The yellow circle represents the approximate circular shape of the Spock galaxy and its position in relation to the caustic. The length of the caustic that intersects the galaxy is estimated to be ~530 pc. At a distance $d$ pc from the caustic, the magnification in the interior region of the caustic drops as $\mu \approx 270/\sqrt{d}$ (pc).</figcaption></figure>

model from Bergamini et al. (2023) at $z_s = 1.0054^5$. The three models show in general good agreement, but also make it obvious that a relatively large degree of uncertainty exists for the precise location of the critical curve around the position of this arc, and at this redshift.

The magnification seen by the Spock galaxy is better shown in the source plane. From the lens model we derive a map of the magnification in the source plane around the predicted position of the Spock galaxy. At the redshift of this galaxy, $1'' = 8.13$ kpc. In order to increase the spatial resolution, we interpolate the smooth deflection field before applying a standard inverse ray shooting algorithm to compute the caustics. The interpolated region is sufficiently large to guarantee that the inverse ray shooting method is not missing pixels in the image plane that land near the predicted position of the Spock galaxy. The resulting caustics are shown in Fig. 8. From this result we find that the Spock galaxy crosses a caustic, as expected from the fact that two of its counterimages form a merging pair of images. When zooming in around the Spock galaxy, we find that several caustics converge in this region. This explains the unusually large elongation of the main arc, despite the small intrinsic size of the galaxy. In the caustic map it is also evident that small changes in the lens configuration can result in multiple caustic crossings for a source placed in the right position. The magnification at the position of Spock scales as $\mu \approx 3/\sqrt{d}$, where $d$ is the distance to the caustic expressed in arcseconds, or $\mu \approx 270/\sqrt{d}$ when $d$ is expressed in parsecs. Using this scaling, we can estimate (later in this work) the magnification at different distances from the caustic.

## 6.1. Luminosity and intrinsic size of Spock galaxy

The third image of the Spock galaxy has an apparent AB magnitude of 26.6 in the F110W filter. When demagnified by a factor 3.5 this corresponds to an AB magnitude of 27.96. Since the galaxy is at $z = 1.0054$, the F110W filter corresponds to the visible part of the spectrum.

To obtain a constraint on the intrinsic size of the galaxy, we use this third image, which has a much lower magnification of $\mu = 3.5$. Unlike the elongated arc which shows only a portion of the source galaxy, the third image includes the full extent of the galaxy so it gives us the full picture of the demagnified source. At the position of the third image we obtain for the tangential and radial components of the magnification $\mu_t = 2.2$ and $\mu_r = 1.6$, respectively. The third image can be described, to first order, as an ellipse with semimajor axis $a \approx 0.2''$ and $b \approx 0.14''$. When demagnified into the source plane this corresponds to an almost circular shape with radius $r \approx 360$ pc.

Now we derive a separate constraint on the portion of the galaxy that is multiply imaged, but using the elongated arc instead, which contains only a portion of the Spock galaxy. By matching the magnification predicted by the caustics with the magnification of the arc observed in the image plane at the two extremes of the longest arc ($\mu \approx 14$ at each end and half the total magnification when one considers a double pair of counterimages), and using the scaling law for the total magnification, $\mu = 270/\sqrt{d}$ (pc), we find that the edge of the Spock galaxy must be ~100 pc away from the caustic. The right panel of Fig. 8 shows our best possible interpretation of the geometry of the source in the caustic region. As noted above, a significant portion of the galaxy must remain outside the high-magnification region in order to not be multiply lensed. Only the small portion inside the caustic region gets imaged into the elongated arc. For

<footnote n="5">⁵ This model is available at http://bazinga.fe.infn.it:6005/SLOT</footnote>
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Portion of the caustic region with the magnification shown as gray scale. Left: small region of 4'' × 4'' showing the caustics for a source at $z = 1.0054$. The yellow circle marks the predicted position of the Spock galaxy in relation to the caustics. The green square marks the region that is enlarged in the right panel. Right: zoom-in of a small portion in the caustics obtained after interpolating the deflection field to a resolution of 2 milliarcseconds per pixel and applying the inverse ray shooting method. The yellow circle represents the approximate circular shape of the Spock galaxy and its position in relation to the caustic. The length of the caustic that intersects the galaxy is estimated to be ~530 pc. At a distance $d$ pc from the caustic, the magnification in the interior region of the caustic drops as $\mu \approx 270/\sqrt{d}$ (pc).</figcaption></figure>

model from Bergamini et al. (2023) at $z_s = 1.0054^5$. The three models show in general good agreement, but also make it obvious that a relatively large degree of uncertainty exists for the precise location of the critical curve around the position of this arc, and at this redshift.

The magnification seen by the Spock galaxy is better shown in the source plane. From the lens model we derive a map of the magnification in the source plane around the predicted position of the Spock galaxy. At the redshift of this galaxy, $1'' = 8.13$ kpc. In order to increase the spatial resolution, we interpolate the smooth deflection field before applying a standard inverse ray shooting algorithm to compute the caustics. The interpolated region is sufficiently large to guarantee that the inverse ray shooting method is not missing pixels in the image plane that land near the predicted position of the Spock galaxy. The resulting caustics are shown in Fig. 8. From this result we find that the Spock galaxy crosses a caustic, as expected from the fact that two of its counterimages form a merging pair of images. When zooming in around the Spock galaxy, we find that several caustics converge in this region. This explains the unusually large elongation of the main arc, despite the small intrinsic size of the galaxy. In the caustic map it is also evident that small changes in the lens configuration can result in multiple caustic crossings for a source placed in the right position. The magnification at the position of Spock scales as $\mu \approx 3/\sqrt{d}$, where $d$ is the distance to the caustic expressed in arcseconds, or $\mu \approx 270/\sqrt{d}$ when $d$ is expressed in parsecs. Using this scaling, we can estimate (later in this work) the magnification at different distances from the caustic.

## 6.1. Luminosity and intrinsic size of Spock galaxy

The third image of the Spock galaxy has an apparent AB magnitude of 26.6 in the F110W filter. When demagnified by a factor 3.5 this corresponds to an AB magnitude of 27.96. Since the galaxy is at $z = 1.0054$, the F110W filter corresponds to the visible part of the spectrum.

To obtain a constraint on the intrinsic size of the galaxy, we use this third image, which has a much lower magnification of $\mu = 3.5$. Unlike the elongated arc which shows only a portion of the source galaxy, the third image includes the full extent of the galaxy so it gives us the full picture of the demagnified source. At the position of the third image we obtain for the tangential and radial components of the magnification $\mu_t = 2.2$ and $\mu_r = 1.6$, respectively. The third image can be described, to first order, as an ellipse with semimajor axis $a \approx 0.2''$ and $b \approx 0.14''$. When demagnified into the source plane this corresponds to an almost circular shape with radius $r \approx 360$ pc.

Now we derive a separate constraint on the portion of the galaxy that is multiply imaged, but using the elongated arc instead, which contains only a portion of the Spock galaxy. By matching the magnification predicted by the caustics with the magnification of the arc observed in the image plane at the two extremes of the longest arc ($\mu \approx 14$ at each end and half the total magnification when one considers a double pair of counterimages), and using the scaling law for the total magnification, $\mu = 270/\sqrt{d}$ (pc), we find that the edge of the Spock galaxy must be ~100 pc away from the caustic. The right panel of Fig. 8 shows our best possible interpretation of the geometry of the source in the caustic region. As noted above, a significant portion of the galaxy must remain outside the high-magnification region in order to not be multiply lensed. Only the small portion inside the caustic region gets imaged into the elongated arc. For

<footnote n="5">⁵ This model is available at http://bazinga.fe.infn.it:6005/SLOT</footnote>
our discussion, the most interesting constraint is the length of the caustic that intersects the Spock galaxy, which we estimate as ~530 pc (see right panel of Fig. 8).

## 6.2. Luminosity of lensed stars

After deriving the size of the source and how the magnification scales with distance to the caustic, we can now estimate the probability of observing microlensing events, similar to the two events of Rodney et al. (2018) and the two events of Kelly et al. (2022). First, we focus on S1 and F1, which are expected to have the largest magnifications from our lens model. In what follows, we assume that S1 and F1 are two different stars, but it is also possible that they are the same star with the critical curve passing through the midpoint between them. This distinction will not have a significant impact in the discussion below. Since the stars allegedly responsible for these events are only observed during short microlensing events, lasting days to weeks, they must fall below the detection limit when only the macromodel magnification is magnifying these stars. For the two lens models shown in Fig. 7, the predicted magnification (without microlensing) is at least 350 for S1 and 180 for F1.

For the time being, let us adopt a fiducial value $\mu = 270$ for the macromodel magnification of each star, and assume F1 is a counterimage of S1; other values of $\mu$ are discussed later. Based on the scaling for the magnification, this corresponds to a distance of 0.25 pc from the caustic (the total magnification is 540, which corresponds to a separation of 0.25 pc from the caustic). If the magnification at each of the two positions is 270, the star must be fainter than AB magnitude 34.6, such that in the absence of a microlensing boost in the magnification, the star remains undetected below the threshold of 28.5 mag for a point source in the HFF. This constraint is easily satisfied by the vast majority of stars at $z = 1.0054$. During a microlensing event, these stars can get a momentary boost of ~1–3 mag (e.g., Kelly et al. 2018; Diego et al. 2018; Venumadhav et al. 2017; Meena et al. 2022), making stars as faint as AB magnitude 37 at $z = 1$ detectable. The distance modulus at this redshift is 42.63 mag, hence the star needs to be fainter than absolute magnitude -8 in order to remain undetected without the aid of microlensing (ignoring specific K-corrections, that we address in more detail below). Absolute magnitude $M_B = -8$ is the realm of bright SG and especially hypergiant stars that can approach absolute magnitudes $M_V \approx -10$ (Humphreys 1978). Then, the star responsible for S1 and F1 is unlikely a hypergiant, since otherwise it would be detectable at all times in the Spock galaxy, even without the occasional microlensing boost. More specifically, the star needs to have bolometric luminosity $L < 3 \times 10^7 L_\odot$ to remain undetected in the F814W filter, or $L < 1.2 \times 10^7 L_\odot$ to be undetected in the F200LP Flashlights observations. On the other hand, since microlensing can only boost the observed flux by 1–3 mag, the stars cannot be fainter than absolute magnitude $M_V \approx -5$; otherwise, they would be undetectable even with microlensing. Stars with this luminosity are in the red and blue SG class, so this is possibly our best candidate for S1 (and F1). SG stars can have a wide range of temperatures, but can reach luminosities between ~$10^4 L_\odot$ and ~$10^5 L_\odot$.

Similar arguments can be followed for S2 and F2, except for this case the magnification from our lens model is almost an order of magnitude smaller ($\mu \approx 40$). Hence, without a second caustic crossing through S2–F2, and assuming S2 is a counterimage of F2, the star responsible for these two events must be 2 mag brighter, having absolute magnitude $M_V \approx -7$. This type of star is very rare, making this possibility less likely than a sec-

ond caustic crossing between S2 and F2 that would increase the magnification and consequently lower the required luminosity of the background star.

For the remainder of this work, we assume that the stars detected in the Spock arc are SGs, with temperatures ranging from $T \approx 3500$ K to $T \approx 30000$ K. This is in good agreement with the photometry observed for S1 and S2, and shown in Fig. 1. For the cooler red SG stars, and according to the Stefan Boltzmann law, the lower temperature is compensated by their much larger radius in such a way that $L \propto R^2 T^4$ can remain very high and relatively uniform. The hotter stars are smaller in radius, but emit significantly more energy at the shortest wavelengths. Related to the wide range of possible temperatures, we need to also consider the detection band, since stars at $z \approx 1$ with $T \approx 3500$ K will have their peak emission in the IR bands (assuming a blackbody spectrum and Wien's law), while a hot (mostly B-type) SG at the same redshift with $T \approx 30000$ K will peak in the UV part of the spectrum. We then consider two HST filters, F814W in the near-IR (used for the detection of the original Spock events S1 and S2) and the bluer F200LP (used in the Flashlights program to detect F1 and F2). F814W is well suited for stars showing a strong Balmer break ($10000 K < T < 20000$ K), while F200LP is better suited for hotter stars ($T > 20000$ K).

To compute the apparent magnitudes of the lensed stars in these filters, we need to adopt a model for their spectrum, $f_\lambda$. We assume the spectrum of these stars is well described (to first order) by a blackbody model with temperature $T$. We are hence ignoring spectral features such as absorption or emission lines, as well as discontinuities such as the Balmer break, but this simple approximation will still produce reasonable predictions. Setting the detection limit in the band to a given magnitude, $m_{\min}$, we can, for a given temperature, compute the minimum luminosity for a star at redshift $z = 1.0054$,

$$L_{\min} = \frac{L_f}{10^{0.4\Delta_{AB}}}, \qquad (2)$$

where $\Delta_{AB} = m_{\min} - m_f$ and $m_f$ is the apparent magnitude of a star at redshift $z = 1.0054$ in the chosen filter and with a fiducial bolometric luminosity $L_f$ (that without loss of generality we fix to $L_f = 10^7 L_\odot$):

$$m_f = -2.5 \log_{10}(F_f) - 48.6, \qquad (3)$$

with $F_f$ the usual flux received in the chosen filter and corrected by the filter response (and redshift dependence),

$$F_f(\nu) = (1 + z) \frac{\int f_f(\frac{\lambda}{1+z})S(\lambda)\lambda d\lambda}{\int S(\lambda)(c/\lambda)d\lambda}, \qquad (4)$$

where $f_f(\lambda)$ is the flux density and $S(\lambda)$ is the filter and detector response (Bessell & Murphy 2012). The multiplicative factor $(1+z)$ accounts for the redshift correction to the flux density. This flux is given simply by

$$f_f(\lambda) = \frac{L_f}{4\pi D_l(z)^2} \frac{BB(\lambda, T)}{\int BB(\lambda, T)d\lambda}, \qquad (5)$$

where $D_l(z)$ is the luminosity distance and $BB(\lambda, T)$ is the black-body spectrum for temperature $T$. Hence, the value of $L_{\min}$ represents the bolometric luminosity of the star (after accounting for magnification) needed to match the chosen magnitude limit, so stars with bolometric luminosity $L$ and amplified by a factor $\mu$ are not detected if $(\mu \times L) < L_{\min}$.
our discussion, the most interesting constraint is the length of the caustic that intersects the Spock galaxy, which we estimate as ~530 pc (see right panel of Fig. 8).

## 6.2. Luminosity of lensed stars

After deriving the size of the source and how the magnification scales with distance to the caustic, we can now estimate the probability of observing microlensing events, similar to the two events of Rodney et al. (2018) and the two events of Kelly et al. (2022). First, we focus on S1 and F1, which are expected to have the largest magnifications from our lens model. In what follows, we assume that S1 and F1 are two different stars, but it is also possible that they are the same star with the critical curve passing through the midpoint between them. This distinction will not have a significant impact in the discussion below. Since the stars allegedly responsible for these events are only observed during short microlensing events, lasting days to weeks, they must fall below the detection limit when only the macromodel magnification is magnifying these stars. For the two lens models shown in Fig. 7, the predicted magnification (without microlensing) is at least 350 for S1 and 180 for F1.

For the time being, let us adopt a fiducial value $\mu = 270$ for the macromodel magnification of each star, and assume F1 is a counterimage of S1; other values of $\mu$ are discussed later. Based on the scaling for the magnification, this corresponds to a distance of 0.25 pc from the caustic (the total magnification is 540, which corresponds to a separation of 0.25 pc from the caustic). If the magnification at each of the two positions is 270, the star must be fainter than AB magnitude 34.6, such that in the absence of a microlensing boost in the magnification, the star remains undetected below the threshold of 28.5 mag for a point source in the HFF. This constraint is easily satisfied by the vast majority of stars at $z = 1.0054$. During a microlensing event, these stars can get a momentary boost of ~1–3 mag (e.g., Kelly et al. 2018; Diego et al. 2018; Venumadhav et al. 2017; Meena et al. 2022), making stars as faint as AB magnitude 37 at $z = 1$ detectable. The distance modulus at this redshift is 42.63 mag, hence the star needs to be fainter than absolute magnitude -8 in order to remain undetected without the aid of microlensing (ignoring specific K-corrections, that we address in more detail below). Absolute magnitude $M_B = -8$ is the realm of bright SG and especially hypergiant stars that can approach absolute magnitudes $M_V \approx -10$ (Humphreys 1978). Then, the star responsible for S1 and F1 is unlikely a hypergiant, since otherwise it would be detectable at all times in the Spock galaxy, even without the occasional microlensing boost. More specifically, the star needs to have bolometric luminosity $L < 3 \times 10^7 L_\odot$ to remain undetected in the F814W filter, or $L < 1.2 \times 10^7 L_\odot$ to be undetected in the F200LP Flashlights observations. On the other hand, since microlensing can only boost the observed flux by 1–3 mag, the stars cannot be fainter than absolute magnitude $M_V \approx -5$; otherwise, they would be undetectable even with microlensing. Stars with this luminosity are in the red and blue SG class, so this is possibly our best candidate for S1 (and F1). SG stars can have a wide range of temperatures, but can reach luminosities between ~$10^4 L_\odot$ and ~$10^5 L_\odot$.

Similar arguments can be followed for S2 and F2, except for this case the magnification from our lens model is almost an order of magnitude smaller ($\mu \approx 40$). Hence, without a second caustic crossing through S2–F2, and assuming S2 is a counterimage of F2, the star responsible for these two events must be 2 mag brighter, having absolute magnitude $M_V \approx -7$. This type of star is very rare, making this possibility less likely than a sec-

ond caustic crossing between S2 and F2 that would increase the magnification and consequently lower the required luminosity of the background star.

For the remainder of this work, we assume that the stars detected in the Spock arc are SGs, with temperatures ranging from $T \approx 3500$ K to $T \approx 30000$ K. This is in good agreement with the photometry observed for S1 and S2, and shown in Fig. 1. For the cooler red SG stars, and according to the Stefan Boltzmann law, the lower temperature is compensated by their much larger radius in such a way that $L \propto R^2 T^4$ can remain very high and relatively uniform. The hotter stars are smaller in radius, but emit significantly more energy at the shortest wavelengths. Related to the wide range of possible temperatures, we need to also consider the detection band, since stars at $z \approx 1$ with $T \approx 3500$ K will have their peak emission in the IR bands (assuming a blackbody spectrum and Wien's law), while a hot (mostly B-type) SG at the same redshift with $T \approx 30000$ K will peak in the UV part of the spectrum. We then consider two HST filters, F814W in the near-IR (used for the detection of the original Spock events S1 and S2) and the bluer F200LP (used in the Flashlights program to detect F1 and F2). F814W is well suited for stars showing a strong Balmer break ($10000 K < T < 20000$ K), while F200LP is better suited for hotter stars ($T > 20000$ K).

To compute the apparent magnitudes of the lensed stars in these filters, we need to adopt a model for their spectrum, $f_\lambda$. We assume the spectrum of these stars is well described (to first order) by a blackbody model with temperature $T$. We are hence ignoring spectral features such as absorption or emission lines, as well as discontinuities such as the Balmer break, but this simple approximation will still produce reasonable predictions. Setting the detection limit in the band to a given magnitude, $m_{\min}$, we can, for a given temperature, compute the minimum luminosity for a star at redshift $z = 1.0054$,

$$L_{\min} = \frac{L_f}{10^{0.4\Delta_{AB}}}, \qquad (2)$$

where $\Delta_{AB} = m_{\min} - m_f$ and $m_f$ is the apparent magnitude of a star at redshift $z = 1.0054$ in the chosen filter and with a fiducial bolometric luminosity $L_f$ (that without loss of generality we fix to $L_f = 10^7 L_\odot$):

$$m_f = -2.5 \log_{10}(F_f) - 48.6, \qquad (3)$$

with $F_f$ the usual flux received in the chosen filter and corrected by the filter response (and redshift dependence),

$$F_f(\nu) = (1 + z) \frac{\int f_f(\frac{\lambda}{1+z})S(\lambda)\lambda d\lambda}{\int S(\lambda)(c/\lambda)d\lambda}, \qquad (4)$$

where $f_f(\lambda)$ is the flux density and $S(\lambda)$ is the filter and detector response (Bessell & Murphy 2012). The multiplicative factor $(1+z)$ accounts for the redshift correction to the flux density. This flux is given simply by

$$f_f(\lambda) = \frac{L_f}{4\pi D_l(z)^2} \frac{BB(\lambda, T)}{\int BB(\lambda, T)d\lambda}, \qquad (5)$$

where $D_l(z)$ is the luminosity distance and $BB(\lambda, T)$ is the black-body spectrum for temperature $T$. Hence, the value of $L_{\min}$ represents the bolometric luminosity of the star (after accounting for magnification) needed to match the chosen magnitude limit, so stars with bolometric luminosity $L$ and amplified by a factor $\mu$ are not detected if $(\mu \times L) < L_{\min}$.
Next, we need to find the maximum possible luminosity in the portion of the Spock galaxy that is being lensed. For this, we can use the fact that during the majority of the observations of the Spock arc, we observe no lensed star. We interpret this as the star being magnified by the most likely value during these observations, which is given by the mode of the curves shown in the left panel of Fig. 10. If a star is not observed when magnified by the mode, $\tilde{\mu}$, then its flux must be smaller than

$$L_{\max} = 1.3 \frac{L_{\min}}{\tilde{\mu}}. \qquad (6)$$

The factor 1.3 in the above expression is to account for the possibility that the star can also have magnifications 30% smaller than the mode. Setting this value to 1 has a small impact on our results, but we adopt the value of 1.3 as conservative. If $L_{\min}$ were interpreted as the minimum bolometric luminosity (after accounting for magnification) that a star has to have in order to be detected, $L_{\max}$ must be interpreted as the maximum luminosity (before accounting for magnification) that a star can have so it is not observed during the majority of the observations. The two constraints can be combined into the final condition for the bolometric luminosity, $L$, of the lensed stars:

$$\frac{L_{\min}}{\mu} < L < L_{\max} = 1.3 \frac{L_{\min}}{\tilde{\mu}}. \qquad (7)$$

Finally, we relate the bolometric luminosity with the mass of the star and the mass of the star with the IMF of the stellar population in the Spock galaxy. For the mass–luminosity relation ($M–L$) of main-sequence (MS) stars, we follow Duric (2004) with a steep $M–L$ for the least massive stars and an Eddington-limit type of $M–L$ relation for stars more massive than 55 $M_{\odot}$:

$$\frac{L}{L_{\odot}} = A \left( \frac{M}{M_{\odot}} \right)^{\beta}, \qquad (8)$$

where $A = (0.23, 1.0, 1.4$, and $3.2) \times 10^4$ for the mass intervals $M < 0.43 M_{\odot}$, $0.43 M_{\odot} < M < 2 M_{\odot}$, $2 M_{\odot} < M < 55 M_{\odot}$, and $55 M_{\odot} < M$ (respectively), and $\beta = 2.3, 4, 3.5$, and $1$ (respectively) for the same mass intervals. In this work, and since we require intrinsically high luminosities, only masses above ~20 $M_{\odot}$ are relevant for our calculations. This relation between mass and luminosity is valid for stars on the MS. Red giants and SGs or white dwarfs do not follow this relation. White dwarfs are not relevant for our study, but red giants and SGs are. In order to relate masses with the IMF, we make the simplification that the abundance of red SGs is comparable to the abundance of hot-ter stars with similar luminosity. That is, given a luminosity, we derive the mass from the above relation, and from that mass we derive the abundance of stars with that mass. Finally, we assume that stars with this mass can have a range of temperatures where the only parameter is the slope ($\delta$) of the distribution (related to the blue SG to red SG ratio, or simply $B/R$ ratio):

$$dN/dT \propto T^{\delta}. \qquad (9)$$

Ideally, we would instead use a luminosity function for the red SGs, but this is uncertain and dependent on variables such as mass loss or binarity (Neugent et al. 2020b; Massey et al. 2023). In the expression above, $\delta = 0$ corresponds to a uniform mix of red and blue SG stars as a function of temperature. A value $\delta = -1$ results in ~10 times more red SGs than blue SGs, and consistent with some local measurements, while $\delta = 1$ would be more representative of a very young population of luminous

<figure><img src="image_11.png" /><figcaption>Fig. 9. Stellar luminosity function inferred from the SED fitting featured in the lower-left corner. We adopt as a fiducial model the one corresponding to the black curve in the large plot.</figcaption></figure>

stars. As described below, our results depend mostly on a relatively narrow range of stars with luminosities close to (but below) the HD limit, and through the exponent $\delta$ in Eq. (9) we can control the relative abundances of blue and red SGs.

As a simple parameterization, we model the HD limit as a function of mass and temperature and then connect the mass to the luminosity through Eq. (8). In particular, we define the HD limit equivalent mass as $M_{HD} = 40 + T/2500$. With this relation and for $T = 3500$ K the corresponding luminosity would be $L_{HD} = 6.4 \times 10^5 L_{\odot}$, while for $T = 30000$ K we get $L_{HD} = 1.2 \times 10^6 L_{\odot}$. This dependence of the maximum luminosity on the temperature of the SG star is roughly consistent with the loci defined by the observations of nearby galaxies Humphreys (1983), Davies et al. (2018), McDonald et al. (2022).

Hence, the relevant variables for us are the HD limit, given by the approximation above, and the relative ratio of blue SG and red SG stars, $B/R$, which for simplicity we also parame-terize as a function of temperature with Eq. (9). The net num-ber of stars with luminosity $L$ is then given by the number of stars with the corresponding mass $M$, according to the $M–L$ relation above and as predicted by the IMF. This is the strongest assumption in our model, but one that is easily rescalable (our results can be rescaled by just varying the number of SG stars in the Spock arc and controlling the $B/R$ ratio with the parame-ter $\delta$). This simplification avoids the complexity of modeling an evolved population of intermediate-mass stars (such as red SGs) and reduces the number of free parameters in the model such as age, star rotation, and especially metallicity (Langer & Maeder 1995; Wagle et al. 2020), substituting these parameters with one single degree of freedom ($\delta$). For simplicity, we begin by pre-senting results with $\delta = 0$, which is close to the $B/R$ super-giant ratio found for some mass range and age (Schaller et al. 1992; Ekström et al. 2012). The range $-1 < \delta < 1$ reproduces well the $B/R$ ratio found in nearby young open clusters and galaxies, especially in their outskirts (Eggenberger et al. 2002; Dohm-Palmer & Skillman 2002), as in the case of the Spock galaxy where only the outer region of the galaxy is being mul-tiply lensed. We discuss these alternative scenarios with varying $\delta$ values in Sect. 9. Our next step is to discuss how many of these luminous stars we expect in the small section of the Spock galaxy that is only parsecs away from the caustic.
Next, we need to find the maximum possible luminosity in the portion of the Spock galaxy that is being lensed. For this, we can use the fact that during the majority of the observations of the Spock arc, we observe no lensed star. We interpret this as the star being magnified by the most likely value during these observations, which is given by the mode of the curves shown in the left panel of Fig. 10. If a star is not observed when magnified by the mode, $\tilde{\mu}$, then its flux must be smaller than

$$L_{\max} = 1.3 \frac{L_{\min}}{\tilde{\mu}}. \qquad (6)$$

The factor 1.3 in the above expression is to account for the possibility that the star can also have magnifications 30% smaller than the mode. Setting this value to 1 has a small impact on our results, but we adopt the value of 1.3 as conservative. If $L_{\min}$ were interpreted as the minimum bolometric luminosity (after accounting for magnification) that a star has to have in order to be detected, $L_{\max}$ must be interpreted as the maximum luminosity (before accounting for magnification) that a star can have so it is not observed during the majority of the observations. The two constraints can be combined into the final condition for the bolometric luminosity, $L$, of the lensed stars:

$$\frac{L_{\min}}{\mu} < L < L_{\max} = 1.3 \frac{L_{\min}}{\tilde{\mu}}. \qquad (7)$$

Finally, we relate the bolometric luminosity with the mass of the star and the mass of the star with the IMF of the stellar population in the Spock galaxy. For the mass–luminosity relation ($M–L$) of main-sequence (MS) stars, we follow Duric (2004) with a steep $M–L$ for the least massive stars and an Eddington-limit type of $M–L$ relation for stars more massive than 55 $M_{\odot}$:

$$\frac{L}{L_{\odot}} = A \left(\frac{M}{M_{\odot}}\right)^{\beta}, \qquad (8)$$

where $A = (0.23, 1.0, 1.4$, and $3.2) \times 10^4$ for the mass intervals $M < 0.43 M_{\odot}$, $0.43 M_{\odot} < M < 2 M_{\odot}$, $2 M_{\odot} < M < 55 M_{\odot}$, and $55 M_{\odot} < M$ (respectively), and $\beta = 2.3, 4, 3.5$, and $1$ (respectively) for the same mass intervals. In this work, and since we require intrinsically high luminosities, only masses above ~20 $M_{\odot}$ are relevant for our calculations. This relation between mass and luminosity is valid for stars on the MS. Red giants and SGs or white dwarfs do not follow this relation. White dwarfs are not relevant for our study, but red giants and SGs are. In order to relate masses with the IMF, we make the simplification that the abundance of red SGs is comparable to the abundance of hot-ter stars with similar luminosity. That is, given a luminosity, we derive the mass from the above relation, and from that mass we derive the abundance of stars with that mass. Finally, we assume that stars with this mass can have a range of temperatures where the only parameter is the slope ($\delta$) of the distribution (related to the blue SG to red SG ratio, or simply $B/R$ ratio):

$$dN/dT \propto T^{\delta}. \qquad (9)$$

Ideally, we would instead use a luminosity function for the red SGs, but this is uncertain and dependent on variables such as mass loss or binarity (Neugent et al. 2020b; Massey et al. 2023). In the expression above, $\delta = 0$ corresponds to a uniform mix of red and blue SG stars as a function of temperature. A value $\delta = -1$ results in ~10 times more red SGs than blue SGs, and consistent with some local measurements, while $\delta = 1$ would be more representative of a very young population of luminous

<figure><img src="image_10.png" /><figcaption>Fig. 9. Stellar luminosity function inferred from the SED fitting featured in the lower-left corner. We adopt as a fiducial model the one corresponding to the black curve in the large plot.</figcaption></figure>

stars. As described below, our results depend mostly on a relatively narrow range of stars with luminosities close to (but below) the HD limit, and through the exponent $\delta$ in Eq. (9) we can control the relative abundances of blue and red SGs.

As a simple parameterization, we model the HD limit as a function of mass and temperature and then connect the mass to the luminosity through Eq. (8). In particular, we define the HD limit equivalent mass as $M_{HD} = 40 + T/2500$. With this relation and for $T = 3500$ K the corresponding luminosity would be $L_{HD} = 6.4 \times 10^5 L_{\odot}$, while for $T = 30000$ K we get $L_{HD} = 1.2 \times 10^6 L_{\odot}$. This dependence of the maximum luminosity on the temperature of the SG star is roughly consistent with the loci defined by the observations of nearby galaxies Humphreys (1983), Davies et al. (2018), McDonald et al. (2022).

Hence, the relevant variables for us are the HD limit, given by the approximation above, and the relative ratio of blue SG and red SG stars, $B/R$, which for simplicity we also parame-terize as a function of temperature with Eq. (9). The net num-ber of stars with luminosity $L$ is then given by the number of stars with the corresponding mass $M$, according to the $M–L$ relation above and as predicted by the IMF. This is the strongest assumption in our model, but one that is easily rescalable (our results can be rescaled by just varying the number of SG stars in the Spock arc and controlling the $B/R$ ratio with the parame-ter $\delta$). This simplification avoids the complexity of modeling an evolved population of intermediate-mass stars (such as red SGs) and reduces the number of free parameters in the model such as age, star rotation, and especially metallicity (Langer & Maeder 1995; Wagle et al. 2020), substituting these parameters with one single degree of freedom ($\delta$). For simplicity, we begin by pre-senting results with $\delta = 0$, which is close to the $B/R$ super-giant ratio found for some mass range and age (Schaller et al. 1992; Ekström et al. 2012). The range $-1 < \delta < 1$ reproduces well the $B/R$ ratio found in nearby young open clusters and galaxies, especially in their outskirts (Eggenberger et al. 2002; Dohm-Palmer & Skillman 2002), as in the case of the Spock galaxy where only the outer region of the galaxy is being mul-tiply lensed. We discuss these alternative scenarios with varying $\delta$ values in Sect. 9. Our next step is to discuss how many of these luminous stars we expect in the small section of the Spock galaxy that is only parsecs away from the caustic.
### 6.3. Abundance of giant stars near the caustic in the Spock galaxy

We use the length of the Spock galaxy that intersects the caustic ($l = 530$ pc) and the scaling of magnification with distance found in a previous section, $\mu = 270/\sqrt{d}$ (pc), to compute the expected number of stars with a given magnification factor. For example, the area that is being magnified with factors greater than $\mu \approx 270$ (or twice 135) is then ~530 pc$^2$. We can find the number of luminous stars by fitting an IMF plus a scaling relation between mass and luminosity to the delensed AB magnitude of the Spock galaxy (27.96 in the V band after correction for magnification; see Sect. 6.1).

To study the stellar population in the Spock galaxy, we carry out annulus photometry on the third image to obtain the SED, as it is the only counterpart showing the entire Spock galaxy without any contamination. The circular aperture chosen is 3 times the PSF of HST, with an annulus of 1.5 times of aperture. The SED is then corrected by the lens-model magnification of $\mu = 3.5$.

We use Bagpipes (Carnall et al. 2018) to carry out SED fitting, adopting a single starburst with allowed nebular emission for the fit. Bagpipes returns a young stellar population without any line emission, where the age of the Spock galaxy is 67 ± 12 Myr, consistent with the blue color observed in the images. The best-fit current stellar mass is $9.1 \pm 1.2 \times 10^6 M_{\odot}$ with a supersolar metallicity of $Z = 1.55 Z_{\odot}$. These results, derived from photometry of the third counterimage, are consistent with those derived on the Spock arc by Rodney et al. (2018), who find a stellar mass of $13.8 \pm 4.5 \times 10^6 M_{\odot}$ and age 290 ± 500 Myr.

After having the best-fit SED, we sample stellar populations following the Kroupa IMF, similar to what Bagpipes used. We find that the Spock galaxy has $2.8 \times 10^4$ stars brighter than $5 \times 10^4 L_{\odot}$ (see Fig. 9). Among these, we expect ~36 stars in the 530 pc$^2$ area considered at the beginning of this subsection, and that is magnified by factors $\mu > 270$.

For the age derived above, we should not expect very massive (O-type) stars to still be present; however, given the fact that the Spock arc is magnifying a very small fraction of the entire galaxy, it is possible that a pocket of more recent star formation in that portion of the galaxy still harbors one of these massive and very hot stars. Hence, for the remainder of this work we assume that these massive stars can still exit in this portion of the galaxy.

## 7. Microlensing simulations for the Spock galaxy

Spock harbors several alleged bright stars$^6$, some of which are close enough to the cluster caustic and can be detected individually when undergoing microlensing events. Two of these microlensing events were reported by Rodney et al. (2018) and two additional events by Kelly et al. (2022). With JWST observations of this cluster becoming public soon, revealing more details about this interesting arc, here we pay special attention to this galaxy and the bright stars in its interior.

The macromodel magnification represents only the expected value of the magnification in a given region of the source plane. In reality, the ubiquitous presence of microlenses in the image plane creates a network of microcaustics in the source plane. A new magnification pattern emerges from reshuffling the macromodel magnification, and transforming its distribution from a delta function (the macromodel value) to a broad probability distribution.

Microcaustics from microlenses in the cluster can temporarily intersect a star in the background galaxy, providing an extra boost in magnification and making these stars detectable. In order to test the likelihood of detecting the star through a microlensing event, we have to calculate the probability of intersecting microcaustics.

The simulation with microlenses involves two main ingredients: (i) the macromodel, which is described by the tangential and radial magnifications (or similarly by $\kappa$ and $\gamma$), and (ii) the micromodel, which is basically given by the surface mass density of microlenses (we assume a standard Salpeter IMF for the distribution of masses with a lower cutoff in mass at $0.02 M_{\odot}$). For the macromodel we use the values of $\kappa = 0.691$ and $\gamma = 0.303$ obtained at the position S1 from the macromodel. The exact values of $\kappa$ and $\gamma$ are not relevant as long as they are representative of the values expected around this position. This results in magnification $\mu = \mu_t\mu_r = 166.7 \times 1.64 = 273.4$, which is in the range of magnifications expected for S1. To simulate different values of the magnification we simply change $\kappa$ very slightly until the desired magnification is reached.

For the microlens model, we assume that the contribution from the intracluster light (ICL) to the total mass at the distance of S1 from the main BCG (~50 kpc) is ~2%, or $\kappa_* = 0.01 \kappa = 0.01382$. This number is in agreement with recent results (Montes 2022; Diego et al. 2023). However, given the bimodal morphology of MACS0416, and the possible recent merger activity, the ICL could be even higher, so one should consider the 2% contribution as a conservative limit. This possibility is reinforced by the fact that the ICL is still clearly visible in HST images at the position of S1. Since the critical surface mass density for the redshifts of the cluster and the Spock galaxy is $\Sigma_{crit} = 2818 M_{\odot} \text{pc}^{-2}$, the surface mass density of microlenses corresponds to $\Sigma_* = 19.45 M_{\odot} \text{pc}^{-2}$. This estimate is in good agreement with the upper limit for the ICL derived by Rodney et al. (2018). They find $3.2 M_{\odot} \text{pc}^{-2} < \Sigma_* < 19.4 M_{\odot} \text{pc}^{-2}$ for the position of the Spock elongated arc. In this work we do not consider exotic forms of dark matter that could also act as microlenses, such as primordial black holes (PBHs). Earlier work based on microlensing of distant stars has constrained the abundance of PBHs to a small fraction (<10%) of the total projected mass (Oguri et al. 2018).

We make a set of simulations varying the macromodel magnification, $\mu = \mu_T\mu_R$, and the surface mass density of microlenses. The resulting probability of magnification for each simulation is presented in Fig. 10. The left panel shows the distribution of the magnification when microlenses are present. Each curve corresponds to a different scenario where the macromodel tangential magnification $\mu_T$ (which is the one that varies most along the arc) and/or the surface mass density of microlenses is varied. In all cases, the radial magnification is fixed to $\mu_R = 1.63$. To first order, the most likely magnification is close to the macromodel value, while the width of the probability distribution scales with the effective surface mass density of microlenses, $\Sigma_{eff} = \mu\Sigma_*$. This behavior was observed in earlier work (Diego 2019), and it is studied in more detail by Palencia et al. (2023). For the smaller values of the magnification, the probability scales as the expected $\mu^3$ at larger magnifications (or $\mu^2$ when expressed in log($\mu$) intervals). For the most extreme values of the macromodel magnification, the probability resembles a log-normal, and the median of the probability falls below the macromodel value as already described in earlier work (Diego et al. 2018; Diego 2019; Palencia et al. 2023). In the right panel of Fig. 10 we show the cumulative distribution for the magnification, and for the same scenarios shown in the left panel.

<footnote>$^6$ See Rodney et al. (2018) for possible alternative interpretations.</footnote>
### 6.3. Abundance of giant stars near the caustic in the Spock galaxy

We use the length of the Spock galaxy that intersects the caustic ($l = 530$ pc) and the scaling of magnification with distance found in a previous section, $\mu = 270/\sqrt{d}$ (pc), to compute the expected number of stars with a given magnification factor. For example, the area that is being magnified with factors greater than $\mu \approx 270$ (or twice 135) is then ~530 pc$^2$. We can find the number of luminous stars by fitting an IMF plus a scaling relation between mass and luminosity to the delensed AB magnitude of the Spock galaxy (27.96 in the V band after correction for magnification; see Sect. 6.1).

To study the stellar population in the Spock galaxy, we carry out annulus photometry on the third image to obtain the SED, as it is the only counterpart showing the entire Spock galaxy without any contamination. The circular aperture chosen is 3 times the PSF of HST, with an annulus of 1.5 times of aperture. The SED is then corrected by the lens-model magnification of $\mu = 3.5$.

We use Bagpipes (Carnall et al. 2018) to carry out SED fitting, adopting a single starburst with allowed nebular emission for the fit. Bagpipes returns a young stellar population without any line emission, where the age of the Spock galaxy is 67 ± 12 Myr, consistent with the blue color observed in the images. The best-fit current stellar mass is $9.1 \pm 1.2 \times 10^6$ M$_\odot$ with a supersolar metallicity of $Z = 1.55 Z_\odot$. These results, derived from photometry of the third counterimage, are consistent with those derived on the Spock arc by Rodney et al. (2018), who find a stellar mass of $13.8 \pm 4.5 \times 10^6$ M$_\odot$ and age 290 ± 500 Myr.

After having the best-fit SED, we sample stellar populations following the Kroupa IMF, similar to what Bagpipes used. We find that the Spock galaxy has $2.8 \times 10^4$ stars brighter than $5 \times 10^4 L_\odot$ (see Fig. 9). Among these, we expect ~36 stars in the 530 pc$^2$ area considered at the beginning of this subsection, and that is magnified by factors $\mu > 270$.

For the age derived above, we should not expect very massive (O-type) stars to still be present; however, given the fact that the Spock arc is magnifying a very small fraction of the entire galaxy, it is possible that a pocket of more recent star formation in that portion of the galaxy still harbors one of these massive and very hot stars. Hence, for the remainder of this work we assume that these massive stars can still exit in this portion of the galaxy.

## 7. Microlensing simulations for the Spock galaxy

Spock harbors several alleged bright stars$^6$, some of which are close enough to the cluster caustic and can be detected individually when undergoing microlensing events. Two of these microlensing events were reported by Rodney et al. (2018) and two additional events by Kelly et al. (2022). With JWST observations of this cluster becoming public soon, revealing more details about this interesting arc, here we pay special attention to this galaxy and the bright stars in its interior.

The macromodel magnification represents only the expected value of the magnification in a given region of the source plane. In reality, the ubiquitous presence of microlenses in the image plane creates a network of microcaustics in the source plane. A new magnification pattern emerges from reshuffling the macromodel magnification, and transforming its distribution from a delta function (the macromodel value) to a broad probability distribution.

Microcaustics from microlenses in the cluster can temporarily intersect a star in the background galaxy, providing an extra boost in magnification and making these stars detectable. In order to test the likelihood of detecting the star through a microlensing event, we have to calculate the probability of intersecting microcaustics.

The simulation with microlenses involves two main ingredients: (i) the macromodel, which is described by the tangential and radial magnifications (or similarly by $\kappa$ and $\gamma$), and (ii) the micromodel, which is basically given by the surface mass density of microlenses (we assume a standard Salpeter IMF for the distribution of masses with a lower cutoff in mass at $0.02 M_\odot$). For the macromodel we use the values of $\kappa = 0.691$ and $\gamma = 0.303$ obtained at the position S1 from the macromodel. The exact values of $\kappa$ and $\gamma$ are not relevant as long as they are representative of the values expected around this position. This results in magnification $\mu = \mu_t\mu_r = 166.7 \times 1.64 = 273.4$, which is in the range of magnifications expected for S1. To simulate different values of the magnification we simply change $\kappa$ very slightly until the desired magnification is reached.

For the microlens model, we assume that the contribution from the intracluster light (ICL) to the total mass at the distance of S1 from the main BCG (~50 kpc) is ~2%, or $\kappa_* = 0.01 \kappa = 0.01382$. This number is in agreement with recent results (Montes 2022; Diego et al. 2023). However, given the bimodal morphology of MACS0416, and the possible recent merger activity, the ICL could be even higher, so one should consider the 2% contribution as a conservative limit. This possibility is reinforced by the fact that the ICL is still clearly visible in HST images at the position of S1. Since the critical surface mass density for the redshifts of the cluster and the Spock galaxy is $\Sigma_{crit} = 2818 M_\odot$ pc$^{-2}$, the surface mass density of microlenses corresponds to $\Sigma_* = 19.45 M_\odot$ pc$^{-2}$. This estimate is in good agreement with the upper limit for the ICL derived by Rodney et al. (2018). They find $3.2 M_\odot$ pc$^{-2} < \Sigma_* < 19.4 M_\odot$ pc$^{-2}$ for the position of the Spock elongated arc. In this work we do not consider exotic forms of dark matter that could also act as microlenses, such as primordial black holes (PBHs). Earlier work based on microlensing of distant stars has constrained the abundance of PBHs to a small fraction (<10%) of the total projected mass (Oguri et al. 2018).

We make a set of simulations varying the macromodel magnification, $\mu = \mu_T\mu_R$, and the surface mass density of microlenses. The resulting probability of magnification for each simulation is presented in Fig. 10. The left panel shows the distribution of the magnification when microlenses are present. Each curve corresponds to a different scenario where the macromodel tangential magnification $\mu_T$ (which is the one that varies most along the arc) and/or the surface mass density of microlenses is varied. In all cases, the radial magnification is fixed to $\mu_R = 1.63$. To first order, the most likely magnification is close to the macromodel value, while the width of the probability distribution scales with the effective surface mass density of microlenses, $\Sigma_{eff} = \mu\Sigma_*$. This behavior was observed in earlier work (Diego 2019), and it is studied in more detail by Palencia et al. (2023). For the smaller values of the magnification, the probability scales as the expected $\mu^3$ at larger magnifications (or $\mu^2$ when expressed in log($\mu$) intervals). For the most extreme values of the macromodel magnification, the probability resembles a log-normal, and the median of the probability falls below the macromodel value as already described in earlier work (Diego et al. 2018; Diego 2019; Palencia et al. 2023). In the right panel of Fig. 10 we show the cumulative distribution for the magnification, and for the same scenarios shown in the left panel.

<footnote>$^6$ See Rodney et al. (2018) for possible alternative interpretations.</footnote>
<figure><img src="image_12.png" /><figcaption>**Fig. 10.** Microlensing probability of magnification for the Spock galaxy. Left: relative probability of magnification in different scenarios where the macromodel tangential magnification ($\mu_T$) or the surface mass density of microlenses ($\Sigma_*$ expressed in units of $M_{\odot} pc^{-2}$) is varied. In all cases, the radial component of the macromodel magnification is fixed to $\mu_R = 1.63$. Different values of the tangential magnification can be interpreted as different distances to the critical curve. At the lowest magnification values considered, the probability of large magnification follows the expected $\mu^{-2}$ law. For larger values of $\mu_T$, the probability resembles a log-normal. Right: fraction of area in the source plane with magnification factors above a certain value. The models are the same as in the left panel.</figcaption></figure>

If we consider the probability of magnification to be over a factor 1000, we see from this figure that it is above ~2% for any star that is within 1 pc from the cluster caustic (or macromodel magnification greater than $166 \times 1.63 = 270$). Surprisingly, the case where the surface mass density of microlenses is $\Sigma_* = 5 M_{\odot} pc^{-2}$ has the same probability as the case where $\Sigma_* = 20 M_{\odot} pc^{-2}$. The probability of magnification normally scales with $\Sigma_*$ (Palencia et al. 2023), so one would naively expect the case with higher $\Sigma_*$ to have larger magnification. This unexpected behavior is explained by the tail of the magnification of the case with $\Sigma_* = 5 M_{\odot} pc^{-2}$, which maintains the $\mu^{-2}$ scaling longer (i.e., toward higher magnifications) than the case with $\Sigma_* = 20 M_{\odot} pc^{-2}$, that for $\mu \approx 1000$ is already falling faster than $\mu^{-2}$, as shown in the left panel. As discussed by Palencia et al. (2023), the breaking of the scaling takes place when the effective surface mass density $\Sigma_{\text{eff}} = \mu\Sigma_* \approx \Sigma_{\text{crit}}$. Similar results were found in earlier work (Schechter & Wambsganss 2002; Schechter et al. 2004; Vernardos et al. 2014). For macromodel magnification $\mu = 666 \times 1.63 = 1086$, the probability of having magnification greater than $\mu = 1000$ is already 20%, but the star spends most of its time with magnification values below $\mu = 800$.

One final ingredient relevant for the discussion below is the typical timescale of a microlensing event of a distant star near a critical curve. Here we are interested in the timescale for crossing one microcaustic, much smaller than the time required to cross the entire microlens, which can span decades (Mosquera & Kochanek 2011). The duration of a high-magnification event produced by a microlens depends on the size of the source, the relative velocity between the source and the web of microcaustics (a redshift-dependent variable; Kayser et al. 1986; Kochanek 2004), the direction of motion, the mass of the microlenses, and the macromodel magnification. For the redshifts involved in this work, a typical velocity is $v = 500 \text{ km s}^{-1}$ (Neira et al. 2020)<fnref n="7" />. At this velocity, the relative position between a source at $z = 1$ and the web of microcaustics changes by 0.17 nanoarcseconds per day. This is approximately equivalent to $100 R_{\odot}$ at $z = 1$. During a star-caustic crossing, the magnification plateaus with a timescale that is proportional to the

radius of the lensed star ($R_*$) and a maximum magnification that is inversely proportional to $\sqrt{R_*}$ (Miralda-Escude 1991). That is, a star with radius ~100 $R_{\odot}$ can cross the microcaustic in ~1 day, during which time it will be at its maximum possible magnification. A typical microlens near a critical curve with macromodel magnification of order 100 produces high magnification factors ($\mu > 1000$) at distances between 1 and 5 nanoarcseconds from the microcaustic (depending on the mass of the microlens and where the crossing takes place near the cusps or folds; Diego 2022). This results in timescales between ~1 week and ~1 month (both modulated by the factor 500/$v$), during which luminous stars can be detected as their brightness increases by ~2–3 mag. This estimated timescale is consistent with alleged microlensing events such as Icarus (Kelly et al. 2018) or the two events in the Spock galaxy reported by Rodney et al. (2018). This timescale agrees also with previous estimates derived for QSO microlensing (Neira et al. 2020).

## 8. Probability of detecting a microlensing event in the Spock galaxy

Now we have all of the ingredients to do our final calculation and compute the expected number of stars that are detected above a given detection threshold and in a specific filter. For simplicity, we consider five regions of fixed magnification, each at a given distance from the caustic. The width of each region, $\Delta_{\text{pc}}$, is determined from the scaling law $\mu = 270/\sqrt{d}$ (pc), by imposing that at the boundaries of the region the magnification changes by 30% with respect to its central value. We fix $\Sigma_* = 20 M_{\odot} pc^{-2}$ and $\mu_R = 1.63$ in all regions, and vary only the tangential magnification which takes the values $\mu_T = 41, 83, 166, 333$, and 666 (all shown in Fig. 10, except $\mu_T = 83$, which for simplicity is omitted). Then we take the model derived from the SED fitting and compute the expected number of stars in the region of length 530 pc and width $\Delta_{\text{pc}}$ (this is an expected number, so it can be less than 1). Finally, all of the stars in this region are magnified by the corresponding constant magnification factor in the region (half the total magnification). We consider different stellar masses (or luminosities) and temperatures. For each combination of mass and temperature, we compute its bolometric luminosity using Eq. (8) and its apparent magnitude in the

<footnote n="7">⁷ This velocity involves the relative transverse motion of three planes: observer, lens, and source; see Appendix B of Kayser et al. (1986).</footnote>
<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Microlensing probability of magnification for the Spock galaxy. Left: relative probability of magnification in different scenarios where the macromodel tangential magnification ($\mu_T$) or the surface mass density of microlenses ($\Sigma_*$ expressed in units of $M_{\odot} pc^{-2}$) is varied. In all cases, the radial component of the macromodel magnification is fixed to $\mu_R = 1.63$. Different values of the tangential magnification can be interpreted as different distances to the critical curve. At the lowest magnification values considered, the probability of large magnification follows the expected $\mu^{-2}$ law. For larger values of $\mu_T$, the probability resembles a log-normal. Right: fraction of area in the source plane with magnification factors above a certain value. The models are the same as in the left panel.</figcaption></figure>

If we consider the probability of magnification to be over a factor 1000, we see from this figure that it is above ~2% for any star that is within 1 pc from the cluster caustic (or macromodel magnification greater than $166 \times 1.63 = 270$). Surprisingly, the case where the surface mass density of microlenses is $\Sigma_* = 5 M_{\odot} pc^{-2}$ has the same probability as the case where $\Sigma_* = 20 M_{\odot} pc^{-2}$. The probability of magnification normally scales with $\Sigma_*$ (Palencia et al. 2023), so one would naively expect the case with higher $\Sigma_*$ to have larger magnification. This unexpected behavior is explained by the tail of the magnification of the case with $\Sigma_* = 5 M_{\odot} pc^{-2}$, which maintains the $\mu^{-2}$ scaling longer (i.e., toward higher magnifications) than the case with $\Sigma_* = 20 M_{\odot} pc^{-2}$, that for $\mu \approx 1000$ is already falling faster than $\mu^{-2}$, as shown in the left panel. As discussed by Palencia et al. (2023), the breaking of the scaling takes place when the effective surface mass density $\Sigma_{\text{eff}} = \mu\Sigma_* \approx \Sigma_{\text{crit}}$. Similar results were found in earlier work (Schechter & Wambsganss 2002; Schechter et al. 2004; Vernardos et al. 2014). For macromodel magnification $\mu = 666 \times 1.63 = 1086$, the probability of having magnification greater than $\mu = 1000$ is already 20%, but the star spends most of its time with magnification values below $\mu = 800$.

One final ingredient relevant for the discussion below is the typical timescale of a microlensing event of a distant star near a critical curve. Here we are interested in the timescale for crossing one microcaustic, much smaller than the time required to cross the entire microlens, which can span decades (Mosquera & Kochanek 2011). The duration of a high-magnification event produced by a microlens depends on the size of the source, the relative velocity between the source and the web of microcaustics (a redshift-dependent variable; Kayser et al. 1986; Kochanek 2004), the direction of motion, the mass of the microlenses, and the macromodel magnification. For the redshifts involved in this work, a typical velocity is $v = 500 \text{ km s}^{-1}$ (Neira et al. 2020)<fnref n="7" />. At this velocity, the relative position between a source at $z=1$ and the web of microcaustics changes by 0.17 nanoarcseconds per day. This is approximately equivalent to $100 R_{\odot}$ at $z=1$. During a star-caustic crossing, the magnification plateaus with a timescale that is proportional to the

radius of the lensed star ($R_*$) and a maximum magnification that is inversely proportional to $\sqrt{R_*}$ (Miralda-Escude 1991). That is, a star with radius ~100 $R_{\odot}$ can cross the microcaustic in ~1 day, during which time it will be at its maximum possible magnification. A typical microlens near a critical curve with macromodel magnification of order 100 produces high magnification factors ($\mu > 1000$) at distances between 1 and 5 nanoarcseconds from the microcaustic (depending on the mass of the microlens and where the crossing takes place near the cusps or folds; Diego 2022). This results in timescales between ~1 week and ~1 month (both modulated by the factor 500/$v$), during which luminous stars can be detected as their brightness increases by ~2–3 mag. This estimated timescale is consistent with alleged microlensing events such as Icarus (Kelly et al. 2018) or the two events in the Spock galaxy reported by Rodney et al. (2018). This timescale agrees also with previous estimates derived for QSO microlensing (Neira et al. 2020).

## 8. Probability of detecting a microlensing event in the Spock galaxy

Now we have all of the ingredients to do our final calculation and compute the expected number of stars that are detected above a given detection threshold and in a specific filter. For simplicity, we consider five regions of fixed magnification, each at a given distance from the caustic. The width of each region, $\Delta_{\text{pc}}$, is determined from the scaling law $\mu = 270/\sqrt{d}$ (pc), by imposing that at the boundaries of the region the magnification changes by 30% with respect to its central value. We fix $\Sigma_* = 20 M_{\odot} pc^{-2}$ and $\mu_R = 1.63$ in all regions, and vary only the tangential magnification which takes the values $\mu_T = 41, 83, 166, 333$, and 666 (all shown in Fig. 10, except $\mu_T = 83$, which for simplicity is omitted). Then we take the model derived from the SED fitting and compute the expected number of stars in the region of length 530 pc and width $\Delta_{\text{pc}}$ (this is an expected number, so it can be less than 1). Finally, all of the stars in this region are magnified by the corresponding constant magnification factor in the region (half the total magnification). We consider different stellar masses (or luminosities) and temperatures. For each combination of mass and temperature, we compute its bolometric luminosity using Eq. (8) and its apparent magnitude in the

<footnote n="7">⁷ This velocity involves the relative transverse motion of three planes: observer, lens, and source; see Appendix B of Kayser et al. (1986).</footnote>
<figure><img src="image_13.png" /></figure>

**Fig. 11.** Number of stars vs. distance. The large plot shows the cumulative number of stars observed at different distances from the caustic, as a function of the star mass. In this and other plots, the mass has to be interpreted as the equivalent mass of a MS star with the same luminosity. SG stars normally have smaller masses than its equivalent MS star. We assume observations reach a limiting magnitude of 28.5 in HST's F814W filter. Stars from the MS with equivalent masses below ~20 $M_{\odot}$ are too faint to be detected, even at the largest magnification factors considered here. Each curve stops when $L_{\text{max}}$ is reached. The total number of stars expected at each distance is shown in the inset at the bottom right. Each bin in the inset is centered at a given magnification (or distance to the caustic) and the width is computed from the $270/\sqrt{d(\text{pc})}$ relation taking 30% above and 30% below the magnification of the bin. colors in the inset correspond to colors in the larger plot. In all cases except for the dashed curve, the assumed surface mass density of microlenses is $\Sigma_* = 20 M_{\odot} \text{pc}^{-2}$, and all stars are assumed to have the same temperature of 10 000 K. The dashed curve corresponds to separations of ~1 pc from the caustic, but for a surface mass density of microlenses of $\Sigma_* = 5 M_{\odot} \text{pc}^{-2}$. The maximum probability of seeing a lensed star is found at ~16 pc from the caustic (yellow bin), if one ignores the Humphreys-Davidson (HD) limit (vertical dotted line). If one considers this limit, the most likely distance to observe stars is at 1 pc (dark blue bin corresponding to total magnification $\mu \approx 270$, or $\mu = 135$ for each counterimage), in good agreement with the observations and lens model.

selected filter from Eq. (3). If the star is brighter than $m_{\text{min}}$, then it is counted as a detection.

The result is shown in Fig. 11 for the five different regions, and for stars with temperature $T = 10\,000$ K being observed with the HST F814W filter and above limiting magnitude $m_{\text{min}} = 28.5$ in this filter. In the larger plot, we show as solid lines the cumulative number of detected stars in each one of the five distance bins, and as a function of the mass of the star (or luminosity as shown in the top scale). The curves stop where the luminosity reaches $L_{\text{max}}$ (see Eq. (6)). At short distances from the caustic (black curve), the macromodel magnification is so large that faint stars can be more easily detected. The black curve also assumes that the double image of the lensed star forms an unresolved source, so its flux doubles. This situation is similar to the case of “Earendel” or “Godzilla”, where given the large inferred magnification factors, the double image is expected to appear as a single unresolved source (Welch et al. 2022; Diego et al. 2022). In this situation, the cumulative function grows faster for less luminous stars. On the opposite extreme, at large distances from the caustic (yellowish curve), where the macromodel magnification is smaller, only the brightest stars can be detected. These are rarer than less luminous stars, but since one is integrating over larger areas, the smaller probability of seeing these stars is partially compensated by the larger area. The prediction at

the smaller macromodel magnification values assumes stars can exist beyond the HD limit (marked by a vertical dotted line). If the HD limit still holds at this redshift, the number of stars detected in the yellow bin would drop by a factor of ~50. In a scenario where there are no stars above $L = 6 \times 10^5 L_{\odot}$ for temperatures $T = 10\,000$ K, we expect the maximum probability of detecting a star at distances of ~1 pc from the caustic (dark-blue curve), or at magnification factors $\mu \approx 270/2$ (the factor 2 accounting for the double image, each carrying half the magnification), in excellent agreement with S1 and F1. Accounting for the contribution from neighboring distances, and after imposing the HD limit, we find that the expected number of lensed stars at the position of S1 and F1 is of order 1, in very good agreement with the observations.

In the same plot, we also show as a dashed line the case where the surface mass density of microlenses is reduced by a factor of four – that is, $\Sigma_* = 5 M_{\odot}$. Despite having fewer microlenses, the probability of large magnification is similar to the scenario where we have four times more microlenses, as discussed earlier, resulting in a similar number of predicted lensed stars. Hence, our results are relatively insensitive to the particular choice of $\Sigma_*$, as long as it is within the constrained range.

The final number of stars in each region is the result of two competing effects. On one side, larger magnification factors make the detection of the more abundant fainter stars easier, but on the other side, regions with larger magnification factors are considerably smaller. There is a sweet spot where the combination of the two effects maximizes the probability of seeing a star. In the inset plot of Fig. 11, we show the total number of detected stars as a function of the distance to the caustic. The width of the bins in the inset corresponds to the width $\Delta_{\text{pc}}$ used in the calculation. Again, in this plot we have assumed that stars can exist beyond the HD limit. If one sets the HD limit as the maximum luminosity, the light blue and green bins would have ~5 and ~50 times fewer lensed stars, respectively. In that case, we deduce that the best chance of finding lensed stars is when the total magnification of the double image is $100 < \mu < 300$, corresponding to a physical distance of one to a few parsecs from the caustic (dark-blue ad light-blue curves, respectively). The fact that we do not observe numerous lensed stars at larger distances from the critical curve is a direct indication that the HD limit must also be applicable at this redshift.

The dependence with the mass (or luminosity) of the star is better appreciated in Fig. 12, where we show the expected number of detected stars above AB magnitude 28.5 in HST's filter F814W, and per mass interval. At short distances from the main caustic (large magnification factors), the fraction of stars that is detected is proportional to the adopted IMF. The black solid curve corresponds to the case with maximum magnification and assuming both counterimages of the lensed star form an unresolved pair (so the observed flux doubles, making the detection of fainter stars possible). At larger distances from the main caustic, as the magnification drops, lensed events naturally select the brightest stars. The most extreme case is given by the yellow curve (macromodel magnification $\mu \approx 67$), where the most likely stars to be detected are the most luminous ones. This dependency with the macromodel magnification offers a unique opportunity to probe the most luminous (and rare) stars in distant galaxies, since these would be the only ones visible at relatively large distances from the critical curves of powerful lenses. In the particular case of the Spock galaxy, since we do not see the predicted abundance of lensed stars at these large distances from the critical curve, we can conclude that these stars cannot exist, as discussed above.
<figure><img src="image_12.png" /></figure>

**Fig. 11.** Number of stars vs. distance. The large plot shows the cumulative number of stars observed at different distances from the caustic, as a function of the star mass. In this and other plots, the mass has to be interpreted as the equivalent mass of a MS star with the same luminosity. SG stars normally have smaller masses than its equivalent MS star. We assume observations reach a limiting magnitude of 28.5 in HST's F814W filter. Stars from the MS with equivalent masses below ~20 $M_\odot$ are too faint to be detected, even at the largest magnification factors considered here. Each curve stops when $L_{\text{max}}$ is reached. The total number of stars expected at each distance is shown in the inset at the bottom right. Each bin in the inset is centered at a given magnification (or distance to the caustic) and the width is computed from the $270/\sqrt{d(\text{pc})}$ relation taking 30% above and 30% below the magnification of the bin. colors in the inset correspond to colors in the larger plot. In all cases except for the dashed curve, the assumed surface mass density of microlenses is $\Sigma_* = 20 M_\odot \text{pc}^{-2}$, and all stars are assumed to have the same temperature of 10 000 K. The dashed curve corresponds to separations of ~1 pc from the caustic, but for a surface mass density of microlenses of $\Sigma_* = 5 M_\odot \text{pc}^{-2}$. The maximum probability of seeing a lensed star is found at ~16 pc from the caustic (yellow bin), if one ignores the Humphreys-Davidson (HD) limit (vertical dotted line). If one considers this limit, the most likely distance to observe stars is at 1 pc (dark blue bin corresponding to total magnification $\mu \approx 270$, or $\mu = 135$ for each counterimage), in good agreement with the observations and lens model.

selected filter from Eq. (3). If the star is brighter than $m_{\text{min}}$, then it is counted as a detection.

The result is shown in Fig. 11 for the five different regions, and for stars with temperature $T = 10\,000$ K being observed with the HST F814W filter and above limiting magnitude $m_{\text{min}} = 28.5$ in this filter. In the larger plot, we show as solid lines the cumulative number of detected stars in each one of the five distance bins, and as a function of the mass of the star (or luminosity as shown in the top scale). The curves stop where the luminosity reaches $L_{\text{max}}$ (see Eq. (6)). At short distances from the caustic (black curve), the macromodel magnification is so large that faint stars can be more easily detected. The black curve also assumes that the double image of the lensed star forms an unresolved source, so its flux doubles. This situation is similar to the case of “Earendel” or “Godzilla”, where given the large inferred magnification factors, the double image is expected to appear as a single unresolved source (Welch et al. 2022; Diego et al. 2022). In this situation, the cumulative function grows faster for less luminous stars. On the opposite extreme, at large distances from the caustic (yellowish curve), where the macromodel magnification is smaller, only the brightest stars can be detected. These are rarer than less luminous stars, but since one is integrating over larger areas, the smaller probability of seeing these stars is partially compensated by the larger area. The prediction at

the smaller macromodel magnification values assumes stars can exist beyond the HD limit (marked by a vertical dotted line). If the HD limit still holds at this redshift, the number of stars detected in the yellow bin would drop by a factor of ~50. In a scenario where there are no stars above $L = 6 \times 10^5 L_\odot$ for temperatures $T = 10\,000$ K, we expect the maximum probability of detecting a star at distances of ~1 pc from the caustic (dark-blue curve), or at magnification factors $\mu \approx 270/2$ (the factor 2 accounting for the double image, each carrying half the magnification), in excellent agreement with S1 and F1. Accounting for the contribution from neighboring distances, and after imposing the HD limit, we find that the expected number of lensed stars at the position of S1 and F1 is of order 1, in very good agreement with the observations.

In the same plot, we also show as a dashed line the case where the surface mass density of microlenses is reduced by a factor of four – that is, $\Sigma_* = 5 M_\odot$. Despite having fewer microlenses, the probability of large magnification is similar to the scenario where we have four times more microlenses, as discussed earlier, resulting in a similar number of predicted lensed stars. Hence, our results are relatively insensitive to the particular choice of $\Sigma_*$, as long as it is within the constrained range.

The final number of stars in each region is the result of two competing effects. On one side, larger magnification factors make the detection of the more abundant fainter stars easier, but on the other side, regions with larger magnification factors are considerably smaller. There is a sweet spot where the combination of the two effects maximizes the probability of seeing a star. In the inset plot of Fig. 11, we show the total number of detected stars as a function of the distance to the caustic. The width of the bins in the inset corresponds to the width $\Delta_{\text{pc}}$ used in the calculation. Again, in this plot we have assumed that stars can exist beyond the HD limit. If one sets the HD limit as the maximum luminosity, the light blue and green bins would have ~5 and ~50 times fewer lensed stars, respectively. In that case, we deduce that the best chance of finding lensed stars is when the total magnification of the double image is $100 < \mu < 300$, corresponding to a physical distance of one to a few parsecs from the caustic (dark-blue ad light-blue curves, respectively). The fact that we do not observe numerous lensed stars at larger distances from the critical curve is a direct indication that the HD limit must also be applicable at this redshift.

The dependence with the mass (or luminosity) of the star is better appreciated in Fig. 12, where we show the expected number of detected stars above AB magnitude 28.5 in HST's filter F814W, and per mass interval. At short distances from the main caustic (large magnification factors), the fraction of stars that is detected is proportional to the adopted IMF. The black solid curve corresponds to the case with maximum magnification and assuming both counterimages of the lensed star form an unresolved pair (so the observed flux doubles, making the detection of fainter stars possible). At larger distances from the main caustic, as the magnification drops, lensed events naturally select the brightest stars. The most extreme case is given by the yellow curve (macromodel magnification $\mu \approx 67$), where the most likely stars to be detected are the most luminous ones. This dependency with the macromodel magnification offers a unique opportunity to probe the most luminous (and rare) stars in distant galaxies, since these would be the only ones visible at relatively large distances from the critical curves of powerful lenses. In the particular case of the Spock galaxy, since we do not see the predicted abundance of lensed stars at these large distances from the critical curve, we can conclude that these stars cannot exist, as discussed above.
<figure><img src="image_14.png" /><figcaption>**Fig. 12.** Mass dependence vs. distance in HST's F814W filter at depth 28.5 mag. This figure is the differential version of Fig. 11. The colors and linestyles of the curves have the same interpretation as in the previous figure (that is, yellow is for the farthest distance to the caustic and black for the closest). The dashed line again represents the case where Σ* = 5 M⊙ pc⁻² as in Fig. 11. The plateau in the yellow curve marks the transition in the M–L relation at M ≈ 55 M⊙.</figcaption></figure>

We observe a general trend with the depth of the observa-
tion that can be appreciated when comparing the previous result
with a shallower observation using the same filter, and assum-
ing the same temperature for the lensed star, but reaching only
27.5 mag. In Fig. 13, we have allowed Lmax to be larger than the
HD limit. This explains the counterintuitive result shown, where
shallower observations predict more events in the light-blue bin
than in the same bin, but for the deeper observations of Fig. 11.
This is due to the fact that Lmax is larger in the shallower obser-
vation, hence the total number of lensed stars is integrated over
a wider range of luminosities. That is, luminosities as high as
5 × 10⁶ L⊙ are permitted in the calculation of the shallow case,
but they are not permitted in the deeper case, which is better con-
strained and already excludes these luminosities because they
are not observed at magnification factors μ = μ⊙. As in the pre-
vious case, if the HD limit is adopted, the number of detections
drops significantly and all bins contain <1 star; the only excep-
tion is the black curve (closest distance to the caustic), where we
expect ~0.2 lensed stars. Hence, we can conclude that in shal-
lower observations of the Spock arc, we expect to detect lensed
stars at the critical curve, or very close to it, but with a relatively
low probability. In particular, we expect that one in five obser-
vations at depth F814W = 27.5 mag may show a lensed star.
But this estimate is made assuming all stars more massive than
M ≈ 20 M⊙ have T = 10 000 K. Cooler or much hotter stars are
less likely to be detected (for the same bolometric luminosity)
since the peak of their emission does not fall in the F814W fil-
ter; hence, this estimate needs to be corrected by the fraction of
stars with T ≈ 10 000 K, making the detection of such stars even
less likely.

In contrast, deep observations (as in Fig. 11) naturally
increase the number of detections, with the bulk of them mov-
ing farther away from the critical curve. The deeper the observa-
tion, the farther one can reach down the luminosity function, and
also the farther away from the critical curve one expects to see
lensed stars. In this case, the larger area of smaller-magnification
regions dominates over the luminosity of the brightest stars, and
we predict a larger number of detections in regions where the
magnification of a lensed star is μ ≤ 150 than in regions with
larger magnification factors. For galaxies at larger redshifts, stars

<figure><img src="image_15.png" /><figcaption>**Fig. 13.** Same as Fig. 11 but for a shallower observation where the limiting magnitude is 27.5 (T = 10 000 K). Most lensed stars are predicted at larger separations from the critical curve, but if one imposes as maximum luminosity the HD limit, lensed stars are observed mostly at the largest magnification factors – that is, closer to the caustic and critical curve.</figcaption></figure>

are fainter and one requires larger magnification factors, making
the situation similar to the shallow observations at lower redshift
discussed above.

For an even deeper observation of the Spock arc, we con-
sidered a depth of 29 mag similar to the one reached by the
Flashlights program with the wide filter F200LP. The result
is presented in Fig. 14. For easier comparison with the previ-
ous results, we required that all stars have the same temperature
(T = 10 000 K). This is a bluer filter that is more sensitive to hot-
ter stars, so despite the increased depth, we observe a reduction
in the expected number of lensed stars when compared with the
shallower observations of F814W at 28.5 AB mag depth. Ignor-
ing events above the HD limit (at M ≈ 41 M⊙), we expect about
40% fewer lensed stars with T = 10 000 K in the Flashlights
observation. However, as we see later, this filter performs much
better for hotter stars, and for which the HD limit is also higher.
The range of magnifications (or distances from the critical curve
or caustic) at which we expect to see lensed stars in this filter is
similar to the case of the F814W filter at 28.5 mag depth. After
integrating over the entire range of distances, the number of pre-
dicted stars is ~1, similar to the case of F814W, and in good
agreement with the observations.

We conclude that in general, each arc, redshift, filter, and
depth needs to be studied individually since there is no general
rule that applies to all. Next we study how these results depend
on the effective temperature of the star.

8.1. Dependence on star temperature

The results presented above are obtained assuming all lensed
stars have the same temperature. Obviously, this is an unrealistic
assumption since SG stars can have a wide range of tempera-
tures from ~3500 K to several tens of thousands of degrees. For
a given redshift and temperature of the lensed stars, there is an
optimal filter for which the star is more likely to be detected.
We present the expected number of detections above AB mag-
nitude 28.5 as a function of the star temperature for the HST
filter F814W in Fig. 15. For this particular filter, star redshift,
and depth, the most likely star to be detected has T ≈ 10 000 K,
and as discussed above, we should expect to find these stars with
magnifications in the range 100 < μ < 300. We emphasise again
<figure><img src="image_13.png" /><figcaption>**Fig. 12.** Mass dependence vs. distance in HST's F814W filter at depth 28.5 mag. This figure is the differential version of Fig. 11. The colors and linestyles of the curves have the same interpretation as in the previous figure (that is, yellow is for the farthest distance to the caustic and black for the closest). The dashed line again represents the case where $Σ_*= 5 M_⊙ pc^{-2}$ as in Fig. 11. The plateau in the yellow curve marks the transition in the $M-L$ relation at $M ≈ 55 M_⊙$.</figcaption></figure>

We observe a general trend with the depth of the observation that can be appreciated when comparing the previous result with a shallower observation using the same filter, and assuming the same temperature for the lensed star, but reaching only 27.5 mag. In Fig. 13, we have allowed $L_{max}$ to be larger than the HD limit. This explains the counterintuitive result shown, where shallower observations predict more events in the light-blue bin than in the same bin, but for the deeper observations of Fig. 11. This is due to the fact that $L_{max}$ is larger in the shallower observation, hence the total number of lensed stars is integrated over a wider range of luminosities. That is, luminosities as high as $5 × 10^6 L_⊙$ are permitted in the calculation of the shallow case, but they are not permitted in the deeper case, which is better constrained and already excludes these luminosities because they are not observed at magnification factors $\mu = \bar{\mu}$. As in the previous case, if the HD limit is adopted, the number of detections drops significantly and all bins contain <1 star; the only exception is the black curve (closest distance to the caustic), where we expect ~0.2 lensed stars. Hence, we can conclude that in shallower observations of the Spock arc, we expect to detect lensed stars at the critical curve, or very close to it, but with a relatively low probability. In particular, we expect that one in five observations at depth F814W = 27.5 mag may show a lensed star. But this estimate is made assuming all stars more massive than $M ≈ 20 M_⊙$ have $T = 10 000$ K. Cooler or much hotter stars are less likely to be detected (for the same bolometric luminosity) since the peak of their emission does not fall in the F814W filter; hence, this estimate needs to be corrected by the fraction of stars with $T ≈ 10 000$ K, making the detection of such stars even less likely.

In contrast, deep observations (as in Fig. 11) naturally increase the number of detections, with the bulk of them moving farther away from the critical curve. The deeper the observation, the farther one can reach down the luminosity function, and also the farther away from the critical curve one expects to see lensed stars. In this case, the larger area of smaller-magnification regions dominates over the luminosity of the brightest stars, and we predict a larger number of detections in regions where the magnification of a lensed star is $\mu \le 150$ than in regions with larger magnification factors. For galaxies at larger redshifts, stars

<figure><img src="image_14.png" /><figcaption>**Fig. 13.** Same as Fig. 11 but for a shallower observation where the limiting magnitude is 27.5 ($T = 10 000$ K). Most lensed stars are predicted at larger separations from the critical curve, but if one imposes as maximum luminosity the HD limit, lensed stars are observed mostly at the largest magnification factors – that is, closer to the caustic and critical curve.</figcaption></figure>

are fainter and one requires larger magnification factors, making the situation similar to the shallow observations at lower redshift discussed above.

For an even deeper observation of the Spock arc, we considered a depth of 29 mag similar to the one reached by the Flashlights program with the wide filter F200LP. The result is presented in Fig. 14. For easier comparison with the previous results, we required that all stars have the same temperature ($T = 10 000$ K). This is a bluer filter that is more sensitive to hotter stars, so despite the increased depth, we observe a reduction in the expected number of lensed stars when compared with the shallower observations of F814W at 28.5 AB mag depth. Ignoring events above the HD limit (at $M ≈ 41 M_⊙$), we expect about 40% fewer lensed stars with $T = 10 000$ K in the Flashlights observation. However, as we see later, this filter performs much better for hotter stars, and for which the HD limit is also higher. The range of magnifications (or distances from the critical curve or caustic) at which we expect to see lensed stars in this filter is similar to the case of the F814W filter at 28.5 mag depth. After integrating over the entire range of distances, the number of predicted stars is ~1, similar to the case of F814W, and in good agreement with the observations.

We conclude that in general, each arc, redshift, filter, and depth needs to be studied individually since there is no general rule that applies to all. Next we study how these results depend on the effective temperature of the star.

## 8.1. Dependence on star temperature

The results presented above are obtained assuming all lensed stars have the same temperature. Obviously, this is an unrealistic assumption since SG stars can have a wide range of temperatures from ~3500 K to several tens of thousands of degrees. For a given redshift and temperature of the lensed stars, there is an optimal filter for which the star is more likely to be detected. We present the expected number of detections above AB magnitude 28.5 as a function of the star temperature for the HST filter F814W in Fig. 15. For this particular filter, star redshift, and depth, the most likely star to be detected has $T ≈ 10 000$ K, and as discussed above, we should expect to find these stars with magnifications in the range $100 < \mu < 300$. We emphasise again
<figure><img src="image_16.png" /><figcaption>**Fig. 14.** Number of stars detected per mass interval as a function of distance to the critical curve (or magnification) in filter F200LP at depth 29 mag (Flashlights). All stars are assumed to have the same temperature *T* = 10 000 K. The colors and styles of the lines are the same as in previous figures (yellow for the farthest distance to the caustic and black for the closest one).</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>**Fig. 15.** Number of stars vs. star temperature and macromodel magnification (or distance to the critical curve) for HST's F814W filter and at limiting AB magnitude of 28.5. For each bin we have assumed that all lensed stars have the temperature of that bin, so one needs to correct the total number in each bin by the fraction of stars with that particular temperature.</figcaption></figure>

that at the lowest magnification values considered in this figure, only stars above the HD limit can be detected, so if these stars do not exist with the abundances we have assumed, the number of stars in these bins should be considerably reduced.

Regarding the bluer F200LP filter, the number of detections as a function of temperature and magnification is shown in Fig. 16). The number of lensed stars in the temperature range 10 000 K < *T* < 20 000 K is very similar to the shallower F814 = 28.5 mag observation, but we observe an increase in the predicted number of events for stars with *T* > 20 000 K. The reduc-

<figure><img src="image_18.png" /><figcaption>**Fig. 16.** Similar to Fig. 15 but for HST's F200LP filter and limiting AB magnitude 29.</figcaption></figure>

tion in number of events, when one considers the HD limit, is still applicable for the bins with the smallest magnification and *T* < 15 000 K, but hotter stars can be more luminous than $6 \times 10^5 L_\odot$. Hence, we expect to see more lensed stars in F200LP at 29 AB mag depth than in F814W at 28.5 mag depth. This matches well the fact that in each of the two epochs of the Flashlights program a lensed star was spotted, while the observations with F814W had a smaller success rate.

When we account for all of the temperatures and possible magnifications, the expected number of detected stars is summarized in Table 1. The numbers in the columns with $\delta = 0$ represent the sum of all the numbers in Figs. 15 and 16 divided by 5, which is the number of temperature bins. This factor 5 is the fraction of stars that fall in each temperature bin when one takes $\delta = 0$ in Eq. (9). That is, we are assuming that the number of stars per temperature bin (before magnification) is homogeneous. In the table we show the case where the HD limit is imposed and when stars above this limit are allowed to exist. The numbers in the case where the HD limit is imposed are relatively consistent with the observations (as discussed later), while the numbers when the HD limit is not imposed are not consistent. Hence, our results in Table 1 clearly favor a scenario where the HD limit is still valid at *z* = 1. We also include the case of JWST observations at 2 µm and reaching 29 mag. The gain provided by JWST is obvious, providing a factor ~8 more detections. If our model is correct, in every single pointing of JWST we should observe ~5 lensed stars. Different observations will show new stars, or variations in the flux of stars observed in previous observation, or a combination of both (new stars and flux fluctuations). The high sensitivity of JWST to the most luminous red SGs at wavelengths $\lambda > 2$ µm will allow a constraint on their maximum luminosity and establish the validity of the HD limit at this redshift. The last two columns in Table 1 show the predictions for alternatives values of the *B/R* ratio (or $\delta$). These are discussed in more detail in Sect. 9.
<figure><img src="image_15.png" /><figcaption>**Fig. 14.** Number of stars detected per mass interval as a function of distance to the critical curve (or magnification) in filter F200LP at depth 29 mag (Flashlights). All stars are assumed to have the same temperature *T* = 10 000 K. The colors and styles of the lines are the same as in previous figures (yellow for the farthest distance to the caustic and black for the closest one).</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>**Fig. 15.** Number of stars vs. star temperature and macromodel magnification (or distance to the critical curve) for HST's F814W filter and at limiting AB magnitude of 28.5. For each bin we have assumed that all lensed stars have the temperature of that bin, so one needs to correct the total number in each bin by the fraction of stars with that particular temperature.</figcaption></figure>

that at the lowest magnification values considered in this figure, only stars above the HD limit can be detected, so if these stars do not exist with the abundances we have assumed, the number of stars in these bins should be considerably reduced.

Regarding the bluer F200LP filter, the number of detections as a function of temperature and magnification is shown in Fig. 16). The number of lensed stars in the temperature range 10 000 K < *T* < 20 000 K is very similar to the shallower F814 = 28.5 mag observation, but we observe an increase in the predicted number of events for stars with *T* > 20 000 K. The reduc-

<figure><img src="image_17.png" /><figcaption>**Fig. 16.** Similar to Fig. 15 but for HST's F200LP filter and limiting AB magnitude 29.</figcaption></figure>

tion in number of events, when one considers the HD limit, is still applicable for the bins with the smallest magnification and *T* < 15 000 K, but hotter stars can be more luminous than $6 \times 10^5 L_\odot$. Hence, we expect to see more lensed stars in F200LP at 29 AB mag depth than in F814W at 28.5 mag depth. This matches well the fact that in each of the two epochs of the Flashlights program a lensed star was spotted, while the observations with F814W had a smaller success rate.

When we account for all of the temperatures and possible magnifications, the expected number of detected stars is summarized in Table 1. The numbers in the columns with $\delta = 0$ represent the sum of all the numbers in Figs. 15 and 16 divided by 5, which is the number of temperature bins. This factor 5 is the fraction of stars that fall in each temperature bin when one takes $\delta = 0$ in Eq. (9). That is, we are assuming that the number of stars per temperature bin (before magnification) is homogeneous. In the table we show the case where the HD limit is imposed and when stars above this limit are allowed to exist. The numbers in the case where the HD limit is imposed are relatively consistent with the observations (as discussed later), while the numbers when the HD limit is not imposed are not consistent. Hence, our results in Table 1 clearly favor a scenario where the HD limit is still valid at *z* = 1. We also include the case of JWST observations at 2 µm and reaching 29 mag. The gain provided by JWST is obvious, providing a factor ~8 more detections. If our model is correct, in every single pointing of JWST we should observe ~5 lensed stars. Different observations will show new stars, or variations in the flux of stars observed in previous observation, or a combination of both (new stars and flux fluctuations). The high sensitivity of JWST to the most luminous red SGs at wavelengths $\lambda > 2$ µm will allow a constraint on their maximum luminosity and establish the validity of the HD limit at this redshift. The last two columns in Table 1 show the predictions for alternatives values of the *B/R* ratio (or $\delta$). These are discussed in more detail in Sect. 9.