TY - JOUR
AU1 - Troxel, M, A
AU2 - Long,, H
AU3 - Hirata, C, M
AU4 - Choi,, A
AU5 - Jarvis,, M
AU6 - Mandelbaum,, R
AU7 - Wang,, K
AU8 - Yamamoto,, M
AU9 - Hemmati,, S
AU1 - Capak,, P
AB - ABSTRACT The Nancy Grace Roman Space Telescope (Roman) mission is expected to launch in the mid-2020s. Its weak lensing program is designed to enable unprecedented systematics control in photometric measurements, including shear recovery, point spread function (PSF) correction, and photometric calibration. This will enable exquisite weak lensing science and allow us to adjust to and reliably contribute to the cosmological landscape after the initial years of observations from other concurrent Stage IV dark energy experiments. This potential requires equally careful planning and requirements validation as the mission prepares to enter its construction phase. We present a suite of image simulations based on galsim that are used to construct a complex, synthetic Roman weak lensing survey that incorporates realistic input galaxies and stars, relevant detector non-idealities, and the current reference 5-yr Roman survey strategy. We present a first study to empirically validate the existing Roman weak lensing requirements flowdown using a suite of 12 matched image simulations, each representing a different perturbation to the wavefront or image motion model. These are chosen to induce a range of potential static and low- and high-frequency time-dependent PSF model errors. We analyse the measured shapes of galaxies from each of these simulations and compare them to a reference, fiducial simulation to infer the response of the shape measurement to each of these modes in the wavefront model. We then compare this to existing analytic flowdown requirements, and find general agreement between the empirically derived response and that predicted by the analytic model. gravitational lensing: weak, techniques: image processing, large-scale structure of Universe 1 INTRODUCTION The nature of dark energy, which drives the accelerated expansion of the Universe, remains one of the most fundamental mysteries in physics twenty years after its discovery (Riess et al. 1998; Perlmutter et al. 1999; Albrecht et al. 2006; Frieman, Turner & Huterer 2008; Weinberg et al. 2013). A number of new experiments have been undertaken to probe dark energy using a variety of physical phenomena, including baryon acoustic oscillations, numbers and masses of galaxy clusters, galaxy clustering, redshift-space distortions, Type Ia supernovae, and weak gravitational lensing. Current-generation experiments are limited to some subset of these probes, but have already begun to expose interesting questions about the soundness of our standard cosmological model, Lambda-Cold Dark Matter, which will require more and better data in all of these probes to resolve. The Nancy Grace Roman Space Telescope (Roman)1|$^,$| 2 has been designed to take advantage of all of these probes to study dark energy and test general relativity with unprecedented systematic control (Spergel et al. 2015; Akeson et al. 2019; Dore et al. 2019). Weak gravitational lensing is a particularly powerful cosmological probe that is sensitive to both the expansion of the Universe and the growth of large-scale structure (Bartelmann & Schneider 2001; Mandelbaum 2018). In the past few years, the current generation of ground-based weak lensing experiments like the Dark Energy Survey (DES),3 Hyper-Suprime Cam (HSC) survey,4 and Kilo-Degree Survey (KiDS)5 have reached levels of precision that rival the previously best possible cosmological constraints when including a free dark energy equation of state (Troxel et al. 2018; DES Collaboration 2019; Hikage et al. 2019; Hildebrandt et al. 2020). These surveys have spurred the development of new algorithms and methods for galaxy shape measurement and weak lensing analysis (e.g. Huff & Mandelbaum 2017; Sheldon & Huff 2017; Zuntz et al. 2018), enhancing the potential power of weak lensing to unravel the fundamental mysteries we face in cosmology today. By the planned launch of Roman in the mid-2020s, we will have final results from the ongoing generation of weak lensing experiments (DES, HSC, and KiDS) and preliminary results from the Dark Energy Spectroscopic Instrument (DESI),6 the Large Synoptic Survey Telescope (LSST),7 and the Euclid mission.8 Faced with the unknown discovery potential of these experiments in the early 2020s, it is vital to maintain the agility of the Roman mission to respond with the best possible science, particularly in what is likely to be a systematics-dominated weak lensing field. The process of quantifying empirically the robustness of the design requirements of the Roman mission for weak lensing in the current phase of mission development is a critical task that this paper will partly address. Precise control of these systematics at the statistical precision offered by current Roman mission forecasts (Eifler et al. 2020a, b) will enable Roman to make crucial contributions to the study of new discoveries made in the early years of LSST and Euclid, and to the resolution of any remaining disagreements between surveys. Towards this goal, we describe in this paper a simulation framework designed to enable the empirical study of requirements flowing down from the Roman wide-field imaging survey, in particular for weak lensing. This simulation pipeline can incorporate a realistic simulated survey strategy, galaxy properties, and instrument effects to create a synthetic Roman wide-field imaging survey. We present in this paper a set of synthetic Roman imaging surveys covering approximately 6 deg2 to full 5-yr survey depth in one filter: a fiducial survey and 12 variations incorporating ways in which the point spread function (PSF) could be mis-estimated. The simulation incorporates realistic distributions of photometric properties for galaxies and stars; complex analytic galaxy models; a simulated observing strategy for a reference 5 yr, 2000 deg2 survey; and realistic detector effects, PSF models, and world coordinate system (WCS) solutions that match current Roman design specifications. We use a blending-free version of this simulation to test the impact on weak lensing science of these simulated wavefront modelling errors, including static, low-, and high-frequency biases. We discuss the Roman weak lensing survey, the current Reference Survey structure, and the weak lensing requirements process in Section 2. The synthetic survey simulation suite is described generally in Section 3, where we outline the simulated survey strategy, input galaxy, and star catalogues, and the galsim implementation of the Roman instrument used to simulate images. We discuss the specific simulation runs produced for this work to study wavefront error propagation in Section 4 and discuss the resulting biases and how these compare to the analytic requirements flowdown in Section 5. We discuss future plans for using this simulation suite in validating Roman requirements and algorithm design in Section 6 and conclude in Section 7. 2 ROMAN BACKGROUND We now proceed to describe requirements and the role of this suite of image simulations in verifying that the requirements flowdown is correct. We begin with a description of weak lensing in Roman that emphasizes the issues most closely tied to the image simulations (Section 2.1), and a high-level review of the requirements process in a cosmology project (Section 2.2). There we describe where in this process we need the mapping between the wavefront error and galaxy ellipticities, ∂ei/∂ψj (where ψj denotes a Zernike mode of the wavefront error). This mapping was obtained using a simplified analytic model, calibrated by toy simulations, in the Phase A requirements flowdown; this approach is described at a high level in Section 2.3, with technical details placed in the appendices. In the rest of this paper, we will use much more advanced image simulations, based on the galsim package, to estimate ∂ei/∂ψj. 2.1 Roman weak lensing The Roman weak lensing program has undergone significant evolution over the past decade (Green et al. 2011, 2012; Spergel et al. 2013, 2015; Doré et al. 2018; Akeson et al. 2019), but the basic philosophy has not changed. The next major advance in cosmology from weak lensing will require unprecedented control of systematic errors in photometric measurements (this includes, but is not limited to, shape measurement and PSF corrections). Roman will make this measurement with a thermally controlled telescope from beyond low Earth orbit, where the PSF can be made both stable and small. The imaging observations will be carried out in multiple filters and will have a cross-linked observing strategy within each filter to enable multiple internal cross-checks in the weak lensing signal. The current Reference Survey in the Roman Science Requirements Document (SRD)9 envisions shape measurements in three filters (J129, H158, and F184), where the PSF is at least half-Nyquist sampled (i.e. pixel size <λ/D; Nyquist sampling would be λ/2D). Here, F184 is the reddest filter on Roman, spanning 1.68–2.00 μm; it is between ground-based H and K, and was chosen based on the thermal constraints from the previously existing telescope hardware that was transferred to the program. Photometric redshift determination requires bluer filters as well. Roman itself will do a photometric survey in the Y106 filter since there was no ground-based option that would reach the required depth. Ground-based observations will be required for the z and bluer filters; the primary option for collecting these data will be LSST (LSST Science Collaboration 2009; Ivezić et al. 2019). The expected imaging depth is 26.9/26.95/26.9/26.25 mag AB in Y106/J129/H158/F184 (5σ point source; the limiting magnitude for the weak lensing samples is typically ∼2 mag shallower and depends on source size). The expected galaxy number density is 35 galaxies/arcmin2 (H158-band, the best for shape measurement) or 50 galaxies/arcmin2 (co-added bands). The Reference Survey also includes 10 per cent of the time devoted to medium-deep fields, which have 10× the exposure time over 1 per cent of the overall survey area, to calibrate the properties of the source galaxies. The Reference Survey area is limited to 2000 deg2 due to the need for internal redundancy (e.g. two passes over the sky in each of four filters means each region must be observed eight times) and the medium-deep fields, and the need to carry out many other observing programs as well in a 5-yr prime mission. This area is less than considered in some previous studies. Options for larger survey area have been considered, and could include a wider layer with less redundancy (e.g. an H158-band survey overlaid with LSST data), an extended mission (Roman has no consumable cryogens, and carries propellant for at least 10 yr), or both (Eifler et al. 2019). The actual survey – which may look different from the Reference Survey and be informed by developments in the coming years – will be chosen closer to launch. However, from the requirements point of view, we focus on enabling the Reference Survey. 2.2 The requirements process Every precision cosmology project has a requirements process to control both its statistical and systematic errors and ensure that the overall mission can achieve its science objectives. In the case of Roman, requirements on the Project (e.g. flight hardware and software or ground system support) were baselined early in the mission (the SRD was placed under configuration control in 2018), but requirements on science analyses are more flexible and will be fixed at a later date. The statistical error requirements are usually formulated in terms of survey area, depth and image quality in each filter, cadence (for time-domain programs), etc.; their relation to the science reach of the mission is handled by forecasting tools to be described (Eifler et al. 2020a, b). Systematic error control is much more difficult, and the approach may differ depending on whether a source of systematic error is observational or astrophysical. Usually, observational systematics (e.g. PSF calibration for weak lensing) can be budgeted within the systems requirements framework of a large project, whereas astrophysical systematics (e.g. baryonic feedback) are addressed through a combination of nuisance parameters, additional observations, and theory/simulation. These astrophysical systematics are important science team responsibilities but are not part of Project requirements and engineering reviews. In general, it is important to distinguish between known systematic biases in both categories, which can be calibrated and removed from the data, and uncertainties on that calibration, which cannot be removed and must either be small enough to ignore or marginalized over in an analysis. In the description of systematic errors here, we are referring to this residual uncertainty. We focus now on the approach to observational systematic errors; our focus is on the Roman process, but note that something similar has been done for other large weak lensing programs such as LSST and Euclid (Euclid Study Scientist & the Science Advisory Team 2010; Vavrek et al. 2016; Ivezic & LSST Science Collaboration 2018; LSST Dark Energy Science Collaboration 2018; Claver & the Systems Engineering Integrated Project Team 2019). First, one identifies a data vector that will contain the cosmological information. For setting Roman weak lensing requirements, the data vector is the concatenated list of shear power spectra and cross-power spectra C across tomographic bins. Other choices, such as including higher-order statistics, using all 3 × 2-point information, or working in correlation function space are possible, but given the tools available at the time of Project start these would have required additional tool development that did not fit in the schedule. Secondly, one identifies an error metric that summarizes the impact of a systematic error on the data vector. We have chosen the error metric |$Z^2 = \Delta {\bf C} {\boldsymbol{\Sigma }}^{-1}\Delta {\bf C}$|⁠, where ΔC is the bias on the data vector and |${\boldsymbol{\Sigma}}$| is the statistics-only covariance matrix. The metric Z is essentially a metric for the ratio of the systematic to the statistical error, and this depends on the solid angle Ω covered by the survey (⁠|$Z\propto \sqrt{\Omega }$|⁠). One also sets a limit on the maximum allowed error Z; in our case, we set Z = 0.5 at 2500 deg2 (or Z = 1 at 10 000 deg2), which means that the observational systematic errors are required to be below 50 per cent of the statistical errors in a 2500 deg2 survey and below 100 per cent of the statistical errors if the survey were to be extended to 10 000 deg2. Thirdly, we note that each category of observational systematic error contributes to Z. In cases where the errors are presumed independent, the Z2 values can be summed (i.e. Z obeys root-sum-square or RSS addition), and the ‘top-level’ budget for Z can be broken down into contributions from different sources. If a source of observational systematic error is parametrized by a parameter p (e.g. overall shear calibration), then a requirement on knowledge of p (parametrized by the 1σ uncertainty Δp) can be obtained by computing the sensitivity dC/dp and setting $$\begin{eqnarray*} \sqrt{\frac{\mathrm{ d}{\bf C}}{\mathrm{ d}p} \Delta p {\boldsymbol{\Sigma }}^{-1} \frac{\mathrm{ d}{\bf C}}{\mathrm{ d}p} \Delta p} \end{eqnarray*}$$(1) equal to the allocation for Z from that contribution. An important aspect of this budgeting is that, like a requirements flowdown, it is hierarchical – a top-level requirement on observational systematics may contain an allocation for shear calibration (one of several contributions), which itself may contain a branch for PSF size (one of several contributions), which itself may contain a branch for detector non-linearity, etc. In the life cycle of a cosmology project, more detail will be filled in first on the branches that have hardware impacts, and then branches related to algorithms or simulations later on. The details of our data vector, covariance, and systematics models are described in Appendix A. The systematic errors in the shear γ are broken down into additive biases (c) and multiplicative biases (m) in accordance with $$\begin{eqnarray*} \gamma ({\boldsymbol{\theta }},z;{\rm obs}) = [1+m(z)] \gamma ({\boldsymbol{\theta }},z;{\rm true}) + c({\boldsymbol{\theta }},z). \end{eqnarray*}$$(2) Appendix A then allocates the systematic budget for Z to residual uncertainties Δc (in different angular bins) and Δm. One challenge is that the shear biases may be redshift-dependent. Fortunately, when we start assigning portions of the shear systematic error budget to underlying root causes, we usually know something about the redshift dependence (for example, most PSF-related errors grow with redshift because the galaxies get smaller). Therefore, we have assigned each possible redshift dependence a weighting factor S, which represents the ratio of what fraction (in an RSS sense) of the error budget is taken up by a systematic with a given redshift dependence, relative to a systematic that is redshift-independent with the same maximum amplitude. A redshift-independent systematic has S = 1; due to covariance between redshift bins, it is possible to have S > 1. 2.3 Mapping from wavefront error to galaxy ellipticities – analytic approach The requirements based on Z are described in terms of shear systematics, but in order to be useful for engineering, we need to write a requirement in terms of wavefront errors. The key step to doing this is to write the derivative of the observed shear γobs with respect to the wavefront error ψj (where j denotes a Zernike mode). Because the PSF size and ellipticity are quadratic rather than linear in the wavefront error, it is necessary to take a quadratic expansion; then ∂γobs,i/∂ψj is linear in the static wavefront error ψj (Noecker 2010). Our approach is to use symmetries to categorize the possible quadratic terms, which gives four independent coefficients if we take the first 11 Zernike modes (up through spherical aberration). We can then use a suite of simple simulations to determine the four coefficients. Then – given a limit on the static wavefront error |$\Vert {\boldsymbol{\psi }}\Vert \le 92$| nm rms, set to achieve diffraction-limited imaging in J-band – we can analytically search the space of possible static wavefront errors and find the maximum possible |∂γobs,i/∂ψj| (units: nm−1). This worst-case sensitivity can be used to set requirements on knowledge of the Roman wavefront. A similar process can be used for changes in the PSF within an exposure (either line-of-sight motion, or wavefront jitter – i.e. beyond the tip-tilt modes – due to vibrations). The baseline plan for Roman will be to independently fit the line-of-sight motion contribution to the PSF in each exposure (Jurling & Content 2012), but not the wavefront jitter. This implies a requirement on the wavefront jitter to make its contribution to the PSF negligible. This requires a computation of the derivative of γobs with respect to the second moments of the PSF (units: mas−2); with respect to the variance or covariance of the wavefront jitter (units: nm−2); or with respect to the covariance of wavefront jitter and line-of-sight motion (units: nm−1 mas−1). In both of these cases, the source of bias in the shear measurement is in practice due to errors in the wavefront model leading to mis-estimation of the PSF model that is used for convolution of the galaxy model when fitting the model shape. All of these calculations are described in detail in Appendix B. We emphasize that while the time-dependent wavefront error and line-of-sight motion are two of the most difficult aspects of weak lensing, they are only a portion of the overall shear measurement error budget. Some other contributions related to the wavefront could come from small-scale field dependence of the wavefront due to figure errors on the fold mirrors (which are closer to an intermediate focus than a pupil) or flatness of the detectors; chromatic dependence of the wavefront; and polarization dependence of the wavefront (Lin et al. 2020). There are also sources associated with the calibration of the Roman detectors (Mosby et al. 2020), including but not limited to inter-pixel capacitance (IPC), persistence, count-rate dependent non-linearity, flat-field and dark current uncertainties, and the brighter-fatter effect. In practice, as we gain additional knowledge of as-built components, new terms are added (e.g. we are currently working on adding the vertical trailing pixel effect; e.g. Freudenburg et al. 2020), so there must be margin to cover these new developments in the top-level error budget for Z (the full budget is being updated and is beyond the scope of this paper). It is possible to group some of these terms together and form an intermediate level requirement on knowledge of the PSF moments, and then have, e.g. time-dependent wavefront errors as a sub-allocation, as was done for Euclid by Cropper et al. (2013). Given the Roman working group structure, with different groups focused on specific elements (e.g. detectors, filters, stability of the optical chain) and with the science team representatives in these groups ultimately looking after the shear bias requirements, we chose instead to treat all instrument-related systematics as sub-allocations of the shear bias requirements. 2.4 Limitations of the analytic approach The analytic approach to estimating the sensitivity to wavefront errors has some advantages: it is simple, maintains a close link to underlying physical principles, enables rapid exploration of the parameter space, and was available at an earlier stage of the project than the image simulations. However, it has some drawbacks: The analytic approach deals with single images, so it does not represent what happens when images are combined. This is especially relevant when the input images are undersampled at the native resolution of the Roman pixels. (The pixel scale is 110 mas, whereas Nyquist sampling would be λ/2D = 56, 68, or 80 mas at the average wavelength of J129, H158, or F184 bands, respectively.) The analytic approach computes the derivatives ∂γobs,i/∂ψj at one point in the focal plane. Therefore, it does not capture the correlations across the focal plane or tiling patterns; the distribution of systematic shear in two-point correlation function space or in power spectrum space is not captured. The analytic approach cannot be extended to include interaction of PSF errors with other aspects of the data, such as noise, detector systematics, blending/selection, etc., in the way that is possible with image simulations. For these reasons, we have also estimated the mapping from wavefront error to galaxy ellipticities using pixel-level image simulations with the galsim package. The version of the simulations used here is highly idealized – for example, the matching to the ‘truth catalogue’ means that selection/blending effects are not realistically implemented, some detector effects were not implemented, and the input galaxies have artificially prescribed shears and do not come from a realistic large-scale structure distribution. This is useful in the current study to enable us to uniquely isolate the impacts of wavefront errors on shear recover. Nevertheless, galsim as a tool is extensible and could be configured to use a realistic Roman input catalogue for future systematics studies. 3 SIMULATION SUITE To empirically test weak lensing requirements, methods, and algorithms in Roman, we have designed a synthetic survey suite that, while not entirely realistic in all object properties, contains sufficiently complex and representative objects so as to enable informative tests and preliminary algorithm development. This synthetic survey utilizes several external simulation and data sources, and generates Roman-like imaging using the galsim framework and its Roman module. The simulation framework is generally capable of producing a full Roman High-Latitude Imaging Survey (HLIS) in all filters matching Cycle 7 specifications.10 The code is publicly available.11 An example SCA image is shown in Fig. 1. The fiducial simulation run is available for download – this public data set is described in Appendix C. Figure 1. Open in new tabDownload slide A simulated ∼140 s exposure Sensor Chip Array (SCA) image, chosen for the presence of several large, bright galaxies and stars. Each SCA (HgCdTe H4RG) has a useable pixel grid of 4088 × 4088, with a pixel scale of 0.11 arcsec. A total of 18 SCAs make up the Roman camera. For comparison, the size of the Hubble Space Telescope Wide Field Camera 3 is shown as the blue outline. Six diffraction spikes due to the Roman secondary mirror support struts are clearly visible for most stars. Figure 1. Open in new tabDownload slide A simulated ∼140 s exposure Sensor Chip Array (SCA) image, chosen for the presence of several large, bright galaxies and stars. Each SCA (HgCdTe H4RG) has a useable pixel grid of 4088 × 4088, with a pixel scale of 0.11 arcsec. A total of 18 SCAs make up the Roman camera. For comparison, the size of the Hubble Space Telescope Wide Field Camera 3 is shown as the blue outline. Six diffraction spikes due to the Roman secondary mirror support struts are clearly visible for most stars. This approach to producing (to varying degrees) realistic, synthetic survey realizations is a common approach for weak lensing experiments, both at the catalogue level (MacCrann et al. 2018; Korytov et al. 2019) and the image level (Suchyta et al. 2016; Fenech Conti et al. 2017; Mandelbaum et al. 2018; Samuroff et al. 2018). These synthetic surveys can serve as sources of calibration or characterization, validation, or increasingly as end-to-end integration tests for measurement and analysis algorithms and pipelines. Our approach here is similar to the approach being implemented in parallel by the LSST Dark Energy Science Collaboration (DESC; Korytov et al. 2019; LSST Dark Energy Survey Collaboration 2020), with comparable levels of morphological complexity for weak lensing algorithm testing, but less complex true object properties. This approach is described in detail in the following subsections. 3.1 Simulation stages The simulation is broken into several stages: 3.1.1 Truth catalogue generation A truth catalogue is generated from the simulated input galaxy distribution, photometric galaxy catalog, and Milky Way simulation. The following true object properties are assigned to each simulated galaxy: (1) The sky position in right ascension (RA) and declination (Dec) from the simulated galaxy distribution; (2) Photometric properties (consistent Y106/J129/H158/F184 magnitudes, size, and redshift) drawn from a random object in the photometric galaxy catalogue; (3) Intrinsic ellipticity components drawn from a Gaussian distribution of width 0.27 (truncated at ±0.7); (4) A random rotation angle; (5) The ratio of fluxes in each of the three galaxy components: (a) de Vaucouleurs bulge, (b) exponential disc, and (c) random-walk star-forming knots (a maximum of 25 per cent of the flux assigned to the disc component can exist in the knots); (6) The gravitational lensing shear applied to the object, drawn from a discrete list of (e1, e2) ∈ {0, ±0.1}. Further details on the provenance of the galaxy catalogues and Milky Way simulation can be found in Sections 3.3 and 3.4, respectively. The true properties for all objects are saved in a single FITS table that is accessed by the following stages. 3.1.2 Image generation In this stage, an empty SCA image is initialized (4088 × 4088 pixels), and a model is built for each galaxy and star in turn, then drawn into the image. The galaxy models are built chromatically from the truth parameters for the object, with each component being assigned a different representative SED of types: S0 (bulge), SBa (disc), and Im (knots), respectively. The assigned SED is the same for all objects, since after redshifting the spectrum and applying the appropriate flux and size in each component, the model is converted to be achromatic in each passband to speed up the drawing (this is discussed further in Section 3.2.2). The intrinsic ellipticity, random rotation, and gravitational shear are then applied. We model stars as point sources with the SED of Alpha Lyra. Stars are also converted to be achromatic before drawing. Both stars and galaxies are then convolved with the appropriate PSF for the SCA (constant across the SCA in the fiducial simulation). An example of the PSF model for an object is shown in Fig. 2, and the PSF model is discussed in more detail in Section 3.2.2. We save images of the true PSF model both at native pixel scale and oversampled by a factor of 8, in stamps of native pixel size 8 × 8 at the position of each galaxy. Figure 2. Open in new tabDownload slide PSF model for SCA #1. The top row shows the model in native pixel scale, while the bottom row is oversampled by a factor of 8. From left to right: a comparison of the high-resolution (‘true’) model, the low-resolution model used in the simulation, and the difference of the two models. The colour bars are defined by the range of the high-resolution model. More detail about the PSF model is given in Section 3.2.2. The difference between the low- and high-resolution PSF models is negligible at the level required for this study. Figure 2. Open in new tabDownload slide PSF model for SCA #1. The top row shows the model in native pixel scale, while the bottom row is oversampled by a factor of 8. From left to right: a comparison of the high-resolution (‘true’) model, the low-resolution model used in the simulation, and the difference of the two models. The colour bars are defined by the range of the high-resolution model. More detail about the PSF model is given in Section 3.2.2. The difference between the low- and high-resolution PSF models is negligible at the level required for this study. The models are drawn in dynamically sized square stamps, the sizes of which are chosen to include at least 99.5 per cent of the flux. These stamps are then added to the SCA image and saved separately (if drawing a galaxy) to provide an isolated image of each simulated galaxy to allow for tests of the impact of blending. Objects that would have a postage stamp that overlaps the SCA image are drawn, such that light from objects in chip gaps are appropriately drawn on to the SCA, but we only save postage stamps for objects that have a centroid that falls on the SCA. We do not save isolated postage stamps of objects that have a stamp size of greater than 288 × 288 pixels, but they are drawn into the images. Finally, each isolated postage stamp is processed through the steps described in Section 3.2.3 to simulate the WFIRST observatory and detectors and written to disc. This means that blended objects will be modelled differently in the isolated postage stamp and full SCA images, since some detector effects are sensitive to the total flux in nearby pixel. When all objects are added to the full SCA image, it is also processed through these steps and written to a FITS image file. 3.1.3 MEDS creation We then compile the output across pointings of the isolated object stamps into MEDS (Multi-Epoch Data Structure) files.12 These files concatenate all exposures of unique objects to allow for fast access for object-by-object data processing (like shape measurement). Each MEDS file also stores for each object (and stamp) its original SCA, the object position and the stamp position within the SCA, the WCS for each stamp, the PSF model for each object, and other ancillary information and metadata. Each MEDS file contains all objects within a nside = 512 Healpixel13 (Górski et al. 2005; Zonca et al. 2019). 3.1.4 Shape measurement The galaxy shape is measured by jointly fitting a two-component model, de Vaucouleurs bulge and exponential disc, across all suitable exposures. Exposures where more than 20 per cent of the pixels are masked (i.e. the centroid falls too close to the edge of the SCA) are rejected. The model fit has seven parameters: e1, 2, px, y, half-light radius, flux, and bulge flux fraction, where e1, 2 is the component of the ellipticity and px, y is the pixel centroid offset. Both model components are constrained to have the same centroid, half-light radius, and shape. The minimization is performed using the ngmix 14 and mof 15 packages (Sheldon 2014). We also measure the PSF size and shape in the oversampled PSF model images using an adaptive moments method (Hirata & Seljak 2003). This stage writes a set of FITS files containing the galaxy and PSF measurement results and relevant truth catalogue information. 3.2 GALSIM The images in the simulation are rendered using the galsim software package (Rowe et al. 2015). This package has been extensively tested and has been shown to yield very accurate rendered images of galaxies and stars. Notably, the image rendering process has been shown to impart biases in the shapes of galaxies at a level much less than 10−3 for the kinds of objects we are simulating here. The galsim package is mostly generic with respect to the telescope and observational strategy, allowing for a wide variety of options in performing the simulation. However, it does have a sub-module (galsim.wfirst) that has a number of WFIRST-specific implementation details. Some of the code in this module pre-dates this work (e.g. Kannawadi et al. 2016), but some of it was developed specifically for this project, especially updating some of the details to match Cycle 7 information, and to reflect new information from laboratory tests of persistence in Roman sensors. The values used for this project correspond to the galsim.wfirst module in galsim release version 2.2.0. 3.2.1 World coordinate system The galsim.wfirst module has code to provide an estimate of the Roman WCS for each SCA given a rotation angle, date, and pointing direction. The WCS gives the 2D mapping from (x, y) coordinates on the image to RA and Dec on the sky. The specific orientations and gaps between the sensors were updated to match Cycle 7 specifications as part of the development work for this project. We create our scene of objects, including their surface brightness profiles, in sky coordinates (RA, Dec). galsim automatically accounts for the Jacobian of the WCS transformation when rendering the surface brightness profiles on each sensor’s pixels. Details such as the telescope distortion and variable pixel area are correctly accounted for in this process. 3.2.2 Point spread function For the PSF, we use a model of the Roman PSF from the galsim.wfirst module. While this module includes a high-resolution Cycle 7 estimate of the Roman spider pattern (i.e. the obscuration of the struts and camera in the pupil plane), we use a faster, low-resolution approximation, which gets the qualitative features correct, but has a slightly different detailed diffraction pattern. For the purposes of this study, we are insensitive to the differences between the two spider patterns, so we did not enable the slower, more accurate option. The PSF uses position-dependent (Zernike) aberration polynomials (Noll 1976), based on an investigation of the field-dependent wavefront errors used in the original Cycle 7 documentation. Aberrations between the tabulated positions are estimated using bilinear interpolation of the tabulated values. More details of the PSF model approximation and its implementation are described in Appendix E. The wavelength-dependent features of the PSF, such as the width of the Airy diffraction pattern, and the wavelength-dependence of the aberrations, are taken at the effective wavelength of the observation bandpass. This is an approximation, which leads to an enormous speed up in the rendering time. However, it does omit some interesting and subtle chromatic effects as different parts of a galaxy, with different effective SEDs, would be convolved by slightly different effective PSFs. There are plans to improve the implementation of this aspect of galsim, but it cannot currently simulate such effects efficiently enough for our needs. There are also plans to enable the use of webbpsf 16 in galsim.wfirst to leverage the work being done on that project to simulate the Roman PSF. The webbpsf model is qualitatively similar to what we are using from galsim.wfirst, but there are slight differences. We expect that the webbpsf model is probably more accurate, but this will be explored in future work. 3.2.3 Implemented detector effects Most of the development of galsim has been driven by the need to render simulations of CCD images. The HgCdTe detectors used by Roman are qualitatively similar, but there are significant differences in the physics, which lead to differences in some of the simulation steps. We discuss the implementation of some of these effects in detail below. For this work, each image is processed through the following stages, simulating what physically happens in the detector: (1) the Poisson background of stray light and thermal emission from the telescope is generated and a ‘sky’ background image is created that also undergoes stages 2–9, (2) the impact of reciprocity failure is added, (3) the electron counts are quantized, (4) dark current is added to the image, (5) non-linear response to flux is applied, (6) the effect of IPC is applied, (7) instrument read noise is applied, (8) electron counts are converted to ADU, and (9) the ADU value is quantized. In this work, we subtract the final background image from the SCA image, simulating a perfect background subtraction algorithm. Reciprocity failure (Biesiadzinski et al. 2011) is a non-linear relationship between the voltage response in the detector to the incident flux of photons at low light levels. The exact mechanism of this effect is unknown and hence we lack a good theoretical model. galsim uses a power law $$\begin{eqnarray*} \frac{p}{p_\mathrm{nominal}} = \left(\frac{p_\mathrm{nominal}}{f_0 t_\mathrm{exp}}\right)^{\frac{\alpha }{\log (10)}} , \end{eqnarray*}$$(3) where pnominal is the pixel response (in electrons) that would have occurred in the absence of reciprocity failure, p is the actual observed response due to reciprocity failure, f0 is the base flux rate (in electrons/sec) at which the nominal gain was calibrated, texp is the exposure time, and α is taken to be 6.5 × 10−3 for the Roman sensors. A particularly pernicious effect present in the HgCdTe detectors is known as ‘persistence’ (Smith et al. 2008; Anderson et al. 2014; McLeod & Smith 2016). In a series of images taken sequentially, some small fraction of the charge accumulated in earlier exposures apparently remains in the sensor and appears in later exposures. The effect lasts for many minutes across multiple reset cycles. Therefore, for simulating the effect, we need to keep track of the precise order and time of each observation, and the electron-level (i.e. pre-read-out) images of multiple prior exposures. The exact functional form of this effect is not very well understood, although some progress is being made in laboratory tests. The functional form for this effect was updated during the Cycle 7 updates, and galsim now uses a Fermi profile when the deposited flux is above the half-well level, and linear when below. Above the half-well level, the functional form is $$\begin{eqnarray*} n_\mathrm{persist} = \frac{A \left(n/n_0\right)^a \left(\frac{t}{1000 \mathrm\,{\mathrm{ sec}}}\right)^{-r}}{ \exp (- \frac{n-n0}{dn})+ 1} \end{eqnarray*}$$(4) where A, n0, a, r, and dn are constants estimated from laboratory measurements (and stored in the galsim.wfirst module). The persistence modelling was not available when this project started, and so is not implemented in the current simulations used in this paper. In addition to the non-linear pixel response, known as reciprocity failure, there is also a non-linearity in the conversion of accumulated charge to the measured voltage (Biesiadzinski et al. 2011; Plazas et al. 2017; Etienne et al. 2018). This is a different effect, which occurs at a different point in the simulation – namely, after the application of dark current (Beletic et al. 2008; Piquette et al. 2014; Zandian et al. 2016) and persistence. galsim treats this as a modification in the effective number of electrons: $$\begin{eqnarray*} n_\mathrm{ e}^\prime = n_\mathrm{ e} - 6 \times 10^{-7} n_\mathrm{ e}^2 \end{eqnarray*}$$(5) where ne is the actual number of electrons accumulated and |$n_\mathrm{ e}^\prime$| is the effective number to account for the voltage response non-linearity. The coefficient 6 × 10−7 is appropriate for one of the WFIRST development detectors measured in the lab (Choi & Hirata 2020). IPC (Kannawadi et al. 2016) essentially amounts to a convolution of the image by a 3 × 3 kernel in pixel coordinates. However, the timing of the convolution is during the readout process, which means that some (but not all) of the noise has already occurred. Thus it cannot be treated as part of the PSF for the purpose of the simulation. It needs to be applied separately after the dark current and Poisson shot noise have been applied, but before the read noise. The IPC coefficients have been measured in the lab for Roman detectors; the values used in the galsim.wfirst module come from the Cycle 5 estimates. 3.3 Galaxy catalogues The input galaxy catalogue is created using a simulated galaxy distribution on the sky taken from one realization of the Buzzard simulation (DeRose et al. 2019; Wechsler et al., in preparation), to introduce realistic galaxy clustering. Each galaxy is then assigned a random set of photometric properties matching a galaxy from a sample based on the Cosmic Assembly Near-infrared Deep Extragalactic Legacy Survey (CANDELS) that simulates the fiducial Roman weak lensing selection (Hemmati et al. 2019). We imposed selection cuts on the lensing source galaxies based on the Exposure Time Calculator (Hirata et al. 2012). The cuts require matched filter S/N ratio >18 in combined J + H, ellipticity error per component <0.2 (in the Bernstein & Jarvis 2002 convention), and resolution factor >0.4 (again in the Bernstein & Jarvis 2002 convention); note that this results in a limiting magnitude that depends on galaxy size. These cuts are also discussed in Hemmati et al. (2019). These selections are made on the input catalogue properties, which improves the efficiency of the simulation. This prevents us from exploring the impact of selection effects, but this is not important to this work and we can use different input galaxy property distributions in future simulation runs. The galaxy distribution, which has a mean galaxy density of approximately 40 arcmin−2, is shown in Fig. 3. In Fig. 4, we show the distributions of size, redshift, and H158 magnitude in the CANDELS sample. We discard less than 1 per cent of the largest objects in the shape measurement stage, however, due to a maximum postage stamp size restriction. In general, the input distribution and properties of galaxies can be easily modified by configuration (i.e. specifying a different input galaxy catalogue or a realistic shear field). Figure 3. Open in new tabDownload slide The distribution of simulated galaxies. The mean galaxy density is 40 arcmin−2. Figure 3. Open in new tabDownload slide The distribution of simulated galaxies. The mean galaxy density is 40 arcmin−2. Figure 4. Open in new tabDownload slide The true (blue) and recovered (orange) distributions of galaxy half-light radius, redshift, and H158 magnitude for galaxies, with a comparison of the magnitude distribution of stars (green). The recovered galaxy magnitude and half-light radius are the distributions inferred from the shape measurement process, while the distribution in redshift simply shows where in redshift objects do not have a valid shape fit – mostly at low redshift, where some large objects are not used. In general, the measured size and magnitude agree well with the true values. Star magnitudes are currently capped at 14 to avoid visual artefacts in the drawn images, which has no impact on the current or most plausible weak lensing studies. Figure 4. Open in new tabDownload slide The true (blue) and recovered (orange) distributions of galaxy half-light radius, redshift, and H158 magnitude for galaxies, with a comparison of the magnitude distribution of stars (green). The recovered galaxy magnitude and half-light radius are the distributions inferred from the shape measurement process, while the distribution in redshift simply shows where in redshift objects do not have a valid shape fit – mostly at low redshift, where some large objects are not used. In general, the measured size and magnitude agree well with the true values. Star magnitudes are currently capped at 14 to avoid visual artefacts in the drawn images, which has no impact on the current or most plausible weak lensing studies. 3.4 Star catalogue We simulate the positions and magnitudes of input stars in Roman passbands using the galaxy simulation Galaxia17 (Sharma et al. 2011). Galaxia uses an analytic model (Robin et al. 2003) to simulate stars in the galaxy that includes a thin and thick disc with warp and flaring, bulge, and halo components. Stars are simulated to 27th magnitude in V band, extinction is added, and they are uniformly translated to Roman bandpasses using the stellar SED of Alpha Lyra derived from HST CALSPEC as packaged with galsim. The star distribution, which has a mean stellar density of approximately 2.5 arcmin−2, is shown in Fig. 5. Figure 5. Open in new tabDownload slide The distribution of simulated stars. The mean stellar density is 2.5 arcmin−2. Figure 5. Open in new tabDownload slide The distribution of simulated stars. The mean stellar density is 2.5 arcmin−2. 3.5 Survey strategy We considered a reference HLS observing strategy consisting of two passes in each of the four HLS imaging filters (plus four passes for the grism). To construct each pass, we take a sequence of n exposures (2 ≤ n ≤ 4, depending on the filter/grism choice), with a small diagonal step after each exposure to cover gaps between the SCAs. These steps also ensure that in each of the n exposures, an image of a star or galaxy does not land on the same small chip defect, nor in the same readout channel of an SCA, and does not interact with its persistence image from the previous exposure.18 After these exposures, we do a ‘field step’ along the short axis of the field19 (∼0.4 degrees), and repeat the exposures. This produces a strip of observed sky; strips are tiled to cover a region on the sky. Subsequent strips are observed in opposite directions (i.e. we alternate ‘up’ versus ‘down’). The HLS is broken down into eight such regions (plus deep fields), each with its own tiling. The H158 filter exposure sequence that overlaps the patch of sky simulated for this work is shown in Fig. 6 and the total number of exposures that overlap each simulated galaxy is shown in Fig. 7. Figure 6. Open in new tabDownload slide A visualization of the individual pointings of the telescope using the H158 filter in the 5-yr simulated Reference Survey that overlap the region of the sky we are simulating images for (the non-shaded region). There are a total of 189 pointings in H158 that overlap this region. Each marker is an individual pointing, whose colour represents the focal plane position angle. Each cluster of pointings typically contains 3–4 very small translational dithers to cover chip gaps. The dither pattern in other filters overlaps in other directions to produce a more homogeneous coverage than is indicated in this figure. Figure 6. Open in new tabDownload slide A visualization of the individual pointings of the telescope using the H158 filter in the 5-yr simulated Reference Survey that overlap the region of the sky we are simulating images for (the non-shaded region). There are a total of 189 pointings in H158 that overlap this region. Each marker is an individual pointing, whose colour represents the focal plane position angle. Each cluster of pointings typically contains 3–4 very small translational dithers to cover chip gaps. The dither pattern in other filters overlaps in other directions to produce a more homogeneous coverage than is indicated in this figure. Figure 7. Open in new tabDownload slide The number of simulated exposures per object (blue) compared to the number eventually used in the shape measurement (orange). The median number of exposures per object both simulated and used is six over this part of the 5-yr simulated Reference Survey. Figure 7. Open in new tabDownload slide The number of simulated exposures per object (blue) compared to the number eventually used in the shape measurement (orange). The median number of exposures per object both simulated and used is six over this part of the 5-yr simulated Reference Survey. The two passes over each region of the HLS are on grids that are rolled relative to each other. This strategy increases the number of exposures, and more importantly ensures that astronomical sources observed on one SCA have repeated observations on other SCAs. This is needed for ‘ubercalibration’ internal to the HLS (e.g. Padmanabhan et al. 2008), and will be helpful in developing a correction in the event that a few SCAs exhibit unusual behaviours (e.g. larger than normal hysteresis). The overall survey strategy has to schedule each pass over each region, while being consistent with the needs of the other surveys and the observing constraints as well. We developed tools to do this early in Roman planning, especially since both L2 and geosynchronous orbits were under consideration, with the latter having complex Earth and Moon avoidance constraints (Spergel et al. 2015). The constraints are much more slowly varying at L2, but we still have the slowly varying Sun avoidance constraint (Roman observes between 54–126° degrees from the Sun), a roll angle constraint (the observatory can roll up to ±15° from the optimal orientation on the solar array; this is very important when attempting to tile a large region of the sky). Moreover, there are cutouts for the microlensing seasons and – during the middle of the reference mission – a 30-h supernova observing session every 5 d. This results in the need to cut each pass into shorter segments that can be observed all at once, in sequence. The strategy described here is an output from an update of the code used in section 3.10 of Spergel et al. (2015). 3.6 Simulation implementation for this study In this work, we study the impact of how a variety of biases in the PSF model propagate to shape measurement and the weak lensing signal. To study this, we produce a set of 13 image simulations that are identical, including noise, modulo a single PSF model change relative to the fiducial simulation in each case. The details of these changes and their impacts are described in more detail in Sections 4 and 5. Shape measurement is then performed on the images with some PSF model bias, but using the fiducial PSF model for convolution in the galaxy shape fitter, to simulate an unknown wavefront error. Several simplifications are employed relative to the generic synthetic survey generation described in Section 3.1 to accommodate the computational load of the many realizations of the survey we are producing. We simulate objects in a 2.5 × 2.5 deg2 patch of the sky. We only simulate pointings targeted for the H158 filter. Since we are not simulating chromatic effects, the specific filter choice does not make a large difference in our results. We use a lower-resolution version of the PSF, which significantly speeds up the convolution. The impact of this approximation on the PSF model, in both native and oversampled pixels, can be seen in Fig. 2, but is not important for this work. To better isolate the effects of PSF errors, we only utilize the isolated object postage stamps in shape measurement. We do not simulate objects with photometry that would fall outside the fiducial weak lensing selection criteria. We do not implement a shear calibration scheme like metacalibration (Sheldon & Huff 2017), since we only care about changes to the recovered shape between simulation runs. Work on applying a method like metacalibration to these simulations is ongoing. We simulate a total of 907 170 unique galaxies and 56 128 unique stars across 189 pointings in each of the runs. The number of exposures per galaxy and the distribution of PSF properties are shown in Figs 7 and 8. Figure 8. Open in new tabDownload slide The measured PSF ellipticity and size, and the signal-to-noise of the galaxy measurements. The nominal signal-to-noise cut for the Roman weak lensing sample is 18, which agrees with the peak of the measured distribution. The PSF properties are measured using the oversampled version of the PSF model images (as in Fig. 2), which we determined to be the minimum required oversampling in the model to achieve unbiased size and ellipticity measurements. The typical Roman PSF tends to be more strongly elliptical in the e2 direction, and the median (and mean) recovered FWHM is 0.178 arcsec. Figure 8. Open in new tabDownload slide The measured PSF ellipticity and size, and the signal-to-noise of the galaxy measurements. The nominal signal-to-noise cut for the Roman weak lensing sample is 18, which agrees with the peak of the measured distribution. The PSF properties are measured using the oversampled version of the PSF model images (as in Fig. 2), which we determined to be the minimum required oversampling in the model to achieve unbiased size and ellipticity measurements. The typical Roman PSF tends to be more strongly elliptical in the e2 direction, and the median (and mean) recovered FWHM is 0.178 arcsec. 4 WAVEFRONT MODEL ERRORS In this paper, we focus on empirical tests of weak lensing requirements for wavefront model control (i.e. the PSF) in Roman. These are used to empirically derive the relationship between recovered shear and wavefront error modes ∂ei/∂ψj, which allows us to validate earlier Phase A analytic estimates of the requirements flowdown for Roman. In the absence of shear biases, or when comparing between runs that should have identical intrinsic shear biases, ∂ei/∂ψj is equivalent to ∂γobs,i/∂ψj. We simulate 13 identical 2.5 × 2.5 deg2 Reference Survey cutouts: a fiducial survey that represents perfect knowledge of the PSF and 12 iterations to simulate various types of errors in the PSF reconstruction. These are split into three types of errors in the wavefront model: (1) static biases in the model, which are constant as a function of time, (2) high-frequency biases in the model, which correspond to rapidly changing conditions compared to the time-scale of a single exposure, and (3) low-frequency biases in the model, which change over the lifetime of the mission, but can be considered static over the time-scale of a single exposure. In each static and low-frequency mode, the (rms) amplitude of the wavefront bias corresponds to 0.005 wavelengths (a fiducial wavelength is taken to be 1293 nm), which is equivalent to approximately 6.5 nm. These PSF changes are summarized in Table 1. Table 1. A summary of the 13 simulation runs. Run name . PSF change . Mode . Notes . fiducial – – – focus ψ4 Static – astig ψ5 Static – coma ψ7 Static – gradz4 ψ4 Static Gradient in focal plane gradz6 ψ6 Static Gradient in focal plane piston ψ4 Static Random per SCA tilt ψ4 Static Rand. gradient per SCA isojitter Gaussian High-Freq. Isotropic anijitter Gaussian High-Freq. Anisotropic ranjitter Gaussian High-Freq. 15 per cent of pointings oscz4 ψ4 Low-Freq. Time-dependent oscz7 ψ7 Low-Freq. Time-dependent Run name . PSF change . Mode . Notes . fiducial – – – focus ψ4 Static – astig ψ5 Static – coma ψ7 Static – gradz4 ψ4 Static Gradient in focal plane gradz6 ψ6 Static Gradient in focal plane piston ψ4 Static Random per SCA tilt ψ4 Static Rand. gradient per SCA isojitter Gaussian High-Freq. Isotropic anijitter Gaussian High-Freq. Anisotropic ranjitter Gaussian High-Freq. 15 per cent of pointings oscz4 ψ4 Low-Freq. Time-dependent oscz7 ψ7 Low-Freq. Time-dependent Open in new tab Table 1. A summary of the 13 simulation runs. Run name . PSF change . Mode . Notes . fiducial – – – focus ψ4 Static – astig ψ5 Static – coma ψ7 Static – gradz4 ψ4 Static Gradient in focal plane gradz6 ψ6 Static Gradient in focal plane piston ψ4 Static Random per SCA tilt ψ4 Static Rand. gradient per SCA isojitter Gaussian High-Freq. Isotropic anijitter Gaussian High-Freq. Anisotropic ranjitter Gaussian High-Freq. 15 per cent of pointings oscz4 ψ4 Low-Freq. Time-dependent oscz7 ψ7 Low-Freq. Time-dependent Run name . PSF change . Mode . Notes . fiducial – – – focus ψ4 Static – astig ψ5 Static – coma ψ7 Static – gradz4 ψ4 Static Gradient in focal plane gradz6 ψ6 Static Gradient in focal plane piston ψ4 Static Random per SCA tilt ψ4 Static Rand. gradient per SCA isojitter Gaussian High-Freq. Isotropic anijitter Gaussian High-Freq. Anisotropic ranjitter Gaussian High-Freq. 15 per cent of pointings oscz4 ψ4 Low-Freq. Time-dependent oscz7 ψ7 Low-Freq. Time-dependent Open in new tab We emphasize that the purpose of these simulations is to measure the sensitivity (i.e. partial derivatives) of the shear biases mi and ci with respect to the PSF parameters. We want to do this with an area that is much less than the Reference Survey area; we choose a change in wavefront Δψ that is considerably greater than the expected requirement so that the partial derivative is not swamped by noise. Similar considerations apply to the pointing jitter. 4.1 Static biases We simulate seven static sources of bias in the PSF model. Three of these simulations include a coherent change in the PSF model Zernike coefficients, where the fiducial value is changed by 0.005 wavelengths in each of defocus ψ4 (focus), oblique astigmatism ψ5 (astig), and vertical coma ψ7 (coma). Two simulations include a coherent gradient in the defocus ψ4 (gradz4) and vertical astigmatism ψ6 (gradz6) across the focal plane with equivalent rms of 0.005 wavelengths. For speed, these are simulated such that the PSF is constant within a single SCA. Finally, two simulations approximate errors in the mounting of the SCAs: (1) a random vertical mounting offset of up to 0.005 wavelengths is assigned to each SCA (piston), and (2) a random tilt in the x or y direction is assigned to each SCA (tilt), with equivalent rms of up to 0.005 wavelengths. These are modelled as changes in the ψ4 coefficient, with the PSF being evaluated based on the object x–y position within the SCA (i.e. each object is assigned a different PSF consistent with this random tilt of the SCA). Potential correlated biases in the WCS model due to these changes are ignored in this work, but should be considered in future studies of the WCS model recovery. 4.2 High-frequency biases Three high-frequency resonant modes are simulated to represent residual vibrations of the telescope after orienting to a new pointing. These are represented by an additional convolution of the image with a Gaussian PSF. We simulate three cases: (1) an isotropic (about the pointing axis) vibration (isojitter), (2) an anisotropic vibration (anijitter), and (3) only applying this anisotropic vibration to a random 15 per cent of pointings (ranjitter). The additional second moments are conserved, which means |$\theta _x\theta _y=\theta _{\mathrm{orig}}^2=15^2$| mas2. For anijitter and ranjitter, we applied a shear |$e_1=\frac{\theta _x-\theta _y}{\theta _x+\theta _y}=0.3$|⁠, with Zernike amplitude change in this case |$\mathrm{ d}\psi =\theta _x^2-\theta _y^2=297$| mas2. In the case of isojitter, θx = θy, which leads to |$\mathrm{ d}\psi =\theta _x^2+\theta _y^2=450$| mas2. 4.3 Low-frequency biases Two low-frequency biases are simulated to represent thermal drift throughout the lifetime of the mission. Thermal perturbations propagate into ψ4 (oscz4) and ψ7 (oscz7) modes. We generate a random time-dependent function f(t) with rms amplitude 0.005 wavelengths following a given power spectrum to quantify the perturbation of the Zernike coefficients over time. The power spectrum of thermal drift noise is taken to be a Lorenzian function $$\begin{eqnarray*} P(\nu)=\frac{A}{1+(\nu /\nu _0)^2}, \end{eqnarray*}$$(6) with normalization factor A. The rms variance of f(t) can be expressed as $$\begin{eqnarray*} \sigma ^2 = \int _{-\infty }^{+\infty }P(\nu)\mathrm{ d}\nu =\pi A\nu _0, \end{eqnarray*}$$(7) which leads to |$A={\sigma ^2}/{\pi \nu _0}$|⁠. |$\nu _0=\frac{1}{2\pi \tau }\approx 3.14\times 10^{-4}$| Hz, with a time constant τ = 1 h. This time-scale is typical of thermal variations that have been seen in integrated modelling (e.g. Spergel et al. 2015); we plan to use actual integrated modelling outputs for the reference observing scenario in a future version of this study. 5 RESULTS Each simulation is analysed in an identical way, except that shape measurement for each simulation assumes the fiducial PSF model is the true model, which simulates the impact of mis-estimating the PSF model and introduces varying levels of bias. All estimates of the multiplicative and additive bias will be explored relative to the fiducial simulation run, since we have not employed an absolute calibration scheme. This is justified to first order, since we are only interested in the relative impacts of the PSF model biases. We find that 3 per cent of objects are not included in the shape measurement stage in the fiducial simulation, due to being too large/bright or because too large a fraction of all cutouts are masked (fall off the edge of an SCA) – see Section 3.1 for more details on these selections. Since the simulated objects have already been pre-selected as objects that should pass the fiducial Roman weak lensing selection, we are able to successfully recover a shape fit for more than 99 per cent of the remaining objects – a total of 871 841 galaxies. We do not make an additional selection on objects that would pass the fiducial Roman shape selection based on measured properties, since we expect all objects to be within this selection if we were to simulate all remaining pointings in other bandpasses. The recovered multiplicative shear bias is only approximately 2 per cent smaller and the mean shear is unchanged if we make this selection, which removes an additional 35 per cent of objects, almost exclusively due to the signal-to-noise cut. We present results for the non-fiducial simulations only for objects that lie in the intersection of successful shape measurement between each simulation and the fiducial simulation, to allow for 1–1 comparison of the shapes and cancellation of shape noise and sources of photon noise, which are identical in each simulation. We neglect the impact of selection biases here, since the intersection criteria exclude on average only 0.3 per cent of objects. 5.1 Summary statistics The bias in an ensemble shear measurement is typically characterized in the weak limit by $$\begin{eqnarray*} e_i^{\mathrm{obs}} = (1+m_i) e_i^{\mathrm{true}} + c_i. \end{eqnarray*}$$(8) We find the following multiplicative and additive biases in the fiducial simulation: $$\begin{eqnarray*} m_1 &=& (-7.56 \pm 0.19)\times 10^{-2}\\ m_2 &=& (-9.40 \pm 0.19)\times 10^{-2}\\ c_1 &=& (1.20 \pm 0.17)\times 10^{-3}\\ c_2 &=& (-1.57 \pm 0.16)\times 10^{-3}. \end{eqnarray*}$$ In some cases, biases are instead parametrized in terms of the PSF leakage as |$e_i^{\mathrm{obs}} = (1+m_i) e_i^{\mathrm{true}} + \alpha _i e^{\mathrm{PSF}}_i + c_i$|⁠. The constraints on m used to interpret requirements in this paper are unchanged in either parametrization. The difference in measured shape versus true input shape (intrinsic shape and shear) is shown in Fig. 9. Figure 9. Open in new tabDownload slide The difference in the measured (obs) shape relative to the true sheared intrinsic shape of the galaxies in the fiducial simulation. The inferred multiplicative and additive shear bias is discussed in Section 5. Figure 9. Open in new tabDownload slide The difference in the measured (obs) shape relative to the true sheared intrinsic shape of the galaxies in the fiducial simulation. The inferred multiplicative and additive shear bias is discussed in Section 5. For each simulation, we compare the recovered shear to the fiducial simulation in several ways. First, we calculate how the inferred values of m and c change from the fiducial result, which is shown in Table 2. More importantly, we are interested in how the inferred shear changes as a function of the induced wavefront error. This allows us to draw a direct connection to the analytic requirements predictions. We show this shear response relative to the wavefront error in Table 3. Finally, we are ultimately interested in how these biases will propagate to the shear correlation function – that is, how any coherent scale dependence of the effects will impact cosmology. Table 2. Additive and multiplicative bias parameter changes in each version of the simulation relative to the fiducial simulation. Note that the wavefront errors injected into these simulations (6.5 nm rms) are much greater than the stability requirements; the variation simulations are designed to measure the partial derivatives of the shear biases with respect to each contribution to the PSF. Run name . △m1 × 103 . △m2 × 103 . △c1 × 104 . △c2 × 104 . focus 4.65 ± 0.24 4.36 ± 0.23 −5.57 ± 0.24 −15.87 ± 0.30 astig −0.23 ± 0.24 0.90 ± 0.21 −29.14 ± 0.30 −0.55 ± 0.20 coma −1.33 ± 0.26 −1.47 ± 0.21 −6.66 ± 0.26 0.22 ± 0.34 gradz4 −2.63 ± 0.27 −2.52 ± 0.23 −0.51 ± 0.21 6.71 ± 0.42 gradz6 0.22 ± 0.23 −0.22 ± 0.23 14.72 ± 0.62 0.44 ± 0.23 piston −13.6 ± 5.0 −15.0 ± 5.3 −3.70 ± 0.31 −7.51 ± 0.34 tilt −38.0 ± 8.0 −36.6 ± 8.0 −0.40 ± 0.38 0.73 ± 0.38 isojitter −40.8 ± 8.6 −39.5 ± 8.5 0.24 ± 0.84 0.43 ± 0.91 anijitter −11.4 ± 1.1 −11.5 ± 1.0 43.27 ± 0.85 −0.13 ± 0.86 ranjitter −12.8 ± 4.4 −12.6 ± 4.1 7.50 ± 0.64 0.33 ± 0.43 oscz4 −7.2 ± 3.1 −6.5 ± 3.1 0.89 ± 0.31 2.15 ± 0.49 oscz7 −12.1 ± 6.8 −11.4 ± 6.3 −0.08 ± 0.31 −0.17 ± 0.44 Run name . △m1 × 103 . △m2 × 103 . △c1 × 104 . △c2 × 104 . focus 4.65 ± 0.24 4.36 ± 0.23 −5.57 ± 0.24 −15.87 ± 0.30 astig −0.23 ± 0.24 0.90 ± 0.21 −29.14 ± 0.30 −0.55 ± 0.20 coma −1.33 ± 0.26 −1.47 ± 0.21 −6.66 ± 0.26 0.22 ± 0.34 gradz4 −2.63 ± 0.27 −2.52 ± 0.23 −0.51 ± 0.21 6.71 ± 0.42 gradz6 0.22 ± 0.23 −0.22 ± 0.23 14.72 ± 0.62 0.44 ± 0.23 piston −13.6 ± 5.0 −15.0 ± 5.3 −3.70 ± 0.31 −7.51 ± 0.34 tilt −38.0 ± 8.0 −36.6 ± 8.0 −0.40 ± 0.38 0.73 ± 0.38 isojitter −40.8 ± 8.6 −39.5 ± 8.5 0.24 ± 0.84 0.43 ± 0.91 anijitter −11.4 ± 1.1 −11.5 ± 1.0 43.27 ± 0.85 −0.13 ± 0.86 ranjitter −12.8 ± 4.4 −12.6 ± 4.1 7.50 ± 0.64 0.33 ± 0.43 oscz4 −7.2 ± 3.1 −6.5 ± 3.1 0.89 ± 0.31 2.15 ± 0.49 oscz7 −12.1 ± 6.8 −11.4 ± 6.3 −0.08 ± 0.31 −0.17 ± 0.44 Open in new tab Table 2. Additive and multiplicative bias parameter changes in each version of the simulation relative to the fiducial simulation. Note that the wavefront errors injected into these simulations (6.5 nm rms) are much greater than the stability requirements; the variation simulations are designed to measure the partial derivatives of the shear biases with respect to each contribution to the PSF. Run name . △m1 × 103 . △m2 × 103 . △c1 × 104 . △c2 × 104 . focus 4.65 ± 0.24 4.36 ± 0.23 −5.57 ± 0.24 −15.87 ± 0.30 astig −0.23 ± 0.24 0.90 ± 0.21 −29.14 ± 0.30 −0.55 ± 0.20 coma −1.33 ± 0.26 −1.47 ± 0.21 −6.66 ± 0.26 0.22 ± 0.34 gradz4 −2.63 ± 0.27 −2.52 ± 0.23 −0.51 ± 0.21 6.71 ± 0.42 gradz6 0.22 ± 0.23 −0.22 ± 0.23 14.72 ± 0.62 0.44 ± 0.23 piston −13.6 ± 5.0 −15.0 ± 5.3 −3.70 ± 0.31 −7.51 ± 0.34 tilt −38.0 ± 8.0 −36.6 ± 8.0 −0.40 ± 0.38 0.73 ± 0.38 isojitter −40.8 ± 8.6 −39.5 ± 8.5 0.24 ± 0.84 0.43 ± 0.91 anijitter −11.4 ± 1.1 −11.5 ± 1.0 43.27 ± 0.85 −0.13 ± 0.86 ranjitter −12.8 ± 4.4 −12.6 ± 4.1 7.50 ± 0.64 0.33 ± 0.43 oscz4 −7.2 ± 3.1 −6.5 ± 3.1 0.89 ± 0.31 2.15 ± 0.49 oscz7 −12.1 ± 6.8 −11.4 ± 6.3 −0.08 ± 0.31 −0.17 ± 0.44 Run name . △m1 × 103 . △m2 × 103 . △c1 × 104 . △c2 × 104 . focus 4.65 ± 0.24 4.36 ± 0.23 −5.57 ± 0.24 −15.87 ± 0.30 astig −0.23 ± 0.24 0.90 ± 0.21 −29.14 ± 0.30 −0.55 ± 0.20 coma −1.33 ± 0.26 −1.47 ± 0.21 −6.66 ± 0.26 0.22 ± 0.34 gradz4 −2.63 ± 0.27 −2.52 ± 0.23 −0.51 ± 0.21 6.71 ± 0.42 gradz6 0.22 ± 0.23 −0.22 ± 0.23 14.72 ± 0.62 0.44 ± 0.23 piston −13.6 ± 5.0 −15.0 ± 5.3 −3.70 ± 0.31 −7.51 ± 0.34 tilt −38.0 ± 8.0 −36.6 ± 8.0 −0.40 ± 0.38 0.73 ± 0.38 isojitter −40.8 ± 8.6 −39.5 ± 8.5 0.24 ± 0.84 0.43 ± 0.91 anijitter −11.4 ± 1.1 −11.5 ± 1.0 43.27 ± 0.85 −0.13 ± 0.86 ranjitter −12.8 ± 4.4 −12.6 ± 4.1 7.50 ± 0.64 0.33 ± 0.43 oscz4 −7.2 ± 3.1 −6.5 ± 3.1 0.89 ± 0.31 2.15 ± 0.49 oscz7 −12.1 ± 6.8 −11.4 ± 6.3 −0.08 ± 0.31 −0.17 ± 0.44 Open in new tab Table 3. The changes of ellipticity with respect to changes in line-of-sight motion for jitter cases and wavefront error for the other modes. Run name . ∂e1/∂ψ × 10−4 . ∂e2/∂ψ × 10−4 . Units . focus −0.87 ± 0.032 −2.5 ± 0.031 nm−1 astig −4.5 ± 0.031 −0.10 ± 0.030 nm−1 coma −1.0 ± 0.032 0.030 ± 0.032 nm−1 gradz4 −0.089 ± 0.030 1.0 ± 0.029 nm−1 gradz6 2.2 ± 0.029 0.044 ± 0.028 nm−1 piston −0.560 ± 0.048 −1.2 ± 0.048 nm−1 tilt −0.036 ± 0.063 0.10 ± 0.063 nm−1 isojitter 0.0004 ± 0.0024 0.0006 ± 0.0023 mas−2 anijitter 0.145 ± 0.003 −0.001 ± 0.003 mas−2 ranjitter 0.0246 ± 0.0016 0.00087 ± 0.00159 mas−2 oscz4 0.13 ± 0.039 0.30 ± 0.039 nm−1 oscz7 −0.019 ± 0.045 −0.044 ± 0.044 nm−1 Run name . ∂e1/∂ψ × 10−4 . ∂e2/∂ψ × 10−4 . Units . focus −0.87 ± 0.032 −2.5 ± 0.031 nm−1 astig −4.5 ± 0.031 −0.10 ± 0.030 nm−1 coma −1.0 ± 0.032 0.030 ± 0.032 nm−1 gradz4 −0.089 ± 0.030 1.0 ± 0.029 nm−1 gradz6 2.2 ± 0.029 0.044 ± 0.028 nm−1 piston −0.560 ± 0.048 −1.2 ± 0.048 nm−1 tilt −0.036 ± 0.063 0.10 ± 0.063 nm−1 isojitter 0.0004 ± 0.0024 0.0006 ± 0.0023 mas−2 anijitter 0.145 ± 0.003 −0.001 ± 0.003 mas−2 ranjitter 0.0246 ± 0.0016 0.00087 ± 0.00159 mas−2 oscz4 0.13 ± 0.039 0.30 ± 0.039 nm−1 oscz7 −0.019 ± 0.045 −0.044 ± 0.044 nm−1 Open in new tab Table 3. The changes of ellipticity with respect to changes in line-of-sight motion for jitter cases and wavefront error for the other modes. Run name . ∂e1/∂ψ × 10−4 . ∂e2/∂ψ × 10−4 . Units . focus −0.87 ± 0.032 −2.5 ± 0.031 nm−1 astig −4.5 ± 0.031 −0.10 ± 0.030 nm−1 coma −1.0 ± 0.032 0.030 ± 0.032 nm−1 gradz4 −0.089 ± 0.030 1.0 ± 0.029 nm−1 gradz6 2.2 ± 0.029 0.044 ± 0.028 nm−1 piston −0.560 ± 0.048 −1.2 ± 0.048 nm−1 tilt −0.036 ± 0.063 0.10 ± 0.063 nm−1 isojitter 0.0004 ± 0.0024 0.0006 ± 0.0023 mas−2 anijitter 0.145 ± 0.003 −0.001 ± 0.003 mas−2 ranjitter 0.0246 ± 0.0016 0.00087 ± 0.00159 mas−2 oscz4 0.13 ± 0.039 0.30 ± 0.039 nm−1 oscz7 −0.019 ± 0.045 −0.044 ± 0.044 nm−1 Run name . ∂e1/∂ψ × 10−4 . ∂e2/∂ψ × 10−4 . Units . focus −0.87 ± 0.032 −2.5 ± 0.031 nm−1 astig −4.5 ± 0.031 −0.10 ± 0.030 nm−1 coma −1.0 ± 0.032 0.030 ± 0.032 nm−1 gradz4 −0.089 ± 0.030 1.0 ± 0.029 nm−1 gradz6 2.2 ± 0.029 0.044 ± 0.028 nm−1 piston −0.560 ± 0.048 −1.2 ± 0.048 nm−1 tilt −0.036 ± 0.063 0.10 ± 0.063 nm−1 isojitter 0.0004 ± 0.0024 0.0006 ± 0.0023 mas−2 anijitter 0.145 ± 0.003 −0.001 ± 0.003 mas−2 ranjitter 0.0246 ± 0.0016 0.00087 ± 0.00159 mas−2 oscz4 0.13 ± 0.039 0.30 ± 0.039 nm−1 oscz7 −0.019 ± 0.045 −0.044 ± 0.044 nm−1 Open in new tab We can study the difference in the measured ellipticity relative to the fiducial simulation in both focal plane and sky coordinates. The mean ellipticity difference binned in the focal plane is shown in Figs 10 and 11, for e1 and e2, respectively. We observe coherent, and sometimes large, biases in the mean ellipticity across the focal plane or individual chips for all static wavefront errors. The time-dependent wavefront errors are generally less pronounced, except in the case of a non-random anisotropic jitter that produces a very strong, coherent bias in e1, the direction of the anisotropy in the smearing. Figure 10. Open in new tabDownload slide The binned mean difference of measured e1 compared to the fiducial run in the focal plane. Each colour bar is in units of 1 × 10−4. Gradients across the focal plane or chips are visible for all static PSF model biases, while the mean difference is particularly large for the anisotropic jitter case, where the large mean e1 difference corresponds to the direction of anisotropy in the Gaussian smearing. Figure 10. Open in new tabDownload slide The binned mean difference of measured e1 compared to the fiducial run in the focal plane. Each colour bar is in units of 1 × 10−4. Gradients across the focal plane or chips are visible for all static PSF model biases, while the mean difference is particularly large for the anisotropic jitter case, where the large mean e1 difference corresponds to the direction of anisotropy in the Gaussian smearing. Figure 11. Open in new tabDownload slide The binned mean difference of measured e2 compared to the fiducial run in the focal plane. Each colour bar is in units of 1 × 10−4. There are visible gradients for most cases of static PSF model biases. Figure 11. Open in new tabDownload slide The binned mean difference of measured e2 compared to the fiducial run in the focal plane. Each colour bar is in units of 1 × 10−4. There are visible gradients for most cases of static PSF model biases. Fig. 12 shows the two-point correlation function ξ+ of the ellipticity difference in sky coordinates, where $$\begin{eqnarray*} \xi _{\pm } = \langle \Delta e_{+}\Delta e_{+} \rangle \pm \langle \Delta e_{\times }\Delta e_{\times } \rangle . \end{eqnarray*}$$(9) Δe+ and Δe× are the tangential and cross components, respectively, of the ellipticity difference relative to the fiducial simulation along the projected separation vector between each pair of galaxies on the sky. As expected for a wavefront error, the ξ− correlation of the differences are all consistent with zero. Like the mean shear in the focal plane, all static wavefront error cases lead to significantly non-zero ξ+ at varying magnitudes. The time-dependent errors have ξ+ consistent with zero, except for the non-random anisotropic jitter case, which shows the largest impact in ξ+ of any case as additional smear only applies to e1 component. The non-random anisotropic jitter, and the static focus and astigmatism errors, all produce a nearly constant ξ+ correlation with angular scale, showing that the results are dominated by uniformly distributed e1, e2. Figure 12. Open in new tabDownload slide The correlation function ξ+ of the ellipticity difference between each run, indicated by the label for each panel, and the fiducial run. Negative points are shown as crosses. The ξ− values are all consistent with zero, and are not shown. Figure 12. Open in new tabDownload slide The correlation function ξ+ of the ellipticity difference between each run, indicated by the label for each panel, and the fiducial run. Negative points are shown as crosses. The ξ− values are all consistent with zero, and are not shown. 5.2 Comparison of analytic to numerical results The partial derivatives |∂e/∂ψ| of the ellipticity with respect to wavefront should be less than |$\Vert {\boldsymbol{\Lambda }} \Vert \Vert {\boldsymbol{\psi }} \Vert$| (where |${\boldsymbol{\Lambda}}$| is the analytically derived matrix defined in equation (B11) of Appendix B), which is |$8.98\times 10^{-4}\,$| nm−1 for the H158-band (this is an RSS of the two ellipticity components). For the partial derivatives of the ellipticity with respect to the second moments of the jitter pattern, |∂e/∂ψ| should be less than Kθθ (where Kθθ is the analytically derived sensitivity to jitter; see equations B19 andB20 of Appendix B), which is 1.38 × 10−5 mas−2 for the H158-band. In the numerical results presented in Table 3, the largest partial derivatives are |$4.5\times 10^{-4}\,$| nm−1 (for wavefront errors) and 1.45 × 10−5 mas−2 (for jitter). For the wavefront drift case, this is consistent with the analytic expectations. For the anisotropic jitter case, the sensitivity |∂e/∂ψ| determined from the simulations is 5 ± 2 per cent larger than the analytic bound (we would expect a sensitivity equal to the analytic bound since the anijitter run is the worst-case jitter pattern: the anisotropy is in the e1 component for every exposure). The 5 ± 2 per cent difference may simply represent the approximations made in the analytic calculation (e.g. no treatment of undersampling/image combination, a different shape measurement algorithm, etc.). A similar comparison is possible for the two-point correlation functions ξ±(θ) of the ellipticity changes Δe. Since we put in a wavefront change of 6.5 nm rms, the analytic prediction is that these correlation functions should be $$\begin{eqnarray*} \xi _\pm (\theta) \le (8.98\times 10^{-4}\, {\rm nm}^{-1} \times 6.5\, {\rm nm})^2 = 3.4\times 10^{-5}. \end{eqnarray*}$$(10) As seen in Fig. 12, this inequality is indeed satisfied. A similar result can be written for the jitter cases. The most stressing case is the anijitter case, which has |$\langle \theta _x^2 - \theta _y^2\rangle = 297\,$| mas2 and hence should satisfy $$\begin{eqnarray*} \xi _\pm (\theta) \le (1.38\times 10^{-5}\, {\rm mas}^{-2}\times 297\, {\rm mas}^2)^2 = 1.68\times 10^{-5}. \end{eqnarray*}$$(11) The anijitter panel in Fig. 12 shows a numerical result that is very close to this. The central values of ξ+(θ) range from (1.79–1.95) × 10−5, which are slightly larger than the analytical estimate, although the 1σ error bars include 1.68 × 10−5. If this is not a statistical fluctuation, it is likely due to the same simplifying approximations in the analytic calculation as described above for |∂e/∂ψ|. In either case, the numerical calculation gives a result that is near the analytically estimated upper bound. 6 FUTURE DEVELOPMENT PLANS The simulation framework described here is a substantial step forward in the development of a significant synthetic Roman imaging survey, which incorporates a realistic set of photometric properties and distributions of galaxies and stars, complex morphological properties for galaxies, and most known detector non-idealities present in HgCdTe H4RG detectors. There are many advances that are still necessary, however, many of which are currently in progress. These include increased simulation volume for more precise end-to-end tests, increased fidelity in the simulation of the data accumulation processes within the detector simulation, and more realistic input galaxy and star information. Together, they will enable more precise and advanced tests of algorithm development and pipeline integration testing for the Roman weak lensing program. In terms of advances in detector physics simulation within galsim.wfirst, galsim.wfirst now includes a model for the persistence effect in the detectors based on measurements from preliminary engineering detectors that will be implemented in future versions of the survey simulation. The current simulations do not include an implementation for the brighter-fatter effect. galsim has an implementation of this for silicon CCDs, but not for the HgCdTe detectors used by Roman. It is possible that using the CCD implementation would be sufficiently accurate for future simulation runs, but this needs to be further investigated. We now also have the engineering data to implement measured realizations of the correlated noise fields derived from detector flats and darks. One significant difference relative to how data will be taken with the WFIRST SCAs is the lack of ‘up the ramp’ information as charge accumulates within the SCA. In practice, we will have access to several linear combinations of intermediate read-outs from the SCAs, which is not currently implemented. Other plans for galsim.wfirst are also discussed in Section 3. These improvements will enable us to use galsim to update our knowledge requirements for detector effects (beyond the analytic estimates used during the formulation phase of the mission). On the mock galaxy catalogue side, we have produced test runs where we interface our weak lensing survey simulation pipeline for WFIRST with the existing LSST DESC Data Challenge 2 (DC2) mock galaxy catalogue (Korytov et al. 2019), CosmoDC2, to produce WFIRST imaging over the same simulated universe as is currently being used for DESC image simulations (LSST Dark Energy Survey Collaboration 2020). The result of this work will be described in a future paper. CosmoDC2 provides a deeper mock catalogue than is currently being used, with synthetic spectra provided for each object to enable fully chromatic studies of weak lensing shape recovery. With planned improvements to the recovery of the near-infrared colours for objects in CosmoDC2, this will also enable a powerful joint-simulation with matched imaging as expected for both LSST and Roman. These matched simulations will enable a range of joint-processing tests at the pixel level to test combinations of ground-based imaging from LSST with space-based imaging from Roman. While the current synthetic survey volume used in this paper is relatively small, due to the necessity of simulating it many times, we plan the production of much larger public simulations in the near future. This will include many of the improvements described above, including multiband imaging across tens of square degrees at full Roman 5-yr Reference Survey depth matched to LSST imaging from DESC. 7 CONCLUSION The Roman observatory will be an exquisite tool for the study of cosmology using weak gravitational lensing. Launching in the mid-2020s, it can harness a unique combination of agility in potential survey design, coupled with a unique range of capabilities and power, to clarify new discoveries and resolve disagreements between the Stage IV surveys that precede it in the early 2020s. To ensure we are able to take full advantage of the potential of Roman, we must develop the necessary tools to both validate instrument requirements and their flowdown from weak lensing cosmology and to enable pixel-level algorithm development and ultimate integration testing of our measurement pipelines. In this paper, we have described a simulation framework to produce a realistic, synthetic Roman imaging survey populated with suitably complex objects that can serve these functions at the current level of necessary realism. This framework combines a simulated 5-yr Reference Survey, an appropriate mock galaxy and star population that would be observed by Roman, and a simulation of most relevant properties of the HgCdTe H4RG detectors to be integrated into the Roman camera. We present a set of 13 matched 2.5 × 2.5 deg2 image simulations to full depth of the reference 5-yr survey, each with the wavefront model perturbed in some way. These perturbations can be classified into three broad categories: static, high-, and low-frequency. We study the galaxy shape recovery in these simulations to empirically measure the relative bias in weak gravitational lensing shear estimates due to these errors in wavefront reconstruction, in order to compare to what is anticipated from the analytical requirements flowdown that was previously developed for Roman. We present quantitative comparisons of the change in the recovered ellipticity due to these various errors in the wavefront model relative to the fiducial simulation. These are presented in terms of both mean shear as a function of focal plane position and the correlation function ξ± of the ellipticity difference as a function of angular separation on the sky. Finally, we derive the response of the change in ellipticity relative to the wavefront model error mode, which we use to evaluate differences relative to previous analytical requirements forecasts. We find general agreement with the analytic requirements flowdown, though note that the empirical measurement of the bias induced in the non-random anisotropic jitter case is typically larger than predicted by the analytic flowdown. We do not consider this to be a significant concern for continued reference to the baseline, analytic requirements flowdown used by the mission, as these differences are at the 1–2σ level, depending on the type of comparison, and thus generally consistent with the analytically predicted upper bound of the effect. We have outlined in Section 6 several future expansions to the validation framework described in this paper for the Roman weak lensing analysis. These include updates to methodology, the incorporation of new flight-candidate detector measurements, and improvements in the fidelity of the image simulations to represent the full range of both properties of objects that will be observed by Roman and the full range of non-idealities in the detector systems. As the Roman mission approaches its construction phase, we expect these simulations to also begin to play a substantial role as the basis for integration tests of measurement pipeline development over the next several years. ACKNOWLEDGEMENTS We thank Dave Content, Jeff Kruk, Alice Liu, Hui Kong, and Erin Sheldon for many useful conversations, and the anonymous referee for insightful suggestions. This work supported by NASA Grant 15-WFIRST15-0008 as part of the Roman Cosmology with the High-Latitude Survey Science Investigation Team (https://www.wfirst-hls-cosmology.org). It used resources on the CCAPP condo of the Ruby Cluster at the Ohio Supercomputing Center (OSC 1987). This research was also done using resources provided by the Open Science Grid (Pordes et al. 2007; Sfiligoi et al. 2009), which is supported by the National Science Foundation and the U.S. Department of Energy’s Office of Science. Plots in this manuscript were produced partly with matplotlib (Hunter 2007), and it has been prepared using NASA’s Astrophysics Data System Bibliographic Services. DATA AVAILABILITY The data underlying this article and further supporting documentation will be shared upon request to the corresponding author. Details are provided in Appendix C. Footnotes 1 Roman was formerly named the Wide-Field Infrared Survey Telescope (WFIRST). 2 http://roman.gsfc.nasa.gov 3 http://www.darkenergysurvey.org/ 4 http://hsc.mtk.nao.ac.jp/ssp/ 5 http://kids.strw.leidenuniv.nl 6 https://www.desi.lbl.gov/ 7 http://www.lsst.org 8 http://sci.esa.int/euclid 9 Document reference number WFIRST-SYS-REQ-0020. 10 These can be found at https://wfirst.gsfc.nasa.gov/science/WFIRST_Reference_Information.html. Note that updates to match Phase B payload design have not been incorporated into the simulation described in this paper, but this is not expected to impact the results of this paper. 11 https://github.com/matroxel/wfirst_imsim 12 https://github.com/esheldon/meds 13 https://healpix.jpl.nasa.gov 14 https://github.com/esheldon/ngmix 15 https://github.com/esheldon/mof 16 https://webbpsf.readthedocs.io/en/stable/ 17 http://galaxia.sourceforge.net 18 Because the same set of offsets is used each time we do a field step, it is possible for all n images of galaxy G to land on the n persistence artefacts from a previous star S. We intend to solve this problem by introducing some pseudo-randomness in the diagonal step sizes, but this has not yet been incorporated in survey simulations. 19 We choose the short axis for two reasons. First, the slew times are shorter, resulting in a more efficient survey. Secondly, the ‘arced’ layout of the focal plane means that we can make a strip with smoother edges by stepping on the short than the long axis; the resulting strips fit together much better when tiling a curved sky. 20 To see why, let |$\xi _\ell = (2\ell +1)C_\ell ^{mm}/(4\pi \sigma _m^2)$| so that ξℓ ≥ 0 and ∑ℓξℓ = 1. Then since Z2 is a positive semidefinite function of the |$C_\ell ^{mm}$|⁠, it is a convex function of the |$C_\ell ^{mm}$|⁠, and Jensen’s inequality (section 3.1.8) shows that |$Z^2 \le \sum _\ell \xi _\ell Z^2_\ell$|⁠, where |$Z^2_{\ell _0}$| is the value of Z2 with the Kronecker delta scale dependence, |$\xi _\ell = \delta _{\ell \ell _0}$|⁠. It follows that Z2 is less than or equal to the maximum of the |$Z^2_{\ell _0}$|⁠. 21 We thank the anonymous referee for encouraging us to think about this more carefully. 22 This was measured on an H2RG. At the time we had to fix this for Phase A requirements flow-down, we did not have a measurement on the H4RG. We now know the charge diffusion for Roman detectors is smaller than this number (Mosby et al. 2020), but it makes a small enough difference that we have not redone the requirements flowdown. 23 It is known that an odd-order mode in the phase can mix with other asymmetric phase modes to produce PSF ellipticity, e.g. if one introduces a large trefoil t then the ellipticity develops a linear term in coma, proportional to tc* (Noecker 2010). However, an amplitude feature with 3-fold or other odd symmetry, such as the spider, does not lead to such an effect. 24 https://roman.gsfc.nasa.gov/science/Roman_Reference_Information.html REFERENCES Ade P. A. R. et al. , 2016 , A&A , 594 , A13 10.1051/0004-6361/201525830 Crossref Search ADS Akeson R. et al. , 2019 , preprint (arXiv:1902.05569) Albrecht A. et al. , 2006 , preprint (arXiv:0609591) Albrecht A. et al. , 2009 , preprint (arXiv:0901.0721) Anderson R. E. , Regan M., Valenti J., Bergeron E., 2014 , preprint (arXiv:1402.4181) Barron N. et al. , 2007 , PASP , 119 , 466 10.1086/517620 Crossref Search ADS Bartelmann M. , Schneider P., 2001 , Phys. Rep. , 340 , 291 10.1016/S0370-1573(00)00082-X Crossref Beletic J. W. et al. , 2008 , in Dorn D. A., Holland A. D., eds, Proc. SPIE Conf. Ser. Vol. 7021, High Energy, Optical, and Infrared Detectors for Astronomy III . SPIE , Bellingham , p. 161 10.1117/12.790382 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Bernstein G. M. , Jarvis M., 2002 , AJ , 123 , 583 10.1086/338085 Crossref Search ADS Biesiadzinski T. , Lorenzon W., Newman R., Schubnell M., Tarlé G., Weaverdyck C., 2011 , PASP , 123 , 179 10.1086/658282 Crossref Search ADS Blas D. , Lesgourgues J., Tram T., 2011 , J. Cosmol. Astropart. Phys. , 7 , 034 10.1088/1475-7516/2011/07/034 Crossref Search ADS Choi A. , Hirata C. M., 2020 , PASP , 132 , 014502 10.1088/1538-3873/ab4504 Crossref Search ADS Claver C. F. , 2019 , The Systems Engineering Integrated Project Team 2019, LSST System Requirements . Available at https://docushare.lsstcorp.org/docushare/dsweb/Get/LSE-29 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Cropper M. et al. , 2013 , MNRAS , 431 , 3103 10.1093/mnras/stt384 Crossref Search ADS DeRose J. et al. , 2019 , preprint (arXiv:1901.02401) DES Collaboration , 2019 , Phys. Rev. Lett. , 122 , 171301 10.1103/PhysRevLett.122.171301 Crossref Search ADS PubMed Doré O. et al. , 2018 , preprint (arXiv:1804.03628) Dore O. et al. , 2019 , BAAS , 51 , 341 Eifler T. et al. , 2019 , BAAS , 51 , 418 Eifler T. et al. , 2020a , preprint (arXiv:2004.04702) Eifler T. et al. , 2020b , preprint (arXiv:2004.05271) Etienne A. et al. , 2018 , in Holland A. D., Beletic J., eds, Proc. SPIE Conf. Ser. Vol. 10709, High Energy, Optical, and Infrared Detectors for Astronomy VIII . SPIE , Bellingham , p. 437 10.1117/12.2314475 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Euclid Collaboration , 2020 , A&A , 635 , A139 10.1051/0004-6361/201936980 Crossref Search ADS Euclid Study Scientist the Science Advisory Team , 2010 , Euclid Science Requirements Document . Available at https://sci.esa.int/documents/33220/36137/1567257215944-Euclid_SciRD_DEM-SA-DC-0001_4_0_2010-03-22.pdf Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Fenech Conti I. , Herbonnet R., Hoekstra H., Merten J., Miller L., Viola M., 2017 , MNRAS , 467 , 1627 10.1093/mnras/stx200 Crossref Freudenburg J. K. C. et al. , 2020 , PASP , 132 , 074504 10.1088/1538-3873/ab9503 Crossref Search ADS Frieman J. , Turner M., Huterer D., 2008 , ARA&A , 46 , 385 10.1146/annurev.astro.46.060407.145243 Crossref Search ADS Górski K. M. , Hivon E., Banday A. J., Wandelt B. D., Hansen F. K., Reinecke M., Bartelmann M., 2005 , ApJ , 622 , 759 10.1086/427976 Crossref Search ADS Green J. et al. , 2011 , preprint (arXiv:1108.1374) Green J. et al. , 2012 , preprint (arXiv:1208.4012) Hemmati S. et al. , 2019 , ApJ , 877 , 117 10.3847/1538-4357/ab1be5 Crossref Search ADS Hikage C. et al. , 2019 , PASJ , 71 , 43 10.1093/pasj/psz010 Crossref Search ADS Hildebrandt H. et al. , 2020 , A&A , 633 , A69 10.1051/0004-6361/201834878 Crossref Search ADS Hirata C. , Seljak U., 2003 , MNRAS , 343 , 459 10.1046/j.1365-8711.2003.06683.x Crossref Search ADS Hirata C. M. , Gehrels N., Kneib J.-P., Kruk J., Rhodes J., Wang Y., Zoubian J., 2012 , preprint (arXiv:1204.5151) Huff E. , Mandelbaum R., 2017 , preprint (arXiv:1702.02600) Hunter J. D. , 2007 , Comput. Sci. Eng. , 9 , 90 10.1109/MCSE.2007.55 Crossref Search ADS Ivezic Z. , LSST Science Collaboration , 2018 , LSST System Science Requirements Document . Available at https://docushare.lsstcorp.org/docushare/dsweb/Get/LPM-17 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Ivezić Ž. et al. , 2019 , ApJ , 873 , 111 10.3847/1538-4357/ab042c Crossref Search ADS Jarvis M. , Jain B., 2004 , preprint (arXiv:astro-ph/0412234) Jee M. J. , Blakeslee J. P., Sirianni M., Martel A. R., White R. L., Ford H. C., 2007 , PASP , 119 , 1403 10.1086/524849 Crossref Search ADS Jurling A. S. , Content D. A., 2012 , in Clampin M. C., Fazio G. G., MacEwen H. A., Oschmann J. M., eds, Proc. SPIE Conf. Ser. Vol. 8442, Space Telescopes and Instrumentation 2012: Optical, Infrared, and Millimeter Wave . SPIE , Bellingham , p. 844210 10.1117/12.925089 Google Scholar Crossref Search ADS Google Preview WorldCat COPAC Kannawadi A. , Shapiro C. A., Mandelbaum R., Hirata C. M., Kruk J. W., Rhodes J. D., 2016 , PASP , 128 , 095001 10.1088/1538-3873/128/967/095001 Crossref Search ADS Kitching T. D. et al. , 2012 , MNRAS , 423 , 3163 10.1111/j.1365-2966.2012.21095.x Crossref Search ADS Kitching T. D. , Taylor A. N., Cropper M., Hoekstra H., Hood R. K. E., Massey R., Niemi S., 2016 , MNRAS , 455 , 3319 10.1093/mnras/stv2523 Crossref Search ADS Korytov D. et al. , 2019 , ApJS , 245 , 26 10.3847/1538-4365/ab510c Crossref Search ADS Lin C.-H. , Tan B., Mandelbaum R., Hirata C. M., 2020 , MNRAS , 496 , 532 10.1093/mnras/staa1298 Crossref Search ADS LSST Science Collaboration , 2009 , preprint (arXiv:0912.0201) LSST Dark Energy Science Collaboration , 2018 , preprint (arXiv:1809.01669) LSST Dark Energy Science Collaboration , 2020 , preprint (arXiv:2010.05926) MacCrann N. et al. , 2018 , MNRAS , 480 , 4614 10.1093/mnras/sty1899 Crossref Search ADS McLeod B. A. , Smith R., 2016 , in Holland A. D., Beletic J., eds, Proc. SPIE Conf. Ser. Vol. 9915, High Energy, Optical, and Infrared Detectors for Astronomy VII . SPIE , Bellingham , p. 99150G 10.1117/12.2233083 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Mandelbaum R. , 2018 , ARA&A , 56 , 393 10.1146/annurev-astro-081817-051928 Crossref Search ADS Mandelbaum R. et al. , 2018 , MNRAS , 481 , 3170 10.1093/mnras/sty2420 Crossref Search ADS Massey R. et al. , 2013 , MNRAS , 429 , 661 10.1093/mnras/sts371 Crossref Search ADS Mosby Gregory J. et al. , 2020 , preprint (arXiv:2005.00505) Neyrinck M. C. , Szapudi I., Szalay A. S., 2009 , ApJ , 698 , L90 10.1088/0004-637X/698/2/L90 Crossref Search ADS Noecker C. , 2010 , in Oschmann J. M., Clampin M. C, MacEwen H. A., eds, Proc. SPIE Conf. Ser. Vol. 7731, Space Telescopes and Instrumentation 2010: Optical, Infrared, and Millimeter Wave . SPIE , Bellingham , p. 77311E 10.1117/12.857744 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Noll R. J. , 1976 , J. Opt. Soc. Am. , 66 , 207 10.1364/JOSA.66.000207 Crossref Search ADS OSC , 1987 , Ohio Supercomputer Center . Available at http://osc.edu/ark:/19495/f5s1ph73 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Padmanabhan N. et al. , 2008 , ApJ , 674 , 1217 10.1086/524677 Crossref Search ADS Paulin-Henriksson S. , Amara A., Voigt L., Refregier A., Bridle S. L., 2008 , A&A , 484 , 67 10.1051/0004-6361:20079150 Crossref Search ADS Perlmutter S. et al. , 1999 , ApJ , 517 , 565 10.1086/307221 Crossref Search ADS Piquette E. C. , McLevige W., Auyeung J., Wong A., 2014 , in Holland A. D., Beletic J., eds, Proc. SPIE Conf. Ser. Vol. 9154, High Energy, Optical, and Infrared Detectors for Astronomy VI . SPIE , Bellingham , p. 776 10.1117/12.2057308 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Plazas A. A. , Shapiro C., Smith R., Rhodes J., Huff E., 2017 , J. Instrum. , 12 , C04009 10.1088/1748-0221/12/04/C04009 Crossref Search ADS Pordes R. et al. , 2007 , J. Phys.: Conf. Ser. , 78 , 012057 10.1088/1742-6596/78/1/012057 Crossref Search ADS Riess A. G. et al. , 1998 , AJ , 116 , 1009 10.1086/300499 Crossref Search ADS Robin A. C. , Reylé C., Derrière S., Picaud S., 2003 , A&A , 409 , 523 10.1051/0004-6361:20031117 Crossref Search ADS Rowe B. T. P. et al. , 2015 , Astron. Comput. , 10 , 121 10.1016/j.ascom.2015.02.002 Crossref Search ADS Samuroff S. et al. , 2018 , MNRAS , 475 , 4524 10.1093/mnras/stx3282 Crossref Search ADS Seo H.-J. , Sato M., Dodelson S., Jain B., Takada M., 2011 , ApJ , 729 , L11 10.1088/2041-8205/729/1/L11 Crossref Search ADS Sfiligoi I. , Bradley D. C., Holzman B., Mhashilkar P., Padhi S., Wurthwein F., 2009 , WRI World Congr. Comput. Sci. Inf. Eng. , 2 , 428 10.1109/CSIE.2009.950 Crossref Sharma S. , Bland-Hawthorn J., Johnston K. V., Binney J., 2011 , ApJ , 730 , 3 10.1088/0004-637X/730/1/3 Crossref Search ADS Sheldon E. S. , 2014 , MNRAS , 444 , L25 10.1093/mnrasl/slu104 Crossref Search ADS Sheldon E. S. , Huff E. M., 2017 , ApJ , 841 , 24 10.3847/1538-4357/aa704b Crossref Search ADS Smith R. M. , Zavodny M., Rahmer G., Bonati M., 2008 , in Dorn D. A., Holland A. D., eds, Proc. SPIE Conf. Ser. Vol. 7021, High Energy, Optical, and Infrared Detectors for Astronomy III . SPIE , Bellingham , p. 192 10.1117/12.789619 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Spergel D. et al. , 2013 , preprint (arXiv:1305.5422) Spergel D. et al. , 2015 , preprint (arXiv:1503.03757) Suchyta E. et al. , 2016 , MNRAS , 457 , 786 10.1093/mnras/stv2953 Crossref Search ADS Troxel M. A. et al. , 2018 , Phys. Rev. D , 98 , 043528 10.1103/PhysRevD.98.043528 Crossref Vavrek R. D. et al. , 2016 , in Angeli G. Z., Dierickx P., eds, Proc. SPIE Conf. Ser. Vol. 9911, Modeling, Systems Engineering, and Project Management for Astronomy VII . SPIE , Bellingham , p. 991105 10.1117/12.2233015 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Weinberg D. H. , Mortonson M. J., Eisenstein D. J., Hirata C., Riess A. G., Rozo E., 2013 , Phys. Rep. , 530 , 87 10.1016/j.physrep.2013.05.001 Crossref Zandian M. et al. , 2016 , in Holland A. D., Beletic J., eds, Proc. SPIE Conf. Ser. Vol. 9915, High Energy, Optical, and Infrared Detectors for Astronomy VII . SPIE , Bellingham , p. 148 10.1117/12.2233664 Google Scholar Google Preview OpenURL Placeholder Text WorldCat COPAC Crossref Zentner A. R. , Rudd D. H., Hu W., 2008 , Phys. Rev. D , 77 , 043507 10.1103/PhysRevD.77.043507 Crossref Zentner A. R. , Semboloni E., Dodelson S., Eifler T., Krause E., Hearin A. P., 2013 , Phys. Rev. D , 87 , 043509 10.1103/PhysRevD.87.043509 Crossref Zonca A. , Singer L., Lenz D., Reinecke M., Rosset C., Hivon E., Gorski K., 2019 , J. Open Source Softw. , 4 , 1298 10.21105/joss.01298 Crossref Search ADS Zuntz J. et al. , 2018 , MNRAS , 481 , 1149 10.1093/mnras/sty2219 Crossref Search ADS APPENDIX A: OVERVIEW OF WEAK LENSING SYSTEMATICS BUDGETING This appendix describes the requirements flowdown and error budgeting for the weak lensing program on the Roman mission, and documents the detailed rationale behind the summary requirements listed in the Roman SRD. This kind of error budgeting has been performed elsewhere in the literature (Paulin-Henriksson et al. 2008; Massey et al. 2013), but this document focuses on the error terms relevant to Roman. For example, the PSFs are based on an obstructed pupil with low-order aberrations rather than using generic formulae involving second moments (some such formulae, including those used in the Joint Dark Energy Mission and WFIRST Interim Design Reference Mission studies, were for Gaussians). We set most systematics requirements for this mission on the basis of having systematic errors sub-dominant to statistical errors in the weak lensing shear power spectra or cross-power spectra (or any linear combinations thereof). Exceptions to this policy are considered in cases where meeting the original systematic budget becomes a cost or complexity driver, or is not possible. Most measurement biases – including those considered in this paper – fall into the ‘additive’ or ‘multiplicative’ forms (see Section A1) and will be treated according to the formalism therein. A1 Additive and multiplicative biases The cosmic shear measurement is sensitive to two major types of measurement errors. Additive bias or ‘spurious shear’ c is a shear signal that is detected even when none is present. Multiplicative bias or ‘calibration bias’ m is an incorrect response to a real shear, e.g. a shear γ is present in the sky but the measurement yields 1.01γ. Normally, we think of additive biases as resulting from mis-estimation of the PSF ellipticity (or its variation across the sky), whereas multiplicative biases result from mis-estimation of the size of the PSF. However, detector non-linearities, approximations used in the data processing/analysis pipelines, and uncertainties about the distribution of galaxy morphologies in the sky can also contribute to both types of biases. The E-mode shear cross-power spectrum between two redshift bins zi and zj is modified in the presence of these biases: $$\begin{eqnarray*} C_\ell ^{z_i,z_j}({\rm obs}) = (1+m_i)(1+m_j)C_\ell ^{z_i,z_j}({\rm true}) + C_\ell ^{c_i,c_j}, \end{eqnarray*}$$(A1) where we write mi ≡ m(zi) as a shorthand for the bias in bin i. To linear order in the biases, the correction to the power spectrum can be written as $$\begin{eqnarray*} \Delta C_\ell ^{z_i,z_j} &=& C_\ell ^{z_i,z_j}({\rm obs}) - C_\ell ^{z_i,z_j}({\rm true}) \nonumber \\ &=& (m_i+m_j)C_\ell ^{z_i,z_j}+C_\ell ^{c_i,c_j}. \end{eqnarray*}$$(A2) If m is spatially variable, there is an additional contribution (e.g. Kitching et al. 2012, 2016): $$\begin{eqnarray*} \Delta C_\ell ^{z_i,z_j} = \int \frac{d^2{\boldsymbol{\ell }}^{\prime }}{(2\pi)^2} \, C_{\ell ^{\prime }}^{z_i,z_j} C^{\delta m_i,\delta m_j}_{{\boldsymbol{\ell }}-{\boldsymbol{\ell }}^{\prime }} \cos ^2 (2\varphi _{{\boldsymbol{\ell }},{\boldsymbol{\ell }}^{\prime }}), \end{eqnarray*}$$(A3) where δmi is the fluctuation in multiplicative bias of bin i, and |$\varphi _{{\boldsymbol{\ell }},{\boldsymbol{\ell }}^{\prime}}$| is the angle between the indicated wave vectors. This is second order in the m-biases, so we expect it to be small compared to the first-order contribution in equation (A2), which comes from the spatially averaged part of the multiplicative bias. We will briefly discuss this spatially variable contribution again in Section A2.3. A2 Setting requirements The power spectra are arranged into a vector C with a covariance matrix Σ. For the weak lensing power spectrum, with Nz redshift bins and Nℓ angular scale bins, there are NℓNz(Nz + 1)/2 power spectra |$C_\ell ^{z_i,z_j}$|⁠; hence C is a vector of length NℓNz(Nz + 1)/2, and Σ is a matrix of size NℓNz(Nz + 1)/2 × NℓNz(Nz + 1)/2. A contaminant that changes the power spectrum by ΔC can have its significance assessed by $$\begin{eqnarray*} Z = \sqrt{\Delta {\bf C} {\bf \Sigma }^{-1}\Delta {\bf C}}, \end{eqnarray*}$$(A4) which is the number of σs at which one could distinguish the correct power spectrum from the contaminated power spectrum. Note that as the survey area Ω is increased, Z will increase as ∝ Ω1/2, and hence contaminants ΔC must be reduced to keep them below statistical errors. If Z = 1, then the power spectrum is biased at the same level as the statistical errors. We use Z as a metric for contaminants, rather than e.g. biases in (w0, wa)-space, for generality: if Z < 1 then the bias due to ΔC in any cosmological parameter from the combination of the Roman weak lensing power spectrum with any other data set(s) from Roman or other experiments is <1σ; whereas if one based the analysis on biases in (w0, wa) then we would need a separate requirement derived from every cosmological analysis planned on Roman weak lensing data. Using Z as a metric also enables us to write requirements that do not depend on other cosmological probes (e.g. the Roman weak lensing systematic error budget does not change if we discover a new way to reduce the scatter in the SN Ia Hubble diagram), which will help to ensure the stability of our requirements going forward. Technically, the above discussion applies only to the E-mode of spurious shear; we have not set a specific requirement on the B-mode, which contains no cosmological information to linear order and is used as a null test. For the latter reason, we set a requirement on the B-mode that is equal to the requirement on the E-mode, so that the B-mode null test will pass if requirements are met. We also note that the weak lensing analysis includes a range of angular scales, ℓmin,tot ≤ ℓ ≤ ℓmax,tot; requirements apply to sources of systematic error that affect these scales, i.e. are ‘in-band’ for the weak lensing measurement. The ‘in-band’ qualifier is critical: as an example, pixelization errors can cause shape measurement errors in galaxies that depend on whether the galaxy lands on a pixel centre, corner, vertical edge, or horizontal edge. For some shape measurement methods, this error may dramatically exceed the additive systematic error budget, but it is concentrated at very small angular scales (multiples of 2π divided by the pixel scale P, or 2π/P = 1.2 × 107). Our requirements are set on the portion of this power that is within (or mixes into) the band limit, ℓ ≤ ℓmax,tot due to e.g. edge effects, selection effects, etc. Equation (A4) still does not completely define a requirement, since we have not described the redshift or scale dependence of the spurious shear in question. Neither dependence is expected to be trivial: errors in PSF models have a greater impact on shape measurements for higher redshift galaxies, since they tend to be smaller; and the angular power spectrum of PSF model errors should be non-white in a survey strategy that ‘marches’ across the sky, even if heavily cross-linked (there may also be a characteristic scale at the size of the field; for example, a repeating error at the ∼0.8 × 0.4° size of the Roman field has reciprocal lattice frequencies at ℓ = 450 and 900, so a large-scale error in the instrument PSF model that is ‘tessellated’ as we tile the sky will appear at these frequencies or multiples thereof). At first, we considered assuming a particular scale and redshift dependence for the errors, but in order to be conservative we would have to assume the worst combination of angular and redshift dependences. Many of our large sources of systematic error, such as PSF ellipticity due to astigmatism, have predictable dependences (e.g. the systematic error induced in galaxy shears is of the same sign in all redshift bins) that are far from the worst case, and this could lead to overconservatism in the requirements. Therefore, we need a more nuanced approach to the requirements, where the allowed amplitude of each term in the error budget is informed by the structure of the correlations it produces. Our approach to this problem is to write a script that accepts a specific angular and redshift dependence (‘template’) for a systematic error, and returns the amplitude A0 of the systematic error at which we would have Z = 1 (i.e. a 1σ bias on the most-contaminated direction in power spectrum space). For cases where the template is not known (or where we have not done the analysis), the script is capable of searching the space of templates and finding the most conservative choice, i.e. the choice that leads to the smallest value of A0. The combined results enable us to build an error tree, where the overall top-level systematics requirement (a limit on Z) can be flowed down to upper limits on each source of systematic error. Finally, some portions of the systematic error budget sum in quadrature (‘root-sum-square’ or RSS addition) and others linearly; in this document, we carefully account for which is which. A2.1 Data vector and covariance model We build our data vector for the shear power spectra and cross-spectra. We recognize that weak lensing analyses have shifted to ‘3 × 2-point’ data vectors containing shear–shear, galaxy–shear, and galaxy–galaxy correlations, and by the time of Roman the list of standard observables may be even longer. However, for setting requirements on shape measurement, shear–shear provides the most demanding use case, and so for simplicity here we only consider shear–shear. We use for our data vector the NℓNz(Nz + 1)/2 power spectra and cross-spectra. Each ℓ is treated separately, so there are Nℓ = ℓmax,tot − ℓmin,tot angular bins; we use ℓmin,tot = 10 and ℓmax,tot = 3161, thereby covering 2.5 orders of magnitude in scale. Roman provides little cosmological constraining power at the larger scales due both to the finite size of its survey and due to the large cosmic variance of the lowest multipoles. The smallest scales are generally not used in cosmic shear analyses because the baryonic effects are severe (e.g. Zentner, Rudd & Hu 2008; Zentner et al. 2013). We use Nz = 15 redshift slices, as shown in Table A1, which are chosen by the Exposure Time Calculator (Hirata et al. 2012) v17, with the Phase B exposure times (5 × 140.25 s in H158-band). In order to ensure that Roman would not become systematics-limited in an extended mission, we set the top-level requirement on systematics to Z = 1 for a survey of area Ω = 104 deg2 (3.05 sr). Table A1. The effective number density in each redshift bin, in units of galaxies/arcmin2, used for setting requirements. These are per bin, i.e. are dneff/dz × Δz. z . neff . z . neff . z . neff . 0.10 ± 0.10 3.62 1.10 ± 0.10 3.75 2.10 ± 0.10 1.21 0.30 ± 0.10 2.12 1.30 ± 0.10 3.17 2.30 ± 0.10 0.95 0.50 ± 0.10 3.05 1.50 ± 0.10 2.52 2.50 ± 0.10 0.82 0.70 ± 0.10 5.90 1.70 ± 0.10 1.45 2.70 ± 0.10 0.68 0.90 ± 0.10 2.79 1.90 ± 0.10 1.68 2.90 ± 0.10 0.19 z . neff . z . neff . z . neff . 0.10 ± 0.10 3.62 1.10 ± 0.10 3.75 2.10 ± 0.10 1.21 0.30 ± 0.10 2.12 1.30 ± 0.10 3.17 2.30 ± 0.10 0.95 0.50 ± 0.10 3.05 1.50 ± 0.10 2.52 2.50 ± 0.10 0.82 0.70 ± 0.10 5.90 1.70 ± 0.10 1.45 2.70 ± 0.10 0.68 0.90 ± 0.10 2.79 1.90 ± 0.10 1.68 2.90 ± 0.10 0.19 Open in new tab Table A1. The effective number density in each redshift bin, in units of galaxies/arcmin2, used for setting requirements. These are per bin, i.e. are dneff/dz × Δz. z . neff . z . neff . z . neff . 0.10 ± 0.10 3.62 1.10 ± 0.10 3.75 2.10 ± 0.10 1.21 0.30 ± 0.10 2.12 1.30 ± 0.10 3.17 2.30 ± 0.10 0.95 0.50 ± 0.10 3.05 1.50 ± 0.10 2.52 2.50 ± 0.10 0.82 0.70 ± 0.10 5.90 1.70 ± 0.10 1.45 2.70 ± 0.10 0.68 0.90 ± 0.10 2.79 1.90 ± 0.10 1.68 2.90 ± 0.10 0.19 z . neff . z . neff . z . neff . 0.10 ± 0.10 3.62 1.10 ± 0.10 3.75 2.10 ± 0.10 1.21 0.30 ± 0.10 2.12 1.30 ± 0.10 3.17 2.30 ± 0.10 0.95 0.50 ± 0.10 3.05 1.50 ± 0.10 2.52 2.50 ± 0.10 0.82 0.70 ± 0.10 5.90 1.70 ± 0.10 1.45 2.70 ± 0.10 0.68 0.90 ± 0.10 2.79 1.90 ± 0.10 1.68 2.90 ± 0.10 0.19 Open in new tab The power spectra were obtained from class (Blas, Lesgourgues & Tram 2011) using the fiducial cosmology from the Planck 2015 ‘TT,TE,EE+lowP+lensing+ext’ results (Ade et al. 2016). The shape noise contribution was added to construct Ctot according to $$\begin{eqnarray*} C^{{\rm tot},z_i,z_j}_\ell = C^{z_i,z_j}_\ell + \frac{\gamma _{\rm rms}^2}{\bar{n}_i}\delta _{ij}, \end{eqnarray*}$$(A5) where |$\bar{n}_i$| is the mean effective number density in galaxies per steradian in redshift slice i, and γrms is the shape noise expressed as an equivalent RMS shear per component; we take γrms = 0.22. We approximate |${\boldsymbol{\Sigma}}$| using the usual Gaussian covariance matrix formula, $$\begin{eqnarray*} &&\Sigma [C^{z_i,z_j}_\ell , C^{z_k,z_m}_{\ell ^{\prime }}]= \nonumber \\ &&\qquad \frac{\delta _{\ell \ell ^{\prime }} \big[ C^{{\rm tot},z_i,z_k}_\ell C^{{\rm tot},z_j,z_m}_\ell + C^{{\rm tot},z_i,z_m}_\ell C^{{\rm tot},z_j,z_k}_\ell \big]}{(2\ell +1)f_{\rm sky}} , \end{eqnarray*}$$(A6) where fsky = Ω/(4π). The non-Gaussian contributions to the error covariance matrix are turned off, because since the FoMSWG (Albrecht et al. 2009) there has been an ongoing program of using non-linear transformations on the data to remove them (e.g. Neyrinck, Szapudi & Szalay 2009; Seo et al. 2011) and we do not want applications of these novel statistics to Roman data to run into systematic error limits. We also do not include astrophysical systematic errors in |${\boldsymbol{\Sigma}}$|⁠; we envision instead that they will be treated with nuisance parameters in the analysis. Another advantage of this is that the covariance matrix |${\boldsymbol{\Sigma}}$| is block diagonal in ℓ-space (it is formally 375 720 × 375 720 without ℓ-binning), which makes computations possible on a machine with limited memory. Indeed, in the Gaussian case one may write $$\begin{eqnarray*} &&\Delta {\bf C} \cdot {\boldsymbol{\Sigma }}^{-1} \Delta {\bf C}=\nonumber \\ &&\qquad \sum _\ell \frac{(2\ell +1)f_{\rm sky}}{2} \sum _{ijkm} \Delta C_\ell ^{ij} [C^{{\rm tot}-1}_\ell ]_{jk} \Delta C_\ell ^{km} [C^{{\rm tot}-1}_\ell ]_{mi}, \end{eqnarray*}$$(A7) where the matrix inverses are Nz × Nz. A2.2 Implementation: additive systematics Each additive systematic error is taken to have an angular dependence given by some template Tℓ, and a redshift dependence given by a set of weights wi = w(zi). That is, there is a reference additive shear cref, with the additive shear in redshift bin i given by c(zi) = wicref. For example, a systematic error independent of redshift bin would be specified with wi = 1 for all i. The reference signal is taken to have a power spectrum proportional to the template: |$C_\ell ^{c_{\rm ref}} = A_0^2T_\ell$|⁠, and the template is normalized so that cref has variance 1 per component (from in-band fluctuations): $$\begin{eqnarray*} \sum _{\ell =\ell _{\rm min,tot}}^{\ell _{\rm max,tot}} \frac{2\ell +1}{4\pi } T_\ell = 1. \end{eqnarray*}$$(A8) The additive cross-power spectrum is then $$\begin{eqnarray*} C_\ell ^{c_i,c_j} = A_0^2 w_iw_jT_\ell , \end{eqnarray*}$$(A9) and the total RMS per component of the spurious shear in bin i is A0|wi|. The additive systematic errors can have various scale dependences. We therefore consider a suite of Nband disjoint angular templates that cover the shape measurement band. Each template satisfies the normalization rule, equation (A8), and has ℓ(ℓ + 1)Tℓ/(2π) = constant: $$\begin{eqnarray*} T_\ell ^{(\alpha)} = &\left[ \sum _{\ell ^{\prime }=\ell _{\rm min}(\alpha)}^{\ell _{\rm max}(\alpha)} \frac{2\ell ^{\prime }+1}{\ell ^{\prime }(\ell ^{\prime }+1)} \right]^{-1} \frac{4\pi }{\ell (\ell +1)}\nonumber \\ &\times \left\lbrace \begin{array}{cc}1 & \ell _{\rm min}(\alpha)\le \ell \le \ell _{\rm max}(\alpha) \\ 0 & {\rm otherwise} \end{array} \right.\!,~\alpha =0,...N_{\rm band}-1. \end{eqnarray*}$$(A10) The current bands are displayed in Table A2. Each band α is allowed a contribution to the total error Z(α). Since there are no statistical correlations between different ℓs in the covariance matrix Σ, the Z(α) can be quadrature-summed (see equation A4). However, additive systematic error is positive in the sense that it adds rather than subtracts power; thus the power spectrum error vectors ΔC from two sources of additive systematic error contributing to the same angular bin are not orthogonal and the Z’s should be added linearly. Another way to think of this is that since Z is proportional to the square of the RMS shear, |$Z\propto A_0^2$|⁠, quadrature-summation of the additive shear is equivalent to linear summation of the Z-values. Table A2. The requirements for additive and multiplicative systematic errors. There are Nband = 4 additive error bands ranging over a total signal band from ℓmin,tot = 10 to ℓmax,tot = 3161. The fraction of the error budget allocated to each band is also indicated, as are the maximum allowed redshift-independent spurious shear (⁠|$A_0^{\rm flat}(\alpha)$|⁠, RMS per component), and the maximum scaling factors for redshift dependence, Smax,±(α) and Smax,+(α). There is only one row for the multiplicative errors, since the implementation does not contain an ℓ dependence; we quote a requirement on the post-calibration shear multiplicative uncertainty |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$|⁠. Band α . ℓmin(α) . ℓmax(α) . Allocation Z(α) . Sys. err. req’t. . Smax,±(α) . Smax,+(α) . . . . . |$A_0^{\rm flat}(\alpha)$| or |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$| . . . Additive errors 0 31 99 0.2596 7.000 × 10−5 8.489 2.782 1 100 315 0.2539 9.900 × 10−5 5.628 2.041 2 316 999 0.2575 1.400 × 10−4 3.569 1.509 3 1000 3161 0.2538 1.900 × 10−4 2.119 1.149 Multiplicative errors mult 0.4600 3.200 × 10−4 2.186 1.140 Band α . ℓmin(α) . ℓmax(α) . Allocation Z(α) . Sys. err. req’t. . Smax,±(α) . Smax,+(α) . . . . . |$A_0^{\rm flat}(\alpha)$| or |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$| . . . Additive errors 0 31 99 0.2596 7.000 × 10−5 8.489 2.782 1 100 315 0.2539 9.900 × 10−5 5.628 2.041 2 316 999 0.2575 1.400 × 10−4 3.569 1.509 3 1000 3161 0.2538 1.900 × 10−4 2.119 1.149 Multiplicative errors mult 0.4600 3.200 × 10−4 2.186 1.140 Open in new tab Table A2. The requirements for additive and multiplicative systematic errors. There are Nband = 4 additive error bands ranging over a total signal band from ℓmin,tot = 10 to ℓmax,tot = 3161. The fraction of the error budget allocated to each band is also indicated, as are the maximum allowed redshift-independent spurious shear (⁠|$A_0^{\rm flat}(\alpha)$|⁠, RMS per component), and the maximum scaling factors for redshift dependence, Smax,±(α) and Smax,+(α). There is only one row for the multiplicative errors, since the implementation does not contain an ℓ dependence; we quote a requirement on the post-calibration shear multiplicative uncertainty |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$|⁠. Band α . ℓmin(α) . ℓmax(α) . Allocation Z(α) . Sys. err. req’t. . Smax,±(α) . Smax,+(α) . . . . . |$A_0^{\rm flat}(\alpha)$| or |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$| . . . Additive errors 0 31 99 0.2596 7.000 × 10−5 8.489 2.782 1 100 315 0.2539 9.900 × 10−5 5.628 2.041 2 316 999 0.2575 1.400 × 10−4 3.569 1.509 3 1000 3161 0.2538 1.900 × 10−4 2.119 1.149 Multiplicative errors mult 0.4600 3.200 × 10−4 2.186 1.140 Band α . ℓmin(α) . ℓmax(α) . Allocation Z(α) . Sys. err. req’t. . Smax,±(α) . Smax,+(α) . . . . . |$A_0^{\rm flat}(\alpha)$| or |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$| . . . Additive errors 0 31 99 0.2596 7.000 × 10−5 8.489 2.782 1 100 315 0.2539 9.900 × 10−5 5.628 2.041 2 316 999 0.2575 1.400 × 10−4 3.569 1.509 3 1000 3161 0.2538 1.900 × 10−4 2.119 1.149 Multiplicative errors mult 0.4600 3.200 × 10−4 2.186 1.140 Open in new tab The allocations for each bin Z(α) were initially set to |$\sqrt{0.25/N_{\rm band}}$|⁠, so that in an RSS sense 25 per cent of the systematic error budget is allocated to additive shear; with four bands this implies Z(α) = 0.25 (i.e. 6.25 per cent of the error budget in an RSS sense) for each band. There have been some updates of the exposure times, throughputs, and number densities since the SRD requirements were set (2016 December); we have kept the requirements the same, and updated the Z-values, so the latter do not exactly equal 0.25. The construction of Z-values for each angular band and each additive systematic is mathematically sufficient to build the error budget. However, they can be difficult to conceptualize. Therefore, we introduce some equivalent notation to describe the weak lensing error budget. For each angular template, we introduce a limiting amplitude |$A_0^{\rm flat}(\alpha)$|⁠, defined to be the amplitude A0 at which we would saturate the requirement on Z(α) for bin α in the case of a redshift-independent systematic |$w_i=1\, \forall i$|⁠. That is, if the additive systematics did not depend on redshift, we could tolerate a total additive systematic shear of |$A_0^{\rm flat}$| (RMS per component) in band α. We also introduce a scaling factor S[w, α] for a systematic error $$\begin{eqnarray*} S[{\bf w},\alpha ] = \frac{Z(\alpha) \, {\rm for\, this\, }w_i}{Z(\alpha)\, {\rm for\, all\, }w_i=1} \end{eqnarray*}$$(A11) that depends on the redshift dependence wi. An additive systematic error that is independent of redshift will have S = 1. A systematic that is ‘made worse’ by its redshift dependence will have S > 1, and a systematic that is ‘made less serious’ by its redshift dependence will have S < 1. The requirement that the (linear) sum of Zs not exceed Z(α) thus translates into $$\begin{eqnarray*} \sum _{\rm systematics} [A(\alpha)]^2\times S[{\bf w},\alpha ] \le [A_0^{\rm flat}(\alpha)]^2, \end{eqnarray*}$$(A12) where A(α) is the RMS additive shear per component due to that systematic. In most cases, we will take the ‘reference’ additive shear to be the additive shear in the most contaminated redshift slice; in this case, wi = 1 for that slice, and |wi| ≤ 1 for the others. Under such circumstances, we can determine a worst-case scaling factorSmax,±(α), which is the largest value of S[w, α] for any weights satisfying the above inequality. We may also determine a worst-case scaling factor Smax,+(α) conditioned on 0 ≤ wi ≤ 1, i.e. for sources of additive shear that have the same sign in all redshift bins. The search within these spaces is simplified by the fact that – according to equation (A7) – the contribution to S[w, α] considering only a single value of ℓ reduces to a semipositive-definite quadratic function of w (it is proportional to |${\bf w}^{\rm T}{\bf C}_\ell ^{{\rm tot}-1}{\bf w}$|⁠, where |${\bf C}_\ell ^{\rm tot}$| is Nz × Nz). Therefore, the worst-case weights |$\lbrace w_i\rbrace _{i=1}^{N_z}$| always occur at the corners of the allowed cube in Nz-dimensional w-space, and we can simply search the |$2^{N_z}$| corners by brute force. A2.3 Implementation: multiplicative systematics The implementation for the spatial mean part of the multiplicative systematic errors is simpler, since one can work directly with the mi. Once again, we may write mi = mwi, where m is the multiplicative error in the worst bin (largest absolute value) and wi = 1 for that bin, and |wi| ≤ 1 for all bins; the wi thus represents the redshift dependence of the multiplicative error. Once again, we may define a scaling factor S[w, mult] for multiplicative biases analogous to equation (A11): $$\begin{eqnarray*} S[{\bf w},{\rm mult}] = \frac{Z^2({\rm mult}) \, {\rm for\, this\, }w_i}{Z^2({\rm mult})\, {\rm for\, all\, }w_i=1}; \end{eqnarray*}$$(A13) this time, we define this with the Z2 rather than Z so that RSS addition will apply to independent multiplicative errors: $$\begin{eqnarray*} \sum _{\rm systematics} \sigma _m^2 \times S[{\bf w},{\rm mult}] \le [\sigma _{m,\rm req^{\prime }t}^{\rm flat}]^2, \end{eqnarray*}$$(A14) where |$\sigma _{m,\rm req^{\prime }t}^{\rm flat}$| is the requirement on knowledge of m. Fundamentally, the square present here but not in equation (A11) arises because multiplicative biases in the power spectrum are proportional to m but additive biases in the power spectrum are proportional to c2. The worst-case scaling factors Smax,+(mult) (conditioned on 0 ≤ wi ≤ 1) and Smax,±(mult) (allowing either sign) can be defined analogously. In this case, since ΔC is linear in the mi and hence wi, it is actually S2[w, mult] that is a semipositive-definite quadratic function of w instead of S, but the technique of searching the corners by brute force still applies. Once again, we initially set Z(mult) = 0.5, allocating 25 per cent of the systematic error in an RSS sense to multiplicative systematics; due to changes in the model since the SRD was first written, the allocation is no longer exactly 25 per cent. The resulting limits are quoted in Table A2. One may also consider the spatially varying multiplicative systematics, which contribute to the power spectrum via equation (A3). As noted in Section A1, we expected this contribution to be small. As a simple test, we tried assessing a redshift-independent multiplicative bias, with the scale dependence |$\ell (\ell +1)C_\ell ^{mm}=$| constant over the range 8 ≤ ℓ < 1442 (i.e. from one wavelength over the survey region to a half-wavelength across an SCA), and with total variance |$\sigma ^2_m = \sum _\ell (2\ell +1)C_\ell ^{mm}/(4\pi)$|⁠. The resulting contribution to the error budget is Z2 = 0.016(σm/0.01)4. We also tried using a Kronecker delta scale dependence, |$C_\ell ^{mm} = 4\pi \sigma _m^2 \delta _{\ell \ell _0}/(2\ell +1)$|⁠; by searching over all the ℓ0 in the range from 8 ≤ ℓ0 < 1442, we found the worst contribution to the error budget to be Z2 = 0.020(σm/0.01)4 if ℓ0 = 16. This represents an upper bound to Z2 for a given σm.20 It is thus clear that the requirements on spatially varying multiplicative systematics will be much looser than those on the spatially varying additive systematics, and hence it is the latter that will drive PSF and wavefront stability requirements. The allocation in terms of acceptable Z2 for spatially varying multiplicative systematics the overall Roman error budget is under discussion.21 A3 Flow-down to PSF requirements In order to translate a requirement on additive bias c or multiplicative bias m into requirements on lower-level quantities, we need to know how a given effect – e.g. an error in the PSF model – affects the shear measurement. We focus here on the additive biases, which are of interest for this paper; the multiplicative biases can be treated in the same formalism and we comment on how to do this at the end. We need to compute ∂γobs(zi)/∂X, where γobs is the measured shear in a region (and in redshift slice i) and X is any quantity on which we want to set a knowledge requirement. The spurious shear in bin i is taken to be $$\begin{eqnarray*} c_i = c(z_i) = \frac{\partial \gamma _{\rm obs}(z_i)}{\partial X} \Delta X, \end{eqnarray*}$$(A15) where ΔX = Xtrue − Xmodel is the error in knowledge of X. In the context of the additive systematic errors, the ratios of the partial derivatives ∂γobs(zi)/∂X set the redshift slice dependence: if i(max) is the redshift bin with the largest derivative (in absolute value) then $$\begin{eqnarray*} w_i = \frac{\partial \gamma _{\rm obs}(z_i)/\partial X}{\partial \gamma _{\rm obs}(z_{i({\rm max})})/\partial X} \end{eqnarray*}$$(A16) and the reference additive shear is cref = ci(max). In principle, the coefficients ∂γobs(zi)/∂X depend on the base model for the PSF, the population of galaxies, and the shape measurement algorithm. Multiple algorithms should be used for Roman, but a final selection has not been made (given how rapidly the field is maturing, such a choice now would be premature). However, all practical methods of measuring shear have some basic properties in common – if e.g. the true PSF has greater e1 than the model (i.e. is elongated in the x-direction), then the inferred shear in that region of the sky will also have greater c1, and this effect will be greater for larger PSFs or smaller galaxies. In setting requirements, we therefore chose a simple, easily understood model. This model is not, and does not need to be, an accurate description of Roman shape measurement at the few ×10−4 accuracy. Rather, it needs to give us estimates of ∂γobs(zi)/∂X early in the development of Roman, with the understanding that we will not update the optical stability requirements every time we have a better model for the distribution of galaxy morphologies. The very simplest choice would be to work with Gaussian PSFs and galaxies; however, our previous experience has been that the non-Gaussian tails of both PSFs and galaxies matter, and furthermore when we discuss the Zernike description of wavefront errors we have predictions for how various combinations of modes affect the ellipticities of different isophotes of the PSF. Therefore, we go one step beyond the Gaussian approximation and include in our analytic flow-down model: Our galaxies are taken to have an exponential profile, |$f_{\rm circ}({\bf x}) \propto \mathrm{ e}^{-1.67834|{\bf x}|/r_{\rm eff}}$|⁠, where rreff is the half-light radius. It can optionally be sheared by applying a finite shear γ to arrive at the galaxy f(x). The PSF is the Fourier transform of an annular pupil with aberrations appearing as contributions to the phase. The resulting ‘optical’ PSF is then convolved with a detector response that includes a tophat and charge diffusion. For HgCdTe detectors, we take the charge diffusion length to be 2.94 μm rms per axis (Barron et al. 2007).22 Other effects that are likely significant for Roman analyses, such as IPC, the brighter-fatter effect, or polarization, are not included. We turned the spider off. The main effect of the spider is the production of 12 diffraction spikes, but it is the core of the PSF that matters most for shape measurement and has the greatest change as one adjusts the Zernike coefficients. The spider further leads to an asymmetric pupil, i.e. with odd-order modes in the decomposition of the amplitude, but this has no appreciable effects on the relation of ellipticity to low-order Zernike modes.23 We use as our measure of ellipticity the two-component ellipticity eI of the observed image I = f⋆P (where f is the galaxy, P is the PSF, and ⋆ denotes convolution). The ellipticity is determined according to the adaptive moment algorithm of Bernstein & Jarvis (2002), section 3.1. The observed two-component ellipticity eI of the galaxy is related to the shear by a 2 × 2 responsivity matrix $$\begin{eqnarray*} R_{ij} = \frac{\partial e_{I,i}}{\partial \gamma _j} = {\cal R}\delta _{ij} + R^{\rm aniso}_{ij}, \end{eqnarray*}$$(A17) which we have decomposed into an isotropic part |${\cal R}$| and a traceless matrix Raniso characterizing the anisotropic part of the responsivity. The inverse of the responsivity matrix relates a bias in the galaxy ellipticities to a bias in the shear: $$\begin{eqnarray*} c_i = \sum _{j=1}^2 [{\bf R}^{-1}]_{ij} \frac{\partial e_{I,j}}{\partial X} \Delta X. \end{eqnarray*}$$(A18) Since the isotropic part of the responsivity dominates except for extreme PSF ellipticity, anisotropic noise correlations, etc., we take the isotropic part and write $$\begin{eqnarray*} \frac{\partial \gamma _{{\rm obs},i}(z_k)}{\partial X} = \left\langle {\cal R}^{-1} \frac{\partial e_{I,i}}{\partial X} \right\rangle , \end{eqnarray*}$$(A19) where the average is taken over the source galaxies in that redshift bin. The various partial derivatives are easily computed as finite differences of the galaxy simulation and ellipticity measurement process. As this model is intended to be simple, the average is taken only over the distribution of source sizes reff – we do not include the intrinsic source ellipticity or a distribution of Sersić indices. The main difference that occurs with the multiplicative systematic errors is that when one changes the PSF size, one must look at the change in responsivity, i.e. |$m = \partial \ln {\cal R}/\partial X$|⁠. APPENDIX B: REQUIREMENTS ON WAVEFRONT STABILITY FOR THE PSF CALIBRATION The determination of the PSF in imaging mode will be based on an empirical (principal components or more advanced version thereof) approach (these methods have a long history in weak lensing – see e.g. Jarvis & Jain 2004; Jee et al. 2007 – but a large amount of work will be required to adapt them to Roman), or on physical fitting of the optical model (Jurling & Content 2012) with empirical corrections. Central to both of these approaches is that we must limit the number of possible principal components in the data by limiting the number of properties of the PSF that vary from one image to another. The Roman approach is to keep the PSF stable during an exposure so that no parameters are needed to describe time dependence of the PSF during an exposure. We make one exception to this policy for image motion, since at the Roman weak lensing level of precision this is unavoidable. Thus, the requirement is for the optics + image motion PSF to be the convolution of the optics PSF with a kernel coming from the image motion, with small residuals. Here, ‘small’ means that the residual error must fit within the overall error budget for PSF (or shear) errors. We note that since Roman detectors can be read non-destructively, and six sub-exposures will be sent to the ground, that one could imagine building a time-dependent PSF from these sub-exposures. We have chosen not to set a looser requirement based on this expectation, since we plan to calibrate detector non-linearity using the consistency of the sub-exposures; this approach does not work if the PSF is varying in an uncontrolled way. Requirements are derived for the two major sources of wavefront change: slow drifts induced by, e.g. thermal variations (Section B1), and jitter induced by e.g. vibrations from the reaction wheels (Section B2). B1 Wavefront drift B1.1 Flowdown methodology In general, we suppose that there is a vector of parameters p that determines the PSF in each exposure (including its field dependence). Some of these are associated with the equilibrium wavefront – this is the subject of this section – whereas others are associated with image motion, jitter, detector properties, etc. The amplitudes |$\psi _i({\boldsymbol{\theta }})$| of each Zernike component of the wavefront error – which depend on field position |${\boldsymbol{\theta}}$| – are functions of these parameters, and will each have their own time dependence |$\psi _i({\boldsymbol{\theta }};t)$|⁠. This induces a time dependence in the PSF |$G({\bf x};{\boldsymbol{\theta }};t)$|⁠, and hence in the observed shear γobs for an object. We may write the amplitudes ψi at a given position as a vector |${\boldsymbol{\psi}}$| of length NZern, where NZern is the number of Zernike coefficients kept. We normalize the Zernike modes to unit RMS, so that |$|{\boldsymbol{\psi }}({\boldsymbol{\theta }})|$| is the RMS wavefront error at position |${\boldsymbol{\theta}}$|⁠. That is, we write the wavefront error at pupil position |${\boldsymbol{\eta}}$| and field position |${\boldsymbol{\theta}}$| as $$\begin{eqnarray*} \psi ({\boldsymbol{\eta }};{\boldsymbol{\theta }}) &=& \sum _{n=2}^\infty \sum _{m} \sqrt{n+1}\, \psi _{nm}({\boldsymbol{\theta }}) R_n^m(\rho)\nonumber \\ && \times \left\lbrace \begin{array}{@{}c@{\hskip6pt}c}1 & m=0 \\ \sqrt{2}\, \cos m\varphi & m\gt 0 \\ \sqrt{2}\, \sin m\varphi & m\lt 0 \end{array} \right., \end{eqnarray*}$$(B1) where ρ is the radius of the pupil position normalized to 1 at the edge, and φ is the polar angle in the pupil plane, m is summed over integers with the same parity as n (both odd or both even) and |m| ≤ n (so that there are n + 1 terms in the m-sum), and |$R_n^m$| is the Zernike polynomial with normalization |$R_n^m(1)=1$|⁠. The factor of |$\sqrt{n+1}$| and (sometimes) |$\sqrt{2}$| guarantee the unit normalization of the RMS over the unit disc. If the wavefront is drifting over time, then to first order in the drift rate we may write $$\begin{eqnarray*} \psi _i({\boldsymbol{\theta }};t) = \psi _i({\boldsymbol{\theta }};t_0) + \dot{\psi }_i({\boldsymbol{\theta }}) (t-t_0), \end{eqnarray*}$$(B2) where t0 is the central epoch chosen and |$-\frac{1}{2}\Delta t \lt t-t_0 \lt \frac{1}{2}\Delta t$|⁠. Again to linear order in t − t0, the PSF that is determined by a least-squares fit with uniform weighting in time will have an expectation value that is |$G({\bf x};{\boldsymbol{\theta }};t_0)$|⁠. There is then a corresponding error in the shear in a given redshift bin zk: $$\begin{eqnarray*} c_{k,i}(t) = \sum _j \frac{\partial \gamma _{{\rm obs},i}(z_k)}{\partial \psi _j} \dot{\psi }_j({\boldsymbol{\theta }}) (t-t_0), \end{eqnarray*}$$(B3) where in this equation, k denotes a redshift bin and i denotes a component. Taking just the most strongly affected (in the sense of |c|) redshift bin to start as the reference, we see that $$\begin{eqnarray*} |c_{\rm ref}(t)| \le \left\Vert \frac{\partial \gamma _{{\rm obs,ref},i}}{\partial \psi _j} \right\Vert |\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| |t-t_0|, \end{eqnarray*}$$(B4) where ‖ ‖ denotes an operator norm (i.e. the maximum singular value of the 2 × NZern matrix). The variance of c per component (i.e. divided by 2) is $$\begin{eqnarray*} A^2 \equiv \frac{1}{2}\langle |c_{\rm ref}|^2 \rangle \le \frac{1}{2} \left[ \left\Vert \frac{\partial \gamma _{{\rm obs,ref},i}}{\partial \psi _j} \right\Vert |\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| \right]^2 \langle (t-t_0)^2\rangle ; \end{eqnarray*}$$(B5) the last expectation value is |$\frac{1}{12}\Delta t^2$| with the average taken over a uniform interval, leading to $$\begin{eqnarray*} A \le \frac{1}{\sqrt{24}} \left\Vert \frac{\partial \gamma _{{\rm obs,ref},i}}{\partial \psi _j} \right\Vert |\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| \Delta t. \end{eqnarray*}$$(B6) Thus from a requirement on A, a determination of the matrix ∂γobs,ref,i/∂ψj, and an interval of time Δt, we can set a requirement on the wavefront drift rate |$|\dot{\boldsymbol{\psi }}|$|⁠. The matrix ∂γobs,ref,i/∂ψj depends on the static aberration pattern and its determination is described below. The interval Δt for PSF fitting is a free parameter, and the wavefront drift rate requirement is tighter if Δt is increased. This must be traded against the statistical error in the PSF solution, where the target precision is easier to achieve if the time baseline Δt used in fitting the model is increased. B1.2 Sensitivity matrix From equation (B6), we see that a key step is to compute the sensitivity matrix ∂γobs,ref,i/∂ψj. Unfortunately, this matrix depends on the specific combination of static wavefront errors, because |${\boldsymbol{\gamma }}_{\rm obs,ref}$| is not a linear function of |${\boldsymbol{\psi}}$|⁠. Indeed, due to symmetries the possible form of |${\boldsymbol{\gamma }}_{\rm obs,ref}$| is restricted, with the result that ∂γobs,ref,i/∂ψj may be suppressed at zero wavefront error (⁠|${\boldsymbol{\psi }}=0$|⁠) and be much larger in the realistic case where |${\boldsymbol{\psi }}\ne 0$| (e.g. Noecker 2010). We do not know a priori what the static wavefront error will be – we have set a requirement of |$|{\boldsymbol{\psi }}|\lt 92\,$| nm, but until we have the as-built observatory we do not know how this will be distributed among the Zernikes. We address this by expanding the sensitivity matrix ∂γobs,ref,i/∂ψj to linear order in |${\boldsymbol{\psi}}$|⁠, which means we need γobs,ref,i Taylor expanded to quadratic order in |${\boldsymbol{\psi}}$| around |${\boldsymbol{\psi }}=0$|⁠. A consequence of this is that if the static wavefront error |${\boldsymbol{\psi}}$| is larger, then the sensitivity matrix is also larger and the stability requirements are tighter. Then, we search the entire space of possible wavefront errors |${\boldsymbol{\psi}}$| – bounded by the top-level requirement that |$|{\boldsymbol{\psi }}|\lt 92\,$| nm – to find the place where the operator norm is maximized. The fact that the PSF inverts (i.e. preserves ellipticity and hence spurious shear) under |${\boldsymbol{\psi }}\rightarrow -{\boldsymbol{\psi}}$| implies that |${\boldsymbol{\gamma }}_{\rm obs,ref}$| is an even function of |${\boldsymbol{\psi}}$| (this statement remains true even for an asymmetric pupil, due e.g. to the spider). For a circularly symmetric pupil (i.e. an annular pupil, so allowing for a secondary obstruction, but not accounting for the offset of the secondary obstruction when using an off-axis portion of the field, nor for the spider arms), we find the further restrictions that $$\begin{eqnarray*} \gamma _{\rm obs,ref\, 1} &=& C_{fa} \psi _{20} \psi _{22} + C_{sa} \psi _{40} \psi _{22} + C_{cc} \big(\psi _{31}^2-\psi _{3-1}^2\big) \nonumber \\ &&+\, C_{ct} (\psi _{31}\psi _{33} + \psi _{3-1}\psi _{3-3}) + ... \end{eqnarray*}$$(B7) and $$\begin{eqnarray*} \gamma _{\rm obs,ref\, 2} &=& C_{fa} \psi _{20} \psi _{2-2} + C_{sa} \psi _{40} \psi _{2-2} + 2C_{cc} \psi _{31} \psi _{3-1} \nonumber \\ &&+\, C_{ct} (\psi _{31}\psi _{3-3} - \psi _{3-1}\psi _{33}) + ...\, , \end{eqnarray*}$$(B8) where we have taken the lowest-order aberrations (focus, astigmatism, coma, trefoil, and spherical) as these dominate the wavefront stability budget. With the wavefront error vector written in this order, |${\boldsymbol{\psi }} = (\psi _{20}; \psi _{22},\psi _{2-2}; \psi _{31},\psi _{3-1}; \psi _{33},\psi _{3-3}; \psi _{40})$|⁠, we find a sensitivity matrix $$\begin{eqnarray*} {\bf M}^{\rm T} &=& \left[ \frac{\partial \gamma _{{\rm obs,ref},i}}{\partial \psi _j} \right]^{\rm T}\nonumber \\ &=& \begin{pmatrix}C_{fa}\psi _{22} & C_{fa}\psi _{2-2} \\ C_{fa}\psi _{20}+C_{sa}\psi _{40} & 0 \\ 0 & C_{fa}\psi _{20}+C_{sa}\psi _{40} \\ 2C_{cc}\psi _{31}+C_{ct}\psi _{33} & 2C_{cc}\psi _{3-1}+C_{ct}\psi _{3-3} \\ -2C_{cc}\psi _{3-1}+C_{ct}\psi _{3-3} & 2C_{cc}\psi _{31}-C_{ct}\psi _{33} \\ C_{ct}\psi _{31} & -C_{ct}\psi _{3-1} \\ C_{ct}\psi _{3-1} & C_{ct}\psi _{31} \\ C_{sa}\psi _{22} & C_{sa}\psi _{2-2} \end{pmatrix} \end{eqnarray*}$$(B9) (we show the transpose here for ease of display; the operator norm is the same). The real pupil is not circularly symmetric, however, as noted above |${\boldsymbol{\gamma }}_{\rm obs,ref}$| remains an even function of |${\boldsymbol{\psi}}$| even for a pupil of general asymmetry; thus M remains an odd function of |${\boldsymbol{\psi}}$| and the leading term is the linear term. The consequence of asymmetry of the pupil is that the coefficients in equation (B9) may be slightly different in different directions, e.g. in J-band in the z = 1.0–1.2 bin, using the full pupil for SCA #1, we find that the top row of MT is |$(8.273\psi _{22}+0.002\psi _{2-2}, 0.002\psi _{22}+8.289\psi _{2-2})\, \mu$|m−2, versus |$8.337(\psi _{22},\psi _{2-2})\, \mu$|m−2 for the circularly symmetric (annular) pupil. This corresponds to a difference in sensitivity of 0.2 per cent (8.273 versus 8.289) between the two astigmatism modes, and a 0.8 per cent (8.273 versus 8.337) change in sensitivity relative to the annular pupil case. This in principle changes our requirements by 0.8 per cent (or 0.9 per cent, which is the maximum over any of the bands and redshift bins); in practice, given the necessary margin factors, it does not make sense to track stability requirements at the |$\lt 1{{\ \rm per\ cent}}$| level. This statement about the sensitivity holds even though the annular pupil is missing some important features of the real PSF (particularly diffraction spikes). We want a limit on the maximum singular value of equation (B9), subject to a limit on |$|{\boldsymbol{\psi }}|$|⁠. To do so, let us first consider writing the singular value decomposition M = UDVT, where U is a 2 × 2 orthogonal matrix, D has two diagonal non-negative entries in non-increasing order (D11 ≥ D22) and is otherwise zeroes (and has dimension 2 × NZern), and V is NZern × NZern. Here, U is simply a rotation of the shear derivative, and due to circular symmetry can be set to the identity by rotating the entire aberration pattern. Thus, without loss of generality we can consider cases where U is the identity, and then $$\begin{eqnarray*} \Vert {\bf M} \Vert = \sqrt{ \sum _j \left(\frac{\partial \gamma _{{\rm obs,ref},1}}{\partial \psi _j} \right)^2 } = \sqrt{ {\boldsymbol{\psi }}^{\rm T} {\boldsymbol{\Lambda }}{\boldsymbol{\psi }} } \le \Vert {\boldsymbol{\Lambda }} \Vert |{\boldsymbol{\psi }}|, \end{eqnarray*}$$(B10) where we used the fact that M is a linear function of |${\boldsymbol{\psi}}$| and defined the matrix |${\boldsymbol{\Lambda}}$| to be the matrix of derivatives of the first row of M: $$\begin{eqnarray*} {\boldsymbol{\Lambda }} = \begin{pmatrix}0 & C_{fa} & 0 & 0 & 0 & 0 & 0 & 0 \\ C_{fa} & 0 & 0 & 0 & 0 & 0 & 0 & C_{sa} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 2C_{cc} & 0 & C_{ct} & 0 & 0 \\ 0 & 0 & 0 & 0 & -2C_{cc} & 0 & C_{ct} & 0 \\ 0 & 0 & 0 & C_{ct} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & C_{ct} & 0 & 0 & 0 \\ 0 & C_{sa} & 0 & 0 & 0 & 0 & 0 & 0 \end{pmatrix}, \end{eqnarray*}$$(B11) which has norm $$\begin{eqnarray*} \Vert {\boldsymbol{\Lambda }} \Vert = \max \left\lbrace \sqrt{C_{fa}^2 + C_{sa}^2}, ~|C_{cc}|+\sqrt{C_{cc}^2+C_{ct}^2} \right\rbrace . \end{eqnarray*}$$(B12) (There are both even-aberration and odd-aberration sectors of this matrix; the operator norm is determined by whichever has greater leverage on the spurious shear. In most cases, the even sector – the first term – is dominant.) We show the derived coefficients in Table B1. Table B1. The coefficients of spurious shear at |${\boldsymbol{\psi }}=0$|⁠, appearing in equations (B7–B12) for the top half of the table (wavefront drift) and equations (B15–B23) for the bottom half (wavefront jitter). Coefficients are shown for the worst (most contaminated) redshift bin. This redshift bin and the S-factor are shown in the two right-most columns for the J129 band (which is always the most sensitive to wavefront jitter because it has the shortest wavelength). The S-factor is shown for the worst angular bin, which is always the smallest scales (3.0 < log10ℓ < 3.5). The C and K coefficients are the same for the ‘even’ aberrations but different for the ‘odd’ aberrations. Band: . J129 . H158 . F184 . Units . Worst z-bin . Worst S-factor . Wavefront drift coefficients Cfa 10.60 9.22 8.06 10−6 mas−1 nm−1 2.8–3.0 0.549 Ccc 2.36 2.03 1.81 10−6 mas−1 nm−1 2.8–3.0 0.612 Cct 6.23 5.41 4.71 10−6 mas−1 nm−1 2.8–3.0 0.558 Csa 4.14 3.21 2.59 10−6 mas−1 nm−1 2.8–3.0 0.672 |$\Vert {\boldsymbol{\Lambda }}\Vert$| 11.38 9.76 8.46 10−6 mas−1 nm−1 Wavefront jitter coefficients Kθc 7.69 6.16 4.83 10−6 mas−1 nm−1 2.8–3.0 0.592 Kθt 3.29 3.78 4.08 10−6 mas−1 nm−1 1.8–2.0 0.438 Kfa 10.60 9.22 8.06 10−6 nm−2 2.8–3.0 0.549 Kcc 3.52 2.76 2.23 10−6 nm−2 2.8–3.0 0.648 Kct 5.52 4.74 4.10 10−6 nm−2 2.8–3.0 0.568 Ksa 4.14 3.21 2.59 10−6 nm−2 2.8–3.0 0.672 Kθθ 14.32 13.77 13.34 10−6 mas−2 2.8–3.0 0.504 Band: . J129 . H158 . F184 . Units . Worst z-bin . Worst S-factor . Wavefront drift coefficients Cfa 10.60 9.22 8.06 10−6 mas−1 nm−1 2.8–3.0 0.549 Ccc 2.36 2.03 1.81 10−6 mas−1 nm−1 2.8–3.0 0.612 Cct 6.23 5.41 4.71 10−6 mas−1 nm−1 2.8–3.0 0.558 Csa 4.14 3.21 2.59 10−6 mas−1 nm−1 2.8–3.0 0.672 |$\Vert {\boldsymbol{\Lambda }}\Vert$| 11.38 9.76 8.46 10−6 mas−1 nm−1 Wavefront jitter coefficients Kθc 7.69 6.16 4.83 10−6 mas−1 nm−1 2.8–3.0 0.592 Kθt 3.29 3.78 4.08 10−6 mas−1 nm−1 1.8–2.0 0.438 Kfa 10.60 9.22 8.06 10−6 nm−2 2.8–3.0 0.549 Kcc 3.52 2.76 2.23 10−6 nm−2 2.8–3.0 0.648 Kct 5.52 4.74 4.10 10−6 nm−2 2.8–3.0 0.568 Ksa 4.14 3.21 2.59 10−6 nm−2 2.8–3.0 0.672 Kθθ 14.32 13.77 13.34 10−6 mas−2 2.8–3.0 0.504 Open in new tab Table B1. The coefficients of spurious shear at |${\boldsymbol{\psi }}=0$|⁠, appearing in equations (B7–B12) for the top half of the table (wavefront drift) and equations (B15–B23) for the bottom half (wavefront jitter). Coefficients are shown for the worst (most contaminated) redshift bin. This redshift bin and the S-factor are shown in the two right-most columns for the J129 band (which is always the most sensitive to wavefront jitter because it has the shortest wavelength). The S-factor is shown for the worst angular bin, which is always the smallest scales (3.0 < log10ℓ < 3.5). The C and K coefficients are the same for the ‘even’ aberrations but different for the ‘odd’ aberrations. Band: . J129 . H158 . F184 . Units . Worst z-bin . Worst S-factor . Wavefront drift coefficients Cfa 10.60 9.22 8.06 10−6 mas−1 nm−1 2.8–3.0 0.549 Ccc 2.36 2.03 1.81 10−6 mas−1 nm−1 2.8–3.0 0.612 Cct 6.23 5.41 4.71 10−6 mas−1 nm−1 2.8–3.0 0.558 Csa 4.14 3.21 2.59 10−6 mas−1 nm−1 2.8–3.0 0.672 |$\Vert {\boldsymbol{\Lambda }}\Vert$| 11.38 9.76 8.46 10−6 mas−1 nm−1 Wavefront jitter coefficients Kθc 7.69 6.16 4.83 10−6 mas−1 nm−1 2.8–3.0 0.592 Kθt 3.29 3.78 4.08 10−6 mas−1 nm−1 1.8–2.0 0.438 Kfa 10.60 9.22 8.06 10−6 nm−2 2.8–3.0 0.549 Kcc 3.52 2.76 2.23 10−6 nm−2 2.8–3.0 0.648 Kct 5.52 4.74 4.10 10−6 nm−2 2.8–3.0 0.568 Ksa 4.14 3.21 2.59 10−6 nm−2 2.8–3.0 0.672 Kθθ 14.32 13.77 13.34 10−6 mas−2 2.8–3.0 0.504 Band: . J129 . H158 . F184 . Units . Worst z-bin . Worst S-factor . Wavefront drift coefficients Cfa 10.60 9.22 8.06 10−6 mas−1 nm−1 2.8–3.0 0.549 Ccc 2.36 2.03 1.81 10−6 mas−1 nm−1 2.8–3.0 0.612 Cct 6.23 5.41 4.71 10−6 mas−1 nm−1 2.8–3.0 0.558 Csa 4.14 3.21 2.59 10−6 mas−1 nm−1 2.8–3.0 0.672 |$\Vert {\boldsymbol{\Lambda }}\Vert$| 11.38 9.76 8.46 10−6 mas−1 nm−1 Wavefront jitter coefficients Kθc 7.69 6.16 4.83 10−6 mas−1 nm−1 2.8–3.0 0.592 Kθt 3.29 3.78 4.08 10−6 mas−1 nm−1 1.8–2.0 0.438 Kfa 10.60 9.22 8.06 10−6 nm−2 2.8–3.0 0.549 Kcc 3.52 2.76 2.23 10−6 nm−2 2.8–3.0 0.648 Kct 5.52 4.74 4.10 10−6 nm−2 2.8–3.0 0.568 Ksa 4.14 3.21 2.59 10−6 nm−2 2.8–3.0 0.672 Kθθ 14.32 13.77 13.34 10−6 mas−2 2.8–3.0 0.504 Open in new tab We are not quite done because we have not specified the redshift or scale dependence of this systematic. Since Cfa is usually dominant, we adopt its redshift dependence to determine the weights w(zi), with the last bin as the reference bin because it is the most heavily contaminated – the galaxies are smallest in that bin and the C-coefficients are largest. However, the weights w(zi) obtained from Cct (the next largest coefficient) are only slightly different. We find that in J129, H158, and F184 bands, the operator norms are |$\Vert {\boldsymbol{\Lambda }}\Vert = 1.14\times 10^{-5}$|⁠, 9.76 × 10−6, and 8.46 × 10−6 nm−2, respectively. Using the J129-band limit, which has the worst total contamination, we find a limit on the total wavefront error of |$|{\boldsymbol{\psi }}({\boldsymbol{\theta }})|\lt 92\,$| nm, we find $$\begin{eqnarray*} A \le \frac{1}{\sqrt{24}}\Vert {\boldsymbol{\Lambda }} \Vert |{\boldsymbol{\psi }}| |\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| \Delta t = 2.14\times 10^{-4}\, {\rm nm}^{-1}\times |\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| \Delta t. \nonumber\\ \end{eqnarray*}$$(B13) The current wavefront drift sub-allocation is that Δt = 1 exposure (140 s) and |$|\dot{\boldsymbol{\psi }}({\boldsymbol{\theta }})| \Delta t\lt 0.37\,$| nm, which produces a spurious shear of 7.91 × 10−5, RMS per component. However, the S-factor for the focus × astigmatism mode is 0.5488 in the worst angular bin, so the implied spurious shear is AS1/2 = 5.86 × 10−5. The requirements in Table A2 give a top-level error of 2.65 × 10−4; 4.9 per cent of the additive shear systematic error budget, in an RSS sense, is currently being taken up by wavefront drift. For a sub-allocation of 1 nm (instead of 0.37 nm), this would be 36 per cent of the additive shear systematic error budget. B2 Wavefront jitter The wavefront jitter is handled by a similar calculation to the wavefront drift. The principal difference is that we are now interested in the spurious shear from a PSF that is the superposition of many instantaneous PSFs with different wavefronts. Moreover, the PSFs can have different line-of-sight positions, so instead of simply considering the covariance matrix of the Zernike amplitudes, we must also consider the line-of-sight motion (parametrized by θx and θy). The spurious shear thus depends on the full covariance matrix of the Zernike amplitudes |${\boldsymbol{\psi}}$| and the line-of-sight motion |${\boldsymbol{\theta}}$|⁠. Of this covariance matrix, the ‘line-of-sight block’ |${\rm Cov}({\boldsymbol{\theta }},{\boldsymbol{\theta }})$| corresponds to simple image motion, and is not related to wavefront jitter. The PSF modelling procedure for Roman will explicitly allow for image motion to be fit separately in each exposure (Jurling & Content 2012), so the presence of |${\rm Cov}({\boldsymbol{\theta }},{\boldsymbol{\theta }})$| does not represent a bias in fitting the PSF. On the other hand, the blocks |${\rm Cov}({\boldsymbol{\theta }},{\boldsymbol{\psi }})$| and |${\rm Cov}({\boldsymbol{\psi }},{\boldsymbol{\psi }})$| involve the wavefront jitter, and their effects on |${\boldsymbol{\gamma }}_{\rm obs,ref}$| must be treated here. (Due to the large number of parameters, we cannot fit a full covariance matrix of all the Zernikes in each exposure.) We can then write the matrix of second derivatives: $$\begin{eqnarray*} \gamma _{{\rm obs},i}(z_k) &=& \gamma _{{\rm obs},i}(z_k)|_{\rm no~wf~jitter}+ \sum _{aj} K^{\rm LOS,WFE}_{iaj}(z_k) {\rm Cov}({\theta }_a,{\psi }_j)\nonumber \\ &&+ \frac{1}{2} \sum _{jj^{\prime }} K^{\rm WFE,WFE}_{ija}(z_k) {\rm Cov}({\psi }_j,{\psi }_{j^{\prime }}). \end{eqnarray*}$$(B14) The matrix KWFE, WFE, describing how much small high-frequency vibrations of the wavefront impact the shear, has a dependence on redshift bin zk, shear component i, and the Zernike modes j and j′. The matrix KLOS, WFE describes the effects of correlations between LOS motion and wavefront jitter. Although not strictly necessary for an analysis of wavefront jitter, we do also compute the sensitivity to line-of-sight motion, which implies a term |$K^{\rm LOS,LOS}_{iab}(z_k) {\rm Cov}({\theta }_a,{\theta }_b)$| in equation (B14). The matrix K in principle varies with the wavefront error, but since it is a second derivative it is non-zero even for zero aberrations. One option is to take this leading term (i.e. K evaluated at |${\boldsymbol{\psi }}=0$|⁠) to set requirements. Another would be to also include the linear dependences on |${\boldsymbol{\psi}}$|⁠; this would be necessary if we were to separately write requirements on the individual Zernike modes, since due to symmetries some entries in K are exactly zero in the unaberrated case, but we do not expect this to be necessary when setting overall limits on wavefront jitter (we verify this explicitly below). Following the methodology of Section B1.2, and again exploiting the symmetries of the problem and suppressing the zk index, we find that at |${\boldsymbol{\psi }}=0$|⁠, the terms involving the covariance of line-of-sight motion and wavefront jitter are $$\begin{eqnarray*} K_{1aj}^{\rm LOS,WFE} = \begin{pmatrix}0 & 0 & 0 & K_{\theta c} & 0 & K_{\theta t} & 0 & 0 \\ 0 & 0 & 0 & 0 & -K_{\theta c} & 0 & K_{\theta t} & 0 \end{pmatrix} \end{eqnarray*}$$(B15) and $$\begin{eqnarray*} K_{2aj}^{\rm LOS,WFE} = \begin{pmatrix}0 & 0 & 0 & 0 & K_{\theta c} & 0 & K_{\theta t} & 0 \\ 0 & 0 & 0 & K_{\theta c} & 0 & -K_{\theta t} & 0 & 0 \end{pmatrix}, \end{eqnarray*}$$(B16) where the two rows are a = 1, 2 and the eight columns are the low-order Zernikes. Similarly, for the wavefront jitter variance, we have $$\begin{eqnarray*} K_{1jj^{\prime }}^{\rm WFE,WFE} = \begin{pmatrix}0 & K_{fa} & 0 & 0 & 0 & 0 & 0 & 0 \\ K_{fa} & 0 & 0 & 0 & 0 & 0 & 0 & K_{sa} \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 2K_{cc} & 0 & K_{ct} & 0 & 0 \\ 0 & 0 & 0 & 0 & -2K_{cc} & 0 & K_{ct} & 0 \\ 0 & 0 & 0 & K_{ct} & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & K_{ct} & 0 & 0 & 0 \\ 0 & K_{sa} & 0 & 0 & 0 & 0 & 0 & 0 \\ \end{pmatrix} \end{eqnarray*}$$(B17) and $$\begin{eqnarray*} K_{2jj^{\prime }}^{\rm WFE,WFE} = \begin{pmatrix}0 & 0 & K_{fa} & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ K_{fa} & 0 & 0 & 0 & 0 & 0 & 0 & K_{sa} \\ 0 & 0 & 0 & 0 & 2K_{cc} & 0 & K_{ct} & 0 \\ 0 & 0 & 0 & 2K_{cc} & 0 & -K_{ct} & 0 & 0 \\ 0 & 0 & 0 & 0 & -K_{ct} & 0 & 0 & 0 \\ 0 & 0 & 0 & K_{ct} & 0 & 0 & 0 & 0 \\ 0 & 0 & K_{sa} & 0 & 0 & 0 & 0 & 0 \\ \end{pmatrix}. \end{eqnarray*}$$(B18) Finally, for the line-of-sight-motion only, we have $$\begin{eqnarray*} K_{1ab}^{\rm LOS,LOS} = \begin{pmatrix}K_{\theta \theta } & 0 \\ 0 & -K_{\theta \theta } \end{pmatrix} \end{eqnarray*}$$(B19) and $$\begin{eqnarray*} K_{2ab}^{\rm LOS,LOS} = \begin{pmatrix}0 & K_{\theta \theta } \\ K_{\theta \theta } & 0 \end{pmatrix}. \end{eqnarray*}$$(B20) Image simulations are required to determine the specific values of Kθc, Kθt, Kfa, Ksa, Kcc, and Kct. These depend on the galaxy sizes, and hence indirectly on redshift slice zk. The coefficients in the ‘worst’ redshift slice are shown in Table B1. Once again, the maximum value of the apparent shear induced by wavefront error can be determined from the eigenvalues of the K matrices, the RMS wavefront jitter, and the line-of-sight motion per axis. We note that the RMS wavefront jitter is $$\begin{eqnarray*} \sigma _{\rm wfe-jitter} = \sqrt{\sum _j {\rm Var} \psi _j}, \end{eqnarray*}$$(B21) and that covariances between WFE jitter and LOS jitter are limited by the rule that the covariance matrix be positive-definite (in particular, the correlation coefficients cannot exceed 1). This implies the limits $$\begin{eqnarray*} &&\left|\sum _{aj} K^{\rm LOS,WFE}_{1aj}(z_k) {\rm Cov}({\theta }_a,{\psi }_j)\right| \nonumber \\ &&\le \sqrt{K_{\theta c}^2 + K_{\theta t}^2}~\, \sigma _{\rm los-jitter} \sigma _{\rm wfe-jitter} \end{eqnarray*}$$(B22) and $$\begin{eqnarray*} &&\left|\frac{1}{2} \sum _{jj^{\prime }} K^{\rm WFE,WFE}_{ij{j^{\prime }}}(z_k) {\rm Cov}({\psi }_j,{\psi }_{j^{\prime }})\right| \nonumber \\ &&\le \frac{1}{2} \max \left(\sqrt{K_{fa}^2+K_{sa}^2}, |K_{cc}|+\sqrt{K_{cc}^2+K_{ct}^2} \right)\, \sigma _{\rm wfe-jitter}^2, \nonumber\\ \end{eqnarray*}$$(B23) where σlos-jitter is the RMS line-of-sight jitter per axis, i.e. we set |${\rm Cov}(\theta _a,\theta _b) = \sigma _{\rm los-jitter}^2\delta _{ab}$|⁠. (Note that only the jitter contributes here: the controlled motion of the line of sight does not correlate with the wavefront jitter since it is not in the same frequency band.) The sum of equations (B22) and (B23) represents a bound on the RMS spurious shear in the γ1 component (a similar bound applies to γ2). We also computed K for the case of the full pupil (including spider) and one realization of the static wavefront error at the centre of SCA #1.24 In this case, the sparseness pattern of K changes, and equation (B22) changes by replacing |$\sqrt{K_{\theta c}^2 + K_{\theta t}^2}$| with the maximum singular value of |${\bf K}^{\rm LOS,WFE}_{1aj}$|⁠. Over the three bands (J, H, and F184), and the 15 redshift bins, we found a maximum change of 1.9 per cent in the coefficient of equation (B22) for the γ1 component, and 3.1 per cent for the γ2 component; in all cases, the sensitivity actually went down. We therefore conclude that for the purposes of setting requirements, computing K by expansion around the unaberrated annular pupil was a sufficient approximation. The RMS spurious shear per component in J-band, weighted by the worst scaling factor S = 0.672 (which accounts for redshift dependence), is then $$\begin{eqnarray*} \gamma _{\rm rms}\sqrt{S N_{\rm ind}} &\le & 6.85\times 10^{-6}\, {\rm nm}^{-1}\, {\rm mas}^{-1}\, \sigma _{\rm los-jitter}\sigma _{\rm wfe-jitter} \nonumber \\ &&+ 4.67\times 10^{-6} \, {\rm nm}^{-2}\, \sigma _{\rm wfe-jitter}^2. \end{eqnarray*}$$(B24) This should be compared to the requirement of 2.65 × 10−4. The line-of-sight jitter is required to be σlos-jitter = 12 mas rms per axis, and the observing strategy has two passes at epochs separated by many months so we take Nind = 2. From equation (B24), we find that the entire error budget would be taken up for σwfe-jitter = 3.76 nm rms. If σwfe-jitter = 1 nm rms, then γrmsS1/2 = 6.14 × 10−5, and 5.4 per cent of the error budget is used (in an RSS sense). APPENDIX C: SIMULATION DATA ACCESS The simulated data for the fiducial run used in this paper form a 2.5 × 2.5 deg2, full-depth synthetic Roman Reference Survey in the H158-band, which is suitable for a variety of uses in testing algorithms to apply to Roman weak lensing data. The simulated data set will be available for download via a public shared Globus endpoint following publication of this paper. This endpoint directory includes the following sub-directories: images: A set of FITS images for each SCA in each pointing. meds: A set of MEDS files that contain cutouts of each object in each exposure, which do not include any neighbouring galaxy light. These MEDS files were used for the analysis in this paper. truth: A FITS catalogue of true object properties and a FITS catalogue containing information about where each object appeared in each SCA in each pointing. The true centroid of the object in SCA pixel coordinates is offset by 0.5 pixels in x and y relative to the positions recorded in the second FITS catalogue, which must be corrected when doing precision operations with the images like shape measurement. This correction is not needed if using the pre-made MEDS files. APPENDIX D: PERFORMANCE STATISTICS Each stage of the image simulation suite can be trivially parallelized. Disc I/O is not a practical issue in the image generation when running at scale across 4–5 k jobs. Within each image generation job (simulation of a complete SCA), the drawing of each object is currently also trivially parallelized across available threads. This will likely change as the simulation of the detectors becomes more realistic, though. In both modes, the image generation stage achieves about 95 per cent CPU utilization. Stages that process the image output into different data formats are less efficient, and the generation of the MEDS files is generally limited by disc I/O and remote data transfer. Typical timing (per thread), memory usage, and resulting data sizes are provided in Table D1. The total data volume is 66 GB of FITS images, 412 GB of MEDS files, and <1 GB of catalogue data, for a total of about 478 GB. These values will scale approximately linearly with the area of sky simulated and the density of objects. Table D1. Performance information for the major image simulation suite stages. Each tile is 0.013 deg2 on the sky. Each image contains about 2255 galaxies and 140 stars. For the image simulation created in this paper, with an area of 6.25 deg2, these numbers correspond to a total CPU time cost per simulation realization of about 7500 CPU hours for image generation and 9000 CPU hours total. The time required for shape measurement is expected to decrease by at least an order of magnitude, since the current measurement algorithm employed is very slow by current standards, but most of the CPU cost is still in generating the images, which is unlikely to be reduced in the future. Benchmark . Image generation . MEDS creation . Shape . . (per SCA image) . (per tile) . measurement . CPU time 180 min 10–15 min 5 s Memory 2–4 GB 1–2 GB <1 GB Data size 25 MB 1 GB <1 MB Benchmark . Image generation . MEDS creation . Shape . . (per SCA image) . (per tile) . measurement . CPU time 180 min 10–15 min 5 s Memory 2–4 GB 1–2 GB <1 GB Data size 25 MB 1 GB <1 MB Open in new tab Table D1. Performance information for the major image simulation suite stages. Each tile is 0.013 deg2 on the sky. Each image contains about 2255 galaxies and 140 stars. For the image simulation created in this paper, with an area of 6.25 deg2, these numbers correspond to a total CPU time cost per simulation realization of about 7500 CPU hours for image generation and 9000 CPU hours total. The time required for shape measurement is expected to decrease by at least an order of magnitude, since the current measurement algorithm employed is very slow by current standards, but most of the CPU cost is still in generating the images, which is unlikely to be reduced in the future. Benchmark . Image generation . MEDS creation . Shape . . (per SCA image) . (per tile) . measurement . CPU time 180 min 10–15 min 5 s Memory 2–4 GB 1–2 GB <1 GB Data size 25 MB 1 GB <1 MB Benchmark . Image generation . MEDS creation . Shape . . (per SCA image) . (per tile) . measurement . CPU time 180 min 10–15 min 5 s Memory 2–4 GB 1–2 GB <1 GB Data size 25 MB 1 GB <1 MB Open in new tab APPENDIX E: PSF MODEL APPROXIMATION For this set of image simulations, we have saved the PSF model at the position of galaxies in two resolutions: the native pixel scale and at a pixel scale that is smaller by a factor of 8 to enable unbiased measurements of the PSF model size and ellipticity, which is undersampled in the native pixel resolution. Measurements of PSF properties use this oversampled PSF model image. The motivation for choosing a pixel grid of 8 × 8 pixels and an oversampling factor of 8 was to recover the true model ellipticity (e1 and e2) and size (T) to better than 0.1 per cent. We show in Fig. E1 the fractional difference in the PSF size and shape measured with various oversampling factors relative to a ‘true’, high-resolution PSF image. This choice has no impact on the measurement of galaxy shapes as implemented in this paper. Measurements of PSF ellipticity and size are performed using an adaptive moments method (e.g. Hirata & Seljak 2003). In all cases, we compare results from a fast approximation of the PSF model used in these simulations, which is shown in Fig. 2, and not the PSF derived from the real pupil plane image. Future studies do include PSF models inferred from down-sampled versions of the full pupil plane image, which are more accurate. Figure E1. Open in new tabDownload slide The fractional error in the recovered PSF ellipticity (e1 and e2) and size (T) for various factors of PSF model pixel scale oversampling factors of 1, 2, 4, 8, and 16 relative to the native Roman pixel scale, in cutouts of size 8 × 8 native pixels. The fractional error is measured relative to a ‘true’ PSF model, which is represented by a cutout of native pixel size 64 × 64 pixels with a resolution that is oversampled by a factor of 32. Results both with (square) and without (circle) the pixel convolution are included. The points are artificially offset horizontally for clarity. Figure E1. Open in new tabDownload slide The fractional error in the recovered PSF ellipticity (e1 and e2) and size (T) for various factors of PSF model pixel scale oversampling factors of 1, 2, 4, 8, and 16 relative to the native Roman pixel scale, in cutouts of size 8 × 8 native pixels. The fractional error is measured relative to a ‘true’ PSF model, which is represented by a cutout of native pixel size 64 × 64 pixels with a resolution that is oversampled by a factor of 32. Results both with (square) and without (circle) the pixel convolution are included. The points are artificially offset horizontally for clarity. This fast approximation to the true PSF model uses six radially oriented struts to create a generic pupil plane image with the struts and a central obscuration. The resolution of this image is set such that the pixelization effects are significantly smaller than the strut width, so should be negligible. Since the six struts in the Roman pupil plane are all at slightly different angles, this leads to 12 diffraction spikes, rather than six, so visually the ‘approximate struts’ PSF (with only six spikes) is noticeably different from the correct appearance outside of the PSF core. In addition, a bug in how the PSF model rotates with the observatory roll angle was discovered prior to publication, which will be corrected in future simulations. This acts to enhance any biases observed in this suite of simulations, due to the PSF model not averaging as much as it should across multiple exposures, but otherwise does not change the primary conclusions. APPENDIX F: ERROR METRICS FOR SETTING REQUIREMENTS In this paper, we have set requirements using the error metric Z, which is based on the error in the data vector (length Ndata) relative to its covariance matrix (see Section 2.2). This is one of several possible choices. There have also been suggestions to choose an error metric more directly related to the cosmological parameters, since these rather than the data vector are the ultimate science result from the mission. There are several ways to implement this idea. Examples are: A similar error metric, r, defined in the space of cosmological parameters |${\boldsymbol{\theta}}$| (length Ncosmo): |$r = \sqrt{\Delta {\boldsymbol{\theta }} \cdot {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{-1} \Delta {\boldsymbol{\theta }}}$|⁠, where |${\boldsymbol{\Sigma }}_{\boldsymbol{\theta}}$| is the cosmological parameter covariance matrix, and |$\Delta {\boldsymbol{\theta}}$| is the bias in the cosmological parameters. This metric asks not about the data vector, but whether the bias in the cosmological parameters is within the 1σ, 2σ, etc. error ellipsoid. If this is done, one presumably fits not just the cosmological parameters to the data vector, but also a set of nuisance parameters |${\boldsymbol{\nu}}$| (length Nnuis). One might use a similar error metric, but only care about some of the cosmological parameters – e.g. the LSST DESC SRD uses the space (w0, wa) (The LSST Dark Energy Science Collaboration 2018), and Massey et al. (2013) consider the bias on the single parameter w (in which case the error metric is ‘bias divided by sigma’ or |$\sqrt{r}$|⁠). This is mathematically the same problem as using Zcosmo, but one treats other parameters such as Ωm or σ8 as nuisance parameters. Euclid Collaboration (2020) present ‘bias over sigma’ for w0 and wa separately, which can be thought of as two versions of the r metric but each in a different 1D space. One could instead treat a systematic error in the data vector as an effective increase in its covariance matrix |${\boldsymbol{\Sigma}}$|⁠: |${\boldsymbol{\Sigma }} \rightarrow {\boldsymbol{\Sigma }} + \Delta {\bf C}\, \Delta {\bf C}^{\rm T}$|⁠, where ΔC is the bias in the data vector. (One thinks of this as a ‘systematics contribution to the data covariance.’) Then one can ask by what factor the ellipsoid volume in cosmological parameter space is increased: |$f_{\rm deg} = \sqrt{|{\boldsymbol{\Sigma }}_{{\boldsymbol{\theta }},\rm new}|/|{\boldsymbol{\Sigma }}_{{\boldsymbol{\theta }},\rm old}|}$|⁠. Note that if the cosmological parameter space is (w0, wa) (and other parameters are taken as nuisance parameters), this is equivalent to the degradation factor in Dark Energy Task Force Figure of Merit (Albrecht et al. 2006). The error metrics Z and r satisfy Z, r ≥ 0, whereas fdeg ≥ 1, with equality holding for no degradation. Our objective in this appendix is to explore the mathematical properties of these different metrics (Z, r, and fdeg), and their robustness under different circumstances. Do these error metrics add linearly, by RSS addition, or in some other way? What happens when we introduce a new nuisance parameter (which may be astrophysical)? What if we shorten the data vector with a scale cut (which is really a special case of adding nuisance parameters)? What if we combine Roman weak lensing with an external data set? The Z metric – while it is the most conservative and leads to the tightest requirements – has the advantage of obeying straightforward rules under these operations, which make it well suited to error budgeting when these tighter requirements can be met. The other metrics might still be considered as a fall-back option under some circumstances if the Z metric cannot reasonably be met, or at the requirements verification stage if one is focused on a particular set of analyses (e.g. Euclid Collaboration 2020). We note that the fdeg error metric treats the systematic error in a probabilistic sense, whereas the others as written treat the error as deterministic but unknown. If we take the probabilistic point of view instead for Z or r, we should instead consider the RMS error metrics, Zrms = 〈Z2〉1/2 or rrms = 〈r2〉1/2. This approach is particularly useful in error budgeting when adding systematic contributions that are ‘independent’ (indeed, only in the probabilistic point of view does considering systematics to be ‘independent’ make sense). The major results from this appendix are shown in Table F1. Some aspects of these error metrics have been considered before. Section 4.1 of Massey et al. (2013) discusses at length the deterministic versus probabilistic interpretation of r. Appendix B2 of LSST Dark Energy Science Collaboration (2018) has a discussion of the comparison of Z and r. Table F1. The rules governing the systematic error metrics considered here: Z, which is the size of the systematic error relative to the statistical error in the data vector; r, which is the size of the systematic error relative to the statistical error in the cosmological parameter vector; and fdeg, which is the increase in volume in parameter space when the systematic error is included in the covariance matrix. For Z and r, one may consider the systematic error in either a deterministic sense or a probabilistic sense (in which case we take the RMS). Table entries include ‘No rule’ if there is no rule for addition, ‘N/A’ for not applicable, ‘Triangle ineq.’ for triangle inequality addition, ‘RSS’ for RSS addition, and ‘≤RSS’ when RSS addition provides an upper bound. Error metric . Data vector . Cosmological parameter . Volume ratio . . space norms . space norms . . . Z . Zrms . r . rrms . fdeg . Relations Z ≤ r Zrms ≤ rrms |$f_{\rm deg} \le \left(1 + \frac{r^2}{N_{\rm cosmo}} \right)^{N_{\rm cosmo}/2}$| Addition of errors (deterministic) Triangle ineq. N/A Triangle ineq. N/A No rule Addition of errors (probabilistic, N/A RSS N/A RSS No rule independent) Adding nuisance parameter No effect No effect No rule No rule No rule Combining independent data sets (A + B): – general Triangle ineq. ≤RSS No rule No rule No rule – if no common nuisance pars. Triangle ineq. ≤RSS Triangle ineq. ≤RSS |$f_{\rm deg}^{\rm (A+B)} \le f_{\rm deg}^{\rm (A)} f_{\rm deg}^{\rm (B)}$| (equality only if |$f_{\rm deg}^{\rm (A)} = f_{\rm deg}^{\rm (B)}=1$|⁠) Error metric . Data vector . Cosmological parameter . Volume ratio . . space norms . space norms . . . Z . Zrms . r . rrms . fdeg . Relations Z ≤ r Zrms ≤ rrms |$f_{\rm deg} \le \left(1 + \frac{r^2}{N_{\rm cosmo}} \right)^{N_{\rm cosmo}/2}$| Addition of errors (deterministic) Triangle ineq. N/A Triangle ineq. N/A No rule Addition of errors (probabilistic, N/A RSS N/A RSS No rule independent) Adding nuisance parameter No effect No effect No rule No rule No rule Combining independent data sets (A + B): – general Triangle ineq. ≤RSS No rule No rule No rule – if no common nuisance pars. Triangle ineq. ≤RSS Triangle ineq. ≤RSS |$f_{\rm deg}^{\rm (A+B)} \le f_{\rm deg}^{\rm (A)} f_{\rm deg}^{\rm (B)}$| (equality only if |$f_{\rm deg}^{\rm (A)} = f_{\rm deg}^{\rm (B)}=1$|⁠) Open in new tab Table F1. The rules governing the systematic error metrics considered here: Z, which is the size of the systematic error relative to the statistical error in the data vector; r, which is the size of the systematic error relative to the statistical error in the cosmological parameter vector; and fdeg, which is the increase in volume in parameter space when the systematic error is included in the covariance matrix. For Z and r, one may consider the systematic error in either a deterministic sense or a probabilistic sense (in which case we take the RMS). Table entries include ‘No rule’ if there is no rule for addition, ‘N/A’ for not applicable, ‘Triangle ineq.’ for triangle inequality addition, ‘RSS’ for RSS addition, and ‘≤RSS’ when RSS addition provides an upper bound. Error metric . Data vector . Cosmological parameter . Volume ratio . . space norms . space norms . . . Z . Zrms . r . rrms . fdeg . Relations Z ≤ r Zrms ≤ rrms |$f_{\rm deg} \le \left(1 + \frac{r^2}{N_{\rm cosmo}} \right)^{N_{\rm cosmo}/2}$| Addition of errors (deterministic) Triangle ineq. N/A Triangle ineq. N/A No rule Addition of errors (probabilistic, N/A RSS N/A RSS No rule independent) Adding nuisance parameter No effect No effect No rule No rule No rule Combining independent data sets (A + B): – general Triangle ineq. ≤RSS No rule No rule No rule – if no common nuisance pars. Triangle ineq. ≤RSS Triangle ineq. ≤RSS |$f_{\rm deg}^{\rm (A+B)} \le f_{\rm deg}^{\rm (A)} f_{\rm deg}^{\rm (B)}$| (equality only if |$f_{\rm deg}^{\rm (A)} = f_{\rm deg}^{\rm (B)}=1$|⁠) Error metric . Data vector . Cosmological parameter . Volume ratio . . space norms . space norms . . . Z . Zrms . r . rrms . fdeg . Relations Z ≤ r Zrms ≤ rrms |$f_{\rm deg} \le \left(1 + \frac{r^2}{N_{\rm cosmo}} \right)^{N_{\rm cosmo}/2}$| Addition of errors (deterministic) Triangle ineq. N/A Triangle ineq. N/A No rule Addition of errors (probabilistic, N/A RSS N/A RSS No rule independent) Adding nuisance parameter No effect No effect No rule No rule No rule Combining independent data sets (A + B): – general Triangle ineq. ≤RSS No rule No rule No rule – if no common nuisance pars. Triangle ineq. ≤RSS Triangle ineq. ≤RSS |$f_{\rm deg}^{\rm (A+B)} \le f_{\rm deg}^{\rm (A)} f_{\rm deg}^{\rm (B)}$| (equality only if |$f_{\rm deg}^{\rm (A)} = f_{\rm deg}^{\rm (B)}=1$|⁠) Open in new tab In this appendix, we use the symbol A≻0 to denote that the symmetric matrix A is positive definite (or A⪰0 to indicate semipositive definite) and A≻B to denote that A − B is positive definite. We limit our analysis to the Fisher matrix approximation, where the covariance matrix does not depend on the data vector, and the derivatives of the theory data vector with respect to the parameters are constants. We make extensive use of the fact that all of these error metrics are invariant under general invertible linear transformations of the data vector, GL(Ndata), and of the cosmological parameter space, GL(Ncosmo). Note that these are general linear transformations, not just rotations. Many results here are easiest to understand and prove using particular choices of basis. For example, general linear transformations allow one to use a basis where |${\boldsymbol{\Sigma}}$| is the Ndata × Ndata identity matrix |${\mathbb {I}}_{N_{\rm data}}$|⁠. F1 Relations among the error metrics Here, we show that Z is in some sense the ‘most conservative’ of the error metrics, followed by r, and then fdeg. This is essentially because Z counts all systematic errors in the data vector, whereas r counts only those that have a projection on to the cosmological parameters. Finally, fdeg treats systematic errors as a contribution to the covariance and hence allows them to be marginalized out. F1.1 Comparison of Z and r Let’s suppose that the cosmological parameters that we fit, |${\boldsymbol{\theta}}$|⁠, have some dependence on the data vector C: there is a Jacobian, |${\bf R} = \partial {\boldsymbol{\theta }}/\partial {\bf C}$| (matrix dimension: Ncosmo × Ndata, with Ncosmo ≤ Ndata). Then $$\begin{eqnarray*} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }} = {\bf R}{\boldsymbol{\Sigma }}{\bf R}^{\rm T}, \end{eqnarray*}$$(F1) and so $$\begin{eqnarray*} r^2 = \Delta {\bf C}^{\rm T}{\bf R}^{\rm T} ({\bf R}{\boldsymbol{\Sigma }}{\bf R}^{\rm T})^{-1} {\bf R}\Delta {\bf C}. \end{eqnarray*}$$(F2) If we go to a basis where |${\boldsymbol{\Sigma }} = {\mathbb {I}}_{N_{\rm data}}$|⁠, and do a singular value decomposition of R = UDVT, where U and V are orthogonal and D has Ncosmo diagonal entries, then we find |$r^2 = \sum _{i=1}^{N_{\rm cosmo}} ({\bf V}^{\rm T}\Delta {\bf C})_i^2$|⁠, whereas |$Z^2 = \sum _{i=1}^{N_{\rm data}} ({\bf V}^{\rm T}\Delta {\bf C})_i^2$|⁠. This shows that $$\begin{eqnarray*} Z \le r. \end{eqnarray*}$$(F3) It follows that this relation holds in the probabilistic sense as well, Zrms ≤ rrms. F1.2 Comparison of Zcosmo and fdeg We now consider fdeg. First, let’s consider what happens without nuisance parameters. If J is the Ndata × Ncosmo matrix of partial derivatives of the theory data vector with respect to the parameters |$\partial {\bf C}_{\rm th}/\partial {\boldsymbol{\theta}}$|⁠, with a partition, then standard Fisher matrix formulae tell us that the covariance matrix of the cosmological + nuisance parameters, |${\boldsymbol{\Sigma }}_{\boldsymbol{\theta}}$| is given by $$\begin{eqnarray*} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }} = ({\bf J}^{\rm T}{\boldsymbol{\Sigma }}^{-1}{\bf J})^{-1}. \end{eqnarray*}$$(F4) If |${\boldsymbol{\Sigma}}$| is increased by the addition of a systematics contribution, |${\boldsymbol{\Sigma }} \rightarrow {\boldsymbol{\Sigma }} +\Delta {\boldsymbol{\Sigma}}$|⁠, then $$\begin{eqnarray*} f_{\rm deg} = \sqrt{ \frac{ |{\bf J}^{\rm T}{\boldsymbol{\Sigma }}^{-1}{\bf J} |}{| {\bf J}^{\rm T}({\boldsymbol{\Sigma }}+\Delta {\boldsymbol{\Sigma }})^{-1}{\bf J} | } }. \end{eqnarray*}$$(F5) This equation simplifies if we do a general linear transformation on the data vector to make |${\boldsymbol{\Sigma }} = {\mathbb {I}}_{N_{\rm data}}$|⁠. This does not uniquely define the choice of basis for the data vector space; we may further do a singular value decomposition of J = UDVT, and then do a rotation in data vector space to set |${\bf V} = {\mathbb {I}}_{N_{\rm data}}$|⁠, a rotation in cosmological parameter space to set |${\bf U} = {\mathbb {I}}_{N_{\rm data}}$|⁠, and a rescaling of the ith basis vector in cosmological parameter space by Dii to set D to 1’s down the diagonal. This makes J equal to 1’s down the diagonal and 0’s elsewhere. In this basis, we have $$\begin{eqnarray*} f_{\rm deg} = \left\lbrace \det [({\mathbb {I}}_{N_{\rm data}} + \Delta {\boldsymbol{\Sigma }})^{-1}]_{N_{\rm cosmo}\, \rm block} \right\rbrace ^{-1/2}, \end{eqnarray*}$$(F6) where the subscript ‘|$N_{\rm cosmo}\, \rm block$|’ means that we take the upper-left Ncosmo × Ncosmo block of the Ndata × Ndata matrix. Now we want a bound on fdeg in terms of r. In the basis we have chosen here, the error on the parameters is related to the error in the data vector by $$\begin{eqnarray*} \Delta {\boldsymbol{\theta }} = ({\bf J}^{\rm T}{\boldsymbol{\Sigma }}^{-1}{\bf J})^{-1} {\bf J}^{\rm T}{\boldsymbol{\Sigma }}^{-1}\Delta {\bf C}, \end{eqnarray*}$$(F7) so Δθi = ΔCi for 1 ≤ i ≤ Ncosmo and then $$\begin{eqnarray*} r^2 = \sum _{i=1}^{N_{\rm cosmo}} \Delta \theta _i^2 = \sum _{i=1}^{N_{\rm cosmo}} \Delta C_i^2 = {\rm Tr} [\Delta {\boldsymbol{\Sigma }}]_{N_{\rm cosmo}\, \rm block}. \end{eqnarray*}$$(F8) (Every step in this equality is valid in either a deterministic or probabilistic sense; in the latter case, the left-hand side becomes |$r_{\rm rms}^2$|⁠.) It follows trivially from the block inversion formula that for a positive definite matrix, the determinant of the block of an inverse is greater than or equal to the determinant of the inverse of the block. Therefore, if |$\lbrace \lambda _i\rbrace _{i=1}^{N_{\rm cosmo}}$| are the eigenvalues of |$[\Delta {\boldsymbol{\Sigma }}]_{N_{\rm cosmo}\, \rm block}$|⁠, then $$\begin{eqnarray*} f_{\rm deg} \le \sqrt{ \prod _{i=1}^{N_{\rm cosmo}} (1+\lambda _i)} ~~~{\rm and} ~~ r^2 = \sum _{i=1}^{N_{\rm cosmo}} \lambda _i. \end{eqnarray*}$$(F9) We may write the formula for fdeg in the form $$\begin{eqnarray*} 2\ln f_{\rm deg} = \sum _{i=1}^{N_{\rm cosmo}} \ln (1+\lambda _i), \end{eqnarray*}$$(F10) where the last term has a strictly negative second derivative. This means that for a fixed |$\sum _{i=1}^{N_{\rm cosmo}} \lambda _i = r^2$|⁠, the maximum value of 2ln fdeg is obtained when all of the λi are the same: λi = r2/Ncosmo. Therefore, we have $$\begin{eqnarray*} f_{\rm deg} \le \left(1 + \frac{r^2}{N_{\rm cosmo}} \right)^{N_{\rm cosmo}/2}. \end{eqnarray*}$$(F11) This relation still applies if we have nuisance parameters, since inclusion of a nuisance parameter can be thought of as a modification to |${\boldsymbol{\Sigma}}$|⁠: |${\boldsymbol{\Sigma }} \rightarrow {\boldsymbol{\Sigma }} + \sigma _p^2 (\partial {\bf C}_{\rm th}/\partial p)(\partial {\bf C}_{\rm th}/\partial p)^{\rm T}$|⁠. F2 Addition of error terms A common problem in error budgeting is addition of error terms – when there are two or more sources of error that must be included in the budget, how should they be added? And what can we say when those two sources of error are independent? This problem is simplest for Z and r, because they are the ordinary vector lengths of ΔC in |${\mathbb {R}}^{N_{\rm data}}$| and |$\Delta {\boldsymbol{\theta}}$| in |${\mathbb {R}}^{N_{\rm cosmo}}$|⁠, respectively, if we use the bases where |${\boldsymbol{\Sigma}}$| and |${\boldsymbol{\Sigma }}_{\boldsymbol{\theta}}$| are the identity. Therefore, they satisfy the ‘usual’ rules of error addition: If we take the deterministic point of view, then Z and r obey the triangle inequality: if there is one source of error (‘A’) that produces a bias in the data vector |$\Delta {\bf C}^{(\rm A)}$| and another that produces a bias |$\Delta {\bf C}^{(\rm B)}$|⁠, then always the error metrics satisfy |$Z^{(\rm A+B)} \le Z^{(\rm A)} + Z^{(\rm B)}$|⁠. If we take the probabilistic point of view, then Zrms and rrms obey RSS addition with two independent contributions A and B: |$Z^{(\rm A+B)\, 2}_{\rm rms} = Z^{(\rm A)\, 2}_{\rm rms} + Z^{(\rm B)\, 2}_{\rm rms}$|⁠. In contrast, fdeg does not obey any simple error addition rule. A simple counterexample to any proposed inequality can be found with Ncosmo = 1, Ndata = 2, and matrices in the language of Section F1.2: $$\begin{eqnarray*} {\bf J} = \begin{pmatrix}1 \\ 0 \end{pmatrix} ~~{\rm and}~~ {\boldsymbol{\Sigma }} = \begin{pmatrix}1 & 0 \\ 0 & 1 \end{pmatrix}, \end{eqnarray*}$$(F12) with systematic error contributions A and B: $$\begin{eqnarray*} \Delta {\bf C}^{(\rm A)} = \begin{pmatrix}\alpha \\ \beta \end{pmatrix} ~~{\rm and}~~ \Delta {\bf C}^{(\rm B)} = \begin{pmatrix}\alpha \\ -\beta \end{pmatrix}. \end{eqnarray*}$$(F13) We take these contributions to be additions to the data covariance matrix: |${\boldsymbol{\Sigma}}$| gets an extra contribution of |$\Delta {\bf C}^{(\rm A)}\Delta {\bf C}^{(\rm A)\, T}$|⁠, |$\Delta {\bf C}^{(\rm B)}\Delta {\bf C}^{(\rm B)\, T}$|⁠, or |$\Delta {\bf C}^{(\rm A)}\Delta {\bf C}^{(\rm A)\, T} + \Delta {\bf C}^{(\rm B)}\Delta {\bf C}^{(\rm B)\, T}$|⁠. Treating the combination of A and B by adding their contributions to the covariance matrix is a direct consequence of supposing them to be independent. However, in this case one can evaluate equation (F5) and find the degradation factors: $$\begin{eqnarray*} f_{\rm deg}^{(\rm A)} = f_{\rm deg}^{(\rm B)} = \sqrt{\frac{ 1 + \alpha ^2 + \beta ^2 }{ 1 + \beta ^2 }}, \end{eqnarray*}$$(F14) whereas $$\begin{eqnarray*} f_{\rm deg}^{(\rm A+B)} = \sqrt{1 + 2\alpha ^2}. \end{eqnarray*}$$(F15) It is clear that in this case, by making α large enough, we can make |$f_{\rm deg}^{(\rm A+B)}$| as large as we want; but then at fixed α, by making β large enough, we can make |$f_{\rm deg}^{(\rm A)}$| and |$f_{\rm deg}^{(\rm B)}$| as close to 1 as we want. In other words, using the fdeg metric, two systematic error contributions A and B, which are individually arbitrarily small can combine to make an arbitrarily large contribution A + B. This means that if one builds an error budget using fdeg as a metric, one must repeat the full Fisher matrix analysis each time a new term is added and cannot rely on each contribution to the error budget staying within an ‘allocation’ (e.g. stating fdeg < 1.02 or |$\lt 2{{\ \rm per\ cent}}$| loss of figure of merit for a particular systematic). Keeping this level of coordination across different contributions to error budgets is a challenge in a large interdisciplinary team, but it can be considered particularly in cases where error budgets based on r or Z become cost, schedule, or engineering risk drivers. F3 Introduction of new nuisance parameters We now consider what happens when new nuisance parameters are added. This is simplest in the case of Z: since it is built entirely in data vector space, Z does not depend on how the data vector is used and is unaffected when nuisance parameters are added. It turns out that no simple rule occurs for r. To see this, let’s consider the simple case of 1 cosmological parameter θ and 1 nuisance parameter ν. We also approximate the theory data vector as a linear function of the parameters over the range of interest, i.e. |$\partial {\bf C}_{\rm th}/\partial ({\boldsymbol{\theta }},{\boldsymbol{\nu }})$| constant. The χ2 surface for these parameters will generally be of the form $$\begin{eqnarray*} \chi ^2 = \begin{pmatrix}\theta & \nu \end{pmatrix} {\bf A} \begin{pmatrix}\theta \\ \nu \end{pmatrix} - 2 (b_\theta \theta + b_\nu \nu) + c, \end{eqnarray*}$$(F16) where A is the 2 × 2 inverse covariance matrix, and b is a vector. A change in the data vector will lead to a change in |$\Delta {\boldsymbol{b}}$|⁠, the slope (it does not affect A if the derivatives are constant, and we do not care about Δc). The bias in the parameters if we include the nuisance parameter is $$\begin{eqnarray*} \begin{pmatrix}\theta ^{(+)} \\ \nu ^{(+)} \end{pmatrix} = {\bf A}^{-1}\Delta {\boldsymbol{b}}, \end{eqnarray*}$$(F17) so $$\begin{eqnarray*} \Delta \theta ^{(+)} = \frac{A_{\nu \nu }\, \Delta b_\theta + A_{\theta \nu }\, \Delta b_\nu }{A_{\theta \theta } A_{\nu \nu } - A_{\theta \nu }^2}. \end{eqnarray*}$$(F18) If the nuisance parameter is not added, then we look at the χ2 surface with ν fixed to 0, in which case the bias in θ is $$\begin{eqnarray*} \Delta \theta ^{(-)} = \frac{\Delta b_\theta }{A_{\theta \theta }}. \end{eqnarray*}$$(F19) The corresponding values of r2 = Δθ2/Var(θ) are: $$\begin{eqnarray*} r^{(+)}{^2} = \frac{(A_{\nu \nu }\, \Delta b_\theta + A_{\theta \nu }\, \Delta b_\nu)^2}{A_{\nu \nu }(A_{\theta \theta } A_{\nu \nu } - A_{\theta \nu }^2)} \end{eqnarray*}$$(F20) and $$\begin{eqnarray*} r^{(-)}{^2} = \frac{(\Delta b_\theta)^2}{A_{\theta \theta }}. \end{eqnarray*}$$(F21) It is now apparent that there is no inequality relating r(+) to r(−): as long as Aθν ≠ 0, any value of the ordered pair |$(\Delta b_\theta , A_{\nu \nu }\, \Delta b_\theta + A_{\theta \nu }\, \Delta b_\nu)\in {\mathbb {R}}^2$| is possible, and therefore knowledge of r(−) by itself provides no constraint on r(+) and vice versa. So the r error metric may go up, go down, or stay the same when nuisance parameters are added. Finally, we show that there is also no inequality relating |$f_{\rm deg}^{(+)}$| (with a nuisance parameter) to |$f_{\rm deg}^{(-)}$| (without it). Let us consider the same example above with 1 cosmological parameter and a single nuisance parameter. The degradation factor with and without the extra nuisance parameter is $$\begin{eqnarray*} f_{\rm deg}^{(+)} = \sqrt{ \frac{ (A_{\theta \theta }^{\rm (old)} A_{\nu \nu }^{\rm (old)} - [A_{\theta \nu }^{\rm (old)}]^2)/A_{\nu \nu }^{\rm (old)} }{ (A_{\theta \theta }^{\rm (new)} A_{\nu \nu }^{\rm (new)} - [A_{\theta \nu }^{\rm (new)}]^2)/A_{\nu \nu }^{\rm (new)} }} \end{eqnarray*}$$(F22) and $$\begin{eqnarray*} f_{\rm deg}^{(-)} = \sqrt{ \frac{A_{\theta \theta }^{\rm (old)}}{A_{\theta \theta }^{\rm (new)}} }, \end{eqnarray*}$$(F23) where A(old) or A(new) is the old (without degradation) or new (with degradation) inverse covariance matrix. These matrices are given by $$\begin{eqnarray*} {\bf A}^{\rm (old)} = {\bf J}{\boldsymbol{\Sigma }}^{-1}{\bf J}^{\rm T} ~~{\rm and}~~ {\bf A}^{\rm (new)} = {\bf J}({\boldsymbol{\Sigma }} + \Delta {\bf C}\, \Delta {\bf C}^{\rm T})^{-1}{\bf J}^{\rm T}, \end{eqnarray*}$$(F24) where J is the Ndata × 2 matrix of partial derivatives of the data vector. For the sake of constructing an example that shows there is no inequality relating |$f_{\rm deg}^{(+)}$| to |$f_{\rm deg}^{(-)}$|⁠, we take Ndata = 2 and J to be the identity. If we take |${\boldsymbol{\Sigma}}$| to be the 2 × 2 identity and ΔC = (α β)T, then $$\begin{eqnarray*} f_{\rm deg}^{(+)} = \sqrt{1+\alpha ^2} ~~{\rm and}~~ f_{\rm deg}^{(-)} = \sqrt{\frac{1+\alpha ^2+\beta ^2}{1+\beta ^2}}. \end{eqnarray*}$$(F25) We can see that any configuration with |$1 \le f_{\rm deg}^{(-)} \le f_{\rm deg}^{(+)}$| can be obtained by choosing α (to get the desired |$f_{\rm deg}^{(+)}$|⁠) and then β (to get the desired |$f_{\rm deg}^{(-)}$|⁠). However, if alternatively we take $$\begin{eqnarray*} {\boldsymbol{\Sigma }} = \begin{pmatrix}1 & \rho \\ \rho & 1 \end{pmatrix} ~~~{\rm and}~~~ \Delta {\bf C} = \begin{pmatrix}\gamma \\ 0 \end{pmatrix} \end{eqnarray*}$$(F26) with |ρ| < 1, then $$\begin{eqnarray*} f_{\rm deg}^{(+)} = \sqrt{1+\gamma ^2} ~~{\rm and}~~ f_{\rm deg}^{(-)} = \sqrt{\frac{1-\rho ^2+\gamma ^2}{1-\rho ^2}}. \end{eqnarray*}$$(F27) This time, we can choose any configuration with |$1\le f_{\rm deg}^{(+)}\le f_{\rm deg}^{(-)}$| by choosing γ (to get the desired |$f_{\rm deg}^{(+)}$|⁠) and then ρ (to get the desired |$f_{\rm deg}^{(-)}$|⁠). Between these two examples, we see that there is no relation between |$f_{\rm deg}^{(+)}$| and |$f_{\rm deg}^{(-)}$|⁠: any behaviour is possible when we add a nuisance parameter. F4 Combination with external data sets Finally, we consider combinations with an external data set. In the context of Roman weak lensing, this might be another experiment (e.g. cosmic microwave background) or within Roman (e.g. the supernova survey). We assume that the data vector for this external data set is independent. The Z metric undergoes RSS addition between data sets, i.e. if we have data sets A and B, then Z2(A + B) = Z2(A) + Z2(B), because the covariance matrix |${\boldsymbol{\Sigma}}$| of the combined data set A + B is block diagonal. This is true both for the deterministic version of the metric and the probabilistic version, Zrms. What happens to the r and fdeg metrics is more complicated. If there are shared nuisance parameters between the experiments, then in general there is no rule on how Zcosmo and fdeg change when an external data set is added, because one could consider the limiting case where the external data set effectively fixes a nuisance parameter, and then one has the situation described in Section F3. We note that if fdeg is defined in terms of (w0, wa) space, and other cosmological parameters such as Ωm and H0 are treated as nuisance parameters, then almost all practical cases of combined data sets will have shared nuisance parameters. If there are no shared nuisance parameters, then it is possible to say more about r and fdeg. Let’s consider r first. When two data sets A and B are combined, and we have Gaussian likelihoods, the bias in parameters in the combined data set is $$\begin{eqnarray*} \Delta {\boldsymbol{\theta }}^{(\rm A+B)} = [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} +{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1} [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)} \Delta {\boldsymbol{\theta }}^{(\rm A)} + {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)} \Delta {\boldsymbol{\theta }}^{(\rm B)}] \end{eqnarray*}$$(F28) and the covariance matrix is |$[{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} + {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1}$|⁠. Then $$\begin{eqnarray*} r^{\rm (A+B)}{^2}&=& \Delta {\boldsymbol{\theta }}^{(\rm A)}{^{\rm T}} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)} [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} +{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)} \Delta {\boldsymbol{\theta }}^{(\rm A)} \nonumber \\ && + \,\Delta {\boldsymbol{\theta }}^{(\rm B)}{^{\rm T}} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)} [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} +{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)} \Delta {\boldsymbol{\theta }}^{(\rm B)} \nonumber \\ && + \,2 \Delta {\boldsymbol{\theta }}^{(\rm A)}{^{\rm T}} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)} [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} +{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)} \Delta {\boldsymbol{\theta }}^{(\rm B)}. \hskip-10pt \end{eqnarray*}$$(F29) Since |${\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)} \succeq [{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}} +{\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (B)}{^{-1}}]^{-1}$|⁠, we know that the first term is |$\le \Delta {\boldsymbol{\theta }}^{(\rm A)}{^{\rm T}} {\boldsymbol{\Sigma }}_{\boldsymbol{\theta }}^{\rm (A)}{^{-1}}\Delta {\boldsymbol{\theta }}^{(\rm A)} = r^{\rm (A)}{^2}$|⁠. Similarly, the second term is ≤r(B)2. The Cauchy–Schwarz inequality then shows that the third term is less than or equal to twice the geometric mean of the first two, i.e. ≤2r(A)r(B). So it follows that r obeys linear addition: $$\begin{eqnarray*} r^{\rm (A+B)} \le r^{\rm (A)} + r^{\rm (B)}. \end{eqnarray*}$$(F30) If one treats the systematic errors in the two data sets as independent in the probabilistic sense, then when taking an average the third term drops out and rrms obeys RSS addition (with an inequality): $$\begin{eqnarray*} r_{\rm rms}^{\rm (A+B)}{^2} \le r_{\rm rms}^{\rm (A)}{^2} + r_{\rm rms}^{\rm (B)}{^2}. \end{eqnarray*}$$(F31) We may also consider fdeg; for two independent probes A and B and no common nuisance parameters, we have $$\begin{eqnarray*} f_{\rm deg}^{\rm (A+B)} = \sqrt{ \frac{\det ({\bf F}^{\rm (A,old)} + {\bf F}^{\rm (B,old)})}{\det ({\bf F}^{\rm (A,new)} + {\bf F}^{\rm (B,new)})} }, \end{eqnarray*}$$(F32) where ‘old’ and ‘new’ indicate the Fisher matrices without and with the systematic, respectively (and after the nuisance parameters have been marginalized out). We can understand what happens if we turn this into an integral: $$\begin{eqnarray*} \ln f_{\rm deg}^{\rm (A+B)} = \frac{1}{2} \int _0^1 \frac{d}{\mathrm{ d} x} \ln \det [{\bf F}^{\rm (A)}(x) + {\bf F}^{\rm (B)}(x)]\, \mathrm{ d} x, \end{eqnarray*}$$(F33) with x = 0 corresponding to the ‘new’ Fisher matrix, x = 1 corresponding to the old: $$\begin{eqnarray*} {\bf F}^{\rm (A)}(x) = x {\bf F}^{\rm (A,old)} + (1-x) {\bf F}^{\rm (A,new)}. \end{eqnarray*}$$(F34) Since |${\bf F}^{\rm(A,new)}\preceq {\bf F}^{\rm(A,old)}$|⁠, we have dF(A)(x)/dx⪰0. The integral in equation (F33) then represents the degradation of parameter constraint volume as we introduce the systematic. Then the derivative of a log determinant satisfies $$\begin{eqnarray*} \ln f_{\rm deg}^{\rm (A+B)} &=& \frac{1}{2} \int _0^1 {\rm Tr} \Bigl \lbrace [{\bf F}^{\rm (A)}(x) + {\bf F}^{\rm (B)}(x)]^{-1} \nonumber \\ && \times \left[ \frac{\mathrm{ d}{\bf F}^{\rm (A)}(x)}{\mathrm{ d}x} + \frac{\mathrm{ d}{\bf F}^{\rm (B)}(x)}{\mathrm{ d}x} \right] \Bigr \rbrace \, \mathrm{ d} x. \end{eqnarray*}$$(F35) Now for any A⪰0 and B≻D≻0, we may write |${\bf A} = \sum _{i=1}^{N_{\rm cosmo}} {\bf uu}^{\rm T}$| and then |${\rm Tr}({\bf BA}) = \sum _{i=1}^{N_{\rm cosmo}} {\bf u}^{\rm T}{\bf Bu}$|⁠. It follows that Tr(BA) ≥ Tr(DA) (with equality only for A = 0). Applying this to [F(A)(x) + F(B)(x)]−1≺[F(A)(x)]−1 and similarly for B, we show that $$\begin{eqnarray*} \ln f_{\rm deg}^{\rm (A+B)} &\le & \frac{1}{2} \int _0^1 {\rm Tr} \Bigl \lbrace [{\bf F}^{\rm (A)}(x)]^{-1}\frac{\mathrm{ d}{\bf F}^{\rm (A)}(x)}{\mathrm{ d}x} \Bigr \rbrace \, \mathrm{ d} x \nonumber \\ && + \frac{1}{2} \int _0^1 {\rm Tr} \Bigl \lbrace [{\bf F}^{\rm (B)}(x)]^{-1}\frac{\mathrm{ d}{\bf F}^{\rm (B)}(x)}{\mathrm{ d}x} \Bigr \rbrace \, \mathrm{ d} x~~~~ \nonumber \\ &=& \ln f_{\rm deg}^{\rm (A)} + \ln f_{\rm deg}^{\rm (B)}. \end{eqnarray*}$$(F36) (Equality holds only when dF(A)/dx = dF(B) = 0.) We see that for no shared nuisance parameters, $$\begin{eqnarray*} f_{\rm deg}^{\rm (A+B)} \le f_{\rm deg}^{\rm (A)} f_{\rm deg}^{\rm (B)}, \end{eqnarray*}$$(F37) with equality only in the case where the systematic has no effect on the Fisher matrices, i.e. when |$f_{\rm deg}^{\rm (A+B)} = f_{\rm deg}^{\rm (A)} = f_{\rm deg}^{\rm (B)} = 1$|⁠. A particular example showing that this inequality cannot be strengthened is in the space of two cosmological parameters: $$\begin{eqnarray*} {\bf F}^{\rm (A,old)} = \begin{pmatrix}1 & 0 \\ 0 & \epsilon \end{pmatrix}, ~~ {\bf F}^{\rm (B,old)} = \begin{pmatrix}\epsilon & 0 \\ 0 & 1 \end{pmatrix}, \nonumber \end{eqnarray*}$$ $$\begin{eqnarray*} {\bf F}^{\rm (A,new)} = \begin{pmatrix}\alpha ^{-2} & 0 \\ 0 & \epsilon \end{pmatrix}, ~~{\rm and}~~ {\bf F}^{\rm (B,new)} = \begin{pmatrix}\epsilon & 0 \\ 0 & \beta ^{-2} \end{pmatrix}, \end{eqnarray*}$$(F38) with α, β ≥ 1. By taking the limit of ϵ → 0, we see that |$f_{\rm deg}^{\rm (A)}\rightarrow \alpha$|⁠, |$f_{\rm deg}^{\rm (B)}\rightarrow \beta$|⁠, and |$f_{\rm deg}^{\rm (A+B)}\rightarrow \alpha \beta$|⁠. © 2020 The Author(s) Published by Oxford University Press on behalf of Royal Astronomical Society This article is published and distributed under the terms of the Oxford University Press, Standard Journals Publication Model (https://academic.oup.com/journals/pages/open_access/funder_policies/chorus/standard_publication_model)
TI - A synthetic Roman Space Telescope High-Latitude Imaging Survey: simulation suite and the impact of wavefront errors on weak gravitational lensing
JF - Monthly Notices of the Royal Astronomical Society
DO - 10.1093/mnras/staa3658
DA - 2021-01-05
UR - https://www.deepdyve.com/lp/oxford-university-press/a-synthetic-roman-space-telescope-high-latitude-imaging-survey-rxSG0lBEXm
SP - 2044
EP - 2070
VL - 501
IS - 2
DP - DeepDyve
ER -