GoKawiilConventions
Throughout this work we assume a flat Λ (dark energy) cold-dark-matter cosmology with matter density parameter Ω m = 0.315 and a Hubble constant H 0 = 67.4 km s−1 Mpc−1 (ref. 46). All reported magnitudes are in the AB system. Following the lensing model of ref. 13, we adopt a flux magnification factor μ = 6.2 and a shear factor of 3.52 for our source (image A of QSO1). Hence, 1 arcsec in the image plane corresponds to 1.52 physical kiloparsecs. For robustness tests, we use the Bayesian information criterion (BIC), defined as BIC ≡ χ2 + klnn, where k is the total number of model parameters and n is the number of points fitted; a decrease in BIC, ΔBIC ≥ 5, between two models was required for robust preference of one over the other, although our main conclusions remain unchanged even if a stricter ΔBIC ≥ 10 threshold is adopted.
Data reduction
We use data from the BlackTHUNDER NIRSpec integral field unit survey, focusing on the 7.3-hour exposures with the G395H grating, giving a nominal spectral resolution R = 3,700 at the wavelength of Hα15. The NIRSpec integral field unit was centred on image A of QSO1 (right ascension 00:14:19.161; declination −30:24:05.664)12. A detailed description of the reduction procedures is available in refs. 4,15; however, a summary is provided here for context.
The spectra were extracted following the procedures of ref. 47, but using version 1.17.1 of the JWST pipeline. At z = 7.04, the Hα line falls just outside the nominal wavelength coverage; however, the F290LP filter does not cut off longer wavelengths and the detector efficiency allows to recover Hα emission. We perform this recovery by extrapolating the wavelength solution, flat-field curves and the grating-equation-derived line spread function (LSF) out of the nominal range and towards the detector sensitivity limit of λ = 5.34 μm. The peak of the Hα line of QSO1 falls on λ = 5.278 μm; hence, our modification readily recovers the entirety of Hα emission. Although flux calibrations beyond the nominal range may suffer inaccuracies, the primary interest of this work is a kinematics study; hence, our key kinematics results are insensitive to flux calibrations. The BH mass measurements are more affected. However, the square-root dependence of the BH mass on luminosity means that flux calibrations have to be wrong by an order of magnitude to significantly impact the measurements.
The nominal spaxel scale of the processed data was 0.05″; however, utilizing the large number of dithers, we are able to oversample the cube to a scale of 0.02″ per spaxel without incurring significant sampling artefacts. We choose the 0.02″ cube for the main kinematic and spectroastrometric analysis, with the 0.05″ cube used to perform consistency checks, ensuring that our results are not pixel sampling artefacts.
Spectroastrometry of the rotation curve
To constrain the density profile of QSO1, we combine spectroastrometric measurements with resolved kinematics. The technical details of spaxe-by-spaxel fitting and spectroastrometry are given in Supplementary Information sections 1.1 and 1.2; here we summarize that we subtract the broad Hα emission from the cube and create images of different velocity channels of narrow Hα for which centroids can be obtained at sub-point-spread-function scales (provided a sufficient signal to noise48) and used to map dynamics below the nominal instrumental resolution49,50. The fiducial spectroastrometric analysis utilizing two velocity channels for higher signal to noise is shown in Extended Data Fig. 1. However, as shown in Fig. 2, splitting the line into finer bins does not change our results.
We infer the outer parts of the rotation curve by binning the line-of-sight velocity field (shown in Fig. 1) on scales >60 pc to avoid beam smearing. This procedure resulted in 4 bins covering the negative and positive sides of the rotation curve with ⟨vsini⟩ ≈ 10 km s−1 (where v is the line-of-sight velocity and i is the inclination angle) with nominal uncertainties of order 1 km s−1. However, these uncertainties, estimated through the standard root mean square (rms) weighting scheme, do not take into account the velocity field cross-correlation between spaxels of each bin owing to beam smearing and hence are probably significantly underestimated. An a priori derivation of the covariance matrix is intractable as it would require fitting individual dithers, which have far too low signal to noise. We thus use an empirical approach—scaling the naive rms-derived errors until the optimal model in the family of models fitted has \({\chi }_{{\rm{R}}}^{2}=1\). This yields an upper limit on the possible errors as it assumes that the optimal model is the ground truth. Consequently, using this method, we can establish a lower limit on the significance of the optimal model (which turns out to be a point mass) over other models considered.
For the spectroastrometric data points, we use a flux weighted average of the velocity channels, giving ⟨vsini⟩ = 51 ± 4 km s−1. The factor of sini is written to explicitly state that these are the projected values, uncorrected for inclination.
... continue reading