When, where, and how star formation happens in a galaxy pair at cosmic noon using CANUCS JWST/NIRISS grism spectroscopy

ABSTRACT

Spatially resolved studies are key to understanding when, where, and how stars form within galaxies. Using slitless grism spectra and broad-band imaging from the CAnadian NIRISS Unbiased Cluster Survey (CANUCS), we study the spatially resolved properties of a strongly lensed (μ = 5.4|$\pm$|1.8) z = 0.8718 galaxy pair consisting of a blue face-on galaxy (10.2 |$\pm$| 0.2 log(⁠|$M/M_\odot$|⁠)) with multiple star-forming clumps and a dusty red edge-on galaxy (9.9 |$\pm$| 0.3 log(⁠|$M/M_\odot$|⁠)). We produce accurate H |$\alpha$| maps from JWST/NIRISS grism data using a new methodology that accurately models spatially varying continuum and emission line strengths. With spatially resolved indicators, we probe star formation on time-scales of |$\sim$|10 Myr (NIRISS H |$\alpha$| emission line maps) and |$\sim$|100 Myr (UV imaging and broad-band SED fits). Taking the ratio of the H |$\alpha$| to UV flux (⁠|$\eta$|⁠), we measure spatially resolved star formation burstiness. We find that in the face-on galaxy both H |$\alpha$| and broad-band star formation rates (SFRs) drop at large galactocentric radii by a factor of |$\sim$|4.7 and 3.8, respectively, while SFR over the last |$\sim$|100 Myrs has increased by a factor of 1.6. Additionally, of the 20 clumps identified in the galaxy pair we find that 7 are experiencing bursty star formation, while 10 clumps are quenching, and 3 are in equilibrium (either being in a state of steady star formation or post-burst). Our analysis reveals that the blue face-on galaxy disc is predominantly in a quenching or equilibrium phase. However, the most intense quenching within the galaxy is seen in the quenching clumps. This pilot study demonstrates what JWST/NIRISS data can reveal about spatially varying star formation in galaxies at Cosmic Noon.

galaxies: evolution, galaxies: formation

1 INTRODUCTION

Knowing when, where, and how star formation occurs within galaxies is crucial to understanding how galaxies grow, evolve, and quench. A powerful way to constrain when star formation occurs is to look at star formation rate (SFR) indicators over various time-scales. Observable indicators such as H |$\alpha$| and rest-frame UV flux have been used to constrain SFRs (Kennicutt, Tamblyn & Congdon 1994; Kennicutt & Evans 2012) at look-back times of |$\sim$|100 and |$\sim$|10 Myrs, respectively. The difference in the SFR time-scales is due to where the flux originates from: H |$\alpha$| emission line flux is caused by short-lived O and B-type stars, while rest-frame UV flux is dominated by longer-lived stars. By comparing shorter and longer SFR time-scales, we can piece together the recent star formation history within galaxies.

Several spatially unresolved studies have used ratios of H |$\alpha$| and rest-frame UV flux to estimate how bursty a galaxy’s star formation is (Glazebrook et al. 1999; Iglesias-Páramo et al. 2004; Lee et al. 2009, 2011; Weisz et al. 2012; Domínguez et al. 2015; Mehta et al. 2017; Broussard et al. 2019; Emami et al. 2019; Fetherolf et al. 2021; Broussard, Gawiser & Iyer 2022; Asada et al. 2023; Mehta et al. 2023). The ratio of H |$\alpha$| to rest-frame UV flux (hereafter referred to as |$\eta$|⁠) works as a burstiness indicator as H |$\alpha$| flux responds more quickly to changes in the star formation history (SFH) than does rest-frame UV flux. Using SFHs to calibrate |$\eta$| we can determine whether a galaxy is bursting with star formation, is in an equilibrium state, or its SFR is decreasing.

One way to study where star formation occurs in galaxies is to exploit spatially resolved photometric and spectroscopic data. Studies using integral field units (IFUs) have researched the spatially resolved star formation properties of galaxies at z < 0.2 (Cid Fernandes et al. 2013; González Delgado et al. 2016; Belfiore et al. 2018; Medling et al. 2018; Quai et al. 2019; Avila-Reese et al. 2023). These works have studied topics such as inside-out quenching, spatially resolved star formation, the star-forming main sequence of galaxies, and specific star formation rate (sSFR) profiles. IFUs have also offered spatially resolved insights at redshifts ranging from 1 to 2, utilizing KMOS (Stott et al. 2014; Beifiori et al. 2017) and SINFONI (Förster Schreiber et al. 2011; Molina et al. 2017). These studies have explored various topics, including star formation in clumpy star-forming regions and the impact of feedback, the formation ages of massive galaxies, the correlation between metallicity gradients and sSFR, and H |$\alpha$| SFRs. Our goal is to do a similar analysis at Cosmic Noon (⁠|$z\sim 1-3$|⁠) using JWST/NIRISS slitless grism data. Among the key strengths of the JWST (Gardner et al. 2023) are its exceptional sensitivity and spatial resolution of the spectroscopic modes. These attributes enable the spectroscopic observations of distant galaxies with an unprecedented level of detail, surpassing what was previously possible with other facilities, including the Hubble Space Telescope (HST). Imaging studies conducted with JWST have already revealed the remarkable intricacies present in the structure of distant galaxies (e.g. Mowla and Iyer et al. 2022; Abdurro’uf et al. 2023; Asada et al. 2023; Giménez-Arteaga et al. 2023). However, the capabilities of JWST extend beyond broad-band or medium-band imaging alone: JWST is equipped with the Near-Infrared Imager and Slitless Spectrograph (NIRISS; Doyon et al. 2012, 2023) whose Wide Field Slitless Spectroscopy (WFSS; Willott et al. 2022) mode captures spatially resolved low-resolution (R|$\sim$| 150) spectroscopic information from all resolved objects in a 2.2 arcmin |$\times$| 2.2 arcmin field all at once.

Previous works have used HST (Nelson et al. 2016; Matharu et al. 2020, 2022; Mehta et al. 2023) and JWST/NIRISS slitless grism spectroscopy (Matharu et al. 2023) to study the spatial distribution of star formation. These works used a method of extracting emission line maps that ignores the spatially varying nature of galaxy spectral energy distributions. Specifically, until now, it has been common to assume that the underlying stellar populations do not vary across the extent of the galaxy – an assumption that is understandable in the case of low-SNR data. However, JWST/NIRISS observations call for the development of more advanced data analysis techniques to fully harness the immense potential offered by this new instrument. Particularly pressing is the need for a more comprehensive treatment when investigating target galaxies that show spatial variations in their spectral energy distributions (SEDs). The need for an improved approach was recognized before the launch of JWST through simulations that revealed potential biases in emission line maps and integrated line flux measurements (Sorba 2017). Our goal has thus been to develop such an improved approach to extracting information from spatially resolved slitless grism spectroscopy.

In this paper, we present a demonstration of our proof of concept analysis of the spatially resolved star formation properties of the Question Mark Pair (QMP), a complex system of two z = 0.8718 galaxies that we observed with NIRISS slitless spectroscopy and NIRCam imaging behind the lensing cluster MACS J0417.5−1154 as part of the CAnadian NIRISS Unbiased Cluster Survey (CANUCS). Fig. 1 highlights our target. The right-most panel of Fig. 1 shows a false colour image of MACS J0417.5−1154 using JWST/NIRCam filters, with zoomed-in images of the QMP shown using HST/ACS filters (leftmost) and JWST/NIRCam (centre). These images show four (of a total of five) lensed images of the QMP and highlight the complex structure of the pair of galaxies, with one of the galaxies being a heavily dust-obscured, edge-on galaxy, and the other a face-on disc-like system with strings of star-forming clumps. In Section 2, we give more information about the target galaxy pair and describe the data we used from the CANUCS survey. In Section 3, we discuss spatially resolved modelling of the NIRISS grism data. In Section 4, we discuss our results including the derived emission line map. In Section 5, we analyse the spatially resolved properties of the QMP. In Section 6, we summarize our conclusions. Throughout this work, we use AB magnitudes and assume a cosmology with |$\Omega _{m,0}$| = 0.3, |$\Omega _{\lambda ,0}$| = 0.7, and |$H_{0}$| = 70 km s^{|$^{-1}$|}.

Figure 1.

False colour images of the Question Mark Pair and MACS J0417.5−1154 (right-hand panel). The left two panels are zoom-in images of four of the multiply lensed images of the Question Mark Pair taken with HST and JWST. By comparing the JWST and HST images, we see how dusty the red edge-on galaxy is as it is barely visible in the HST/ACS imaging.

2 DATA AND TARGET

2.1 Target: the question mark galaxy pair

Our objective is to study spatially resolved star formation properties derived from JWST/NIRISS emission line maps. The QMP is an excellent case study as this complex system has several clumpy star-forming regions, is composed of two potentially interacting galaxies with different stellar populations, and has increased SNR due to it being highly lensed.

To analyse spatially resolved properties from JWST slitless grism data, we developed a new methodology for the recovery of emission line maps. To demonstrate, the value and effectiveness of our new method, we decided to demonstrate it on the QMP, which is a particularly complex galaxy system.

Fig. 2 shows the NIRCam images of the target of our study along with its NIRISS grism data. In each of the panels, we have marked several objects of interest, including objects 1.A–D which are 4 of the 5 lensed images of the QMP, the 5th is obscured by the |$z=0.441$| brightest cluster galaxy (BCG) which is labelled as object 2. The top right and bottom left panels show the NIRISS grism data (in all three filters) of the field near the QMP. Here, we can clearly see H |$\alpha$| emission in the F115W data for the QMP (1.A for GR150R and 1.A and 1.B in GR150C marked by an arrow). The bottom right panel shows more of the field along with the reason we refer to this galaxy system as the QMP as the red galaxy (not previously seen in the HST data) outlines a question mark, though it should be noted that the point of the question mark is an unrelated galaxy. The QMP system consists of two galaxies that are magnified and slightly distorted by the gravitational potential of the |$z=0.441$| foreground galaxy cluster MACS J0417.5−1154: a blue face-on galaxy that consists of a disc and several clumpy star-forming regions and a highly dust-obscured red edge-on galaxy. The face-on galaxy in our pair has been previously reported by Jauzac et al. (2019) to be a multiply lensed ‘complex ring galaxy’ at |$z=0.8718$| and called ‘The Doughnut’ by these authors. Our deep JWST/NIRCam imaging reveals that this galaxy is in fact not a classic ring galaxy but a disc sprinkled with star-forming clumps. This face-on disc galaxy is accompanied by a previously unknown edge-on red galaxy that is at a spectroscopically confirmed similar redshift (Section 4.3). We dubbed the galaxy pair ‘The Question Mark Pair’ given their prominence in what looks like a question mark-shaped multiply imaged lens system traced by the red galaxy in our NIRCam imaging (see bottom-right panel of Fig. 2).

Figure 2.

JWST/NIRISS grism data and false colour images. In each of the panels, we have marked objects of interest. The bottom right panel shows the multiply imaged QMP (1.A–D) in the MACS J0417.5−1154 cluster, here we see four of the five images of the galaxy, the fifth being under the BCG (2). The top left-hand panel shows the image of the QMP, we will be examining as it is the most magnified and least contaminated of the multiple images. The top right and bottom left panels show the NIRISS grism data in all orients and filters, focusing on the image of the QMP we are using.

We focus our attention on the top-left image of the QMP (Fig. 2; top left panel), which corresponds to Image 1.1 in Jauzac et al. (2019). We do so because this image is the most magnified among the multiple images of the QMP (μ = 5.4|$\pm$|1.8 in the Jauzac et al. 2019 lens model), while at the same time not suffering from strong shear, as can be seen by comparing its appearance with those of the other images of the galaxy pair in Fig. 2, and the least contaminated in the NIRISS grism data. By focusing on this complex, highly resolved system we highlight the main benefits of high-SNR spatially resolved NIRISS grism spectroscopy and showcase our new line-map extraction method.

2.2 Data

We use data from the MACS J0417.5−1154 cluster field taken as part of the CAnadian NIRISS Unbiased Cluster Survey (CANUCS; Willott et al. 2022). From this data set, we utilize broad-band imaging in JWST NIRCam (F090W, F115W, F150W, F200W, F277W, F356W, F410M, F444W) with exposure times of 6.4 ks each with SNR between 5 and 10 for an AB = 29 point source (Asada et al. 2023; Strait et al. 2023). In addition, we use archival HSTACS data (F435W, F606W, F814W). MACS J0417.5−1154 cluster galaxies were modelled and removed from the imaging (Martis et al. in prep). The point spread functions (PSFs) for broad-band imaging were empirically determined by median-stacking non-saturated bright stars. As a result, our broad-band data consist of PSF-convolved images standardized to the resolution of the JWST NIRCam F444W.

The JWST NIRISS grism data were taken in both the GR150R and GR150C grism utilizing the F115W, F150W, and F200W filters with exposure times of 19.2 ks in all filters (Matharu et al. 2023). The grism data were first processed using the grism modelling analysis software Grizli (Brammer 2022) which performs an end-to-end reduction of the data and modelling of sources. Subsequent processing of the data was performed after the 2nd, 1st, 0th, and 1st orders of the cluster galaxy spectra were modelled and removed (see Section 3.4)

3 METHODS

3.1 Segmentation

When using broad-band imaging, spatially resolved studies can treat each pixel or group of pixels as an independent object and fit galaxies pixel by pixel or region by region (e.g. Abraham et al. 1999; Sawicki 2012; Wuyts et al. 2012; Sorba & Sawicki 2015, 2018; Abdurro’uf et al. 2023). However, due to the self-contaminating nature of grism data (van Dokkum et al. 2011; Estrada-Carpenter et al. 2019) the elements cannot be treated as independent along the angle of dispersion. Therefore a pixel-by-pixel fit of the galaxy with the grism data would require that all elements be fit simultaneously. There are thus two main limitations to a pixel-by-pixel fit with grism data. The first of these is computing resource limitations. When fitting, each element must be forward-modelled in each exposure for each filter, and so, depending on the size of the galaxy and the number of exposures, we can easily exhaust our computing resources. The second main limitation is signal-to-noise ratio (SNR). Pixels in the outskirts of the galaxy have low SNR and fitting them will likely only return the priors. Therefore it would be best to group lower SNR pixels to create larger/higher SNR regions that can be fit to derive more robust fits while simultaneously reducing the computation time and resources used (Wuyts et al. 2012, 2013, 2014; Lang et al. 2014; Tadaki et al. 2014; Chan et al. 2016; Fetherolf et al. 2020).

For our spatially resolved analysis, we segment the galaxies into regions, grouping pixels using colours and fluxes from the broaduband images. To do this, we first take the brightest pixel in the direct NIRISS F150W image of the galaxy and using a nearest neighbours algorithm find the distance to all pixels in colour space using all possible combinations of our 11 broad-band filters. We then start to group the nearest neighbours of our brightest pixel until we reach an SNR limit (based on the NIRISS F150W SNR map); note that for the brightest pixels, this will often be satisfied with the single pixel alone. Here, we do not consider positional information and regions do not have to be contiguous. For this proof of concept study, we decided to fit with as many regions as possible therefore we chose an SNR sufficient to create 675 regions (SNR = 100), which was the number that we could fit with our available computing resources. Once we have reached the SNR limit for our region we mask those pixels from the broad-band images and start the process over again using the newest brightest pixel. We repeat this process until all pixels have joined a region. Note that the last region added will always sit below our SNR limit, but this will likely not affect our fits as this region will consist of low SNR background pixels. In the process of developing our nearest neighbours method, we conducted tests against other segmentation methods that utilize Gaussian mixture models and Voronoi tessellation. Our findings indicated that our nearest neighbours method generated regions with more homogeneous properties, such as lower scatter in sSFR, dust, H |$\alpha$| equivalent width, and so on. Furthermore, contamination modelling performed with our nearest neighbours segmentation resulted in more accurate models.

Fig. 3 shows our segmentation of the QMP. In the left-hand panel of Fig. 3, we see the results of our segmentation. In the centre of the galaxy, we have many single-pixel regions (274, |$\sim$| 40  per cent of our regions) with an overall median region size of 2 pixels and our largest region being 494 pixels. The right-hand panel of Fig. 3 shows the SNR map of our regions. Here, we see that most of the galaxy is shaded blue, which sits at our SNR limit of 100, with several regions sitting well above this limit (most of these regions being our 1 or 2-pixel regions).

Figure 3.

Segmentation map for the QMP. The first panel shows the 675 segments that the system has been broken into. In the centre of the system most of the segments are only 1 pixel and towards the outer regions (where the SNR is the lowest) we have larger regions. On the right, we show the SNR map of the system where the SNR has been calculated for each segment of our segmentation map. Here, we see that in the inner region, many of the single-pixel regions sit above our SNR limit while the outer regions sit at the SNR = 100 limit.

3.2 Spatially resolved priors/broad-band fits

Segmenting the galaxy into several regions results in a large parameter space we need to explore. Therefore to streamline the fitting process we derive a set of bespoke priors for each region by fitting to the spatially resolved broad-band images.

We start by using our PSF-matched NIRCam and ACS images, then apply the segmentation we derive in Section 3.1. The Flexible Stellar Population Synthesis code (FSPS; Conroy & Gunn 2010) is used with a combination of MILeS and BaSeL libraries and a Kroupa initial mass function (Kroupa 2001).

We fit each of the 675 regions using a set of 50 000 pre-generated FSPS models (which include nebular lines) varying sSFR, t|$_{50}$|⁠, both metallicities, ionization, and Av, with redshift set to z|$_{spec}$| = 0.8718 (Mahler et al. ; Jauzac et al. 2019) using the non-parametric SFHs described in Iyer & Gawiser (2017) and Iyer et al. (2019). From these fits, we derive posteriors for each parameter for each region. These provide maps of physical properties (such as stellar mass etc). and will also serve as priors for our fits to the grism data.

3.3 Spatially resolved grism modelling

The goal of our spatially resolved methodology is to generate highly accurate models of the grism spectra. As our end goal is to generate the purest emission line maps possible we must model the continuum and each of the emission lines accurately. To do this, each region will be fit with a linear combination of models, allowing for complex spectra to be fit relatively quickly.

Fig. 4 is a schematic of how our forward-modelling works:

• Step (1) Sample from the priors of each region and create multiple models using FSPS. The number of models generated for each region depends on multiple factors but ultimately is limited by computing memory. For the QMP, we generate 4 models per region. In the leftmost panel of Fig. 4, we show one of the models we may use for just three regions.

• Step (2) We forward-model each model spectrum (Estrada-Carpenter et al. 2023) for each region in all orients using grizli. For the QMP that will result in 2700 model grism spectra. The middle panels of Fig. 4 show three regions being forward-modelled in one filter and one orientation using the spectra from Step 1.

• Step (3) We find the linear combination of the model grism spectra which optimally fit the grism data, as seen in the rightmost panel of Fig. 4 (for a single filter).

• Step (4) We then explore our parameter space by repeating Steps 1–3 thousands of times to sample the priors well, once our library of models is built we find the best-fitting model that minimizes |$\chi ^2$|⁠.

Figure 4.

Schematic of the forward-modelling procedure for the multiregion analysis method. First, we show a subset of the input spectra from FSPS we used to model, here we utilize 4 spectra per region which results in a total of 4|$\times$|675  = 2700 spectra. We then forward-model each region using our segmentation map for each filter in both orients, resulting in 675|$\times 3\times$|2 = 129 600 forward-modelled grism spectra. We then sum the linear combination of the models using coefficients that optimize the fit resulting in the complete forward-modelled grism spectrum. For simplicity, only a single orient in one filter (F115W) is shown here, although our procedure simultaneously models all three filters in both orients.

3.4 Cluster galaxy removal

One of the difficulties in working with grism data is dealing with contamination, i.e. spectra from other objects that overlap the spectrum we are interested in (this may come from any order of dispersion). In the MACS J0417.5−1154 data, this is made even more difficult with the presence of large/extended foreground cluster galaxies. We do gain magnification by observing in the cluster field but the trade-off is heavily contaminated grism data. The left-hand panel in Fig. 5 shows a cutout (26 arcsec |$\times$| 26 arcsec) of the MACS J0417.5−1154 grism data centred on the |$z=0.441$| brightest cluster galaxy (BCG). As the cluster galaxies are so bright, we have to not only worry about the 0th and 1st-order dispersions but also the 2nd and -1st orders as well.

Figure 5.

Cluster galaxy removal from the grism data. In each of the panels objects of interest are labeled, these are the same labels as seen in Fig. 2. In the left-hand panel, we show the raw data in NIRISS F115W GR150C, the middle plot shows the NIRISS data with the cluster galaxies removed using the standard single-region method of modelling contamination with polynomials, and the right-hand panel shows the NIRISS data with the cluster galaxies removed using our spatially resolved fitting. Here, we see that our spatially resolved models do a much better job of removing the contaminating sources while leaving behind sources that were previously contaminated by the cluster galaxies (object 1.C).

In the centre panel of Fig. 5, we show the same field but with the single-region models for the cluster galaxies removed. Here, following the standard approach, single-region cluster galaxy models are generated by forward-modelling the galaxy spectra using a set of polynomials (up to and including third-order polynomials) using the image of the galaxy within the segmentation. The result is that the standard single-region models tend to oversubtract the contaminating spectra and fail to model the extended light, as seen particularly in the field centre. This oversubtraction is likely due to spatially varying stellar populations that are not accounted for in a single-region approach. Given the limitations of the standard procedure, a better approach is needed.

To properly subtract the cluster galaxies in the grism data we use our multiregion modelling approach of Section 3.3. To model the cluster galaxies, we use isophotal models of the BCGs produced by Martis et al. in preparation. These isophotal models replace our direct images and we use them to segment the galaxy and produce the forward-modelled grism spectra. We fit these bright foreground galaxies using our multiregion modelling procedure (Section 3.3) and subtract them from the NIRISS grism data. The results of this can be seen in the right-hand panel of Fig. 5. Our models show a clear visual improvement, including reducing oversubtraction and recovery of underlying spectra. In Fig. 5, we have labelled objects of interest using the same numbering scheme from Fig. 2. We see in the centre of the right-hand panel that we can recover grism spectra of galaxies that were completely contaminated by the BCG. Object 1.C is a highly deformed image of the QMP which sits close to the BCG, as seen in Fig. 2. Using the standard contamination subtraction method Object 1.C’s spectrum is unrecoverable, but using our multiregion methods not only do we recover the spectrum but we can also produce emission line maps from the data. Our cluster galaxy multiregion grism models provide a much better contamination removal than the standard process. From our modelling, we see that the BCG added |$\sim$|12 per cent contaminating flux to the QMP grism spectra which we successfully removed.

3.5 Emission line maps

Our goal is to produce as accurate an emission line map as possible to study the spatially resolved star formation properties of the QMP. To achieve this we must model the continuum well in addition to any nearby emission lines present within the galaxy spectra. Doing so is necessary because the way we produce an emission line map from the grism data is to first forward-model the galaxy using its best-fitting model, excluding the line we are interested in extracting. We then subtract the best-fitting model (excluding the line we are extracting) from the grism data and dither the resulting residuals (which will include the two orthogonal orients) to produce the emission line map.

Our multiregion approach allows for a spatially varying continuum and emission line strength, something that is not done using the standard method of forward-modelling the entire galaxy as a single entity. Our approach is superior since not accounting for a spatially varying continuum/emission line strength may lead to over or undersubtraction of the continuum/nearby emission lines leading to missing/excessive flux in the emission line maps, as we discuss next.

Fig. 6 shows an example of our multiregion fitting versus the standard single-region fitting using the best-fitting models for both methods. For the standard single-region method, grizli uses a set of optimally fit templates (13 FSPS models and 32 emission lines) to model the galaxy with the redshift set to z = 0.8718. panel (a) shows the F115W-GR150C grism data for the QMP, panel (b) shows the best-fitting model from the standard single-region method, and panel (c) shows the best-fitting model from our multiregion method. From these panels, we see that the multiregion models do an excellent job of modelling the real grism spectrum both in the continuum and emission lines, while the single-region model does not do as well at modelling the emission lines and the continuum. This is the most apparent in the H |$\alpha$| emission as since the single region method weighs the forward-modelled spectrum by the entire image of the galaxy, the best-fitting model for the single region method contains a large amount of H |$\alpha$| flux in the bulge of the blue face-on galaxy. This does not match what we see in panel (a). On the other hand, the multiregion best-fitting model shows a more similar distribution of H |$\alpha$| flux to the actual data (panel a). Panels (e) and (f) show the residual flux for the single-region and multiregion approaches, respectively. In panel (f), we see no discernible structure within the residuals while in panel (e), we see a considerable amount of residual in both the emission lines (both H |$\alpha$| + N ii and S ii) and along several regions where the continuum has been undersubtracted (both of which will negatively affect the emission line maps). Panel (d) shows the distribution of the residual flux with the multiregion residuals shown in blue and the single-region residuals shown in orange. We see in panel (d) that the multiregion residuals show a smaller range than those of the single-region model; they are also clustered symmetrically around zero, whereas those for the single-region model are systematically offset in the positive direction. In other words, the multiregion residual distribution is consistent with simple shot noise, whereas residuals from the single-region model show a systematic bias. Numerically, the multiregion residuals show a residual flux of 0.0|$^{+0.02}_{-0.02}$| (10|$^{-17}$| erg⁻¹s⁻¹ Å⁻¹ cm⁻²), or a fractional error of 3.5|$^{+9.5}_{-1.0}$| per cent. In contrast, the single-region residuals show a residual flux of 0.04|$^{+0.03}_{-0.03}$| (10|$^{-17}$| erg⁻¹ s⁻¹ Å⁻¹ cm⁻²), or a fractional error of 12.9|$^{+30.7}_{-3.9}$| per cent. It is clear that the multiregion method does better. This overall improvement means our multiregion models will contain much less contamination and will be truer to the actual spatial distribution of H |$\alpha$|

Figure 6.

Panel (a) shows JWST/NIRISS F115W-GR150R grism data for the QMP. Panel (b) shows the best-fitting model using the standard single-region approach, and panel (c) shows our best-fitting model using our multiregion approach. In the top panels, we see that our multiregion model more accurately reproduces the NIRISS grism data. Panel (d) shows the multiregion residual distributions for the single-region method and multiregion method; here, we see that the multiregion method residuals sit at 0.0|$^{+0.02}_{-0.02}$| (10^{|$^{-17}$|} erg⁻¹s⁻¹ Å⁻¹ cm⁻²) – 3.5|$^{+9.5}_{-1.0}$| per cent error – while the single region method produces offset residuals that sit at 0.04|$^{+0.03}_{-0.03}$| (10^{|$^{-17}$|} erg⁻¹s⁻¹ Å⁻¹ cm⁻²) – 12.9|$^{+30.7}_{-3.9}$| per cent error. Panel (e) shows the single-region residuals and panel (f) shows the multiregion residuals. Here, we see structure in the single-region residuals stemming from emission lines and continuum, while we see almost no structure in the multiregion residuals.

The difference in the two models leads to large differences in the emission line maps. Fig. 7 shows the first few steps to creating an emission line map, and shows the large difference in the two models. The first panel in each row of Fig. 7 shows the NIRISS F115W-GR150C data for the QMP. The second panel in each row shows the best-fitting model for the multiregion (panel b) and single-region (panel f) approach excluding H |$\alpha$| emission. By subtracting the second panel from the first, we should be left with only H |$\alpha$| emission flux, this is what is shown in the third panel on each row. Ideally, the residual flux should appear as a stamp of the QMP in H |$\alpha$| emission, the final panel in each row shows the direct image of the QMP in the orient of the grism data and this is what the residual flux should resemble. In panel (c), we see that the residual H |$\alpha$| emission using the multiregion approach does match our expectations as the flux is distributed like the QMP. However, in panel (g) we see that for the single-region approach, we are left with an excess of continuum flux that was not properly subtracted. This will result in an extended H |$\alpha$| emission profile that is not really there. When we complete the process and derive the final emission line maps, we find that the standard single-region method produces a map with 50 per cent more flux than our multiregion method. This extra flux is likely not from the target H |$\alpha$| emission line, but instead is an artefact of the undersubtracted continuum and poorly modelled S ii emission, all of which we see in Fig. 6 panel (e) and Fig. 7 panel (g). This excess flux will end up in the emission line maps which results in a 40 per cent higher SFR in the blue face-on galaxy and 70 per cent higher SFR in the red edge-on galaxy than our multiregion results. The differences in the SFR and flux may distort any result stemming from the standard single-region method.

Figure 7.

H |$\alpha$| residuals for the multi (top row) and single (bottom row) region models. In each row, the first panel is the NIRISS F115W-GR150C data for the QMP, the second is the best-fitting model excluding H |$\alpha$| emission, the third panel shows the residual flux when subtracting panel 2 from panel 1 in each row, and the last panel shows the NIRISS F150W direct image of the galaxy in the orient of the data for reference. Subtracting panel (2) from panel (1) should result in only H |$\alpha$| flux being left behind. When comparing the multiregion (panel c) to the single-region (panel g) H |$\alpha$| residuals we see that the multiregion residuals are very reminiscent of the QMP while the single-region residuals show extended H |$\alpha$| flux which, as it does not match what we see in the final panels, we can only conclude is due to undersubtraction of the continuum flux.

4 PHYSICAL PROPERTY MAPS

In this section, we employ our multiregion line-extraction method to investigate star formation in the QMP, following the successful demonstration of its superiority over the conventional one-region approach.

4.1 The question mark pair

Fig. 8 shows several property maps derived either from stellar population fits to the broad-band fluxes region by region (see Section 3.2) or the extracted H |$\alpha$| emission line map (panel d). Panel (a) has been corrected for lensing while panels (d, f, g, h, and i) have been corrected for lensing and dust. Note that our maps are not in the source plane, but values such as distance, mass, flux, and SFR have all been corrected for using our lens model described in the Appendix. We see in Fig. 8 panel (e) that the two galaxies have drastically different colours with the red edge-on being highly reddened and the blue face-on galaxy having multiple colours likely indicating varying stellar populations.

Figure 8.

Maps of the stellar mass surface density, segmentation of regions plus the H |$\alpha$| emitting clump in the red galaxy, dust (Av), H |$\alpha$| emission line flux surface density and SFR|$_{H\alpha }$| surface density, false colour image, change in SFR from the broad-band SFRs to the H |$\alpha$| SFRs, UV|$_{2300}$| surface density, |$\log (\eta _{2300})$| + equilibrium offset (1.65), and the SFR|$_{H\alpha }$| to SFR|$_{BB}$| ratio. Panel (a) has been corrected for lensing, while panels (d, f, g, h, and i) have been corrected for lensing and dust.

What makes the QMP so interesting is that these galaxies are possibly at the beginning of an interaction (as their morphologies do not seem to be disturbed). An interaction between the galaxy pair could lead to a burst of star formation, and this may be the reason why the blue face-on galaxy contains so many clumpy star-forming regions. To break down what the different components of the two galaxies are doing, we segment the galaxies into the regions shown in Fig. 8 panel (b), where we have split the galaxy pair into the blue galaxy bulge and disc, red bulge and disc, and bursting, equilibrium, and quenching clumps. We identify clumpy regions in the QMP by first applying a Gaussian smoothing algorithm to the NIRCam F115W image to create a smooth model of the galaxy. We then subtract this model from the original image leaving behind clumpy regions, note that this will also include bulges and patchy regions of the galaxy. We then employ Source Extractor for python (Bertin & Arnouts 1996; Barbary 2016) to segment the residual image. The mass map (see Fig. 8 panel a) is utilized to identify the two galaxy bulges. Then, to eliminate patchy/noisy regions, we impose an H |$\alpha$| SNR limit of 10, this should leave only clumpy star-forming regions in our segmentation. The remaining clumpy regions are then categorized into bursty, equilibrium, and quenching phases using the method detailed in Section 5.1. Finally, we distinguish between the blue disc and red disc using a colour-based classification utilizing a Gaussian mixture model.

We know that the two galaxies are at the same redshift as they both contain clumpy star-forming regions with H |$\alpha$| emission at the same observed wavelength. In addition, these two galaxies appear associated, as in the three unobscured lensed images of the QMP (Fig. 2 objects 1.a, 1.b, and 1.d) the two galaxies have the same alignment; if the two galaxies were just a chance projection from different redshifts then the two would align differently in the 3 lensed images.

Using the stellar mass surface density map and segmentation map (Fig. 8 panels a and b, respectively) we have determined that this would be a major merger with a stellar mass ratio of |$\sim$| 2:1 with the higher mass blue face-one galaxy having a stellar mass of 10.2 |$\pm$| 0.2 log(M/M|$_\odot$|⁠) and the lower mass red edge-on galaxy having 9.9 |$\pm$| 0.3 log(M/M|$_\odot$|⁠).

4.2 H |$\alpha$| map corrections

One of the main objectives of this work is to quantify the recent (10 Myr) SFR from the H |$\alpha$| flux. However, the flux in our H |$\alpha$| emission line maps is not purely from star-forming H ii regions and therefore we cannot simply convert it to SFR. In order to extract the SFRs from our H |$\alpha$| map we need to apply several corrections: we correct for dust attenuation, contamination of the H |$\alpha$| line by the neighbouring [N ii] line, and for the contribution of non-star-forming sources of H |$\alpha$| flux.

Many of the property maps we are going to analyse can be heavily affected by dust, so before we begin our analysis, we correct for dust attenuation using dust measurements from the spatially resolved broad-band data. In Fig. 8 panel (c), we show the attenuation map derived from spatially resolved SED fits to the broad-band images (see Section 3.2). In panel c, we see a region of high attenuation (⁠|$\gt $|1 A|$_V$|⁠), which coincides with the red edge-on galaxy (⁠|$1.2^{+0.4}_{-0.7}$| A|$_V$|⁠). In contrast, the blue face-on galaxy is less attenuated (⁠|$0.6^{+0.5}_{-0.2}$| A|$_V$|⁠). We note that as our A|$_V$| measurements are from the broad-band data they do not account for the possibility that the extinction towards the nebular regions may be higher than the continuum (Charlot & Fall 2000). To better account for the differences in dust attenuation, we use the relationship E_{|$_{stars}$|}(B–V) = (0.44 |$\pm$| 0.03)E_{|$_{gas}$|}(B–V) Calzetti (1997) for our dust corrections to H |$\alpha$|⁠. With the attenuation map for the QMP in hand, we apply a position-dependent dust correction to our images of interest. We do this using the Calzetti dust law (Calzetti et al. 2000), as this is the dust parametrization used in our SED fitting.

The flux in each pixel of our H |$\alpha$| maps is actually H |$\alpha$| + [N ii], as these lines will always be blended in the grism spectra. To account for this, we estimate the [N ii] to H |$\alpha$| ratio by deriving the spatially resolved metallicity. We do this by first deriving the ionization map using S32 ([S iii]|$\lambda \lambda$| 9069, 9532/[S ii]|$\lambda \lambda$| 6717, 6731), and then deriving metallicity using S23 ([S iii]|$\lambda \lambda$| 9069, 9532  + [S ii]|$\lambda \lambda$| 6717, 6731/H |$\alpha$|⁠) via equations from Kewley, Nicholls & Sutherland (2019). Note that [S iii]|$\lambda$| 9069 falls outside of the NIRISS F200W coverage, so to convert the [S iii]|$\lambda$| 9532 line flux to the total doublet line flux we use the [S iii]|$\lambda$| 9532/[S iii]|$\lambda$| 9069 flux ratio of |$\sim$|2.5 (Osterbrock & Ferland 2006; Tayal, Zatsarinny & Sossah 2019). As [N ii]/H |$\alpha$| is a metallicity indicator, we can then use the S23 metallicity to estimate [N ii]/H |$\alpha$| using equations from Kewley et al. (2019). This gives us a way to estimate [N ii]/H |$\alpha$| that accounts for metallicity gradients. Fig. 9 panel (a) shows the spatially resolved fraction of H |$\alpha$| (⁠|$f_{H\alpha }$|⁠) using this metallicity correction. We plan to further analyse the spatially resolved metallicity map of the QMP in a future paper.

Figure 9.

Correction maps applied to our H |$\alpha$| emission line maps. Panel (a) shows the correction which accounts for the H |$\alpha$| + [N ii] blended flux using metallicity estimations. Panel (b) shows the correction accounting for H |$\alpha$| emission due to DIG. By multiplying the H |$\alpha$| map by these correcting factors, we end up with a deblended H |$\alpha$| map composed of emission due to star formation.

Last, we need to account for H |$\alpha$| emission from sources other than star formation. The first alternate source we look at is diffuse ionized gas (DIG; Martin 1997; Sanders et al. 2017; Zhang et al. 2017). DIG can account for a significant amount of H |$\alpha$| emission, accounting for 30–60 per cent of H |$\alpha$| in local spiral galaxies (Zurita, Rozas & Beckman 2000; Oey et al. 2007) and is shown in Sanders et al. (2017) to be correlated with the H |$\alpha$| surface density. Using equation (24) from Sanders et al. (2017), we estimate the fraction of H |$\alpha$| originating from DIG (f|$_{DIG}$|⁠) and derive from it our DIG correction seen in Fig. 9 panel (b; shown here as f|$_{H\alpha }$| or 1—f|$_{DIG}$|⁠), where we see that in the outskirts of the galaxy as much as 60 per cent of H |$\alpha$| is from DIG. DIG also affects metallicity estimations as metallicities are calibrated using H ii region models, which do not match the conditions of DIG. Therefore, we also apply a correction to our [S ii] maps using equations from Sanders et al. (2017). We do not apply any corrections to our [S iii] map as DIG primarily enhances low-ionization emission lines (Sanders et al. 2020). Therefore, our metallicity and DIG corrections are intertwined. To estimate these corrections, we iterated over the derivation until the maps converged.

Other sources of possible H |$\alpha$| emission are active galactic nuclei (AGNs) and shocks. Strong AGNs have a very distinctive look in grism emission line maps. When no AGN templates are included the continuum and emission are not properly modelled creating a cross pattern in the emission line map (it is cross as NIRISS uses orthogonal dispersion angles). We do not see this cross pattern in our maps and therefore do not believe there is a strong AGN present. There may be a weak AGN, but we do not have the necessary emission lines needed to create a BPT diagram to identify AGN emission (Baldwin, Phillips & Terlevich 1981), and so do not correct for this possibility. As for shocks, we currently do not have the data needed to estimate the proportion of H |$\alpha$| due to purely shocks as kinematic information and high-resolution spectra would be necessary to make this distinction (Kewley et al. 2019). For this proof of concept work, we assume no additional H|$\alpha$| emission from AGN or shocks. This assumption may lead to biases in our results, and we plan to add AGN and shock templates in our larger sample studies and will study the spatially resolved effects of AGNs and shocks in greater detail in future work.

4.3 H |$\alpha$| emission line map

Fig. 8 panel (d) shows the (dust, metallicity, DIG, and lensing corrected) H |$\alpha$| emission line surface density map for the QMP, PSF-matched to the JWST/NIRCam F444W filter. The map is shown in log space to emphasize the H |$\alpha$| emission in the face-on disc. We see that the red edge-on galaxy shows less H |$\alpha$| flux (5.8 |$\pm$| 0.1 |$\times 10^{-17}$| erg⁻¹s⁻¹cm^{−|$^2$|}) when compared to the blue face-on galaxy (10.0 |$\pm$| 0.04 |$\times 10^{-17}$| erg⁻¹ s⁻¹ cm^{−|$^2$|}). We also see that a large amount of the H |$\alpha$| flux is concentrated in the clumpy star-forming regions (4.1 |$\pm$| 0.02 |$\times 10^{-17}$| erg⁻¹ s⁻¹ cm^{−|$^2$|}), accounting for 41 per cent of the total H |$\alpha$| flux. One caveat we would like to add is that the red edge-on galaxy is highly obscured and we may be missing some H |$\alpha$| flux. As we have no better dust measurement (H |$\beta$| falls below our wavelength range and therefore we cannot correct using the Balmer decrement), we move forward assuming our dust maps do provide sufficient corrections.

4.4 SFR

With the H |$\alpha$| map in hand, we can begin to study additional properties of the galaxy. Using the relationship from Kennicutt et al. (1994), we calculate the SFR from the H |$\alpha$| flux, which details recent star formation (time-scale of approximately |$\sim$|10 Myr). The resulting SFR|$_{H\alpha }$| map is shown in Fig. 8 panel (d). Here, we see that the regions with the highest SFRs in the blue face-on galaxy are concentrated in the inner chain of clumpy regions.

In total, the SFR of the blue face-on galaxy is 6.0 |$\pm$| 0.2 M_{|$_\odot$|} yr⁻¹ (sSFR −9.3 |$\pm$| 0.1 yr^{|$^{-1}$|}) while the red edge-on galaxy has an SFR of 3.6 |$\pm$| 0.3 M|$_\odot$| yr⁻¹ (sSFR −9.2 |$\pm$| 0.2 yr^{|$^{-1}$|}). Breaking the galaxies into their components we see that the red edge-on galaxy, |$\sim$| 87 per cent of its star formation occurs in its disc and |$\sim$|13 per cent in its bulge. For the blue face-on galaxy, |$\sim$| 54 per cent of its star formation occurs in its disc, while |$\sim$| 41 per cent of its star formation happens in its clumpy star-forming regions, with the bulge accounting for the remaining |$\sim$|5 per cent.

5 DISCUSSION

5.1 Star formation burstiness from H |$\alpha$|-to-UV flux ratios

One of the main properties we are interested in studying for the QMP is the spatially resolved burstiness of star formation. We do this by taking the ratio of the dust-corrected H |$\alpha$| and UV|$_{2300}$| (HST/ACS F435W – rest frame 2300 Å) maps (Fig. 8 panels d and g show the dust and lensing-corrected H |$\alpha$| and UV|$_{2300}$| maps, respectively). The largest potential source of error in our burstiness calculations is our dust correction, given the wavelength separation and the fact that the UV is continuum while H |$\alpha$| is an emission line. Please refer to Section 4.2 for our dust correction methods and any associated caveats. Both H |$\alpha$| and UV fluxes are commonly used to estimate SFR at different time-scales (Kennicutt et al. 1994; Kennicutt & Evans 2012). Stemming from gas around massive short-lived stars, H |$\alpha$| is sensitive to time-scales in the tens of millions of years. In contrast, the UV continuum flux originates from longer-lived stars and is sensitive to star formation on time-scales up to hundreds of millions of years. When a burst of star formation occurs, we see that H |$\alpha$| increases quickly, while the UV continuum flux is slower to respond.

$$\begin{eqnarray} \eta _{2300} = F_{H\alpha }/\lambda F_{\lambda ,UV}. \end{eqnarray}$$

(1)

Fig. 10 shows how the |$\eta _{2300}$| ratio responds to several illustrative star formation histories and suggests that we can use this flux ratio to constrain the time-variability of star formation (Lee et al. 2009; Weisz et al. 2012; Emami et al. 2019; Asada et al. 2023, 2024; Mehta et al. 2023). The subscript 2300 in equation (1) indicates that the measurement of the UV continuum is at rest frame |$\sim$| 2300 Å. This is because it is crucial to specify the wavelength at which the rest-frame UV flux is measured as the numerical values of |$\eta$| depend upon it despite the overall qualitative behaviour remaining unchanged.

Figure 10.

The connection between a galaxy SFH and |$\log (\eta _{2300})$|⁠. In each of the panels, the top sub-panel shows the SFH we are examining and on the bottom sub-panel, we show the evolution of |$\log (\eta _{2300})$| for galaxies with log(Z/Z|$_\odot$|⁠) = −2 (dashed line) and log(Z/Z|$_\odot$|⁠) = 0 (solid line) along with the region of equilibrium (shaded region) derived using FSPS models. When |$\eta$| is above the equilibrium region, the SFH is in a bursting phase, in the equilibrium region the SFH is either flat or slightly rising or transitioning from bursting to quenching, and below the equilibrium region the SFH is in a quenching phase (either long-term or short-term quenching).

All panels in Fig. 10 are structured so that the onset of star formation is at age |$=$| 0 with the top plot of each panel showing the SFH and the bottom panel showing the evolution of |$\log (\eta _{2300})$|⁠. Each of the |$\log (\eta _{2300})$| tracks derived using FSPS with the settings stated in Section 3.2, assuming no dust and an ionization of log(U)  = −3. In all panels of Fig. 10, we show a grey region defined as the |$\log (\eta _{2300})$| equilibrium (−1.75 |$\lt $||$\log (\eta _{2300})$||$\lt $| −1.55), which is roughly where the ratio resides during periods of constant or slightly rising star formation, as can be seen in panels A and B, respectively. The size of this region was determined by varying the metallicity of the models, with the lower bound determined by log(Z/Z|$_\odot$|⁠) |$=$| −2 (dashed line in bottom plots of each panel), and the higher bound determined by log(Z/Z|$_\odot$|⁠) |$=$| 0 (solid line in bottom plots of each panel).

Panels (a) and (b) show that equilibrium can occur with a constant or slightly rising SFR, indicating that to leave the equilibrium region a large sudden change in SFR is required. Panels (c) and (d) show such cases as these are bursty SFHs, with panel (c) showing bursts with rapid rises and declines, while panel (d) show bursts with rapid rises and slow declines. Panels (c) and (d) show that the rapid change in SFR will cause |$\log (\eta _{2300})$| to leave the equilibrium region while panel (d) also shows that during the fall of a burst |$\log (\eta _{2300})$| can spend an extended time in the equilibrium region. Panel (d) also shows that the value of |$\log (\eta _{2300})$| is correlated to how rapid the burst and subsequent fall-off are. Panels (e) and (f) show that being below the equilibrium region (the quenching phase) can occur without a burst and can happen through slow quenching (panel e) and rapid quenching (panel f). These SFHs indicate that being above the equilibrium region (bursty phase) only occurs when a burst happens while being below the equilibrium region can occur during the fall of a burst or long-term quenching. Therefore to understand the full context of |$\log (\eta _{2300})$| we must combine these results with spatially resolved SFHs.

5.2 Spatially resolved burstiness

Fig. 8 panel (h) shows the spatially resolved burstiness (⁠|$\log (\eta _{2300})$|⁠, see Section 5.1) of the QMP. Here, we have added 1.65 to |$\log (\eta _{2300})$| to set the equilibrium region to zero, with red indicating quenching, blue indicating bursting, and black indicating being in equilibrium. As noted in Section 5.1 sitting above the equilibrium is associated with bursting star formation while values below the equilibrium correspond to quenching whether that be long-term quenching (i.e. exponentially decaying star formation) or short-term quenching (downturn of the SFH associated with the end of a burst).

In Fig. 8 panel (h), we see that the majority of the blue face-on galaxy is quenching (65 per cent of the pixels below the equilibrium, 18 per cent in equilibrium, and 17 per cent sitting above – with 59 per cent of the mass below the equilibrium, 19 per cent in equilibrium, and 22 per cent above). In contrast, the red edge-on galaxy is bursting with star formation (10 per cent of the pixels below the equilibrium, 19 per cent in equilibrium, and 71 per cent sitting above – with 8 per cent of the mass below the equilibrium, 24 per cent in equilibrium, and 68 per cent above), and the blue galaxy is therefore either experiencing more long-term quenching or post-burst quenching. Asada et al. (2024) showed that the |$\log (\eta _{2300})$| offset below the equilibrium is correlated with how rapid the quenching is, therefore the quenching clumps are experiencing a more rapid quenching event than the disc as the quenching clumps have an overall lower |$\log (\eta _{2300})$| measurements. It may be that clumps that are in equilibrium/quenching may have entered a post-burst phase and it may be that the reason we see smaller |$\log (\eta _{2300})$| values in the non-bursty regions is that these clumpy regions can more efficiently use up their fuel, while the disc is inefficient and therefore is slowly quenching. Larger samples would be needed to say if that was a common trait, something we will address in future work. In Fig. 8 panel (i), we show another parametrization of burstiness using the SFR|$_{H\alpha }$| to the broad-band SFR (SFR|$_{BB}$|⁠) ratio (SFRr|$_{H\alpha /BB}$|⁠). SFR|$_{BB}$| was calculated for each region by taking the posterior SFH of that region, sampling from it 1000 times and then calculating the mean SFR over the final 100 Myr as this time-scale is similar to the time-scale of the UV emission. Panel (i) is shown using the same colour scheme as panel (h) where red indicates quenching, blue indicates bursting, and black indicates being in equilibrium (i.e. SFR|$_{H\alpha }$| and SFR|$_{BB}$| having similar values). The results using SFRr|$_{H\alpha /BB}$| are the same as with |$\log (\eta _{2300})$| where we find that the blue face-on galaxy mostly sits in or below the equilibrium, and the red edge-on galaxy is mostly bursting.

5.3 Spatially resolved properties

5.3.1 Stellar mass dependence

In Fig. 11, we show how the properties of |$\log (\eta _{2300})$| SFR, and sSFR correlate to mass for the regions shown in Fig. 8 panel (b; defined in Section 4.1) with the blue and red galaxy bulges shown as dots (blue and red, respectively), the blue and red disc shown as rings (blue and red, respectively), and the bursting/equilibrium/quenching clumps shown as cyan/black/yellow stars.

Figure 11.

SFH properties versus stellar mass for each of the regions indicated in Fig. 8 panel (b; defined in Section 4.1, blue and red galaxy bulges and discs as well as the bursting and non-bursting clumps). Here, the data points are plotted as arrows with the beginning of the arrow indicating the broad-band measurement and the end of the arrow indicating the H |$\alpha$| derived property. The span of the arrow indicates how large the change is over the relative time-scales (⁠|$\sim$|100 Myr) while the directions tell us if that value is rising or falling. Panel (a) shows the relationship between |$\log (\eta _{2300})$| and mass with the region of equilibrium shown in grey. Panel (b) shows the SFR–mass relationship with the star-forming main sequence from Nelson et al. (2021). In panel (c), we show the sSFR–mass relationship including the star-forming main sequence from Nelson et al. (2021). This figure shows that regions with bursting star formation do not necessarily need to be above the main sequence.

Panel (a) of Fig. 11 shows |$\log (\eta _{2300})$| versus stellar mass. We see no relationship between burstiness and stellar mass for the multiple regions. What we find is that the bursting clumps and the red disc and bulge sit above the equilibrium region, while the equilibrium clumps and blue disc/bulge sit in the equilibrium region, and only the quenching clumps sit below equilibrium. A question which arises is why the blue disc is in equilibrium when the majority of it is quenching (60 per cent of its mass is below equilibrium). The reason for this is that the bursting regions within the blue disc exhibit significantly higher H |$\alpha$| emission, which elevates the |$\log (\eta _{2300})$| from a quenching state to equilibrium. This highlights a drawback of integrated |$\log (\eta _{2300})$| measurements, as a bursting region can easily outshine the rest of the galaxy. This effect is particularly noticeable in larger regions that contain a mixture of bursting and non-bursting areas. However, it is less problematic in smaller, more uniform regions, such as the star-forming clumps.

Fig. 11 panel (b) shows the SFR–mass relationship for the regions along with the galaxy star-forming main sequence from Nelson et al. (2021). Each of the regions is marked with an arrow, with the arrow’s starting point being the SFR measured from the broad-band fits (correlating to a time-scale of |$\sim$| 100 Myr) while the endpoint is the SFR derived from the H |$\alpha$| emission (correlating to a time-scale of |$\sim$|10 Myr). Therefore the direction of the arrows shows us if the region’s SFR has been increasing or decreasing over the last |$\sim$| 100 Myr. In panel (b), we see that the blue and red discs sit in the main sequence, while the blue bulge sits below. We see a correlation between |$\Delta$|SFR and location on the main sequence. Note that we do not have a fit for the main sequence below log(M/M|$_\odot$|⁠) = 9, though relatively the bursting clumps have larger |$\Delta$|SFR values and sit above the equilibrium and quenching clumps which either have minor rises in SFR, or have a negative |$\Delta$| SFR.

In panel (c) of Fig. 11, we show the sSFR–mass relationship along with the star-forming main sequence from Nelson et al. (2021). We find that the regions of the QMP mostly fall into the category of being star-forming with log(sSFRs) between −10.0 and −8.0 (yr^{|$^{-1}$|}). One interesting takeaway from panel (b) is that the bulge of the red galaxy has a higher sSFR than the disc. This may indicate that the red edge-on galaxy is experiencing a starburst induced by a possible interaction between the galaxies.

When we combine the findings from the integrated values in all panels, it becomes clear that there is no correlation between a region’s position on the main sequence and |$\log (\eta _{2300})$| (therefore burstiness) as we see both quenching and bursting galaxies in the main sequence.

5.3.2 Dependence on distance from galaxy centre

In Fig. 12, we show the change in |$\log (\eta _{2300})$|⁠, SFR, and sSFR with the distance from the galaxy centre (bulge). Here, we focus on the blue face-on galaxy as we have not corrected the red edge-on galaxy for inclination. In each panel, the bursting/equilibrium/quenching regions are shown as cyan/black/yellow arrows, the blue face-on galaxy disc (with the clumps omitted) is shown in blue with the solid line indicating values derived from H |$\alpha$|⁠, and the dashed line for values from the broad-band (corresponding to 10 and 100 Myrs, respectively). The distances of each pixel have been corrected for lensing (see Appendix  A).

Figure 12.

The change in SFH properties in reference to the distance from the centre of mass of the galaxy. Only the blue face-on galaxy is shown here as the red face-on galaxy is highly inclined. In all panels bursting/equilibrium/quenching clumps are shown as arrows (with the direction of the arrow indicating the direction of change). The blue face-on galaxies disc is shown as a solid blue line for values derived from our H |$\alpha$| maps and a dashed line for values derived from the broad-band data with the solid region showing the gap between the two values. Panel (a) shows the evolution of |$\log (\eta _{2300})$|⁠, panel (b) shows the evolution in SFR, and panel (c) shows the evolution of sSFR. By comparing the panels, we see that burstiness is best correlated with sSFR, while not being well correlated with SFR or distance from the centre.

In Fig. 12 of panel (a), we see how burstiness evolves with distance from the centre. In both the disc and clumps, we see larger offsets (both bursting and quenching) closer to the centre of the galaxy. To determine how burstiness correlates with distance from the centre, we will need to utilize a larger sample of galaxies, something we plan to do in future work.

In panel (b) of Fig. 12, we see that the SFR for the blue galaxy rises outwardly in the central 1 kpc and then steadily drops going towards the outskirts of the galaxy. We also see that the clumpy star-forming regions in the outskirts have lower SFRs than those towards the centre of the galaxy. When comparing SFR|$_{H\alpha }$| to SFR|$_{BB}$| for the blue face-on galaxy we see that the regions with the largest increases are associated with bursting clumps with larger |$\log (\eta _{2300})$| values (as seen in panel a).

Fig. 12 of panel (c) shows how sSFR changes with distance from the centre, where we see that in the inner 2 kpc sSFR steadily rises from the core of the galaxy and peaks at |$\sim$|2 kpc (the approximate location of the inner chain of star-forming clumps), then lays mostly flat the rising in the outskirts (⁠|$\gtrsim$|4.4 kpc). As for the star-forming clumps, we see that they are mostly flat in sSFR with distance.

6 CONCLUSIONS

We have presented the JWST/NIRISS case study of the spatial distribution of H |$\alpha$| emission in a z = 0.8718 interacting galaxy pair. Our extraction of H |$\alpha$| emission goes beyond the standard analysis that assumes a constant stellar population at all locations, and – exploiting the high-SNR data that JWST can deliver – allows for the variation of stellar populations throughout the galaxy, resulting in more accurate models of the NIRISS grism data. We find that:

• Multiregion models of BCGs better remove the contaminating spectra than their single-region counterparts. Using our approach, it is now possible to analyse once unrecoverable contaminated spectra.

• Our multiregion models better model the JWST/NIRISS grism spectra, with residuals that are unbiased, unlike the single-region approach.

• We find that the single-region emission line maps are heavily contaminated with poorly subtracted continuum and emission line flux. These maps find 50 per cent more flux than our multiregion emission line maps, though this extra flux is all likely contamination.

This first illustration of our new technique is applied to a complex system of two potentially interacting galaxies. Using H |$\alpha$|/UV ratio (⁠|$\eta$|⁠) maps, which characterize the burstiness of star formation by exploiting the different time-scales of H |$\alpha$| (10 Myr) and UV (100 Myr) flux, we interpret our data as showing the following trends:

• The higher mass blue face-on galaxy is mostly quenching with a few patches or bursty star formation given their H |$\alpha$|/UV ratios.

• Notably among the 20 star-forming clumps, 10 have low |$\log (\eta _{2300})$| ratios that suggest star formation in these clumps may have already peaked and has started to decline.

• In contrast to the blue face-on galaxy, the lower-mass red edge-on galaxy has consistently higher |$\log (\eta _{2300})$| values, suggesting that it is currently bursting with star formation.

• When comparing |$\log (\eta _{2300})$| to our other burstiness indicator (SFRr|$_{H\alpha /BB}$|⁠) we find that the two indicators are in agreement that quenching/equilibrium regions in the galaxy align spatially with clumpy star-forming regions.

• We speculate that the reason that the quenching clumps have lower overall |$\log (\eta _{2300})$| values is that these regions may more efficiently use up their fuel, and therefore have entered a post-burst phase. This process may take longer in the disc.

Additionally, we compare spatially resolved SFRs from H |$\alpha$| emission and broad-band SED fits and find:

• The red edge-on galaxy’s SFR has increased over the last 100 Myr, combined with our |$\log (\eta _{2300})$| results this indicates that the galaxy is possibly experiencing a starburst as its bulge has a higher sSFR than its disc.

We have demonstrated the power of spatially resolved NIRISS spectroscopy to study the spatial distribution of star formation properties in distant galaxies using our new technique for extracting emission line maps in the era of high-SNR JWST slitless grism observations. Several works at low redshift (z|$\lt $| 0.2) have studied when/where star formation happens using large samples of IFU data (González Delgado et al. 2016; Belfiore et al. 2018; Medling et al. 2018; Avila-Reese et al. 2023). Our future work will extend the technique we introduced in this paper to larger, statistically significant, and complete samples of galaxies to determine when, where, and how star formation is happening at Cosmic Noon.

Acknowledgement

This research was enabled by grant 18JWST-GTO1 from the Canadian Space Agency and a Discovery Grant from the Natural Sciences and Engineering Research Council of Canada. MB acknowledges support from the ERC Grant FIRSTLIGHT, Slovenian national research agency ARRS through grants N1-0238 and P1-0188, and the program HST-GO-16667, provided through a grant from the STScI under NASA contract NAS5-26555. The authors acknowledge the Texas A|$\&$|M University Brazos HPC cluster and Texas A|$\&$|M High-Performance Research Computing Resources (HPRC, http://hprc.tamu.edu) that contributed to the research reported here. This research used the Canadian Advanced Network For Astronomy Research (CANFAR) operated in partnership with the Canadian Astronomy Data Centre and The Digital Research Alliance of Canada with support from the National Research Council of Canada the Canadian Space Agency, CANARIE, and the Canadian Foundation for Innovation.

DATA AVAILABILITY

Raw JWST data used in this work will be available from the Mikulski Archive for Space Telescopes (https://archive.stsci.edu), at the end of the 1-yr proprietary time (doi: 10.17909/ph4n-6n76). Processed data products will be available on a similar time-scale at http://canucs-jwst.com.

APPENDIX A: LENS MODELS

A strong lensing model is built using Lenstool (Kneib 1993; Jullo et al. 2007) and will be presented in a separate paper (Desprez et al., in preparation). This model is leveraging the multiple images constraints from Mahler et al. (2019) and Jauzac et al. (2019), as well as new multiple images systems and redshifts obtained from the CANUCS-JWST data. The model includes cluster–size mass haloes and galaxy sizes, described as double Pseudo-Isothermal Elliptical profiles (Elíasdóttir et al. 2007). Additionally, six clumps identified in the different images of the QMP are used as model constraints and allow us to produce a highly accurate lens model in the area of the QMP. The quality of our model can be ascertained from the fact that the average distance between the observed and predicted positions of the different multiple QMP clump images is 0.27 arcsec, whereas it is 0.44 arcsec for all the multiple images in the model. This indicates that the model is well-constrained, and particularly so around the QMP.

Author notes

Canada Research Chair.