# The Dark Energy Spectroscopic Instrument: One-dimensional power spectrum from first Lyman-α forest samples with Fast Fourier Transform

Corentin Ravoux, Marie Lynn Abdul Karim, Eric Armengaud, Michael

Walther, Naim Göksel Karaçaylı, Paul Martini, Julien Guy, Jessica Nicole

Aguilar, Steven Ahlen, Stephen Bailey, et al.

► To cite this version:

Corentin Ravoux, Marie Lynn Abdul Karim, Eric Armengaud, Michael Walther, Naim Göksel Karaçaylı, et al.. The Dark Energy Spectroscopic Instrument: One-dimensional power spectrum from first Lyman-α forest samples with Fast Fourier Transform. Monthly Notices of the Royal Astronomical Society, 2023, 526 (4), pp.5118-5140. 10.1093/mnras/stad3008. hal-04143297

HAL Id: hal-04143297

https://hal.science/hal-04143297

Submitted on 24 May 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# The Dark Energy Spectroscopic Instrument: One-dimensional power spectrum from first Lyman-α forest samples with Fast Fourier Transform

Corentin Ravoux, Marie Lynn Abdul Karim, Eric Armengaud, Michael Walther, Naim Göksel Karaçaylı, Paul Martini, Julien Guy, Jessica Nicole Aguilar, Steven Ahlen, Stephen Bailey, et al.

► To cite this version:

Corentin Ravoux, Marie Lynn Abdul Karim, Eric Armengaud, Michael Walther, Naim Göksel Karaçaylı, et al.. The Dark Energy Spectroscopic Instrument: One-dimensional power spectrum from first Lyman-α forest samples with Fast Fourier Transform. Monthly Notices of the Royal Astronomical Society, 2023, 526 (4), pp.5118-5140. 10.1093/mnras/stad3008. hal-04143297

HAL Id: hal-04143297

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

Submitted on 24 May 2024

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

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

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

Distributed under a Creative Commons Attribution 4.0 International License
# The Dark Energy Spectroscopic Instrument: one-dimensional power spectrum from first Ly α forest samples with Fast Fourier Transform

Corentin Ravoux<sup>1,2★</sup> Marie Lynn Abdul Karim,<sup>2</sup> Eric Armengaud<sup>1</sup>,<sup>2</sup> Michael Walther<sup>1</sup>,<sup>3,4</sup>

Naim Göksel Karaçaylı<sup>1</sup>,<sup>5,6,7,8</sup> Paul Martini,<sup>5,7,8</sup> Julien Guy,<sup>9</sup> Jessica Nicole Aguilar,<sup>9</sup> Steven Ahlen<sup>1</sup>,<sup>10</sup>

Stephen Bailey<sup>1</sup>,<sup>9</sup> Julian Bautista<sup>1</sup>,<sup>1</sup> Sergio Felipe Beltran<sup>1</sup>,<sup>11</sup> David Brooks,<sup>12</sup>

Laura Cabayol-Garcia,<sup>13</sup> Solène Chabanier<sup>1</sup>,<sup>9</sup> Edmond Chaussidon<sup>1</sup>,<sup>2</sup> Jonás Chaves-Montero<sup>1</sup>,<sup>13</sup>

Kyle Dawson,<sup>14</sup> Rodrigo de la Cruz<sup>1</sup>,<sup>11</sup> Axel de la Macorra<sup>1</sup>,<sup>15</sup> Peter Doel,<sup>12</sup> Kevin Fanning<sup>1</sup>,<sup>8</sup>

Andreu Font-Ribera<sup>1</sup>,<sup>13</sup> Jaime Forero-Romero<sup>1</sup>,<sup>16,17</sup> Satya Gontcho A Gontcho<sup>1</sup>,<sup>9</sup>

Alma X. Gonzalez-Morales<sup>1</sup>,<sup>11,18</sup> Calum Gordon,<sup>13</sup> Hiram K. Herrera-Alcantar<sup>1</sup>,<sup>11</sup> Klaus Honscheid,<sup>5,6,8</sup>

Vid Iršič<sup>1</sup>,<sup>19</sup> Mustapha Ishak<sup>1</sup>,<sup>20</sup> Robert Kehoe,<sup>20</sup> Theodore Kisner<sup>1</sup>,<sup>9</sup> Anthony Kremin<sup>1</sup>,<sup>9</sup>

Martin Landriau<sup>1</sup>,<sup>9</sup> Laurent Le Guillou<sup>1</sup>,<sup>21</sup> Michael Levi<sup>1</sup>,<sup>9</sup> Zarija Lukić,<sup>9</sup> Christophe Magneville,<sup>2</sup>

Aaron Meisner<sup>1</sup>,<sup>22</sup> Ramon Miquel,<sup>13,23</sup> John Moustakas<sup>1</sup>,<sup>24</sup> Eva-Maria Mueller,<sup>25</sup>

Andrea Muñoz-Gutiérrez,<sup>15</sup> Lucas Napolitano<sup>1</sup>,<sup>26</sup> Jundan Nie<sup>1</sup>,<sup>27</sup> Gustavo Niz<sup>1</sup>,<sup>11,28</sup>

Nathalie Palanque-Delabrouille<sup>1</sup>,<sup>2,9</sup> Will Percival<sup>1</sup>,<sup>29,30,31</sup> Ignasi Pérez-Ràfols<sup>1</sup>,<sup>13</sup> Matthew Pieri,<sup>32</sup>

Claire Poppett,<sup>9,33,34</sup> Francisco Prada<sup>1</sup>,<sup>35</sup> César Ramírez Pérez,<sup>13</sup> Graziano Rossi,<sup>36</sup> Eusebio Sanchez<sup>1</sup>,<sup>37</sup>

David Schlegel,<sup>9</sup> Michael Schubnell,<sup>38,39</sup> Hee-Jong Seo<sup>1</sup>,<sup>40</sup> Francesco Sinigaglia,<sup>41,42</sup> Ting Tan,<sup>21</sup>

Gregory Tarlé<sup>1</sup>,<sup>39</sup> Ben Wang<sup>1</sup>,<sup>43</sup> Benjamin Weaver,<sup>22</sup> Christophe Yèche<sup>1</sup>,<sup>2</sup> and Zhimin Zhou<sup>1</sup><sup>27</sup>

*Affiliations are listed at the end of the paper*

Accepted 2023 October 2. Received 2023 September 13; in original form 2023 June 13

## ABSTRACT

We present the one-dimensional Ly α forest power spectrum measurement using the first data provided by the Dark Energy Spectroscopic Instrument (DESI). The data sample comprises 26 330 quasar spectra, at redshift $z > 2.1$, contained in the DESI Early Data Release and the first 2 months of the main survey. We employ a Fast Fourier Transform (FFT) estimator and compare the resulting power spectrum to an alternative likelihood-based method in a companion paper. We investigate methodological and instrumental contaminants associated with the new DESI instrument, applying techniques similar to previous Sloan Digital Sky Survey (SDSS) measurements. We use synthetic data based on lognormal approximation to validate and correct our measurement. We compare our resulting power spectrum with previous SDSS and high-resolution measurements. With relatively small number statistics, we successfully perform the FFT measurement, which is already competitive in terms of the scale range. At the end of the DESI survey, we expect a five times larger Ly α forest sample than SDSS, providing an unprecedented precise one-dimensional power spectrum measurement.

**Key words:** intergalactic medium – large-scale structure of Universe – cosmology: observations.

## 1 INTRODUCTION

The Ly α (Ly α) forest can be observed from the ground in the optical spectra of distant quasars at redshift between the end phase of reionization ($z \sim 6$) and the peak of galaxy formation $z \sim 2$. The Ly α forest consists of a series of Ly α absorption lines caused by intervening neutral hydrogen located at various redshifts between the quasar and the observer. The Ly α forest is a powerful probe of the

underlying matter density field at redshift $z > 2$, together with the astrophysical state of the intergalactic medium (Gunn & Peterson 1965; Lynds 1971; Meiksin 2009; McQuinn 2016).

In particular, the small-scale distribution of neutral hydrogen ($\sim$Mpc) is imprinted in the fluctuations of the Ly α forest along the line of sight that can be accessed by measuring the one-dimensional (1D) Ly α forest power spectrum (denoted $P_{1D,α}$). This measurement is sensitive to the amplitude and slope of the matter power spectrum at redshift $z > 2$. The impact of cosmological parameters on $P_{1D,α}$ can only be accurately predicted using hydrodynamical simulations (Borde et al. 2014; Bolton et al. 2017; Walther et al. 2021; Puchwein

<footnote>*E-mail: ravoux@cppm.in2p3.fr</footnote>
<figure><img src="image_2.png" /></figure>

# The Dark Energy Spectroscopic Instrument: one-dimensional power spectrum from first Ly α forest samples with Fast Fourier Transform

Corentin Ravoux<sup>1,2★</sup> Marie Lynn Abdul Karim,<sup>2</sup> Eric Armengaud<sup>1</sup>,<sup>2</sup> Michael Walther<sup>1,3,4</sup>

Naim Göksel Karaçaşlı<sup>1</sup>,<sup>5,6,7,8</sup> Paul Martini,<sup>5,7,8</sup> Julien Guy,<sup>9</sup> Jessica Nicole Aguilar,<sup>9</sup> Steven Ahlen<sup>1</sup>,<sup>10</sup>

Stephen Bailey<sup>1</sup>,<sup>9</sup> Julian Bautista<sup>1</sup>,<sup>1</sup> Sergio Felipe Beltran<sup>1</sup>,<sup>11</sup> David Brooks,<sup>12</sup>

Laura Cabayol-Garcia,<sup>13</sup> Solène Chabanier<sup>1</sup>,<sup>9</sup> Edmond Chaussidon<sup>1</sup>,<sup>2</sup> Jonás Chaves-Montero<sup>1</sup>,<sup>13</sup>

Kyle Dawson,<sup>14</sup> Rodrigo de la Cruz<sup>1</sup>,<sup>11</sup> Axel de la Macorra<sup>1</sup>,<sup>15</sup> Peter Doel,<sup>12</sup> Kevin Fanning<sup>1</sup>,<sup>8</sup>

Andreu Font-Ribera<sup>1</sup>,<sup>13</sup> Jaime Forero-Romero<sup>1</sup>,<sup>16,17</sup> Satya Gontcho A Gontcho<sup>1</sup>,<sup>9</sup>

Alma X. Gonzalez-Morales<sup>1</sup>,<sup>11,18</sup> Calum Gordon,<sup>13</sup> Hiram K. Herrera-Alcantar<sup>1</sup>,<sup>11</sup> Klaus Honscheid,<sup>5,6,8</sup>

Vid Iršič<sup>1</sup>,<sup>19</sup> Mustapha Ishak<sup>1</sup>,<sup>20</sup> Robert Kehoe,<sup>20</sup> Theodore Kisner<sup>1</sup>,<sup>9</sup> Anthony Kremin<sup>1</sup>,<sup>9</sup>

Martin Landriau<sup>1</sup>,<sup>9</sup> Laurent Le Guillou<sup>1</sup>,<sup>21</sup> Michael Levi<sup>1</sup>,<sup>9</sup> Zarija Lukić,<sup>9</sup> Christophe Magneville,<sup>2</sup>

Aaron Meisner<sup>1</sup>,<sup>22</sup> Ramon Miquel,<sup>13,23</sup> John Moustakas<sup>1</sup>,<sup>24</sup> Eva-Maria Mueller,<sup>25</sup>

Andrea Muñoz-Gutiérrez,<sup>15</sup> Lucas Napolitano<sup>1</sup>,<sup>26</sup> Jundan Nie<sup>1</sup>,<sup>27</sup> Gustavo Niz<sup>1</sup>,<sup>11,28</sup>

Nathalie Palanque-Delabrouille<sup>1</sup>,<sup>2,9</sup> Will Percival<sup>1</sup>,<sup>29,30,31</sup> Ignasi Pérez-Ràfols<sup>1</sup>,<sup>13</sup> Matthew Pieri,<sup>32</sup>

Claire Poppett,<sup>9,33,34</sup> Francisco Prada<sup>1</sup>,<sup>35</sup> César Ramírez Pérez,<sup>13</sup> Graziano Rossi,<sup>36</sup> Eusebio Sanchez<sup>1</sup>,<sup>37</sup>

David Schlegel,<sup>9</sup> Michael Schubnell,<sup>38,39</sup> Hee-Jong Seo<sup>1</sup>,<sup>40</sup> Francesco Sinigaglia,<sup>41,42</sup> Ting Tan,<sup>21</sup>

Gregory Tarlé<sup>1</sup>,<sup>39</sup> Ben Wang<sup>1</sup>,<sup>43</sup> Benjamin Weaver,<sup>22</sup> Christophe Yèche<sup>1</sup>,<sup>2</sup> and Zhimin Zhou<sup>1</sup><sup>27</sup>

*Affiliations are listed at the end of the paper*

Accepted 2023 October 2. Received 2023 September 13; in original form 2023 June 13

## ABSTRACT

We present the one-dimensional Ly α forest power spectrum measurement using the first data provided by the Dark Energy Spectroscopic Instrument (DESI). The data sample comprises 26 330 quasar spectra, at redshift $z > 2.1$, contained in the DESI Early Data Release and the first 2 months of the main survey. We employ a Fast Fourier Transform (FFT) estimator and compare the resulting power spectrum to an alternative likelihood-based method in a companion paper. We investigate methodological and instrumental contaminants associated with the new DESI instrument, applying techniques similar to previous Sloan Digital Sky Survey (SDSS) measurements. We use synthetic data based on lognormal approximation to validate and correct our measurement. We compare our resulting power spectrum with previous SDSS and high-resolution measurements. With relatively small number statistics, we successfully perform the FFT measurement, which is already competitive in terms of the scale range. At the end of the DESI survey, we expect a five times larger Ly α forest sample than SDSS, providing an unprecedented precise one-dimensional power spectrum measurement.

**Key words:** intergalactic medium – large-scale structure of Universe – cosmology: observations.

## 1 INTRODUCTION

The Ly α (Ly α) forest can be observed from the ground in the optical spectra of distant quasars at redshift between the end phase of reionization ($z \sim 6$) and the peak of galaxy formation $z \sim 2$. The Ly α forest consists of a series of Ly α absorption lines caused by intervening neutral hydrogen located at various redshifts between the quasar and the observer. The Ly α forest is a powerful probe of the

underlying matter density field at redshift $z > 2$, together with the astrophysical state of the intergalactic medium (Gunn & Peterson 1965; Lynds 1971; Meiksin 2009; McQuinn 2016).

In particular, the small-scale distribution of neutral hydrogen ($\sim$Mpc) is imprinted in the fluctuations of the Ly α forest along the line of sight that can be accessed by measuring the one-dimensional (1D) Ly α forest power spectrum (denoted $P_{1D,α}$). This measurement is sensitive to the amplitude and slope of the matter power spectrum at redshift $z > 2$. The impact of cosmological parameters on $P_{1D,α}$ can only be accurately predicted using hydrodynamical simulations (Borde et al. 2014; Bolton et al. 2017; Walther et al. 2021; Puchwein

<footnote>*E-mail: ravoux@cppm.in2p3.fr</footnote>
et al. 2023). The realization of those simulations is made arduous by the large dynamic range needed to model the Ly α forest adequately (Lukić et al. 2015; Chabanier et al. 2022). Fitting data measurements with those simulation predictions provides constraints on the cosmological parameters σ<sub>8</sub>, n<sub>s</sub>, and Ω<sub>m</sub>, as well as on parameters describing the thermal properties of the intergalactic medium. In particular, the simulations described in Walther et al. (2021) are able to predict *P*<sub>1D,α</sub> with sufficient accuracy (at the 1 per cent level) when compared to expected uncertainties from the Dark Energy Spectroscopic Instrument (DESI) measurement.

Due to its sensitivity to the matter fluctuations at small scales, measurements of $P_{1D,\alpha}$ can constrain physics beyond the Standard Model, such as the mass of neutrinos, the mass of warm dark matter candidates, or a possible running of the spectral index due to primordial inflation physics. First, $P_{1D,\alpha}$ is well suited to constrain the sum of neutrino masses that damps the matter power spectrum at small scales (Lesgourgues & Pastor 2006, 2012). Stringent constraints are obtained by coupling $P_{1D,\alpha}$ with hydrodynamical simulations, and by combining it with the cosmic microwave background (Seljak, Slosar & McDonald 2006; Palanque-Delabrouille et al. 2015; Yèche et al. 2017; Palanque-Delabrouille et al. 2020). Secondly, several studies combined high- and moderate-resolution $P_{1D,\alpha}$ measurements to obtain constraints on the warm dark matter mass (Viel et al. 2005, 2008, 2013; Baur et al. 2016, 2017; Yèche et al. 2017; Palanque-Delabrouille et al. 2020). The hydrodynamical simulations used in those studies either directly model neutrinos as particles or using a rescaling of the matter power spectrum to account for neutrinos (Pedersen et al. 2020a, b; Pedersen, Font-Ribera & Gnedin 2023). Finally, other exotic dark matter models such as fuzzy dark matter (Armengaud et al. 2017; Iršič et al. 2017) can also be constrained using $P_{1D,\alpha}$ measurement.

Between the first $P_{1D,\alpha}$ measurements (Croft et al. 1998; McDonald et al. 2000; Croft et al. 2002; Kim et al. 2004) and today, the large increase in observation capabilities brought numerous Ly α forest samples that can be split between moderate-resolution ($\lambda/\Delta\lambda \lesssim 5000$) and high-resolution ($\lambda/\Delta\lambda \gtrsim 20\ 000$) observations.

The computation of $P_{1D,\alpha}$ with high-resolution data sets such as SQUAD (Murphy et al. 2018), KODIAQ (O'Meara et al. 2015; O'Meara 2017), or XQ-100 (Lopez et al. 2016) are performed in Viel et al. (2013), Iršič et al. (2016), Walther et al. (2018), Boera et al. (2019), Day, Tytler & Kambalur (2019), Khaire et al. (2019), Gaikwad et al. (2021), and Karaçaylı et al. (2022). Those measurements use high signal-to-noise quasars to probe the intergalactic medium at very small scales (~100 kpc) but does not provide sufficient statistics to accurately measure the large-scale clustering (>5 Mpc), needed for cosmological interpretation.

Moderate-resolution surveys provide large numbers of Ly α forests, which yield smaller statistical uncertainties on $P_{1D,\alpha}$ estimates. However, the resolution of such spectrographs limits the reach of very small scales. The $P_{1D,\alpha}$ measurement with moderate-resolution spectrograph was first performed on a small sample of the Sloan Digital Sky Survey (SDSS) data in McDonald et al. (2006). Subsequently, the increase of the Ly α forest statistic has largely improved this measurement with the Baryon Oscillation Spectroscopic Survey (BOSS) in Palanque-Delabrouille et al. (2013) and the extended Baryon Oscillation Spectroscopic Survey (eBOSS) in Chabanier et al. (2019) using 43 751 quasar spectra.

Several methods can be used to measure $P_{1D,\alpha}$ from Ly α forest samples. The most straightforward relies on the Fast Fourier Transform (FFT) and was applied in BOSS and eBOSS (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019) analyses. The 1D power spectrum can also be measured with configuration space

estimators such as the quadratic maximum likelihood estimator (QMLE). This method has already been applied to moderate resolution observations (McDonald et al. 2006; Palanque-Delabrouille et al. 2013) and more recently on high-resolution data in Karaçaylı et al. (2022). The QMLE method is applied to the same data than used in the present paper and presented in a companion paper (Karaçaylı et al. 2023a). The FFT method yields a straightforward calculation of $P_{1D,\alpha}$ and offers more control over the different calculation steps. Conversely, the more-complex QMLE estimation is not sensitive to gaps in the quasar spectra. The results between the two methods are presented in the companion paper and are in good agreement. FFT and QMLE results agree at 1 per cent level precision up to half the Nyquist frequency.

The purpose of this work is to compute $P_{1D,\alpha}$ from the first DESI data, following the same methodology as in the latest eBOSS measurement in Chabanier et al. (2019). Using the same method facilitate the comparison between eBOSS and DESI. The $P_{1D,\alpha}$ is sensitive to instrumental properties such as noise and spectral resolution. As the telescopes used and the data are very different, it is essential to characterize the DESI instrument. We improve the algorithms and methodology used in Chabanier et al. (2019) to account for systematic and instrumental differences between eBOSS and DESI. In particular, due to the spectral resolution improvement of DESI, our measurement allows accessing smaller scales than eBOSS.

The outline of this paper is as follows: Section 2 describes the DESI instrument and data processing used to perform this $P_{1D,\alpha}$ measurement. The $P_{1D,\alpha}$ pipeline is presented in Section 3, and the characterization of the DESI instrument in Section 4. We generate synthetic data to validate and correct our measurement in Section 5. The treatment of statistical and systematic uncertainties for the $P_{1D,\alpha}$ measurement is given in Section 6. Finally, we present our measurement on DESI data, as well as a comparison to previous measurements in Section 7, and conclude in Section 8.

## 2 INSTRUMENT AND DATA DESCRIPTION

DESI has as objective to measure the spectra of 40 million galaxies and quasars in a footprint of 14 000 deg<sup>2</sup> over 5 yr (Levi et al. 2013; DESI 2016a; Abareshi et al. 2022). This project aims to continue the cosmic mapping efforts started by SDSS, while drastically increasing its constraining power on the Lambda cold dark matter model and its possible extensions.

We first focus on the description of the data used for our measurement. In the following, we describe the data starting from the instrument (DESI), its associated spectroscopic pipeline, the different data acquisition phases, and the input catalogues of our study.

### 2.1 DESI instrument

The DESI instrument is mounted on the Mayall telescope, located on the Kitt Peak National Observatory (KPNO) in the Tohono O'odham Nation. The Mayall telescope is a reflective prime-focus telescope with a 4-m diameter primary mirror. The DESI instrument (DESI 2016b; Miller et al. 2023) receives photons through an optical corrector designed to increase the field of view to 7.5 deg<sup>2</sup> on the focal plane. The focal plane system, composed of 5000 robotically controlled fibres, can quickly modify its configuration to aim at the targeted objects on a specific footprint (Silber et al. 2022). An optical fibre system redirects the light of the observed objects to a separate climate-controlled enclosure containing 10 spectrographs. Each spectrograph comprises three CCD cameras, whose properties are given in the Table 1. In comparison to SDSS spectrographs,
et al. 2023). The realization of those simulations is made arduous by the large dynamic range needed to model the Ly α forest adequately (Lukić et al. 2015; Chabanier et al. 2022). Fitting data measurements with those simulation predictions provides constraints on the cosmological parameters σ<sub>8</sub>, n<sub>s</sub>, and Ω<sub>m</sub>, as well as on parameters describing the thermal properties of the intergalactic medium. In particular, the simulations described in Walther et al. (2021) are able to predict *P*<sub>1D,α</sub> with sufficient accuracy (at the 1 per cent level) when compared to expected uncertainties from the Dark Energy Spectroscopic Instrument (DESI) measurement.

Due to its sensitivity to the matter fluctuations at small scales, measurements of $P_{1D,\alpha}$ can constrain physics beyond the Standard Model, such as the mass of neutrinos, the mass of warm dark matter candidates, or a possible running of the spectral index due to primordial inflation physics. First, $P_{1D,\alpha}$ is well suited to constrain the sum of neutrino masses that damps the matter power spectrum at small scales (Lesgourgues & Pastor 2006, 2012). Stringent constraints are obtained by coupling $P_{1D,\alpha}$ with hydrodynamical simulations, and by combining it with the cosmic microwave background (Seljak, Slosar & McDonald 2006; Palanque-Delabrouille et al. 2015; Yèche et al. 2017; Palanque-Delabrouille et al. 2020). Secondly, several studies combined high- and moderate-resolution $P_{1D,\alpha}$ measurements to obtain constraints on the warm dark matter mass (Viel et al. 2005, 2008, 2013; Baur et al. 2016, 2017; Yèche et al. 2017; Palanque-Delabrouille et al. 2020). The hydrodynamical simulations used in those studies either directly model neutrinos as particles or using a rescaling of the matter power spectrum to account for neutrinos (Pedersen et al. 2020a, b; Pedersen, Font-Ribera & Gnedin 2023). Finally, other exotic dark matter models such as fuzzy dark matter (Armengaud et al. 2017; Iršič et al. 2017) can also be constrained using $P_{1D,\alpha}$ measurement.

Between the first $P_{1D,\alpha}$ measurements (Croft et al. 1998; McDonald et al. 2000; Croft et al. 2002; Kim et al. 2004) and today, the large increase in observation capabilities brought numerous Ly α forest samples that can be split between moderate-resolution ($\lambda/\Delta\lambda \lesssim 5000$) and high-resolution ($\lambda/\Delta\lambda \gtrsim 20\ 000$) observations.

The computation of $P_{1D,\alpha}$ with high-resolution data sets such as SQUAD (Murphy et al. 2018), KODIAQ (O'Meara et al. 2015; O'Meara 2017), or XQ-100 (Lopez et al. 2016) are performed in Viel et al. (2013), Iršič et al. (2016), Walther et al. (2018), Boera et al. (2019), Day, Tytler & Kambalur (2019), Khaire et al. (2019), Gaikwad et al. (2021), and Karaçaylı et al. (2022). Those measurements use high signal-to-noise quasars to probe the intergalactic medium at very small scales (~100 kpc) but does not provide sufficient statistics to accurately measure the large-scale clustering (>5 Mpc), needed for cosmological interpretation.

Moderate-resolution surveys provide large numbers of Ly α forests, which yield smaller statistical uncertainties on $P_{1D,\alpha}$ estimates. However, the resolution of such spectrographs limits the reach of very small scales. The $P_{1D,\alpha}$ measurement with moderate-resolution spectrograph was first performed on a small sample of the Sloan Digital Sky Survey (SDSS) data in McDonald et al. (2006). Subsequently, the increase of the Ly α forest statistic has largely improved this measurement with the Baryon Oscillation Spectroscopic Survey (BOSS) in Palanque-Delabrouille et al. (2013) and the extended Baryon Oscillation Spectroscopic Survey (eBOSS) in Chabanier et al. (2019) using 43 751 quasar spectra.

Several methods can be used to measure $P_{1D,\alpha}$ from Ly α forest samples. The most straightforward relies on the Fast Fourier Transform (FFT) and was applied in BOSS and eBOSS (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019) analyses. The 1D power spectrum can also be measured with configuration space

estimators such as the quadratic maximum likelihood estimator (QMLE). This method has already been applied to moderate resolution observations (McDonald et al. 2006; Palanque-Delabrouille et al. 2013) and more recently on high-resolution data in Karaçaylı et al. (2022). The QMLE method is applied to the same data than used in the present paper and presented in a companion paper (Karaçaylı et al. 2023a). The FFT method yields a straightforward calculation of $P_{1D,\alpha}$ and offers more control over the different calculation steps. Conversely, the more-complex QMLE estimation is not sensitive to gaps in the quasar spectra. The results between the two methods are presented in the companion paper and are in good agreement. FFT and QMLE results agree at 1 per cent level precision up to half the Nyquist frequency.

The purpose of this work is to compute $P_{1D,\alpha}$ from the first DESI data, following the same methodology as in the latest eBOSS measurement in Chabanier et al. (2019). Using the same method facilitate the comparison between eBOSS and DESI. The $P_{1D,\alpha}$ is sensitive to instrumental properties such as noise and spectral resolution. As the telescopes used and the data are very different, it is essential to characterize the DESI instrument. We improve the algorithms and methodology used in Chabanier et al. (2019) to account for systematic and instrumental differences between eBOSS and DESI. In particular, due to the spectral resolution improvement of DESI, our measurement allows accessing smaller scales than eBOSS.

The outline of this paper is as follows: Section 2 describes the DESI instrument and data processing used to perform this $P_{1D,\alpha}$ measurement. The $P_{1D,\alpha}$ pipeline is presented in Section 3, and the characterization of the DESI instrument in Section 4. We generate synthetic data to validate and correct our measurement in Section 5. The treatment of statistical and systematic uncertainties for the $P_{1D,\alpha}$ measurement is given in Section 6. Finally, we present our measurement on DESI data, as well as a comparison to previous measurements in Section 7, and conclude in Section 8.

## 2 INSTRUMENT AND DATA DESCRIPTION

DESI has as objective to measure the spectra of 40 million galaxies and quasars in a footprint of 14 000 deg<sup>2</sup> over 5 yr (Levi et al. 2013; DESI 2016a; Abareshi et al. 2022). This project aims to continue the cosmic mapping efforts started by SDSS, while drastically increasing its constraining power on the Lambda cold dark matter model and its possible extensions.

We first focus on the description of the data used for our measurement. In the following, we describe the data starting from the instrument (DESI), its associated spectroscopic pipeline, the different data acquisition phases, and the input catalogues of our study.

### 2.1 DESI instrument

The DESI instrument is mounted on the Mayall telescope, located on the Kitt Peak National Observatory (KPNO) in the Tohono O'odham Nation. The Mayall telescope is a reflective prime-focus telescope with a 4-m diameter primary mirror. The DESI instrument (DESI 2016b; Miller et al. 2023) receives photons through an optical corrector designed to increase the field of view to 7.5 deg<sup>2</sup> on the focal plane. The focal plane system, composed of 5000 robotically controlled fibres, can quickly modify its configuration to aim at the targeted objects on a specific footprint (Silber et al. 2022). An optical fibre system redirects the light of the observed objects to a separate climate-controlled enclosure containing 10 spectrographs. Each spectrograph comprises three CCD cameras, whose properties are given in the Table 1. In comparison to SDSS spectrographs,
<table><caption>**Table 1.** Spectral range and effective resolving power ($R = \Delta\lambda/\lambda$) for each channel of the DESI spectrographs (Abareshi et al. 2022).</caption><thead><tr><td>Channel</td><td>Spectral range (Å)</td><td>Resolving power</td></tr></thead><tbody><tr><td>Blue (B)</td><td>3600–5930</td><td>2000–3200</td></tr><tr><td>Red (R)</td><td>5600–7720</td><td>3200–4100</td></tr><tr><td>Near infrared (Z)</td><td>7470–9800</td><td>4100–5100</td></tr></tbody></table>

the effective resolving power ($\Delta\lambda/\lambda$) improved by at least a factor of 2.

## 2.2 DESI spectroscopic pipeline

The high complexity of the DESI survey induces the need for advanced software pipelines and products, including the imaging from the DESI Legacy Imaging Surveys (Zou et al. 2017; Dey et al. 2019; Schlegel et al., in preparation), a pipeline to select the targets to observe (Myers et al. 2022), a pipeline to assign fibres (Raichoor et al., in preparation), a pipeline to parse the survey and to optimize the observation strategy (Schlaflay et al. 2023), and an extensive spectroscopic reduction pipeline (Guy et al. 2022).

This spectroscopic pipeline, called desispec https://github.com/desihub/desispec, transforms the raw CCD images into spectra, and is detailed in Guy et al. (2022). Before extracting the spectra, the images are subtracted by dark and bias calibration frames to remove expected background sources, and to estimate the associated readout noise (noise estimation details are given in Appendix C). The non-uniform CCD pixel response is corrected using a dedicated flat-field slit on the spectrograph and the CCD over-scan is removed. A dedicated software detects and flags cosmic rays or defective CCD pixels.

The spectral extraction is performed using the 'spectroperfectionism' method (Bolton & Schlegel 2010), an optimal spectroscopic extraction that correctly models complex 2D point spread functions (PDFs). This method provides the encoding for each fibre and each wavelength of the non-Gaussian instrument resolution into a resolution matrix (noted $R$ in next sections) used to compute $P_{1D,\alpha}$.

All the spectra are defined on the same wavelength grid without additional resampling. Consequently, the extracted spectra are linearly binned in observed wavelength with a constant separation $\Delta\lambda_{pix} = 0.8$ Å. Conversely, the SDSS/BOSS spectrographs were logarithmically binned with a constant spectral pixel size of $\Delta v = 69$ km s$^{-1}$, which corresponds to $\Delta\lambda/\lambda = \Delta v/c = 2.3 \times 10^{-4}$, equivalent to $\Delta[\log(\lambda)] = 10^{-4}$ (Smee et al. 2013).

Once the spectrum of each fibre has been extracted, several post-processing steps allow removing a variety of further observational effects. The non-uniform response of individual fibres as a function of wavelength is corrected with flat-field frames by observing a white screen attached to the telescope dome and illuminated with a LED array. For all exposures, some fibres are dedicated to observing the sky. The so-called sky spectra associated with these fibres provide the sky level and the intensity of atmospheric emission lines and are subtracted to the spectra associated with targets (sky subtraction). The transmission defaults of the atmosphere and telescope as a whole are corrected by the observation of calibrated star spectra. This step converts CCD units (number of electrons) to observed flux units. Finally, the spectrum of an object is obtained by coadding its different exposures. The resulting spectrum is expressed separately into the three spectrographs bands described in Table 1.

All the software of the pipeline employed for the analysis of DESI data are listed in the repository desihub https://github.com/desihub. In particular, the spectra analysed in this

<table><caption>**Table 2.** Summary of the DESI data sets used in this study, their associated acronyms, a subjective description, and the total number of quasar spectra whose redshifts are between 2.0 and 5.0.</caption><thead><tr><td>Data set</td><td>Acronym</td><td>Description</td><td>Quasar number</td></tr></thead><tbody><tr><td>Target selection validation</td><td>SV1</td><td>Small sample with deeper exposures (up to 16) than the end of DESI</td><td>12 355</td></tr><tr><td>One-percent survey</td><td>SV3</td><td>Same density and depth as for the full 5 yr DESI main survey</td><td>12 686</td></tr><tr><td>First 2 months of the main survey</td><td>DESI-M2(M2)</td><td>Sparse data with small number of exposures</td><td>87 373</td></tr></tbody></table>

article have been processed with the fuji https://github.com/desihub/desispec/releases/tag/0.51.13. version of the spectroscopic pipeline.

## 2.3 DESI data

The spectroscopic pipeline described previously was used to analyse data obtained over different periods. In this article, we use three data sets from the first observations of DESI described in Table 2, which includes a total of 112 414 quasar spectra whose redshifts are between 2.0 and 5.0.

The first two data sets 'Target Selection Validation' (SV1) (DESI 2023a) and 'One-Percent Survey' (SV3) are part of the Early Data Release (EDR) whose complete description is given in DESI (2023b). 'Target Selection Validation' was conducted from 2020 December to 2021 March and includes a large number of exposures (up to 16) for the same targets. The objective of this survey was among others, to study extensively the survey performance as a function of instrumental depth and to build visual inspection truth tables. The SV1 data set also includes 'Secondary Tiles' as detailed in DESI (2023b). Following the completion of 'Target Selection Validation' observations, the 'One-Percent Survey' phase was dedicated more specifically to the evaluation of the survey design. The number of exposures is similar to the main survey at its end (between 4 and 5 for each Ly $\alpha$ quasars) and the goal was to determine the best strategy to cover the sky while limiting fibre loss.

The main DESI survey started in June 2021 and in this article we also use the first 2 months of data, named DESI-M2 (and noted M2 in this paper for conciseness), which is not present in the EDR but will be included in the Data Release 1. In the M2 data set, most quasars have only one exposure. While all three data sets are studied in this article, in the end we removed SV1 due to its different noise properties (see Appendix A). The final measurement is computed on SV3 + M2 data set only.

## 2.4 Input catalogues

The input catalogues used to compute $P_{1D,\alpha}$ were obtained by applying specific procedures to the three data sets previously described.

The targeting of the quasars used in our study (Yèche et al. 2020; Chaussidon et al. 2022) was verified with visual inspection of subsets of early observations (Alexander et al. 2022). For quasars, the DESI pipeline categorizes the observed spectra and estimates their redshifts using the redrock spectral template fitting software https://github.com/desihub/redrock (Bailey et al., in preparation; Brodzeller et al. 2023). In order to optimize, the completeness of the quasar catalogue while keeping a high purity, additionally a broad Mg II line finder mgii_afterburner
<table><caption>**Table 1.** Spectral range and effective resolving power ($R = \Delta\lambda/\lambda$) for each channel of the DESI spectrographs (Abareshi et al. 2022).</caption><thead><tr><td>Channel</td><td>Spectral range (Å)</td><td>Resolving power</td></tr></thead><tbody><tr><td>Blue (B)</td><td>3600–5930</td><td>2000–3200</td></tr><tr><td>Red (R)</td><td>5600–7720</td><td>3200–4100</td></tr><tr><td>Near infrared (Z)</td><td>7470–9800</td><td>4100–5100</td></tr></tbody></table>

the effective resolving power ($\Delta\lambda/\lambda$) improved by at least a factor of 2.

## 2.2 DESI spectroscopic pipeline

The high complexity of the DESI survey induces the need for advanced software pipelines and products, including the imaging from the DESI Legacy Imaging Surveys (Zou et al. 2017; Dey et al. 2019; Schlegel et al., in preparation), a pipeline to select the targets to observe (Myers et al. 2022), a pipeline to assign fibres (Raichoor et al., in preparation), a pipeline to parse the survey and to optimize the observation strategy (Schlaflay et al. 2023), and an extensive spectroscopic reduction pipeline (Guy et al. 2022).

This spectroscopic pipeline, called desispec https://github.com/desihub/desispec, transforms the raw CCD images into spectra, and is detailed in Guy et al. (2022). Before extracting the spectra, the images are subtracted by dark and bias calibration frames to remove expected background sources, and to estimate the associated readout noise (noise estimation details are given in Appendix C). The non-uniform CCD pixel response is corrected using a dedicated flat-field slit on the spectrograph and the CCD over-scan is removed. A dedicated software detects and flags cosmic rays or defective CCD pixels.

The spectral extraction is performed using the 'spectroperfectionism' method (Bolton & Schlegel 2010), an optimal spectroscopic extraction that correctly models complex 2D point spread functions (PDFs). This method provides the encoding for each fibre and each wavelength of the non-Gaussian instrument resolution into a resolution matrix (noted $R$ in next sections) used to compute $P_{1D,\alpha}$.

All the spectra are defined on the same wavelength grid without additional resampling. Consequently, the extracted spectra are linearly binned in observed wavelength with a constant separation $\Delta\lambda_{pix} = 0.8$ Å. Conversely, the SDSS/BOSS spectrographs were logarithmically binned with a constant spectral pixel size of $\Delta v = 69$ km s$^{-1}$, which corresponds to $\Delta\lambda/\lambda = \Delta v/c = 2.3 \times 10^{-4}$, equivalent to $\Delta[\log(\lambda)] = 10^{-4}$ (Smee et al. 2013).

Once the spectrum of each fibre has been extracted, several post-processing steps allow removing a variety of further observational effects. The non-uniform response of individual fibres as a function of wavelength is corrected with flat-field frames by observing a white screen attached to the telescope dome and illuminated with a LED array. For all exposures, some fibres are dedicated to observing the sky. The so-called sky spectra associated with these fibres provide the sky level and the intensity of atmospheric emission lines and are subtracted to the spectra associated with targets (sky subtraction). The transmission defaults of the atmosphere and telescope as a whole are corrected by the observation of calibrated star spectra. This step converts CCD units (number of electrons) to observed flux units. Finally, the spectrum of an object is obtained by coadding its different exposures. The resulting spectrum is expressed separately into the three spectrographs bands described in Table 1.

All the software of the pipeline employed for the analysis of DESI data are listed in the repository desihub https://github.com/desihub. In particular, the spectra analysed in this

<table><caption>**Table 2.** Summary of the DESI data sets used in this study, their associated acronyms, a subjective description, and the total number of quasar spectra whose redshifts are between 2.0 and 5.0.</caption><thead><tr><td>Data set</td><td>Acronym</td><td>Description</td><td>Quasar number</td></tr></thead><tbody><tr><td>Target selection validation</td><td>SV1</td><td>Small sample with deeper exposures (up to 16) than the end of DESI</td><td>12 355</td></tr><tr><td>One-percent survey</td><td>SV3</td><td>Same density and depth as for the full 5 yr DESI main survey</td><td>12 686</td></tr><tr><td>First 2 months of the main survey</td><td>DESI-M2(M2)</td><td>Sparse data with small number of exposures</td><td>87 373</td></tr></tbody></table>

article have been processed with the fuji https://github.com/desihub/desispec/releases/tag/0.51.13. version of the spectroscopic pipeline.

## 2.3 DESI data

The spectroscopic pipeline described previously was used to analyse data obtained over different periods. In this article, we use three data sets from the first observations of DESI described in Table 2, which includes a total of 112 414 quasar spectra whose redshifts are between 2.0 and 5.0.

The first two data sets 'Target Selection Validation' (SV1) (DESI 2023a) and 'One-Percent Survey' (SV3) are part of the Early Data Release (EDR) whose complete description is given in DESI (2023b). 'Target Selection Validation' was conducted from 2020 December to 2021 March and includes a large number of exposures (up to 16) for the same targets. The objective of this survey was among others, to study extensively the survey performance as a function of instrumental depth and to build visual inspection truth tables. The SV1 data set also includes 'Secondary Tiles' as detailed in DESI (2023b). Following the completion of 'Target Selection Validation' observations, the 'One-Percent Survey' phase was dedicated more specifically to the evaluation of the survey design. The number of exposures is similar to the main survey at its end (between 4 and 5 for each Ly $\alpha$ quasars) and the goal was to determine the best strategy to cover the sky while limiting fibre loss.

The main DESI survey started in June 2021 and in this article we also use the first 2 months of data, named DESI-M2 (and noted M2 in this paper for conciseness), which is not present in the EDR but will be included in the Data Release 1. In the M2 data set, most quasars have only one exposure. While all three data sets are studied in this article, in the end we removed SV1 due to its different noise properties (see Appendix A). The final measurement is computed on SV3 + M2 data set only.

## 2.4 Input catalogues

The input catalogues used to compute $P_{1D,\alpha}$ were obtained by applying specific procedures to the three data sets previously described.

The targeting of the quasars used in our study (Yèche et al. 2020; Chaussidon et al. 2022) was verified with visual inspection of subsets of early observations (Alexander et al. 2022). For quasars, the DESI pipeline categorizes the observed spectra and estimates their redshifts using the redrock spectral template fitting software https://github.com/desihub/redrock (Bailey et al., in preparation; Brodzeller et al. 2023). In order to optimize, the completeness of the quasar catalogue while keeping a high purity, additionally a broad Mg II line finder mgii_afterburner
<figure><img src="image_2.png" /><figcaption>**Figure 1.** A particularly high-signal spectrum of a quasar located at a redshift $z = 3.42$ measured by DESI with an exposure time of 2300 s. This quasar was observed on 2021 April 12, in the SV3 programme, on DESI tile 221 (TARGETID = 39627746095137037, RA = 217.263°, Dec. = -1.755°). The quasar flux is represented in blue and its noise in orange. The Ly α forest is shown in green. The side-band regions 1 and 2 pictured in red and yellow are used to estimate the forest contamination by metals.</figcaption></figure>

<figure><img src="image_3.png" /><figcaption>**Figure 2.** Histogram of the quasar redshift (with $z > 2.0$) whose spectra are observed in the SV1, SV3, and M2 data sets. The histogram of the sum of data sets is shown in dashed line. Finally, the number of BAL quasars characterized by a balnicity index higher than zero, and not used in the $P_{1\text{D,α}}$ computation pipeline, is shown with a dotted line.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 3.** Normalized histogram of the effective exposure time in the Ly α forest region for the quasar spectra in the SV1, SV3, and M2 data sets. As mentioned in Section 2.3, there is a wide disparity of exposure time for the three data sets. As a reference, the nominal time of one DESI exposure is set to 1200 s.</figcaption></figure>

https://github.com/desihub/desispec/blob/main/py/desispec/mgii_afterburner.py and a machine learning classifier applying deep convolutional neural networks QuasarNP https://github.com/desihub/QuasarNP (Busca & Balland 2018; Farr, Font-Ribera & Pontzen 2020) are run after redrock. Both post-processing programs are run on all objects targeted by DESI as detailed in Chaussidon et al. (2022).

An example of a quasar spectrum at redshift $z = 3.42$ is given in Fig. 1. This spectrum, with a particularly high flux, is part of the SV3 data set that contains spectra with exposure time equivalent to the end of the DESI survey. The redshifts and effective exposure times for each data set considered in this paper are shown in Figs 2 and 3, respectively. The effective exposure time accounts for nightly observing conditions by normalizing the real exposure time to a reference with airmass 1, zero galactic extinction, a 1.1 arcsec seeing

(FWHM), and zenith dark sky (Guy et al. 2022). The nominal exposure time of one exposure is defined to 1000 s. The large differences in term of exposure time emphasize the need to treat the data sets differently, at least for noise properties. SV1 and SV3 contains a small number of forests but with heterogeneous exposure times, in opposition to M2, which contains many quasars only observed once.

Broad absorption line (BAL) quasars are specific quasars whose spectra exhibit consistent blueshifted absorptions associated with many spectral features. They are identified using the baltools https://github.com/paulmartini/baltools software. It consists of a $\chi^2$ minimizer algorithm that looks for blueshifted C IV or Si IV absorptions in an unabsorbed quasar model. The fit is performed for rest-frame wavelengths between 1260 and 2400 Å. A quasar is considered a BAL type if its spectrum exhibits a region between C IV
<figure><img src="image_3.png" /><figcaption>**Figure 1.** A particularly high-signal spectrum of a quasar located at a redshift $z = 3.42$ measured by DESI with an exposure time of 2300 s. This quasar was observed on 2021 April 12, in the SV3 programme, on DESI tile 221 (TARGETID = 39627746095137037, RA = 217.263°, Dec. = -1.755°). The quasar flux is represented in blue and its noise in orange. The Ly α forest is shown in green. The side-band regions 1 and 2 pictured in red and yellow are used to estimate the forest contamination by metals.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Histogram of the quasar redshift (with $z > 2.0$) whose spectra are observed in the SV1, SV3, and M2 data sets. The histogram of the sum of data sets is shown in dashed line. Finally, the number of BAL quasars characterized by a balnicity index higher than zero, and not used in the $P_{1\text{D,α}}$ computation pipeline, is shown with a dotted line.</figcaption></figure>

<figure><img src="image_5.png" /><figcaption>**Figure 3.** Normalized histogram of the effective exposure time in the Ly α forest region for the quasar spectra in the SV1, SV3, and M2 data sets. As mentioned in Section 2.3, there is a wide disparity of exposure time for the three data sets. As a reference, the nominal time of one DESI exposure is set to 1200 s.</figcaption></figure>

https://github.com/desihub/desispec/blob/main/py/desispec/mgii_afterburner.py and a machine learning classifier applying deep convolutional neural networks QuasarNP https://github.com/desihub/QuasarNP (Busca & Balland 2018; Farr, Font-Ribera & Pontzen 2020) are run after redrock. Both post-processing programs are run on all objects targeted by DESI as detailed in Chaussidon et al. (2022).

An example of a quasar spectrum at redshift $z = 3.42$ is given in Fig. 1. This spectrum, with a particularly high flux, is part of the SV3 data set that contains spectra with exposure time equivalent to the end of the DESI survey. The redshifts and effective exposure times for each data set considered in this paper are shown in Figs 2 and 3, respectively. The effective exposure time accounts for nightly observing conditions by normalizing the real exposure time to a reference with airmass 1, zero galactic extinction, a 1.1 arcsec seeing

(FWHM), and zenith dark sky (Guy et al. 2022). The nominal exposure time of one exposure is defined to 1000 s. The large differences in term of exposure time emphasize the need to treat the data sets differently, at least for noise properties. SV1 and SV3 contains a small number of forests but with heterogeneous exposure times, in opposition to M2, which contains many quasars only observed once.

Broad absorption line (BAL) quasars are specific quasars whose spectra exhibit consistent blueshifted absorptions associated with many spectral features. They are identified using the baltools https://github.com/paulmartini/baltools software. It consists of a $\chi^2$ minimizer algorithm that looks for blueshifted C IV or Si IV absorptions in an unabsorbed quasar model. The fit is performed for rest-frame wavelengths between 1260 and 2400 Å. A quasar is considered a BAL type if its spectrum exhibits a region between C IV
and Si IV emission lines, with at least 10 per cent flux decrement below the continuum and a width greater than 2000 km s⁻¹. BAL quasars spectra (4.18 per cent in the total data set) are removed in the analysis performed in this work.

Damped Ly α absorbers (DLAs) are regions within a quasar spectrum that show oversaturated absorption with prominent Lorentzian wings as the quasar flux intersects the dense, circumgalactic medium of an intervening (proto-)galaxy. They are a subclass of high-column density systems and are a significant contaminant of the Ly α forest signal, particularly because of their wings and the additional contamination by circumgalactic metal absorption lines (McDonald et al. 2005). The correct modeling of such systems in simulations has been proven to be particularly complicated (Pontzen et al. 2008).

We use a catalogue resulting from the combination of a convolutional neural network (CNN) algorithm desi-dlas https://github.com/cosmodesi/desi-dlas (Parks et al. 2017; Zou et al., in preparation) and a Gaussian process (GP) algorithm https://github.com/jibanCat/gp.dla.detection_dr16q-public (Ho, Bird & Garnett 2021). desi-dlas is trained with SDSS spectra to identify candidate high-column density objects for rest-frame wavelengths between 900 and 1346 Å. It returns locations of high-column density systems in the spectra, as well as their H I column density and a confidence parameter. The GP finder provides similar output using the same training set and a Bayesian model selection. We only consider the high-column density objects with column density $N_{HI} > 10^{20.3} \text{ cm}^{-2}$ (DLAs). In accordance with recommendations from Parks et al. (2017) and Ho et al. (2021), we consider CNN confidence level higher than 0.2 as valid DLA detections when the ratio between the quasar continuum and the noise is higher than 3. We take a confidence level limit of 0.3 when this ratio is lower than 3. For the GP model, a 0.9 minimal confidence level is applied. In the case when absorbers are detected by the two models, the combined DLA catalogue uses $N_{HI}$ values and DLA redshifts from GP model.

Although DLAs by themselves constitute tracers of the matter distribution, they have an extended impact on the observed spectra. They increase the correlations of neighbour spectrum pixels, thus artificially increasing $P_{1D,α}$ level. Therefore, we choose to mask the core of DLA regions of the spectra by fixing the transmitted flux fraction to its mean value for spectrum pixels, where the DLA-induced absorption is larger than 20 per cent. In addition, the absorption in the Lorentzian damping wings that remain after the cut is corrected with a Voigt profile following Bautista et al. (2017) and Chabanier et al. (2021).

Finally, we use a catalogue of masks to account for atmospheric and Galactic emission lines, which has been adapted to DESI resolution https://github.com/coretinravoux/p1desi/blob/main/etc/skylines/list_mask.p1d_DESI.EDR.txt. The creation of this catalogue is detailed in Section 4.4.

## 3 1D POWER-SPECTRUM ESTIMATION

The $P_{1D,α}$ estimator is build using the data product described previously in two phases. First, the fitting of the continuum of quasars is used to convert the absolute received flux to a normalized quantity $\delta_F$. Secondly, $P_{1D,α}$ is computed by employing an FFT and by averaging the product of this transformation for all the selected Ly α forests.

### 3.1 Continuum fitting

A standard normalized quantity used in the calculation of correlations and power spectra is the flux contrast $\delta_F$ of the Ly α forest, defined

as

$$ \delta_F(\lambda) = \frac{F(\lambda)}{\bar{F}(\lambda)} - 1 = \frac{f(\lambda)}{C_q(\lambda, z_q)\bar{F}(\lambda)} - 1, \quad (1) $$

where F is the transmitted flux fraction, and $\bar{F}(\lambda)$ is its average value, the mean transmission of the intervening IGM. Note that for the purpose of this work, we do not need to know the individual quasar continua $C_q$, but only the product $C_q(\lambda, z_q)\bar{F}(\lambda)$. Thus, similar to previous Ly α studies based on survey data (Chabanier et al. 2019; du Mas des Bourboux et al. 2020; Ramírez-Pérez et al. 2023), we directly measure this product in our continuum fitting process using the picca https://github.com/igmhub/picca. (du Mas des Bourboux et al. 2021) software package. This software also merges the quasar spectra over different bands to obtain f over the all wavelength range. The continuum of each quasar is modelled as the product of a universal continuum C common to all quasars, and a first-order polynomial term in wavelength:

$$ C_q(\lambda, z_q) = (a_q + b_q\lambda) C \left( \lambda_{rf} = \frac{\lambda}{(1 + z_q)} \right), \quad (2) $$

where $a_q$ and $b_q$ are quasar-dependent constants. In previous studies (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), $b_q = 0$ was assumed for all quasars, i.e. only a wavelength-independent normalization factor was taken into account. We add an additional linear wavelength-dependent term to account for the diversity of quasars after verifying that this change does not impact the mean level of our $P_{1D,α}$ measurement.

The $a_q$ and $b_q$ parameters for each quasar are determined along with C by maximizing the following log-likelihood:

$$ \ln \mathcal{L} = -\frac{1}{2} \sum_i \frac{[f_i - \bar{F}(\lambda_i)C_q(\lambda_i, z_q, a_q, b_q)]^2}{\sigma_q^2(\lambda_i)} - \ln[\sigma_q^2(\lambda_i)], \quad (3) $$

where the sum is run over all the spectrum pixels of the quasar q, and $\sigma_q$ is the standard deviation estimator of the flux f.

In contrast to analyses of the large-scale 3D correlation function such as Bourboux et al. (2020), we want all spectrum pixels to contribute equally to the continuum fitting and the $P_{1D,α}$ computation, as the opposite could bias $P_{1D,α}$. Therefore, we impose noise-independent weights in the continuum fitting procedure:

$$ \sigma_q^2(\lambda) = (\bar{F}(\lambda)C_q(\lambda))^2 \quad (4) $$

This procedure is the same as in the previous $P_{1D,α}$ analyses based on BOSS/eBOSS data (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019). The standard deviation associated with $\delta_F$ at the end of the continuum fitting procedure is defined by

$$ \sigma_{\delta_F}(\lambda) = \frac{\sigma_{\text{pip,q}}(\lambda)}{\bar{F}(\lambda)C_q(\lambda)}, \quad (5) $$

where $\sigma_{\text{pip,q}}$ is the noise provided by the DESI spectroscopic pipeline detailed in Appendix C.

The universal continuum C and the $a_q$ and $b_q$ parameters are computed iteratively. In particular, C is estimated from the average of all spectra, i.e. in a non-parametric way. During the entire fitting procedure, spectrum pixels that are masked due to the presence of a DLA or an atmospheric line are not considered in the fit. We use 7 iterations and have verified that the continuum fits are converged at this point.
and Si IV emission lines, with at least 10 per cent flux decrement below the continuum and a width greater than 2000 km s⁻¹. BAL quasars spectra (4.18 per cent in the total data set) are removed in the analysis performed in this work.

Damped Ly α absorbers (DLAs) are regions within a quasar spectrum that show oversaturated absorption with prominent Lorentzian wings as the quasar flux intersects the dense, circumgalactic medium of an intervening (proto-)galaxy. They are a subclass of high-column density systems and are a significant contaminant of the Ly α forest signal, particularly because of their wings and the additional contamination by circumgalactic metal absorption lines (McDonald et al. 2005). The correct modeling of such systems in simulations has been proven to be particularly complicated (Pontzen et al. 2008).

We use a catalogue resulting from the combination of a convolutional neural network (CNN) algorithm desi-dlas https://github.com/cosmodesi/desi-dlas (Parks et al. 2017; Zou et al., in preparation) and a Gaussian process (GP) algorithm https://github.com/jibanCat/gp.dla.detection_dr16q-public (Ho, Bird & Garnett 2021). desi-dlas is trained with SDSS spectra to identify candidate high-column density objects for rest-frame wavelengths between 900 and 1346 Å. It returns locations of high-column density systems in the spectra, as well as their H I column density and a confidence parameter. The GP finder provides similar output using the same training set and a Bayesian model selection. We only consider the high-column density objects with column density $N_{HI} > 10^{20.3} \text{ cm}^{-2}$ (DLAs). In accordance with recommendations from Parks et al. (2017) and Ho et al. (2021), we consider CNN confidence level higher than 0.2 as valid DLA detections when the ratio between the quasar continuum and the noise is higher than 3. We take a confidence level limit of 0.3 when this ratio is lower than 3. For the GP model, a 0.9 minimal confidence level is applied. In the case when absorbers are detected by the two models, the combined DLA catalogue uses $N_{HI}$ values and DLA redshifts from GP model.

Although DLAs by themselves constitute tracers of the matter distribution, they have an extended impact on the observed spectra. They increase the correlations of neighbour spectrum pixels, thus artificially increasing $P_{1D,α}$ level. Therefore, we choose to mask the core of DLA regions of the spectra by fixing the transmitted flux fraction to its mean value for spectrum pixels, where the DLA-induced absorption is larger than 20 per cent. In addition, the absorption in the Lorentzian damping wings that remain after the cut is corrected with a Voigt profile following Bautista et al. (2017) and Chabanier et al. (2021).

Finally, we use a catalogue of masks to account for atmospheric and Galactic emission lines, which has been adapted to DESI resolution https://github.com/coretinravoux/p1desi/blob/main/etc/skylines/list_mask.p1d_DESI.EDR.txt. The creation of this catalogue is detailed in Section 4.4.

## 3 1D POWER-SPECTRUM ESTIMATION

The $P_{1D,α}$ estimator is build using the data product described previously in two phases. First, the fitting of the continuum of quasars is used to convert the absolute received flux to a normalized quantity $\delta_F$. Secondly, $P_{1D,α}$ is computed by employing an FFT and by averaging the product of this transformation for all the selected Ly α forests.

### 3.1 Continuum fitting

A standard normalized quantity used in the calculation of correlations and power spectra is the flux contrast $\delta_F$ of the Ly α forest, defined

as

$$ \delta_F(\lambda) = \frac{F(\lambda)}{\bar{F}(\lambda)} - 1 = \frac{f(\lambda)}{C_q(\lambda, z_q)\bar{F}(\lambda)} - 1, \quad (1) $$

where F is the transmitted flux fraction, and $\bar{F}(\lambda)$ is its average value, the mean transmission of the intervening IGM. Note that for the purpose of this work, we do not need to know the individual quasar continua $C_q$, but only the product $C_q(\lambda, z_q)\bar{F}(\lambda)$. Thus, similar to previous Ly α studies based on survey data (Chabanier et al. 2019; du Mas des Bourboux et al. 2020; Ramírez-Pérez et al. 2023), we directly measure this product in our continuum fitting process using the picca https://github.com/igmhub/picca. (du Mas des Bourboux et al. 2021) software package. This software also merges the quasar spectra over different bands to obtain f over the all wavelength range. The continuum of each quasar is modelled as the product of a universal continuum C common to all quasars, and a first-order polynomial term in wavelength:

$$ C_q(\lambda, z_q) = (a_q + b_q\lambda) C \left( \lambda_{rf} = \frac{\lambda}{1 + z_q} \right), \quad (2) $$

where $a_q$ and $b_q$ are quasar-dependent constants. In previous studies (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), $b_q = 0$ was assumed for all quasars, i.e. only a wavelength-independent normalization factor was taken into account. We add an additional linear wavelength-dependent term to account for the diversity of quasars after verifying that this change does not impact the mean level of our $P_{1D,α}$ measurement.

The $a_q$ and $b_q$ parameters for each quasar are determined along with C by maximizing the following log-likelihood:

$$ \ln \mathcal{L} = -\frac{1}{2} \sum_i \frac{[f_i - \bar{F}(\lambda_i)C_q(\lambda_i, z_q, a_q, b_q)]^2}{\sigma_q^2(\lambda_i)} - \ln[\sigma_q^2(\lambda_i)], \quad (3) $$

where the sum is run over all the spectrum pixels of the quasar q, and $\sigma_q$ is the standard deviation estimator of the flux f.

In contrast to analyses of the large-scale 3D correlation function such as Bourboux et al. (2020), we want all spectrum pixels to contribute equally to the continuum fitting and the $P_{1D,α}$ computation, as the opposite could bias $P_{1D,α}$. Therefore, we impose noise-independent weights in the continuum fitting procedure:

$$ \sigma_q^2(\lambda) = (\bar{F}(\lambda)C_q(\lambda))^2 \quad (4) $$

This procedure is the same as in the previous $P_{1D,α}$ analyses based on BOSS/eBOSS data (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019). The standard deviation associated with $\delta_F$ at the end of the continuum fitting procedure is defined by

$$ \sigma_{\delta_F}(\lambda) = \frac{\sigma_{\text{pip,q}}(\lambda)}{\bar{F}(\lambda)C_q(\lambda)}, \quad (5) $$

where $\sigma_{\text{pip,q}}$ is the noise provided by the DESI spectroscopic pipeline detailed in Appendix C.

The universal continuum C and the $a_q$ and $b_q$ parameters are computed iteratively. In particular, C is estimated from the average of all spectra, i.e. in a non-parametric way. During the entire fitting procedure, spectrum pixels that are masked due to the presence of a DLA or an atmospheric line are not considered in the fit. We use 7 iterations and have verified that the continuum fits are converged at this point.
The noise level of a Ly α forest is characterized by defining the average signal-to-noise ratio SNR in the Ly α forest region:

$$ \overline{\text{SNR}} = \left\langle \frac{f(\lambda)}{\sigma_{\text{pip,q}}(\lambda)} \right\rangle_{\lambda}. \qquad (6) $$

Only Ly α forests with a SNR larger than 1 are used in the continuum fitting procedure. This procedure is also restricted to the observed redshift range 3600 Å < λ < 7600 Å, to avoid the shorter wavelength range where a large fraction of the quasar spectra is absorbed by the atmosphere. We also select the rest-frame wavelength in the range 1050 Å < λ<rf < 1180 Å, so that the measured contrasts are dominated by the Ly α forest. In particular, we try to avoid the Ly β singlet and the OVI doublet emission regions respectively located at λ<Ly β = 1025.72 Å and λ<Ovi = (1031.912, 1037.613) Å in the rest frame.

The cut λ<rf < 1180 Å facilitates the continuum fitting procedure and mitigates most of the proximity effect: close to a quasar, the neutral hydrogen fraction is indeed influenced by the quasar’s UV radiation in addition to the extragalactic UV radiation background (Bajtlik, Duncan & Ostriker 1988).

As detailed in Section 2.2, the quasar spectra are linearly binned in observed wavelength with Δλ<pix = 0.8 Å. Note that when converting to rest-frame wavelength λ<rf the pixel size will be redshift dependent. For the continuum fitting process, we thus need to rebin our spectra to a uniform grid in λ<rf. As the quasar continuum is relatively smooth and to avoid noisy continuum fits for analyses of relatively small data sets, we chose a grid for the common continuum C that is 10 times coarser than the lowest redshift quasar pixels considered, i.e.:

$$ \Delta\lambda_{\text{pix,rf}} = 10 \frac{\Delta\lambda_{\text{pix}}}{1 + z_{\min}}. \qquad (7) $$

By taking z<min = 2.0 as the lowest redshift, we obtain Δλ<pix,rf = 2.67 Å. With increasing size of the data set, such a rebinning could be relaxed in future DESI measurements. At the end of the continuum fitting procedure, the stacking of all the Lyα contrasts is forced to be equal to zero to avoid introducing flux calibration errors.

## 3.2 FFT power spectrum estimator

Conceptually, the δ<sub>F</sub> quantity can be separated into different contributions in the Ly α forest region:

$$ \delta_F(\lambda) = \delta_{\text{astro}}(\lambda) + \delta_{\text{noise}}(\lambda), \qquad (8) $$

where δ<sub>astro</sub> corresponds to the fluctuations caused by all the elements of the intergalactic medium (including Ly α), and δ<sub>noise</sub> corresponds to noise fluctuations.

Considering instrumental effects, the true underlying flux contrast is modified on its way through the instrument in multiple ways. First, photons traverse the spectrograph leading to the output flux being convolved with the spectrograph line spread function W(λ, **R**, Δλ<pix) which depends on the resolution matrix presented in Section 2.2. The W term also account for the signal pixelization as the photons are counted into CCD pixels of size Δλ<pix. Finally, we need to account for noise sourcing from the processes of photon counting and readout. In total, we can write the measured flux contrast δ<sub>F</sub> as

$$ \delta_F(\lambda) = \delta_{\text{astro}}(\lambda) \oplus W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda). \qquad (9) $$

We assume that the impact of the noise and the resolution term W are decorrelated, i.e. that the noise contrast is not affected by the instrumental effects of pixelization and resolution. Verifying this assumption is beyond the scope of this paper.

The δ<sub>astro</sub> contrast can be decomposed between the Ly α signal and the one from all other elements of the intergalactic medium called metals. The contribution of those metal absorptions to δ<sub>F</sub> can be decomposed into two parts. On the one hand, there are absorption lines with rest-frame wavelength λ ≫ λ<Ly α = 1215.67 Å. Those lines can be independently observed redwards of the Ly α forest, using specific rest-frame spectral regions called side-bands SB1 (1270 < λ<rf < 1380 Å), and SB2 (1410 < λ<rf < 1520 Å), as shown in Fig. 1. We group the absorption of all those metals in a contrast noted δ<sub>metals</sub>.

On the other hand, there are absorption lines with rest-frame wavelength λ ≲ λ<Ly α, such as Si II and Si III elements (λ<SiII = 1, 190 and 1193 Å, and λ<SiIII = 1, 206.50 Å). They cannot be observed independently of the Ly α forest, but will show absorption that is correlated with the Ly α absorption and will lead to an oscillatory feature in the estimated P<sub>1D,α</sub>. For those lines, we adopt the same approach developed in McDonald et al. (2006) and subsequently used in other analyses. We leave the features inside our power spectrum estimates to be corrected by fitting an additional oscillation during parameter inference. We note δ<sub>Ly α</sub> the contrast containing the Ly α forest and the effect of those latter metals. Finally, the flux contrast can be expressed:

$$ \delta_F(\lambda) = (\delta_{\text{Ly}\alpha}(\lambda) + \delta_{\text{metals}}(\lambda)) \oplus W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda). \qquad (10) $$

The 1D Ly α power spectrum (P<sub>1D,α</sub>) can be estimated from this decomposition by applying an FFT algorithm on δ<sub>F</sub> of each Ly α forest. This method was applied in previous measurements (Croft et al. 1998; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019).

The FFT estimator is implemented in picca (Bourboux et al. 2021). For each Ly α forest, the raw power spectrum<fnref n="1" /> is defined from the Fourier transform of δ<sub>F</sub>(r):

$$ (2\pi)\delta_D(k-k')P_{\text{raw}}(k) = \delta_F(k)\delta_F(k'), \qquad (11) $$

where δ<sub>D</sub> is the 1D Dirac distribution.

Applying a Fourier transform on equation (10), the raw power spectrum is expressed by

$$ P_{\text{raw}}(k) = (P_{\text{Ly}\alpha}(k) + P_{\text{Ly}\alpha-\text{Si II/Si III}}(k) + P_{\text{metals}}(k)) \cdot W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}}) + P_{\text{noise}}(k). \qquad (12) $$

In this decomposition, P<sub>Ly α</sub>, P<sub>metals</sub>, and P<sub>noise</sub> are the power spectra, respectively, associated with the contrasts δ<sub>Ly α</sub>, δ<sub>metals</sub>, and δ<sub>noise</sub> with the same definition as in equation (11). We assumed here that the SiII and SiIII power spectra are negligible and that all cross-correlation terms between the contrasts are null. The only non-neglected cross-term is P<sub>Ly α-Si II/Si III</sub>(k) = 2|δ<sub>Ly α</sub>(k)δ<sub>Si II/Si III</sub>(k)|, which corresponds to the correlated absorptions of Ly α with either Si II or Si III. The oscillations induced by this term have a wavenumber 2π/(λ<Ly α - λ<Si II/Si III>) when 'k' is expressed in Å⁻<fnref n="1" />.

The FFT estimator for the 1D power spectrum is computed as an average over all available Ly α forests in the measurement sample. It is designed to match the sum of Ly α and Ly α-Si II/Si III power spectra in equation (12), so that we define

$$ P_{1D,\alpha}(k) = \langle P_{\text{Ly}\alpha}(k) + P_{\text{Ly}\alpha-\text{Si II/Si III}}(k) \rangle, \qquad (13) $$

where 〈·〉 denotes the average over all the Ly α forests used for P<sub>1D,α</sub> calculation. From equation (12), the estimator of P<sub>1D,α</sub> is defined by

$$ P_{1D,\alpha}(k) = \left\langle \frac{P_{\text{raw}}(k) - P_{\text{noise}}(k)}{W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}})} \right\rangle - P_{\text{metals}}(k). \qquad (14) $$

<footnote n="1">¹This is not a proper power spectrum, strictly speaking, as it does not correspond to an average.</footnote>
The noise level of a Ly α forest is characterized by defining the average signal-to-noise ratio SNR in the Ly α forest region:

$$ \overline{\text{SNR}} = \left\langle \frac{f(\lambda)}{\sigma_{\text{pip,q}}(\lambda)} \right\rangle_{\lambda}. \qquad (6) $$

Only Ly α forests with a SNR larger than 1 are used in the continuum fitting procedure. This procedure is also restricted to the observed redshift range 3600 Å < λ < 7600 Å, to avoid the shorter wavelength range where a large fraction of the quasar spectra is absorbed by the atmosphere. We also select the rest-frame wavelength in the range 1050 Å < λ<rf < 1180 Å, so that the measured contrasts are dominated by the Ly α forest. In particular, we try to avoid the Ly β singlet and the OVI doublet emission regions respectively located at λ<Ly β = 1025.72 Å and λ<Ovi = (1031.912, 1037.613) Å in the rest frame.

The cut λ<rf < 1180 Å facilitates the continuum fitting procedure and mitigates most of the proximity effect: close to a quasar, the neutral hydrogen fraction is indeed influenced by the quasar’s UV radiation in addition to the extragalactic UV radiation background (Bajtlik, Duncan & Ostriker 1988).

As detailed in Section 2.2, the quasar spectra are linearly binned in observed wavelength with Δλ<pix = 0.8 Å. Note that when converting to rest-frame wavelength λ<rf the pixel size will be redshift dependent. For the continuum fitting process, we thus need to rebin our spectra to a uniform grid in λ<rf. As the quasar continuum is relatively smooth and to avoid noisy continuum fits for analyses of relatively small data sets, we chose a grid for the common continuum C that is 10 times coarser than the lowest redshift quasar pixels considered, i.e.:

$$ \Delta\lambda_{\text{pix,rf}} = 10 \frac{\Delta\lambda_{\text{pix}}}{1 + z_{\min}}. \qquad (7) $$

By taking z<min = 2.0 as the lowest redshift, we obtain Δλ<pix,rf = 2.67 Å. With increasing size of the data set, such a rebinning could be relaxed in future DESI measurements. At the end of the continuum fitting procedure, the stacking of all the Lyα contrasts is forced to be equal to zero to avoid introducing flux calibration errors.

## 3.2 FFT power spectrum estimator

Conceptually, the δ<sub>F</sub> quantity can be separated into different contributions in the Ly α forest region:

$$ \delta_F(\lambda) = \delta_{\text{astro}}(\lambda) + \delta_{\text{noise}}(\lambda), \qquad (8) $$

where δ<sub>astro</sub> corresponds to the fluctuations caused by all the elements of the intergalactic medium (including Ly α), and δ<sub>noise</sub> corresponds to noise fluctuations.

Considering instrumental effects, the true underlying flux contrast is modified on its way through the instrument in multiple ways. First, photons traverse the spectrograph leading to the output flux being convolved with the spectrograph line spread function W(λ, **R**, Δλ<pix) which depends on the resolution matrix presented in Section 2.2. The W term also account for the signal pixelization as the photons are counted into CCD pixels of size Δλ<pix. Finally, we need to account for noise sourcing from the processes of photon counting and readout. In total, we can write the measured flux contrast δ<sub>F</sub> as

$$ \delta_F(\lambda) = \delta_{\text{astro}}(\lambda) \oplus W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda). \qquad (9) $$

We assume that the impact of the noise and the resolution term W are decorrelated, i.e. that the noise contrast is not affected by the instrumental effects of pixelization and resolution. Verifying this assumption is beyond the scope of this paper.

The δ<sub>astro</sub> contrast can be decomposed between the Ly α signal and the one from all other elements of the intergalactic medium called metals. The contribution of those metal absorptions to δ<sub>F</sub> can be decomposed into two parts. On the one hand, there are absorption lines with rest-frame wavelength λ ≫ λ<Ly α = 1215.67 Å. Those lines can be independently observed redwards of the Ly α forest, using specific rest-frame spectral regions called side-bands SB1 (1270 < λ<rf < 1380 Å), and SB2 (1410 < λ<rf < 1520 Å), as shown in Fig. 1. We group the absorption of all those metals in a contrast noted δ<sub>metals</sub>.

On the other hand, there are absorption lines with rest-frame wavelength λ ≲ λ<Ly α, such as Si II and Si III elements (λ<SiII = 1, 190 and 1193 Å, and λ<SiIII = 1, 206.50 Å). They cannot be observed independently of the Ly α forest, but will show absorption that is correlated with the Ly α absorption and will lead to an oscillatory feature in the estimated P<sub>1D,α</sub>. For those lines, we adopt the same approach developed in McDonald et al. (2006) and subsequently used in other analyses. We leave the features inside our power spectrum estimates to be corrected by fitting an additional oscillation during parameter inference. We note δ<sub>Ly α</sub> the contrast containing the Ly α forest and the effect of those latter metals. Finally, the flux contrast can be expressed:

$$ \delta_F(\lambda) = (\delta_{\text{Ly}\alpha}(\lambda) + \delta_{\text{metals}}(\lambda)) \oplus W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda). \qquad (10) $$

The 1D Ly α power spectrum (P<sub>1D,α</sub>) can be estimated from this decomposition by applying an FFT algorithm on δ<sub>F</sub> of each Ly α forest. This method was applied in previous measurements (Croft et al. 1998; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019).

The FFT estimator is implemented in picca (Bourboux et al. 2021). For each Ly α forest, the raw power spectrum<fnref n="1" /> is defined from the Fourier transform of δ<sub>F</sub>(r):

$$ (2\pi)\delta_D(k-k')P_{\text{raw}}(k) = \delta_F(k)\delta_F(k'), \qquad (11) $$

where δ<sub>D</sub> is the 1D Dirac distribution.

Applying a Fourier transform on equation (10), the raw power spectrum is expressed by

$$ P_{\text{raw}}(k) = (P_{\text{Ly}\alpha}(k) + P_{\text{Ly}\alpha-\text{Si II/Si III}}(k) + P_{\text{metals}}(k)) \cdot W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}}) + P_{\text{noise}}(k). \qquad (12) $$

In this decomposition, P<sub>Ly α</sub>, P<sub>metals</sub>, and P<sub>noise</sub> are the power spectra, respectively, associated with the contrasts δ<sub>Ly α</sub>, δ<sub>metals</sub>, and δ<sub>noise</sub> with the same definition as in equation (11). We assumed here that the SiII and SiIII power spectra are negligible and that all cross-correlation terms between the contrasts are null. The only non-neglected cross-term is P<sub>Ly α-Si II/Si III</sub>(k) = 2 |δ<sub>Ly α</sub>(k)δ<sub>Si II/Si III</sub>(k)|, which corresponds to the correlated absorptions of Ly α with either Si II or Si III. The oscillations induced by this term have a wavenumber 2π/(λ<Ly α - λ<Si II/Si III)) when 'k' is expressed in Å⁻<fnref n="1" />.

The FFT estimator for the 1D power spectrum is computed as an average over all available Ly α forests in the measurement sample. It is designed to match the sum of Ly α and Ly α-Si II/Si III power spectra in equation (12), so that we define

$$ P_{1D,\alpha}(k) = \langle P_{\text{Ly}\alpha}(k) + P_{\text{Ly}\alpha-\text{Si II/Si III}}(k) \rangle, \qquad (13) $$

where 〈·〉 denotes the average over all the Ly α forests used for P<sub>1D,α</sub> calculation. From equation (12), the estimator of P<sub>1D,α</sub> is defined by

$$ P_{1D,\alpha}(k) = \left\langle \frac{P_{\text{raw}}(k) - P_{\text{noise}}(k)}{W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}})} \right\rangle - P_{\text{metals}}(k). \qquad (14) $$

<footnote n="1">¹This is not a proper power spectrum, strictly speaking, as it does not correspond to an average.</footnote>
The $P_{1D,\alpha}$ measurement is split in different redshift bins to take into account its evolution. The Ly $\alpha$ forest is split into subforests that correspond to consecutive and non-overlapping subregions of equal length. This procedure also reduces the correlations between the different redshift bins. We chose to cut Ly $\alpha$ forests into three subforests whose rest-frame wavelength boundaries are $\lambda_{rf} = 1093.3$ and 1136.6 Å, so that the length of each subforest is $L_{sub} = 43.3$ Å. With this subforest separation, a single Ly $\alpha$ forest can contribute to up to three different redshift bins in the $P_{1D,\alpha}$ measurement. The subforest splitting constrains the minimal accessible observed wavenumber to $k_{min} = 2\pi/(L_{sub}(1+z_{min})) = 0.0453$ Å$^{-1}$, by taking the minimal redshift used. Each subforest spans at most $\Delta z = 0.2$ and we choose the same $\Delta z$ to define the redshift binning for $P_{1D,\alpha}$.

For observed wavelength $\lambda \lesssim 3,700$ Å, the noise level is high in comparison to the spectra because of atmospheric absorptions. To minimize the impact of this noise, we remove the spectrum pixels for which the observed wavelength is lower than 3750 Å, which corresponds to Ly $\alpha$ absorbers located at $z = 2.085$. In the future, with a dedicated study to control the noise at shorter wavelength, the $P_{1D,\alpha}$ analysis can be extended to redshift $z \sim 2$.

In accordance with the eBOSS study (Chabanier et al. 2019), we remove subforests shorter than 75 spectrum pixels due to a cut in the UV region or to the presence of a large DLA. We also do not consider the Ly $\alpha$ subforests with more than 120 masked spectrum pixels.

Unlike the analysis in Chabanier et al. (2019), we do not apply a second redshift-dependent SNR cut for the averaging of $P_{1D,\alpha}$. Instead, we develop and test a SNR weighting scheme, as detailed in the appendix B. This procedure is used for all the article except in Section 4.2.1 where the impact of the SNR cut is investigated.

# 4 DESI INSTRUMENTAL CHARACTERIZATION

As the DESI instrument is new, we first focus the analysis on characterizing the instrumental effects on our $P_{1D,\alpha}$ measurement. In particular, we describe the impact of spectral resolution, instrumental noise, metal power spectrum, and atmospheric emission lines in the following.

## 4.1 Spectrograph resolution

### 4.1.1 Resolution correction modelling

As mentioned in Section 2.2, DESI spectra are linearly binned in observed wavelength, so that the natural unit for wavenumbers is Å$^{-1}$. The maximal measurable wavenumber follows the Nyquist-Shannon limit: $k_{max} = k_{Nyq} = \pi/\Delta\lambda_{pix} \approx 3.92$ Å$^{-1}$.

For SDSS $P_{1D,\alpha}$ measurements (Croft et al. 1998; McDonald et al. 2006; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), the pixelization is logarithmically binned in observed wavelength, making it suitable to express the power spectrum in Hubble velocity unit $v$, because $\Delta v \propto \Delta\lambda/\lambda = \Delta\log(\lambda)$. In eBOSS (Chabanier et al. 2019), the resolution correction is modelled by a Gaussian function $W(k, \Delta\lambda) = \exp(-0.5(k[\text{km}^{-1}]\Delta v)^2)$, where $\Delta v$ is the spectral resolution in velocity units.

The spectroscopic resolution of DESI is improved with respect to the SDSS spectrographs. On the DESI blue band (see Table 1) where most of the Ly $\alpha$ forest are observed, the effective resolving power $R = \lambda/\Delta\lambda$ ranges from 2000 to 3200 (Abareshi et al. 2022). In comparison, SDSS spectrographs had a $1500 < R < 2300$ in its blue band ($3600 < \lambda < 6350$ Å). This improved resolution brings the opportunity to probe the clustering of matter at smaller scale by

<figure><img src="image_5.png" /><figcaption>**Figure 4.** Illustration of the resolution improvement between eBOSS and DESI. The DESI spectrum of the quasar represented in Fig. 1 is zoomed on a region of the Ly $\alpha$ forest on the top panel. The spectrum of the same quasar obtained by the SDSS-IV (Gunn et al. 2006; Smee et al. 2013; Dawson et al. 2016; Blanton et al. 2017) in the 16th Data Release (DR16; Ahumada et al. 2020) of the eBOSS. The SDSS name associated with this quasar is 142 903.03–014 519.3. The DESI spectrum is obtained after three individual exposures for a total exposure time of 2300 s. The eBOSS spectrum have 11 exposures for a total of 6300 s. The large-scale absorption structures are similar but due to its improved spectroscopic resolution, the DESI spectrum clearly shows more details at small scales.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 5.** Average resolution correction ($W^2(k, \mathbf{R}, \Delta\lambda_{pix}) = \langle R^2(k) \rangle$) (blue points) using Ly $\alpha$ forest from SV3 + M2 data set. This resolution correction is weakly dependent on the redshift range (shown by shaded dashed black lines along the blue points). Only one shaded black line is above the blue points ($z = 3.8$). All the others are at the same level as the blue points. The mean value over all redshift bins is shown with the points. The shaded area represents the regime for which the associated impact of resolution and pixelization removes more than 80 per cent of the power spectrum. This criteria is chosen to define the maximal wavenumber of our $P_{1D,\alpha}$ measurement, shown with a vertical black line.</figcaption></figure>

measuring the small fluctuations in the Ly $\alpha$ forest, as illustrated in Fig. 4.

As described in Section 2.2, the DESI spectrograph resolution is entirely characterized by the resolution matrix $\mathbf{R}$ (Guy et al. 2022). In opposition to SDSS, the resolution matrix also accounts for the pixelization of the signal. Consequently, we choose to express the resolution correction function $W$ in equation (12) directly as the Fourier transform of the resolution matrix that we note $\mathbf{R}(k)$.

Fig. 5 shows the average correction due to resolution and pixelization for the SV3 + M2 data set. This correction indicates that
The $P_{1D,\alpha}$ measurement is split in different redshift bins to take into account its evolution. The Ly $\alpha$ forest is split into subforests that correspond to consecutive and non-overlapping subregions of equal length. This procedure also reduces the correlations between the different redshift bins. We chose to cut Ly $\alpha$ forests into three subforests whose rest-frame wavelength boundaries are $\lambda_{rf} = 1093.3$ and 1136.6 Å, so that the length of each subforest is $L_{sub} = 43.3$ Å. With this subforest separation, a single Ly $\alpha$ forest can contribute to up to three different redshift bins in the $P_{1D,\alpha}$ measurement. The subforest splitting constrains the minimal accessible observed wavenumber to $k_{min} = 2\pi/(L_{sub}(1+z_{min})) = 0.0453$ Å$^{-1}$, by taking the minimal redshift used. Each subforest spans at most $\Delta z = 0.2$ and we choose the same $\Delta z$ to define the redshift binning for $P_{1D,\alpha}$.

For observed wavelength $\lambda \lesssim 3,700$ Å, the noise level is high in comparison to the spectra because of atmospheric absorptions. To minimize the impact of this noise, we remove the spectrum pixels for which the observed wavelength is lower than 3750 Å, which corresponds to Ly $\alpha$ absorbers located at $z = 2.085$. In the future, with a dedicated study to control the noise at shorter wavelength, the $P_{1D,\alpha}$ analysis can be extended to redshift $z \sim 2$.

In accordance with the eBOSS study (Chabanier et al. 2019), we remove subforests shorter than 75 spectrum pixels due to a cut in the UV region or to the presence of a large DLA. We also do not consider the Ly $\alpha$ subforests with more than 120 masked spectrum pixels.

Unlike the analysis in Chabanier et al. (2019), we do not apply a second redshift-dependent SNR cut for the averaging of $P_{1D,\alpha}$. Instead, we develop and test a SNR weighting scheme, as detailed in the appendix B. This procedure is used for all the article except in Section 4.2.1 where the impact of the SNR cut is investigated.

# 4 DESI INSTRUMENTAL CHARACTERIZATION

As the DESI instrument is new, we first focus the analysis on characterizing the instrumental effects on our $P_{1D,\alpha}$ measurement. In particular, we describe the impact of spectral resolution, instrumental noise, metal power spectrum, and atmospheric emission lines in the following.

## 4.1 Spectrograph resolution

### 4.1.1 Resolution correction modelling

As mentioned in Section 2.2, DESI spectra are linearly binned in observed wavelength, so that the natural unit for wavenumbers is Å$^{-1}$. The maximal measurable wavenumber follows the Nyquist-Shannon limit: $k_{max} = k_{Nyq} = \pi/\Delta\lambda_{pix} \approx 3.92$ Å$^{-1}$.

For SDSS $P_{1D,\alpha}$ measurements (Croft et al. 1998; McDonald et al. 2006; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), the pixelization is logarithmically binned in observed wavelength, making it suitable to express the power spectrum in Hubble velocity unit $v$, because $\Delta v \propto \Delta\lambda/\lambda = \Delta\log(\lambda)$. In eBOSS (Chabanier et al. 2019), the resolution correction is modelled by a Gaussian function $W(k, \Delta\lambda) = \exp(-0.5(k[\text{km}^{-1}]\Delta v)^2)$, where $\Delta v$ is the spectral resolution in velocity units.

The spectroscopic resolution of DESI is improved with respect to the SDSS spectrographs. On the DESI blue band (see Table 1) where most of the Ly $\alpha$ forest are observed, the effective resolving power $R = \lambda/\Delta\lambda$ ranges from 2000 to 3200 (Abareshi et al. 2022). In comparison, SDSS spectrographs had a 1500 < $R$ < 2300 in its blue band (3600 < $\lambda$ < 6350 Å). This improved resolution brings the opportunity to probe the clustering of matter at smaller scale by

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Illustration of the resolution improvement between eBOSS and DESI. The DESI spectrum of the quasar represented in Fig. 1 is zoomed on a region of the Ly $\alpha$ forest on the top panel. The spectrum of the same quasar obtained by the SDSS-IV (Gunn et al. 2006; Smee et al. 2013; Dawson et al. 2016; Blanton et al. 2017) in the 16th Data Release (DR16; Ahumada et al. 2020) of the eBOSS. The SDSS name associated with this quasar is 142 903.03–014 519.3. The DESI spectrum is obtained after three individual exposures for a total exposure time of 2300 s. The eBOSS spectrum have 11 exposures for a total of 6300 s. The large-scale absorption structures are similar but due to its improved spectroscopic resolution, the DESI spectrum clearly shows more details at small scales.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Average resolution correction ($W^2(k, \mathbf{R}, \Delta\lambda_{pix}) = \langle R^2(k) \rangle$) (blue points) using Ly $\alpha$ forest from SV3 + M2 data set. This resolution correction is weakly dependent on the redshift range (shown by shaded dashed black lines along the blue points). Only one shaded black line is above the blue points ($z = 3.8$). All the others are at the same level as the blue points. The mean value over all redshift bins is shown with the points. The shaded area represents the regime for which the associated impact of resolution and pixelization removes more than 80 per cent of the power spectrum. This criteria is chosen to define the maximal wavenumber of our $P_{1D,\alpha}$ measurement, shown with a vertical black line.</figcaption></figure>

measuring the small fluctuations in the Ly $\alpha$ forest, as illustrated in Fig. 4.

As described in Section 2.2, the DESI spectrograph resolution is entirely characterized by the resolution matrix $\mathbf{R}$ (Guy et al. 2022). In opposition to SDSS, the resolution matrix also accounts for the pixelization of the signal. Consequently, we choose to express the resolution correction function $W$ in equation (12) directly as the Fourier transform of the resolution matrix that we note $\mathbf{R}(k)$.

Fig. 5 shows the average correction due to resolution and pixelization for the SV3 + M2 data set. This correction indicates that
resolution and pixelization suppress more than 95 per cent of the signal at $k > k_{res,95} = 2.73 \AA^{-1}$, and 98 per cent at $k > k_{res,98} = 3.15 \AA^{-1}$. This observation drives the maximal wavenumber that of the $P_{1D,\alpha}$ measurement. We choose the conservative value $k_{max} = 2 \AA^{-1}$, for which the average resolution corrections is equal to 80 per cent. We will extend this conservative limit in future studies after a full characterization of resolution on CCD pixel-level simulations.

### 4.1.2 Validation with CCD image simulations

We use CCD image simulations of the DESI instrument to verify our resolution modelling described previously. This method is also used for the QMLE estimation of $P_{1D,\alpha}$ in a companion paper (Karaçaylı et al. 2023a). Those simulations are built using the desisim https://github.com/desisub/desisim package developed to model and validate the spectroscopic extraction pipeline presented in Section 2.2 and detailed in Guy et al. (2022). We use the desisim package to produce realistic realizations of two-dimensional spectroscopic images. Those images simulate various instrumental effects such as the different sources of noise detailed in Appendix C, gain and bias of the CCD amplifiers, the throughput of each spectrograph, PSF of each fibre, and sky emission.

Using the desisim package, we transform Ly $\alpha$ transmissions that follows a given input $P_{1D,\alpha}$ into realistic Ly $\alpha$ forests with apparent magnitudes representative of DESI quasar targets and noise that is representative of a single, 1000 s exposure in nominal conditions. We simulate 45 000 Ly $\alpha$ forest spectra located over ten DESI tiles and process this simulated data set with the spectroscopic pipeline. Since we only want to see the impact of spectral resolution, the true imposed noise level is used to reduce the data.

Two sets of Ly $\alpha$ contrasts $\delta_F$ are generated from this simulated quasar sample. The first set, called RAW, is produced directly from Ly $\alpha$ transmissions. This type of realization does not contain the effect of spectral resolution. In parallel, we create Ly $\alpha$ contrasts from the full CCD image simulations, noted CCD, but using the true imposed continuum for each quasar to only see the impact of resolution. Finally, we run the $P_{1D,\alpha}$ FFT pipeline presented in Section 3.2 with our resolution modelling on both Ly $\alpha$ contrast sets.

The ratio between $P_{1D,\alpha}$ obtained from RAW and CCD sets is shown in Fig. 6. As expected, the main difference between those measurements resides in the smallest scales ($k > 1.0 \AA^{-1}$). This ratio is not redshift dependent, as indicated by the light black curves in the background. We only consider the average overall redshifts. We checked that applying an additional pixelization correction sinc2 ($0.5k\Delta\lambda_{pix}$) similarly to eBOSS increases a lot the discrepancy between RAW and CCD power spectrum; thus, confirming that the resolution matrix accounts the pixelization, at least partly. We derive a correction by fitting a second-order polynomial function to the following averaged ratio:

$$A_{\text{res}}(z, k) = \left\langle \frac{P_{\text{1D},\alpha,\text{RAW}}(k, z)}{P_{\text{1D},\alpha,\text{CCD}}(k, z)} \right\rangle_z. \quad (15)$$

This term is directly multiplied to our $P_{1D,\alpha}$ measurement to account for the miss-estimation of the resolution correction.

## 4.2 Noise power spectrum measurement

### 4.2.1 Comparison of noise estimators with high-wavenumber data

The $P_{1D,\alpha}$ measurement is significantly impacted by the noise power spectrum at small scales. Thus, it is necessary to obtain an accurate estimate of this component to correct for it.

<figure><img src="image_7.png" /><figcaption>**Figure 6.** Ratios between the power spectrum obtained directly from Ly $\alpha$ transmissions (RAW) and the one derived from the CCD image simulations (CCD), for different redshift bins (light black) and averaged over all redshifts bins (blue points with error bars). Each power spectrum is re-binned by a factor of 3 to reduce error bars. A second-order polynomial is fitted to the ratio (blue continuous line). It is used to correct the miss estimation of resolution in the $P_{1D,\alpha}$ measurement.</figcaption></figure>

The noise power spectrum is estimated either directly from the pipeline noise, or by using an exposure difference method. A detailed description of the obtained noise power spectrum estimators, respectively, noted $P_{\text{pipeline}}$ and $P_{\text{diff}}$, can be found in Appendix C.

Additionally, the noise power-spectrum level is determined by taking advantage of the combined effect of the resolution and the pixelization shown in Fig. 5, as well as the Ly $\alpha$ thermal broadening. Those effects erase essentially all ‘signal’ power and thus at large wavenumbers, equation (12) simplifies into $P_{\text{raw}}(k) \simeq P_{\text{noise}}(k)$. The difference (or ratio) of $P_{\text{raw}}$ and $P_{\text{noise}}$ on the largest k-bins accessible can be used to validate the noise estimator and correct it empirically. We define the following asymptotic difference and ratio by averaging those quantities at large wavenumbers. We decide to use the criteria $k > k_{\text{res,98}}$ where $k_{\text{res,98}}$ is defined as the wavenumber for which the resolution and pixelization suppress more than 98 per cent of the signal:

$$\begin{align}
\alpha &= \langle P_{\text{raw}}(k) - P_{\text{noise}}(k) \rangle_{k>k_{\text{res,98}}}, \nonumber \\
\beta &= \left\langle \frac{P_{\text{noise}}(k)}{P_{\text{raw}}(k)} \right\rangle_{k>k_{\text{res,98}}}. \tag{16}
\end{align}$$

Fig. 7 shows the measurement of $\alpha$ and $\beta$ on SV3 data set with $\overline{\text{SNR}} > 3$ (with $\overline{\text{SNR}}$ defined by equation 6), using the pipeline noise to compute $P_{\text{noise}}$. Considering the observed statistical fluctuations, we notice that the asymptotic behaviour of power spectra at high wavenumber enables a good measurement of $\alpha$ (respectively $\beta$), whose value is close to 0 (respectively 1). Additionally, we verified that the variation of $\alpha$ and $\beta$ as a function of redshift is small in comparison to the statistical fluctuations of the ratio and difference.

### 4.2.2 Characterization on DESI data sets

We compute the $\alpha$ and $\beta$ coefficients for the pipeline ($P_{\text{pipeline}}$) and difference ($P_{\text{diff}}$) noise estimators on SV1, SV3, and M2 data sets, while varying the applied $\overline{\text{SNR}}$ cut.

Fig. 8 shows the measured $\alpha$ values on the Ly $\alpha$ forest regions for the different data sets and noise estimators. As previously stated, $\alpha$ should be equal to zero when the noise is perfectly estimated. The values of $\alpha$ are small compared to the absolute level of the noise power spectrum shown in Fig. 7. The SV1 data set exhibits a $\overline{\text{SNR}}$ cut
resolution and pixelization suppress more than 95 per cent of the signal at $k > k_{res,95} = 2.73 \AA^{-1}$, and 98 per cent at $k > k_{res,98} = 3.15 \AA^{-1}$. This observation drives the maximal wavenumber that of the $P_{1D,\alpha}$ measurement. We choose the conservative value $k_{max} = 2 \AA^{-1}$, for which the average resolution corrections is equal to 80 per cent. We will extend this conservative limit in future studies after a full characterization of resolution on CCD pixel-level simulations.

### 4.1.2 Validation with CCD image simulations

We use CCD image simulations of the DESI instrument to verify our resolution modelling described previously. This method is also used for the QMLE estimation of $P_{1D,\alpha}$ in a companion paper (Karaçaylı et al. 2023a). Those simulations are built using the desisim https://github.com/desisub/desisim package developed to model and validate the spectroscopic extraction pipeline presented in Section 2.2 and detailed in Guy et al. (2022). We use the desisim package to produce realistic realizations of two-dimensional spectroscopic images. Those images simulate various instrumental effects such as the different sources of noise detailed in Appendix C, gain and bias of the CCD amplifiers, the throughput of each spectrograph, PSF of each fibre, and sky emission.

Using the desisim package, we transform Ly $\alpha$ transmissions that follows a given input $P_{1D,\alpha}$ into realistic Ly $\alpha$ forests with apparent magnitudes representative of DESI quasar targets and noise that is representative of a single, 1000 s exposure in nominal conditions. We simulate 45 000 Ly $\alpha$ forest spectra located over ten DESI tiles and process this simulated data set with the spectroscopic pipeline. Since we only want to see the impact of spectral resolution, the true imposed noise level is used to reduce the data.

Two sets of Ly $\alpha$ contrasts $\delta_F$ are generated from this simulated quasar sample. The first set, called RAW, is produced directly from Ly $\alpha$ transmissions. This type of realization does not contain the effect of spectral resolution. In parallel, we create Ly $\alpha$ contrasts from the full CCD image simulations, noted CCD, but using the true imposed continuum for each quasar to only see the impact of resolution. Finally, we run the $P_{1D,\alpha}$ FFT pipeline presented in Section 3.2 with our resolution modelling on both Ly $\alpha$ contrast sets.

The ratio between $P_{1D,\alpha}$ obtained from RAW and CCD sets is shown in Fig. 6. As expected, the main difference between those measurements resides in the smallest scales ($k > 1.0 \AA^{-1}$). This ratio is not redshift dependent, as indicated by the light black curves in the background. We only consider the average overall redshifts. We checked that applying an additional pixelization correction sinc2 (0.5k$\Delta\lambda_{pix}$) similarly to eBOSS increases a lot the discrepancy between RAW and CCD power spectrum; thus, confirming that the resolution matrix accounts the pixelization, at least partly. We derive a correction by fitting a second-order polynomial function to the following averaged ratio:

$$A_{\text{res}}(z, k) = \left\langle \frac{P_{\text{1D},\alpha,\text{RAW}}(k, z)}{P_{\text{1D},\alpha,\text{CCD}}(k, z)} \right\rangle_z . \qquad (15)$$

This term is directly multiplied to our $P_{1D,\alpha}$ measurement to account for the miss-estimation of the resolution correction.

## 4.2 Noise power spectrum measurement

### 4.2.1 Comparison of noise estimators with high-wavenumber data

The $P_{1D,\alpha}$ measurement is significantly impacted by the noise power spectrum at small scales. Thus, it is necessary to obtain an accurate estimate of this component to correct for it.

<figure><img src="image_8.png" /><figcaption>**Figure 6.** Ratios between the power spectrum obtained directly from Ly $\alpha$ transmissions (RAW) and the one derived from the CCD image simulations (CCD), for different redshift bins (light black) and averaged over all redshifts bins (blue points with error bars). Each power spectrum is re-binned by a factor of 3 to reduce error bars. A second-order polynomial is fitted to the ratio (blue continuous line). It is used to correct the miss estimation of resolution in the $P_{1D,\alpha}$ measurement.</figcaption></figure>

The noise power spectrum is estimated either directly from the pipeline noise, or by using an exposure difference method. A detailed description of the obtained noise power spectrum estimators, respectively, noted $P_{\text{pipeline}}$ and $P_{\text{diff}}$, can be found in Appendix C.

Additionally, the noise power-spectrum level is determined by taking advantage of the combined effect of the resolution and the pixelization shown in Fig. 5, as well as the Ly $\alpha$ thermal broadening. Those effects erase essentially all ‘signal’ power and thus at large wavenumbers, equation (12) simplifies into $P_{\text{raw}}(k) \simeq P_{\text{noise}}(k)$. The difference (or ratio) of $P_{\text{raw}}$ and $P_{\text{noise}}$ on the largest k-bins accessible can be used to validate the noise estimator and correct it empirically. We define the following asymptotic difference and ratio by averaging those quantities at large wavenumbers. We decide to use the criteria $k > k_{\text{res,98}}$ where $k_{\text{res,98}}$ is defined as the wavenumber for which the resolution and pixelization suppress more than 98 per cent of the signal:

$$\begin{align}
\alpha &= \langle P_{\text{raw}}(k) - P_{\text{noise}}(k) \rangle_{k>k_{\text{res,98}}} , \nonumber \\
\beta &= \left\langle \frac{P_{\text{noise}}(k)}{P_{\text{raw}}(k)} \right\rangle_{k>k_{\text{res,98}}}. \tag{16}
\end{align}$$

Fig. 7 shows the measurement of $\alpha$ and $\beta$ on SV3 data set with $\overline{\text{SNR}} > 3$ (with $\overline{\text{SNR}}$ defined by equation 6), using the pipeline noise to compute $P_{\text{noise}}$. Considering the observed statistical fluctuations, we notice that the asymptotic behaviour of power spectra at high wavenumber enables a good measurement of $\alpha$ (respectively $\beta$), whose value is close to 0 (respectively 1). Additionally, we verified that the variation of $\alpha$ and $\beta$ as a function of redshift is small in comparison to the statistical fluctuations of the ratio and difference.

### 4.2.2 Characterization on DESI data sets

We compute the $\alpha$ and $\beta$ coefficients for the pipeline ($P_{\text{pipeline}}$) and difference ($P_{\text{diff}}$) noise estimators on SV1, SV3, and M2 data sets, while varying the applied $\overline{\text{SNR}}$ cut.

Fig. 8 shows the measured $\alpha$ values on the Ly $\alpha$ forest regions for the different data sets and noise estimators. As previously stated, $\alpha$ should be equal to zero when the noise is perfectly estimated. The values of $\alpha$ are small compared to the absolute level of the noise power spectrum shown in Fig. 7. The SV1 data set exhibits a $\overline{\text{SNR}}$ cut
<figure><img src="image_8.png" /><figcaption>**Figure 7.** Measurement of the asymptotic difference ($\alpha$) and ratio ($\beta$) between noise ($P_{\text{noise}}$) and raw ($P_{\text{raw}}$) power spectra. Those terms are defined by equation (16) and corresponds to the average of the difference and ratio of power spectra for large wavenumbers ($k > k_{\text{res,98}}$). Here, $\alpha$ and $\beta$ are measured for the pipeline noise ($P_{\text{pipeline}}$), using the SV3 observations with a minimal $\overline{\text{SNR}}$ cut of 3.</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>**Figure 8.** Asymptotic differences $\alpha$ between the noise and raw power spectra for SV1 (blue), SV3 (yellow), and M2 (green) data sets, as a function of the minimal $\overline{\text{SNR}}$ cut. This difference is measured for both $P_{\text{noise}}$ estimators from the pipeline ($P_{\text{pipeline,points}}$) and from exposure differences ($P_{\text{diff,stars}}$). The continuous lines are fits of the $\alpha$ values for pipeline noise, whose parameters are given in Table 3.</figcaption></figure>

dependence which is not present for M2 and SV3. We choose to make a data set-dependent correction to remove this residual noise. The miss-estimated noise is higher for data sets with a larger number of exposures (such as SV1). We think this originates from unaccounted common sources of noise coming from the statistical uncertainties in the CCD calibration data (dark current, pixel flat field), which explain why this effect increases with the number of exposures. The

<table><caption>**Table 3.** Additive corrections applied to the pipeline noise for different spectral regions and data sets. An $\overline{\text{SNR}}$ dependence is included in the case of SV1 only. The same parameters are used for both SB1 and SB2.</caption><thead><tr><th>Band</th><th>Data</th><th>P<sub>noise,miss</sub> = α (Å)</th></tr></thead><tbody><tr><td rowspan="3">Ly α</td><td>SV1</td><td>0.026 × (SNR)<sup>-1.77</sup> + 0.00076</td></tr><tr><td>SV3</td><td>0.00127</td></tr><tr><td>M2</td><td>0.00109</td></tr><tr><td rowspan="3">Side-bands</td><td>SV1</td><td>0.018 × (SNR)<sup>-1.52</sup> + 0.000032</td></tr><tr><td>SV3</td><td>0.00048</td></tr><tr><td>M2</td><td>0.00019</td></tr></tbody></table>

dependence in $\overline{\text{SNR}}$ can be explained by the fact that this effect is amplified when we consider noisier spectra.

In Fig. 8, the exposure difference noise estimator $P_{\text{diff}}$ is shown for SV1 and SV3. In M2, quasars are observed with one, or a small number of exposures; thus, $P_{\text{diff}}$ cannot be reliably computed. For the SV1 and SV3 observations, the difference noise power spectrum exhibits the same trend as the pipeline noise estimate. However, $P_{\text{diff}}$ consistently underestimates the noise level compared to $P_{\text{pipeline}}$. We think that this underestimation is due to the fact that $P_{\text{diff}}$ is not accounting for all the common sources of noise between exposures. As a consequence, we only consider the $P_{\text{pipeline}}$ estimation from now on.

In Appendix C3, we perform additional studies to characterize the additivity of the miss-estimated noise and its behaviour for different spectra regions. Those regions, called side-bands, are used in the next section for metal power spectrum estimation. Table 3 summarizes the corrections to the pipeline noise we computed for the different data sets. For SV3 and M2, we choose to apply a constant additive correction ($\alpha$). For SV1, given the observed dependence of $\alpha$ as a function of the minimal $\overline{\text{SNR}}$ cut, we fit it with a power law.

## 4.3 Side-band power spectrum

Following previous $P_{\text{1D,α}}$ studies (McDonald et al. 2006; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), special spectrum regions, called side-bands and that are devoid of Ly $\alpha$ absorption, are used to statistically estimate the power spectrum components $P_{\text{metals}}$ caused by metal absorptions in the Ly $\alpha$ forest. The resulting signal from the side-bands contains information about the abundance, temperature and clustering of metals in the intergalactic medium. In our study, we aim at creating a model to closely reproduce the side-band power spectrum ($P_{\text{SB}}$) so that we can statistically subtract it in the measurement of $P_{\text{1D,α}}$ in equation (14).

We define the side-bands SB1 (1270 Å < $\lambda_{\text{rf}}$ < 1380 Å) and SB2 (1410 Å < $\lambda_{\text{rf}}$ < 1520 Å). In both side-bands, the fraction of transmitted flux contrast can be expressed similarly to equation (10):

$$ \delta_F(\lambda)|_{\text{SB}} = \delta_{\text{metals}}(\lambda)|_{\text{SB}} \otimes W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda)|_{\text{SB}}. \quad (17) $$

The $\delta_{\text{metals}}(\lambda)|_{\text{SB}}$ contrast contains all the fluctuations caused by metals with rest-frame absorption wavelength higher than 1380 Å for SB1 and 1520 Å for SB2. Similarly to the calculation of $P_{\text{1D,α}}$, the side-band power spectrum writes:

$$ P_{\text{SB}}(k) = \left\langle \frac{P_{\text{raw}}(k)|_{\text{SB}} - P_{\text{noise}}(k)|_{\text{SB}}}{W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}})} \right\rangle. \quad (18) $$

The main difference between both side-bands is that SB1 contains Si IV absorption, which is not present in SB2. Consequently, we use the side-band power spectrum of SB1 to estimate $P_{\text{metals}}$ in
<figure><img src="image_9.png" /><figcaption>**Figure 7.** Measurement of the asymptotic difference ($\alpha$) and ratio ($\beta$) between noise ($P_{\text{noise}}$) and raw ($P_{\text{raw}}$) power spectra. Those terms are defined by equation (16) and corresponds to the average of the difference and ratio of power spectra for large wavenumbers ($k > k_{\text{res,98}}$). Here, $\alpha$ and $\beta$ are measured for the pipeline noise ($P_{\text{pipeline}}$), using the SV3 observations with a minimal $\overline{\text{SNR}}$ cut of 3.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Asymptotic differences $\alpha$ between the noise and raw power spectra for SV1 (blue), SV3 (yellow), and M2 (green) data sets, as a function of the minimal $\overline{\text{SNR}}$ cut. This difference is measured for both $P_{\text{noise}}$ estimators from the pipeline ($P_{\text{pipeline,points}}$) and from exposure differences ($P_{\text{diff,stars}}$). The continuous lines are fits of the $\alpha$ values for pipeline noise, whose parameters are given in Table 3.</figcaption></figure>

dependence which is not present for M2 and SV3. We choose to make a data set-dependent correction to remove this residual noise. The miss-estimated noise is higher for data sets with a larger number of exposures (such as SV1). We think this originates from unaccounted common sources of noise coming from the statistical uncertainties in the CCD calibration data (dark current, pixel flat field), which explain why this effect increases with the number of exposures. The

<table><caption>**Table 3.** Additive corrections applied to the pipeline noise for different spectral regions and data sets. An $\overline{\text{SNR}}$ dependence is included in the case of SV1 only. The same parameters are used for both SB1 and SB2.</caption><thead><tr><th>Band</th><th>Data</th><th><i>P</i><sub>noise,miss</sub> = α (Å)</th></tr></thead><tbody><tr><td rowspan="3">Ly α</td><td>SV1</td><td>0.026 × (SNR)<sup>-1.77</sup> + 0.00076</td></tr><tr><td>SV3</td><td>0.00127</td></tr><tr><td>M2</td><td>0.00109</td></tr><tr><td rowspan="3">Side-bands</td><td>SV1</td><td>0.018 × (SNR)<sup>-1.52</sup> + 0.000032</td></tr><tr><td>SV3</td><td>0.00048</td></tr><tr><td>M2</td><td>0.00019</td></tr></tbody></table>

dependence in $\overline{\text{SNR}}$ can be explained by the fact that this effect is amplified when we consider noisier spectra.

In Fig. 8, the exposure difference noise estimator $P_{\text{diff}}$ is shown for SV1 and SV3. In M2, quasars are observed with one, or a small number of exposures; thus, $P_{\text{diff}}$ cannot be reliably computed. For the SV1 and SV3 observations, the difference noise power spectrum exhibits the same trend as the pipeline noise estimate. However, $P_{\text{diff}}$ consistently underestimates the noise level compared to $P_{\text{pipeline}}$. We think that this underestimation is due to the fact that $P_{\text{diff}}$ is not accounting for all the common sources of noise between exposures. As a consequence, we only consider the $P_{\text{pipeline}}$ estimation from now on.

In Appendix C3, we perform additional studies to characterize the additivity of the miss-estimated noise and its behaviour for different spectra regions. Those regions, called side-bands, are used in the next section for metal power spectrum estimation. Table 3 summarizes the corrections to the pipeline noise we computed for the different data sets. For SV3 and M2, we choose to apply a constant additive correction ($\alpha$). For SV1, given the observed dependence of $\alpha$ as a function of the minimal $\overline{\text{SNR}}$ cut, we fit it with a power law.

## 4.3 Side-band power spectrum

Following previous $P_{\text{1D,α}}$ studies (McDonald et al. 2006; Palanque-Delabrouille et al. 2013; Chabanier et al. 2019), special spectrum regions, called side-bands and that are devoid of Ly $\alpha$ absorption, are used to statistically estimate the power spectrum components $P_{\text{metals}}$ caused by metal absorptions in the Ly $\alpha$ forest. The resulting signal from the side-bands contains information about the abundance, temperature and clustering of metals in the intergalactic medium. In our study, we aim at creating a model to closely reproduce the side-band power spectrum ($P_{\text{SB}}$) so that we can statistically subtract it in the measurement of $P_{\text{1D,α}}$ in equation (14).

We define the side-bands SB1 (1270 Å < $\lambda_{\text{rf}}$ < 1380 Å) and SB2 (1410 Å < $\lambda_{\text{rf}}$ < 1520 Å). In both side-bands, the fraction of transmitted flux contrast can be expressed similarly to equation (10):

$$ \delta_F(\lambda)|_{\text{SB}} = \delta_{\text{metals}}(\lambda)|_{\text{SB}} \otimes W(\lambda, \mathbf{R}, \Delta\lambda_{\text{pix}}) + \delta_{\text{noise}}(\lambda)|_{\text{SB}}. \quad (17) $$

The $\delta_{\text{metals}}(\lambda)|_{\text{SB}}$ contrast contains all the fluctuations caused by metals with rest-frame absorption wavelength higher than 1380 Å for SB1 and 1520 Å for SB2. Similarly to the calculation of $P_{\text{1D,α}}$, the side-band power spectrum writes:

$$ P_{\text{SB}}(k) = \left\langle \frac{P_{\text{raw}}(k)|_{\text{SB}} - P_{\text{noise}}(k)|_{\text{SB}}}{W^2(k, \mathbf{R}, \Delta\lambda_{\text{pix}})} \right\rangle. \quad (18) $$

The main difference between both side-bands is that SB1 contains Si IV absorption, which is not present in SB2. Consequently, we use the side-band power spectrum of SB1 to estimate $P_{\text{metals}}$ in
<figure><img src="image_10.png" /><figcaption>**Figure 9.** 1D power spectra measured in the side-band regions SB1 and SB2, using the SV3 + M2 data set after applying the data set-dependent noise correction. Top: Average of $P_{SB}$ over all redshift bins in the rest wavenumber frame $k_{\text{rest},i} = k_{\text{obs},i} \times (1+z)$. The fitted model represented in continuous line is given by equation (19). Bottom: $P_{SB}$ on side-band SB1 as a function of redshift and observed wavenumber. Each redshift bin is fitted using the product between equation (19), with parameters fixed from previous fit, and a first-order polynomial function.</figcaption></figure>

equation (14), and the SB2 power spectrum as a consistency check. The side-band power spectrum is computed at the same observed wavelength range as $P_{1D,\alpha}$. Because of the higher rest-frame absorption wavelength, the quasars employed to calculate the side-band power spectrum are at a lower redshift than the sample used for Ly $\alpha$. In particular, quasars at $z < 2.0$ are employed to calculate $P_{SB}$ in the lowest redshift bins of $P_{1D,\alpha}$. However, the metals in the intergalactic medium that produces the absorptions responsible for the side-band power spectrum ($P_{SB}$) are at the same redshift as those which produced the metal power spectrum ($P_{\text{metals}}$) in the $P_{1D,\alpha}$ calculation for Ly $\alpha$. Consequently, the redshift dependence of metal absorptions is correctly taken into account.

We note that the method we use to remove metal contribution is not perfect. In particular, the blending of metals with Ly $\alpha$ emission pointed out in Day et al. (2019) is not fully accounted here. This second-order effect should be included in future studies for which the precision level will significantly improve.

The measurement of the side-band power spectrum using the SV3 + M2 data set is shown in Fig. 9. The top panel shows the stack of $P_{SB}$ for all redshift bins as a function of wavenumber in rest-frame $k_{\text{rest}} = (1+z) \times k_{\text{obs}}$. The bottom panel shows separated redshift bins as a function of observed wavenumber. A total number of 201 849 and 276 279 sub-forests were used for SB1 and SB2 respective measurements.

In top panel of Fig. 9, the average side-band power is lower for the SB2 than SB1, as expected by the addition of Si IV absorptions.

In the eBOSS measurement (Chabanier et al. 2019), a sixth-degree polynomial is used to fit the shape of the side-band power. For our

measurement, we exploit the stacked $P_{SB}$ profile to design a more physically motivated model.

A complete list of metals present in the Ly $\alpha$ forest and that impact $P_{1D,\alpha}$ is given in Pieri et al. (2014), Yang et al. (2022), and the strongest absorptions are from Si III, Si II, Si IV, and C IV.

The emission peaks, and consequently absorption peaks, of Si IV and C IV are actually two doublets. Their rest-frame wavelength given by NIST (Kramida et al. 2021) are $\lambda_{\text{SiIV}}^a = 1, 393.76$ Å, $\lambda_{\text{SiIV}}^b = 1, 402.77$ Å, $\lambda_{\text{CIV}}^a = 1, 548.202$ Å, and $\lambda_{\text{CIV}}^b = 1, 550.774$ Å. The presence of an absorption doublet in the side-band creates a peak in the 1D correlation function, which translates into an oscillatory pattern in the power spectrum, whose periodicity depends on the doublet separation. This effect is studied more in detail on the same data set in Karaçaylı et al. (2023b) to determine cosmic ion abundance. In the top panel of Fig. 9, both side-band power spectra display a large oscillation caused by the C IV doublet. As expected, the SB1 power spectrum shows an additional oscillation induced by Si IV doublet absorptions. These considerations lead us to model $P_{SB}$ as the sum of a power law including all-metal contributions and oscillations due to Si IV and C IV doublets:

$$P_{\text{SB,m}} = A \times k^{-\epsilon} + \sum_{i} P_{\text{doublet},i}(k, A_i, a_i, k_i, \psi_i). \quad (19)$$

Oscillations induced by a doublet have a frequency characterized by the rest-frame wavenumber $k_{\text{rest},i} = 2\pi/\delta\lambda_i$, where $\delta\lambda_i$ is the doublet separation in Å. We choose to use damped sinusoidal functions to model the doublet oscillations as follows:

$$P_{\text{doublet},i}(k, A_i, a_i, k_i, \psi_i) = A_i e^{-a_i k} \sin \left( 2\pi \left( \frac{k}{k_i} \right) + \psi_i \right), \quad (20)$$

where $k_i$ is a free parameter with a uniform prior centred around $k_{\text{rest},i}$.

The result of the fit on the redshift-averaged $P_{SB}$, taking into account the oscillations of C IV and Si IV for SB1, and only C IV for SB2, is shown in the top panel of Fig. 9.

The SB1 fitted function is used to derive the redshift dependence of the side-band power spectrum, shown in the bottom panel of Fig. 9. For each redshift bin, we fit a product between the global SB1 fitted function (expressed in the observational wavelength frame), and a first-order polynomial. As $P_{SB}$ may also include other uncorrelated contaminations besides metals, we do not seek to interpret the fitted values for each power spectrum. In particular, we note that the $k_i$ are systematically shifted in comparison to their doublet oscillation frequency $k_{\text{rest},i}$. For SB1, the fitted values are $k_{\text{CIV}} = 3.32$ Å$^{-1}$ and $k_{\text{SiIV}} = 0.812$ Å$^{-1}$, whereas rest-frame wavenumber are $k_{\text{rest,CIV}} = 2.44$ Å$^{-1}$ and $k_{\text{rest,SiIV}} = 0.697$ Å$^{-1}$. For SB2, we obtain $k_{\text{CIV}} = 3.23$ Å$^{-1}$. We think that this effect might be due to the blended impact of all metals present in the intergalactic medium. Thus, it is necessary to vary $k_i$ parameters to closely fit our data.

This $P_{SB}$ measurement already represents a clear improvement with respect to that of BOSS (Palanque-Delabrouille et al. 2013), and eBOSS (Chabanier et al. 2019): by eye, Si IV- and C IV-induced oscillatory patterns are seen even for individual redshift bins, and $P_{SB}$ is essentially a decreasing function of wavenumber, even at high wavenumber. This indicates an improvement in the noise modelling.

## 4.4 Atmospheric and galactic emission lines

Atmospheric emission lines are corrected from DESI spectra by the spectral extraction pipeline as described in Section 2.2. The average of 15 000 sky spectra on exposures with optimal observing conditions, noted $\langle f_{\text{sky}} \rangle$, is shown in Fig. 10.
<figure><img src="image_11.png" /><figcaption>**Figure 9.** 1D power spectra measured in the side-band regions SB1 and SB2, using the SV3 + M2 data set after applying the data set-dependent noise correction. Top: Average of $P_{SB}$ over all redshift bins in the rest wavenumber frame $k_{\text{rest},i} = k_{\text{obs},i} \times (1+z)$. The fitted model represented in continuous line is given by equation (19). Bottom: $P_{SB}$ on side-band SB1 as a function of redshift and observed wavenumber. Each redshift bin is fitted using the product between equation (19), with parameters fixed from previous fit, and a first-order polynomial function.</figcaption></figure>

equation (14), and the SB2 power spectrum as a consistency check. The side-band power spectrum is computed at the same observed wavelength range as $P_{1D,\alpha}$. Because of the higher rest-frame absorption wavelength, the quasars employed to calculate the side-band power spectrum are at a lower redshift than the sample used for Ly $\alpha$. In particular, quasars at $z < 2.0$ are employed to calculate $P_{SB}$ in the lowest redshift bins of $P_{1D,\alpha}$. However, the metals in the intergalactic medium that produces the absorptions responsible for the side-band power spectrum ($P_{SB}$) are at the same redshift as those which produced the metal power spectrum ($P_{\text{metals}}$) in the $P_{1D,\alpha}$ calculation for Ly $\alpha$. Consequently, the redshift dependence of metal absorptions is correctly taken into account.

We note that the method we use to remove metal contribution is not perfect. In particular, the blending of metals with Ly $\alpha$ emission pointed out in Day et al. (2019) is not fully accounted here. This second-order effect should be included in future studies for which the precision level will significantly improve.

The measurement of the side-band power spectrum using the SV3 + M2 data set is shown in Fig. 9. The top panel shows the stack of $P_{SB}$ for all redshift bins as a function of wavenumber in rest-frame $k_{\text{rest}} = (1+z) \times k_{\text{obs}}$. The bottom panel shows separated redshift bins as a function of observed wavenumber. A total number of 201 849 and 276 279 sub-forests were used for SB1 and SB2 respective measurements.

In top panel of Fig. 9, the average side-band power is lower for the SB2 than SB1, as expected by the addition of Si IV absorptions.

In the eBOSS measurement (Chabanier et al. 2019), a sixth-degree polynomial is used to fit the shape of the side-band power. For our

measurement, we exploit the stacked $P_{SB}$ profile to design a more physically motivated model.

A complete list of metals present in the Ly $\alpha$ forest and that impact $P_{1D,\alpha}$ is given in Pieri et al. (2014), Yang et al. (2022), and the strongest absorptions are from Si III, Si II, Si IV, and C IV.

The emission peaks, and consequently absorption peaks, of Si IV and C IV are actually two doublets. Their rest-frame wavelength given by NIST (Kramida et al. 2021) are $\lambda_{\text{SiIV}}^a = 1, 393.76$ Å, $\lambda_{\text{SiIV}}^b = 1, 402.77$ Å, $\lambda_{\text{CIV}}^a = 1, 548.202$ Å, and $\lambda_{\text{CIV}}^b = 1, 550.774$ Å. The presence of an absorption doublet in the side-band creates a peak in the 1D correlation function, which translates into an oscillatory pattern in the power spectrum, whose periodicity depends on the doublet separation. This effect is studied more in detail on the same data set in Karaçaylı et al. (2023b) to determine cosmic ion abundance. In the top panel of Fig. 9, both side-band power spectra display a large oscillation caused by the C IV doublet. As expected, the SB1 power spectrum shows an additional oscillation induced by Si IV doublet absorptions. These considerations lead us to model $P_{SB}$ as the sum of a power law including all-metal contributions and oscillations due to Si IV and C IV doublets:

$$P_{\text{SB,m}} = A \times k^{-\epsilon} + \sum_{i} P_{\text{doublet},i}(k, A_i, a_i, k_i, \psi_i). \quad (19)$$

Oscillations induced by a doublet have a frequency characterized by the rest-frame wavenumber $k_{\text{rest},i} = 2\pi/\delta\lambda_i$, where $\delta\lambda_i$ is the doublet separation in Å. We choose to use damped sinusoidal functions to model the doublet oscillations as follows:

$$P_{\text{doublet},i}(k, A_i, a_i, k_i, \psi_i) = A_i e^{-a_i k} \sin \left( 2\pi \left( \frac{k}{k_i} \right) + \psi_i \right), \quad (20)$$

where $k_i$ is a free parameter with a uniform prior centred around $k_{\text{rest},i}$.

The result of the fit on the redshift-averaged $P_{SB}$, taking into account the oscillations of C IV and Si IV for SB1, and only C IV for SB2, is shown in the top panel of Fig. 9.

The SB1 fitted function is used to derive the redshift dependence of the side-band power spectrum, shown in the bottom panel of Fig. 9. For each redshift bin, we fit a product between the global SB1 fitted function (expressed in the observational wavelength frame), and a first-order polynomial. As $P_{SB}$ may also include other uncorrelated contaminations besides metals, we do not seek to interpret the fitted values for each power spectrum. In particular, we note that the $k_i$ are systematically shifted in comparison to their doublet oscillation frequency $k_{\text{rest},i}$. For SB1, the fitted values are $k_{\text{CIV}} = 3.32$ Å$^{-1}$ and $k_{\text{SiIV}} = 0.812$ Å$^{-1}$, whereas rest-frame wavenumber are $k_{\text{rest,CIV}} = 2.44$ Å$^{-1}$ and $k_{\text{rest,SiIV}} = 0.697$ Å$^{-1}$. For SB2, we obtain $k_{\text{CIV}} = 3.23$ Å$^{-1}$. We think that this effect might be due to the blended impact of all metals present in the intergalactic medium. Thus, it is necessary to vary $k_i$ parameters to closely fit our data.

This $P_{SB}$ measurement already represents a clear improvement with respect to that of BOSS (Palanque-Delabrouille et al. 2013), and eBOSS (Chabanier et al. 2019): by eye, Si IV- and C IV-induced oscillatory patterns are seen even for individual redshift bins, and $P_{SB}$ is essentially a decreasing function of wavenumber, even at high wavenumber. This indicates an improvement in the noise modelling.

## 4.4 Atmospheric and galactic emission lines

Atmospheric emission lines are corrected from DESI spectra by the spectral extraction pipeline as described in Section 2.2. The average of 15 000 sky spectra on exposures with optimal observing conditions, noted $\langle f_{\text{sky}} \rangle$, is shown in Fig. 10.
<figure><img src="image_11.png" /></figure>

**Figure 10.** Average of the sky spectra of 15 000 sky fibres with optimal observation conditions (speed >105, effective exposure time >1, 100 s, seeing <1.05 deg, and airmass <1.3). These sky spectra originate from three exposures in the SV1 and SV3 data sets. The different spectral bands of DESI are represented (B in blue, R in orange, and Z in green). The median smoothing of this average sky spectrum multiplied by a threshold $\Lambda_1$ = 2.5, shown as a dashed red line, is used to select atmospheric emission lines we want to mask (light purple lines). The line located at 4360 Å was added manually, considering its large impact on noise seen in Fig. 11.

The noise of spectrum pixels associated with intense atmospheric lines is strongly increased. It induces additional oscillations in the Ly $\alpha$ contrasts and increases the level of $P_{\text{ID,}\alpha}$. We need to correct this effect as those atmospheric lines are not linked to IGM physics. We choose to mask the major atmospheric lines as in previous measurements (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019). The catalogue of lines in these studies https://github.com/igmhub/picca/blob/master/etc/list_veto_line_Pk1D.txt was adapted to the spectral resolution of the SDSS instrument. The improved resolution of DESI makes it possible to reduce the masking size for narrow atmospheric lines, decreasing the impact of masking on the $P_{\text{ID,}\alpha}$ measurement.

We develop an algorithm similar to Lee et al. (2013) to compute an atmospheric line catalog adapted to the DESI instrument. A median smoothing $M(d_\lambda)$ of spectral width $d_\lambda$ = 160 Å is applied on the average sky spectrum. Atmospheric lines are selected when the average sky spectrum is larger than the product of the smoothed sky flux, $M(d_\lambda) \otimes \langle f_{\text{sky}} \rangle$, by a threshold $\Lambda_1$. In Fig. 10, the dashed red line represents this product for $\Lambda_1$ = 2.5 and $d_\lambda$ = 160 Å.

A second threshold $\Lambda_w$ = 1.2 defines the width of atmospheric lines. The upper and lower wavelength limits of an atmospheric line are defined as the first wavelengths on each side whose average sky spectrum is lower than $\Lambda_w \times M(d_\lambda) \otimes \langle f_{\text{sky}} \rangle$. To remain conservative and prevent numerical effects potentially caused by masking at a spectrum pixel position, the line widths are increased by 1 Å on each side.

In this atmospheric line catalogue, we also add the galactic absorption lines, which correspond to relatively broad absorptions made by dust in the Milky Way. We take the same lines as eBOSS: Ca II H and K lines at 3968 and 3933 Å, and the Na D doublet at 5893 Å. The DESI atmospheric emission line catalogue built from this procedure is available online.

We verify that the produced atmospheric line catalogue correctly masks the DESI noise. For this purpose, we compute an average noise by stacking the pipeline noise $\sigma_{\text{pip}}$ of 125 477 objects categorized by redrock as quasars or luminous red galaxies on SV3 observations. This average noise with the DESI atmospheric line catalogue such that ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) is shown in the bottom panel of Fig. 11. A zoom on an atmospheric line is shown in the top panels of Fig. 11.

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

**Figure 11.** Comparison between the atmospheric emission line mask catalogues used for eBOSS and the new DESI catalogue we designed in our study. Bottom: Average of the DESI pipeline noise for 125 477 objects categorized by redrock as quasars or luminous red galaxies on SV3 observations (blue curve). The DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) line catalogue is shown on top with purple vertical lines. Top: Zoom on a specific atmospheric emission line on the DESI stacked noise (blue curve). The mask used on this specific line for DESI is shown in the left panel in purple and for eBOSS in orange on the right panel. The DESI mask decreases the masked length in accordance with stacked DESI noise.

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

**Figure 12.** Percentage of spectral length masked by atmospheric lines as a function of redshift, for eBOSS and DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) catalogues.

with eBOSS and DESI masks. The eBOSS mask is too wide for DESI stacked noise, which highlights the benefit of creating a new catalogue. After a visual inspection of most atmospheric lines, we validate that all the spectrum pixels showing an increase in noise are masked by setting the width threshold to $\Lambda_w$ = 1.2.

Comparing the average noise (Fig. 11) to the average sky flux (Fig. 10), there is a consistency between atmospheric emission lines and observed peaks in the pipeline noise. The feature at 4360 Å is an exception, as it appears wide and relatively high in the DESI noise and not in the average sky flux. Its wavelength is inside a known transmission dip around 4400 Å due to an issue with DESI's spectrograph collimator coating (Guy et al. 2022). For this specific line, we take the same value as the eBOSS catalogue and force the algorithm to consider it as an atmospheric line even if it does not pass the $\Lambda_1$ requirement.

A comparison of the percentage of Ly $\alpha$ forest masked for eBOSS and DESI catalogue, as a function of redshift, is given in Fig. 12. To remain conservative, we chose the value of $\Lambda_1$ = 2.5 to obtain a catalogue of atmospheric emission lines with a spectral length
<figure><img src="image_12.png" /></figure>

**Figure 10.** Average of the sky spectra of 15 000 sky fibres with optimal observation conditions (speed >105, effective exposure time >1, 100 s, seeing <1.05 deg, and airmass <1.3). These sky spectra originate from three exposures in the SV1 and SV3 data sets. The different spectral bands of DESI are represented (B in blue, R in orange, and Z in green). The median smoothing of this average sky spectrum multiplied by a threshold $\Lambda_1$ = 2.5, shown as a dashed red line, is used to select atmospheric emission lines we want to mask (light purple lines). The line located at 4360 Å was added manually, considering its large impact on noise seen in Fig. 11.

The noise of spectrum pixels associated with intense atmospheric lines is strongly increased. It induces additional oscillations in the Ly $\alpha$ contrasts and increases the level of $P_{\text{ID,}\alpha}$. We need to correct this effect as those atmospheric lines are not linked to IGM physics. We choose to mask the major atmospheric lines as in previous measurements (Palanque-Delabrouille et al. 2013; Chabanier et al. 2019). The catalogue of lines in these studies https://github.com/igmhub/picca/blob/master/etc/list_veto_line_Pk1D.txt was adapted to the spectral resolution of the SDSS instrument. The improved resolution of DESI makes it possible to reduce the masking size for narrow atmospheric lines, decreasing the impact of masking on the $P_{\text{ID,}\alpha}$ measurement.

We develop an algorithm similar to Lee et al. (2013) to compute an atmospheric line catalog adapted to the DESI instrument. A median smoothing $M(d_\lambda)$ of spectral width $d_\lambda$ = 160 Å is applied on the average sky spectrum. Atmospheric lines are selected when the average sky spectrum is larger than the product of the smoothed sky flux, $M(d_\lambda) \otimes \langle f_{\text{sky}} \rangle$, by a threshold $\Lambda_1$. In Fig. 10, the dashed red line represents this product for $\Lambda_1$ = 2.5 and $d_\lambda$ = 160 Å.

A second threshold $\Lambda_w$ = 1.2 defines the width of atmospheric lines. The upper and lower wavelength limits of an atmospheric line are defined as the first wavelengths on each side whose average sky spectrum is lower than $\Lambda_w \times M(d_\lambda) \otimes \langle f_{\text{sky}} \rangle$. To remain conservative and prevent numerical effects potentially caused by masking at a spectrum pixel position, the line widths are increased by 1 Å on each side.

In this atmospheric line catalogue, we also add the galactic absorption lines, which correspond to relatively broad absorptions made by dust in the Milky Way. We take the same lines as eBOSS: Ca II H and K lines at 3968 and 3933 Å, and the Na D doublet at 5893 Å. The DESI atmospheric emission line catalogue built from this procedure is available online.

We verify that the produced atmospheric line catalogue correctly masks the DESI noise. For this purpose, we compute an average noise by stacking the pipeline noise $\sigma_{\text{pip}}$ of 125 477 objects categorized by redrock as quasars or luminous red galaxies on SV3 observations. This average noise with the DESI atmospheric line catalogue such that ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) is shown in the bottom panel of Fig. 11. A zoom on an atmospheric line is shown in the top panels of Fig. 11.

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

**Figure 11.** Comparison between the atmospheric emission line mask catalogues used for eBOSS and the new DESI catalogue we designed in our study. Bottom: Average of the DESI pipeline noise for 125 477 objects categorized by redrock as quasars or luminous red galaxies on SV3 observations (blue curve). The DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) line catalogue is shown on top with purple vertical lines. Top: Zoom on a specific atmospheric emission line on the DESI stacked noise (blue curve). The mask used on this specific line for DESI is shown in the left panel in purple and for eBOSS in orange on the right panel. The DESI mask decreases the masked length in accordance with stacked DESI noise.

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

**Figure 12.** Percentage of spectral length masked by atmospheric lines as a function of redshift, for eBOSS and DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) catalogues.

with eBOSS and DESI masks. The eBOSS mask is too wide for DESI stacked noise, which highlights the benefit of creating a new catalogue. After a visual inspection of most atmospheric lines, we validate that all the spectrum pixels showing an increase in noise are masked by setting the width threshold to $\Lambda_w$ = 1.2.

Comparing the average noise (Fig. 11) to the average sky flux (Fig. 10), there is a consistency between atmospheric emission lines and observed peaks in the pipeline noise. The feature at 4360 Å is an exception, as it appears wide and relatively high in the DESI noise and not in the average sky flux. Its wavelength is inside a known transmission dip around 4400 Å due to an issue with DESI's spectrograph collimator coating (Guy et al. 2022). For this specific line, we take the same value as the eBOSS catalogue and force the algorithm to consider it as an atmospheric line even if it does not pass the $\Lambda_1$ requirement.

A comparison of the percentage of Ly $\alpha$ forest masked for eBOSS and DESI catalogue, as a function of redshift, is given in Fig. 12. To remain conservative, we chose the value of $\Lambda_1$ = 2.5 to obtain a catalogue of atmospheric emission lines with a spectral length
masked similar to the eBOSS catalogue. A complete study on synthetic data, out of this paper's scope, will be done to decrease the length masked, and consequently the impact on $P_{1D,α}$.

# 5 SYNTHETIC DATA CORRECTIONS

Synthetic data (otherwise called mocks) are generated to characterize the impact on $P_{1D,α}$ of continuum fitting, spectral resolution, noise modelling, and spectrum pixel masking. From this, we derive empirical corrections of these effects and apply them to the data measurement.

## 5.1 Synthetic data sets

We generated a set of DESI-Lite (Karaçaylı, Font-Ribera & Padmanabhan 2020) mocks, specifically designed for $P_{1D,α}$. The full description of these mocks is given in Karaçaylı et al. (2023a). The DESI-LITE software produces uncorrelated Ly α forests that mimic the redshift and noise distribution of the SV3 + M2 data set. Ten independent realizations are generated with different initial conditions. For each realization, a random catalogue of DLAs is created to follow the redshift and column density distribution of the latest eBOSS catalogue (Chabanier et al. 2021).

The quickquasars software https://github.com/desihub/desisim/blob/main/py/desisim/scripts/quickquasars.py (Herrera-Alcantar et al., in preparation), included in the desisim package, transforms Ly α transmission into spectra with observational and astrophysical contaminants. For each realization, two sets of spectra are generated by imprinting DLAs according to the catalog aforementioned or not.

We run five different $P_{1D,α}$ variations to study the impact of different Ly α contaminants:

(i) TRUECONT: The true quasar continuum imposed by quickquasars is applied instead of the continuum fitting procedure described in Section 3.1. In comparison to RAW mocks, this realization is impacted by finite noise and resolution.

(ii) CONT: The Ly α contrasts are calculated using the pipeline detailed in Section 3.1. This type of mocks includes the impact of continuum fitting.

(iii) DLAm: Realization for which the DLAs are not added to forests at the quickquasars stage, though we mask spectrum pixels as if they were present. The objective of this kind of mocks is to characterize the impact of DLA masking.

(iv) LINEm: Similarly to DLAm but masking the atmospheric emission lines catalogue built in Section 4.4.

(v) DLA: Realization for which the DLA are applied to the spectra without masking them. The objective of this mocks is to measure the impact of DLA to compute a DLA completeness systematic error in Section 6.

For all the mocks, we take the same procedure as for the FFT calculation on observational data. In particular, the same SNR-weighting is applied, and the number of sub-forest for each realization is around 81 500 with a small statistical variation between realizations. This is slightly larger than the number of sub-forest of the data sample given in Section 7.

In the next sections, all the results are shown for the combination of ten independent realizations. To decrease the error bars, we also performed a linear rebin that provides a wavenumber binning three times coarser than what we used for observational data. The error bars of the presented ratios are computed in quadrature.

<figure><img src="image_14.png" /><figcaption>Figure 13. Ratios between the power spectrum obtained using true continuum (TRUECONT) and the one derived with our pipeline (CONT) on the combination of 10 mocks. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Fitting functions are represented by continuous lines, and used to correct the $P_{1D,α}$ measurement. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

## 5.2 Continuum-fitting correction

The continuum-fitting procedure defined in Section 3.1 systematically distort the measured $C_q(\lambda, z_q)\bar{F}(\lambda)$ term by suppressing large-scale modes, and may bias the $P_{1D,α}$ measurement. This is a well-known effect in BAO measurements (Bourboux et al. 2020). To create a correction which contains this effect, we compare the mock computed using the true continuum (TRUECONT) with the one which follows the standard continuum fitting procedure (CONT):

$$A_{\text{cont}}(z, k) = \frac{P_{1D,\alpha,\text{TRUECONT}}(k, z)}{P_{1D,\alpha,\text{CONT}}(k, z)}. \quad (21)$$

This correction is shown for the combination of ten DESI-Lite mocks in Fig. 13. We use a second-order polynomial function to fit this correction and apply it to the $P_{1D,α}$ measurement.

This correction differs in amplitude compared to the eBOSS measurement (Chabanier et al. 2019). As for eBOSS, the 1D power spectrum with continuum fitting is higher than that measured with the true continuum. However, in our case, the impact is much smaller than eBOSS, for which this ratio was near 4 per cent (without using a first-order polynomial function in the continuum fitting). Furthermore, we do not have a large-scale impact as significant as eBOSS.

## 5.3 Spectrum pixel masking

For both DLA and atmospheric line masking, we remove some data points from the measured spectra. This does not impact studies performed on real-space spectra, such as the continuum fitting or the QMLE for $P_{1D,α}$ (Karaçaylı et al. 2022). On the other hand, the FFT calculation requires that spectrum pixels are equally spaced. Consequently, when computing the Fourier transform, we impose a value of $\delta_F = 0$ (equivalent to mean transmitted flux fraction value for
masked similar to the eBOSS catalogue. A complete study on synthetic data, out of this paper's scope, will be done to decrease the length masked, and consequently the impact on $P_{1D,α}$.

# 5 SYNTHETIC DATA CORRECTIONS

Synthetic data (otherwise called mocks) are generated to characterize the impact on $P_{1D,α}$ of continuum fitting, spectral resolution, noise modelling, and spectrum pixel masking. From this, we derive empirical corrections of these effects and apply them to the data measurement.

## 5.1 Synthetic data sets

We generated a set of DESI-Lite (Karaçaylı, Font-Ribera & Padmanabhan 2020) mocks, specifically designed for $P_{1D,α}$. The full description of these mocks is given in Karaçaylı et al. (2023a). The DESI-LITE software produces uncorrelated Ly α forests that mimic the redshift and noise distribution of the SV3 + M2 data set. Ten independent realizations are generated with different initial conditions. For each realization, a random catalogue of DLAs is created to follow the redshift and column density distribution of the latest eBOSS catalogue (Chabanier et al. 2021).

The quickquasars software https://github.com/desihub/desisim/blob/main/py/desisim/scripts/quickquasars.py (Herrera-Alcantar et al., in preparation), included in the desisim package, transforms Ly α transmission into spectra with observational and astrophysical contaminants. For each realization, two sets of spectra are generated by imprinting DLAs according to the catalog aforementioned or not.

We run five different $P_{1D,α}$ variations to study the impact of different Ly α contaminants:

(i) TRUECONT: The true quasar continuum imposed by quickquasars is applied instead of the continuum fitting procedure described in Section 3.1. In comparison to RAW mocks, this realization is impacted by finite noise and resolution.

(ii) CONT: The Ly α contrasts are calculated using the pipeline detailed in Section 3.1. This type of mocks includes the impact of continuum fitting.

(iii) DLAm: Realization for which the DLAs are not added to forests at the quickquasars stage, though we mask spectrum pixels as if they were present. The objective of this kind of mocks is to characterize the impact of DLA masking.

(iv) LINEm: Similarly to DLAm but masking the atmospheric emission lines catalogue built in Section 4.4.

(v) DLA: Realization for which the DLA are applied to the spectra without masking them. The objective of this mocks is to measure the impact of DLA to compute a DLA completeness systematic error in Section 6.

For all the mocks, we take the same procedure as for the FFT calculation on observational data. In particular, the same SNR-weighting is applied, and the number of sub-forest for each realization is around 81 500 with a small statistical variation between realizations. This is slightly larger than the number of sub-forest of the data sample given in Section 7.

In the next sections, all the results are shown for the combination of ten independent realizations. To decrease the error bars, we also performed a linear rebin that provides a wavenumber binning three times coarser than what we used for observational data. The error bars of the presented ratios are computed in quadrature.

<figure><img src="image_15.png" /><figcaption>Figure 13. Ratios between the power spectrum obtained using true continuum (TRUECONT) and the one derived with our pipeline (CONT) on the combination of 10 mocks. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Fitting functions are represented by continuous lines, and used to correct the $P_{1D,α}$ measurement. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

## 5.2 Continuum-fitting correction

The continuum-fitting procedure defined in Section 3.1 systematically distort the measured $C_q(\lambda, z_q)\bar{F}(\lambda)$ term by suppressing large-scale modes, and may bias the $P_{1D,α}$ measurement. This is a well-known effect in BAO measurements (Bourboux et al. 2020). To create a correction which contains this effect, we compare the mock computed using the true continuum (TRUECONT) with the one which follows the standard continuum fitting procedure (CONT):

$$A_{\text{cont}}(z, k) = \frac{P_{1D,\alpha,\text{TRUECONT}}(k, z)}{P_{1D,\alpha,\text{CONT}}(k, z)} \quad (21)$$

This correction is shown for the combination of ten DESI-Lite mocks in Fig. 13. We use a second-order polynomial function to fit this correction and apply it to the $P_{1D,α}$ measurement.

This correction differs in amplitude compared to the eBOSS measurement (Chabanier et al. 2019). As for eBOSS, the 1D power spectrum with continuum fitting is higher than that measured with the true continuum. However, in our case, the impact is much smaller than eBOSS, for which this ratio was near 4 per cent (without using a first-order polynomial function in the continuum fitting). Furthermore, we do not have a large-scale impact as significant as eBOSS.

## 5.3 Spectrum pixel masking

For both DLA and atmospheric line masking, we remove some data points from the measured spectra. This does not impact studies performed on real-space spectra, such as the continuum fitting or the QMLE for $P_{1D,α}$ (Karaçaylı et al. 2022). On the other hand, the FFT calculation requires that spectrum pixels are equally spaced. Consequently, when computing the Fourier transform, we impose a value of $\delta_F = 0$ (equivalent to mean transmitted flux fraction value for
<figure><img src="image_15.png" /><figcaption>Figure 14. Ratio between the unmasked (CONT) and masked (LINEm) power spectra (equation 22) for atmospheric line masking on the combination of the 10 mocks. The DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) atmospheric line catalogue is used. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Second-order polynomial functions are employed to fit the corrections in each redshift bin. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

F) and infinite standard deviation to the masked spectrum pixels. This masking introduces a k-dependent bias, which we need to quantify.

In order to determine and correct this bias in our $P_{1D,\alpha}$ measurement, we compare mocks for which DLAs or atmospheric lines are masked (respectively DLAm and LINEm) with mocks where no masking is applied (CONT). On those two mocks, the DLAs and atmospheric emission lines are not imposed on spectra. We want to derive only the corrections of masking in order to apply them on data. The coefficients used for both masking corrections are defined as the ratio between the unmasked and the masked power spectra:

$$A_{\text{line}}(k, z) = \frac{P_{1D,\alpha,\text{CONT}}(k, z)}{P_{1D,\alpha,\text{LINEm}}(k, z)}. \qquad (22)$$

$$A_{\text{dla}}(k, z) = \frac{P_{1D,\alpha,\text{CONT}}(k, z)}{P_{1D,\alpha,\text{DLAm}}(k, z)}. \qquad (23)$$

### 5.3.1 Atmospheric emission lines

The correction induced by the DESI atmospheric line mask, as defined in Section 4.4, is shown in Fig. 14. We verified that, for all redshifts where the masked Ly $\alpha$ forest length of DESI is close to eBOSS, the impact of masking is lower in the DESI case. It indicates that applying thinner masks to our measurement mitigates the impact of masking.

As expected, the correction roughly scales with the number of masked spectrum pixels. The effect of masking is a relatively smooth function of wavenumber, and its main impact is at low wavenumber. The most impacted redshift bins are $z = 2.2$ (Cosgalactic absorption lines), $z = 2.6$ (lines at 4360 Å in the transmission dip), and at high redshift for which many atmospheric lines need to be masked. As shown in Fig. 12, redshifts $z = 2.4, 2.8, 3.0, 3.2$ have no masks

applied, and only a few for $z = 3.4$. It is also in agreement with the level of corrections. The impact of atmospheric line masking is qualitatively in agreement with the eBOSS results in Chabanier et al. (2019). We choose to model $A_{\text{line}}(k, z)$ by a second-order polynomial fit and use this correction in the final calculation of $P_{1D,\alpha}$.

### 5.3.2 DLA masking

DLAs are added at random locations in the Ly $\alpha$ forest during the creation of the mocks. For this study, we do not attempt to characterize the completeness of the DLA finder applied to the data, and we use a ‘truth’ DLA catalogue for masking.

We mask the ‘truth’ catalog with the same parameters as the DLA data catalogue, i.e. for $N_{\text{HI}} > 10^{20.3} \text{ cm}^{-2}$. The correction induced by the masking, $A_{\text{dla}}(k, z)$, is represented in Fig. 15. As it was already seen in the eBOSS measurement (Chabanier et al. 2019), the DLA masking has a small impact compare to atmospheric emission lines. This is due to the random distribution of DLAs and the smaller masking in terms of Ly $\alpha$ forest length. As the impact is very similar for all wavenumbers, we apply a k-independent correction $A_{\text{dla}}(k, z) = A_{\text{dla}}(z)$, whose amplitude is 0.5 per cent on average.

## 6 UNCERTAINTY ESTIMATION

The statistical uncertainty of our averaged $P_{1D,\alpha}$ measurement, noted $\sigma_{\text{stat}}$, is obtained during the SNRweighting scheme presented in Appendix B. For each $(k, z)$ bin, a binned histogram of standard deviation as a function of SNR is derived. Fitting this histogram provides a function $\sigma_{k,z}(\text{SNR})$ that is used to define the statistical uncertainty:

$$\sigma_{\text{stat}}(k, z) = \sqrt{\frac{1}{\sum_{i} (\sigma_{k,z}(\text{SNR}_i))^{-2}}}, \qquad (24)$$

where $\sigma$ is the fitted function, and the i index runs over the SNR bins chosen.

The obtained statistical uncertainties are shown in Fig. 16. Despite using a SNR-dependent weighting in the $P_{1D,\alpha}$ calculation compared to a redshift-dependent SNR cut as eBOSS, we find similar trends as in Chabanier et al. (2019). The statistical uncertainty depends mainly on the number of sub-forests for each redshift bin, so that $\sigma_{\text{stat}}$ is an increasing function of redshift. This error bar depends also on the power spectrum level, as pointed out by the right panel of Fig. 16, for which low-redshift bins are not separated. In the case of $z \sim 2.2-2.4$, at small scales ($k \gtrsim 1.5 \text{ Å}^{-1}$), the statistical uncertainties are crossing each other. It is caused by the large noise increase in the blue spectral band due to atmospheric absorptions.

Looking at wavenumber dependence, $\sigma_{\text{stat}}$ increases as a function of k for small scales ($k \gtrsim 1.5 \text{ Å}^{-1}$). This is due to the resolution correction, which effectively increases the rms of individual $P_{1D,\alpha}$. At large scales ($k \lesssim 1 \text{ Å}^{-1}$), $\sigma_{\text{stat}}$ is a decreasing function of k, mainly due to the decrease in Fourier modes available to compute $P_{1D,\alpha}$.

In our study, we characterized the impact of several instrumental and astrophysical contaminants. From this extensive study, we associate systematic errors, noted $\sigma_{\text{syst}}$, to our $P_{1D,\alpha}$ measurement. Fig. 17 shows the systematic uncertainties for different redshift bins and their relative values with respect to statistical errors. Similarly to Chabanier et al. (2019), we made conservative choices to define those uncertainties:

(i) Noise estimation: As presented in Section 4.2, the pipeline noise is corrected using the $\alpha$ corrective term, which depends on
<figure><img src="image_16.png" /><figcaption>**Figure 14.** Ratio between the unmasked (CONT) and masked (LINEm) power spectra (equation 22) for atmospheric line masking on the combination of the 10 mocks. The DESI ($\Lambda_1$ = 2.5, $\Lambda_w$ = 1.2) atmospheric line catalogue is used. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Second-order polynomial functions are employed to fit the corrections in each redshift bin. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

F) and infinite standard deviation to the masked spectrum pixels. This masking introduces a k-dependent bias, which we need to quantify.

In order to determine and correct this bias in our $P_{1D,\alpha}$ measurement, we compare mocks for which DLAs or atmospheric lines are masked (respectively DLAm and LINEm) with mocks where no masking is applied (CONT). On those two mocks, the DLAs and atmospheric emission lines are not imposed on spectra. We want to derive only the corrections of masking in order to apply them on data. The coefficients used for both masking corrections are defined as the ratio between the unmasked and the masked power spectra:

$$A_{\text{line}}(k, z) = \frac{P_{1D,\alpha,\text{CONT}}(k, z)}{P_{1D,\alpha,\text{LINEm}}(k, z)}. \qquad (22)$$

$$A_{\text{dla}}(k, z) = \frac{P_{1D,\alpha,\text{CONT}}(k, z)}{P_{1D,\alpha,\text{DLAm}}(k, z)}. \qquad (23)$$

### 5.3.1 Atmospheric emission lines

The correction induced by the DESI atmospheric line mask, as defined in Section 4.4, is shown in Fig. 14. We verified that, for all redshifts where the masked Ly $\alpha$ forest length of DESI is close to eBOSS, the impact of masking is lower in the DESI case. It indicates that applying thinner masks to our measurement mitigates the impact of masking.

As expected, the correction roughly scales with the number of masked spectrum pixels. The effect of masking is a relatively smooth function of wavenumber, and its main impact is at low wavenumber. The most impacted redshift bins are $z = 2.2$ (Cosgalactic absorption lines), $z = 2.6$ (lines at 4360 Å in the transmission dip), and at high redshift for which many atmospheric lines need to be masked. As shown in Fig. 12, redshifts $z = 2.4, 2.8, 3.0, 3.2$ have no masks

applied, and only a few for $z = 3.4$. It is also in agreement with the level of corrections. The impact of atmospheric line masking is qualitatively in agreement with the eBOSS results in Chabanier et al. (2019). We choose to model $A_{\text{line}}(k, z)$ by a second-order polynomial fit and use this correction in the final calculation of $P_{1D,\alpha}$.

### 5.3.2 DLA masking

DLAs are added at random locations in the Ly $\alpha$ forest during the creation of the mocks. For this study, we do not attempt to characterize the completeness of the DLA finder applied to the data, and we use a ‘truth’ DLA catalogue for masking.

We mask the ‘truth’ catalog with the same parameters as the DLA data catalogue, i.e. for $N_{\text{HI}} > 10^{20.3} \text{ cm}^{-2}$. The correction induced by the masking, $A_{\text{dla}}(k, z)$, is represented in Fig. 15. As it was already seen in the eBOSS measurement (Chabanier et al. 2019), the DLA masking has a small impact compare to atmospheric emission lines. This is due to the random distribution of DLAs and the smaller masking in terms of Ly $\alpha$ forest length. As the impact is very similar for all wavenumbers, we apply a k-independent correction $A_{\text{dla}}(k, z) = A_{\text{dla}}(z)$, whose amplitude is 0.5 per cent on average.

## 6 UNCERTAINTY ESTIMATION

The statistical uncertainty of our averaged $P_{1D,\alpha}$ measurement, noted $\sigma_{\text{stat}}$, is obtained during the SNRweighting scheme presented in Appendix B. For each $(k, z)$ bin, a binned histogram of standard deviation as a function of SNR is derived. Fitting this histogram provides a function $\sigma_{k,z}(\text{SNR})$ that is used to define the statistical uncertainty:

$$\sigma_{\text{stat}}(k, z) = \sqrt{\frac{1}{\sum_{i} (\sigma_{k,z}(\text{SNR}_i))^{-2}}}, \qquad (24)$$

where $\sigma$ is the fitted function, and the i index runs over the SNR bins chosen.

The obtained statistical uncertainties are shown in Fig. 16. Despite using a SNR-dependent weighting in the $P_{1D,\alpha}$ calculation compared to a redshift-dependent SNR cut as eBOSS, we find similar trends as in Chabanier et al. (2019). The statistical uncertainty depends mainly on the number of sub-forests for each redshift bin, so that $\sigma_{\text{stat}}$ is an increasing function of redshift. This error bar depends also on the power spectrum level, as pointed out by the right panel of Fig. 16, for which low-redshift bins are not separated. In the case of $z \sim 2.2-2.4$, at small scales ($k \gtrsim 1.5 \text{ Å}^{-1}$), the statistical uncertainties are crossing each other. It is caused by the large noise increase in the blue spectral band due to atmospheric absorptions.

Looking at wavenumber dependence, $\sigma_{\text{stat}}$ increases as a function of k for small scales ($k \gtrsim 1.5 \text{ Å}^{-1}$). This is due to the resolution correction, which effectively increases the rms of individual $P_{1D,\alpha}$. At large scales ($k \lesssim 1 \text{ Å}^{-1}$), $\sigma_{\text{stat}}$ is a decreasing function of k, mainly due to the decrease in Fourier modes available to compute $P_{1D,\alpha}$.

In our study, we characterized the impact of several instrumental and astrophysical contaminants. From this extensive study, we associate systematic errors, noted $\sigma_{\text{syst}}$, to our $P_{1D,\alpha}$ measurement. Fig. 17 shows the systematic uncertainties for different redshift bins and their relative values with respect to statistical errors. Similarly to Chabanier et al. (2019), we made conservative choices to define those uncertainties:

(i) Noise estimation: As presented in Section 4.2, the pipeline noise is corrected using the $\alpha$ corrective term, which depends on
<figure><img src="image_16.png" /><figcaption>**Figure 15.** Ratio between the unmasked (CONT) and masked (DLAm) power spectra (equation 23) for DLA masking on the combination of the 10 mocks. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Continuous lines show constant fits are used for correction. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

the data set considered. We assign a systematic uncertainty equal to 30 per cent of the average $\alpha$ for each redshift bin.

(ii) *Resolution*: We fit the average resolution correction given in Fig. 5 by a simplified model $\exp(-0.5(k\Delta\lambda)^2) \cdot \text{sinc}(0.5k\Delta\lambda_{\text{pix}})$ with $\Delta\lambda_{\text{pix}}$ fixed to the DESI spectral pixel separation. This procedure allows to determine the effective spectral resolution $\Delta\lambda$ in Å. Measurements of DESI PSF stability, shown in Abareshi et al. (2022, fig. 30), indicate that its fractional change is less than 1 per cent over all spectrographs. We therefore assign a conservative systematic uncertainty $\sigma_{\Delta\lambda} = 1$ per cent$\Delta\lambda$. Using the above-mentioned

simplified resolution model, this translates into a $P_{1\text{D},\alpha}$ uncertainty equal to $2k^2\Delta\lambda\sigma_{\Delta\lambda} \cdot P_{1\text{D},\alpha}(k)$.

(iii) *Resolution correction:* We apply a correction to the resolution modeling as presented in Section 4.1, by multiplying $A_{\text{res}}$ to $P_{1\text{D},\alpha}$. We add an associated systematic error defined as 30 per cent of this correction.

(iv) *Side-band:* The fitted side-band power spectrum $P_{\text{SB1,m}}$ measured in Section 4.3 is subtracted to $P_{1\text{D},\alpha}$ to account essentially for metal absorptions in the Ly $\alpha$ forest region. We associate with this correction a systematic uncertainty equal to the statistical errors of the measured SB1 power spectrum. This is a conservative choice, as the modelling performed in Section 4.3 closely reproduces $P_{\text{SB}}$.

(v) *Spectrum pixel masking:* The impact of masking DLAs and atmospheric emission lines on the $P_{1\text{D},\alpha}$ measurement was determined with synthetic data in Section 5.3. Spectrum pixel masking is corrected by multiplying $A_{\text{line}}(z, k) \cdot A_{\text{dla}}(z, k)$ to the $P_{1\text{D},\alpha}$ estimator. We define the systematic error associated with each masking as 30 per cent of this correction.

(vi) *Continuum fitting:* Similarly, we assign a systematic error of 30 per cent times the $A_{\text{cont}}(z, k)$ correction computed in Section 5.2.

(vii) *DLA completeness:* Using the synthetic data described in Section 5, we derived the impact of DLA on the 1D power spectrum as the ratio between mock with DLA (DLA) and without (CONT). We fit this ratio with an adapted function provided by Rogers et al. (2017). As detailed previously in Section 2.4, our DLA catalogue of data results from the combination of two finders. The trend of this ratio is reported on the penultimate panel of Fig. 17. In Chabanier et al. (2021), the authors perform a full study on eBOSS data and provide the completeness of the CNN finder. The completeness of this finder is higher than 85 per cent for $\log(N_{\text{H}}) > 20.3$. To be conservative, we choose to associate an uncertainty of 15 per cent of the total impact of DLAs on $P_{1\text{D},\alpha}$ to the incompleteness of our catalogue. We stress that this uncertainty is overestimated, as the CNN finder has a higher completeness for DLAs with higher column density and since we are using an additional GP algorithm.

For all the corrections we applied on our $P_{1\text{D},\alpha}$ measurement, the choice of 30 per cent in the associated systematic uncertainties is motivated by the fact that we consider a shift randomly ranging between no correction and 100 per cent of the correction. It is

<figure><img src="image_17.png" /><figcaption>**Figure 16.** Statistical uncertainties (left) and its relative value with respect to $P_{1\text{D},\alpha}$ (right) of the DESI SV3 + M2 measurement in Å, as a function of wavenumber.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>**Figure 15.** Ratio between the unmasked (CONT) and masked (DLAm) power spectra (equation 23) for DLA masking on the combination of the 10 mocks. Each power spectrum is re-binned by a factor of 3 to reduce error bars. Continuous lines show constant fits are used for correction. For clarity, we artificially offset the points corresponding to different redshifts.</figcaption></figure>

the data set considered. We assign a systematic uncertainty equal to 30 per cent of the average $\alpha$ for each redshift bin.

(ii) *Resolution*: We fit the average resolution correction given in Fig. 5 by a simplified model $\exp(-0.5(k\Delta\lambda)^2) \cdot \text{sinc}(0.5k\Delta\lambda_{\text{pix}})$ with $\Delta\lambda_{\text{pix}}$ fixed to the DESI spectral pixel separation. This procedure allows to determine the effective spectral resolution $\Delta\lambda$ in Å. Measurements of DESI PSF stability, shown in Abareshi et al. (2022, fig. 30), indicate that its fractional change is less than 1 per cent over all spectrographs. We therefore assign a conservative systematic uncertainty $\sigma_{\Delta\lambda} = 1$ per cent$\Delta\lambda$. Using the above-mentioned

simplified resolution model, this translates into a $P_{1\text{D},\alpha}$ uncertainty equal to $2k^2\Delta\lambda\sigma_{\Delta\lambda} \cdot P_{1\text{D},\alpha}(k)$.

(iii) *Resolution correction*: We apply a correction to the resolution modeling as presented in Section 4.1, by multiplying $A_{\text{res}}$ to $P_{1\text{D},\alpha}$. We add an associated systematic error defined as 30 per cent of this correction.

(iv) *Side-band*: The fitted side-band power spectrum $P_{\text{SB1,m}}$ measured in Section 4.3 is subtracted to $P_{1\text{D},\alpha}$ to account essentially for metal absorptions in the Ly $\alpha$ forest region. We associate with this correction a systematic uncertainty equal to the statistical errors of the measured SB1 power spectrum. This is a conservative choice, as the modelling performed in Section 4.3 closely reproduces $P_{\text{SB}}$.

(v) *Spectrum pixel masking*: The impact of masking DLAs and atmospheric emission lines on the $P_{1\text{D},\alpha}$ measurement was determined with synthetic data in Section 5.3. Spectrum pixel masking is corrected by multiplying $A_{\text{line}}(z, k) \cdot A_{\text{dla}}(z, k)$ to the $P_{1\text{D},\alpha}$ estimator. We define the systematic error associated with each masking as 30 per cent of this correction.

(vi) *Continuum fitting*: Similarly, we assign a systematic error of 30 per cent times the $A_{\text{cont}}(z, k)$ correction computed in Section 5.2.

(vii) *DLA completeness*: Using the synthetic data described in Section 5, we derived the impact of DLA on the 1D power spectrum as the ratio between mock with DLA (DLA) and without (CONT). We fit this ratio with an adapted function provided by Rogers et al. (2017). As detailed previously in Section 2.4, our DLA catalogue of data results from the combination of two finders. The trend of this ratio is reported on the penultimate panel of Fig. 17. In Chabanier et al. (2021), the authors perform a full study on eBOSS data and provide the completeness of the CNN finder. The completeness of this finder is higher than 85 per cent for $\log(N_{\text{H}}) > 20.3$. To be conservative, we choose to associate an uncertainty of 15 per cent of the total impact of DLAs on $P_{1\text{D},\alpha}$ to the incompleteness of our catalogue. We stress that this uncertainty is overestimated, as the CNN finder has a higher completeness for DLAs with higher column density and since we are using an additional GP algorithm.

For all the corrections we applied on our $P_{1\text{D},\alpha}$ measurement, the choice of 30 per cent in the associated systematic uncertainties is motivated by the fact that we consider a shift randomly ranging between no correction and 100 per cent of the correction. It is

<figure><img src="image_18.png" /><figcaption>**Figure 16.** Statistical uncertainties (left) and its relative value with respect to $P_{1\text{D},\alpha}$ (right) of the DESI SV3 + M2 measurement in Å, as a function of wavenumber.</figcaption></figure>
<figure><img src="image_18.png" /><figcaption>**Figure 17.** Value of the systematic uncertainties $\sigma_{\text{syst}}$ in Å, for different redshift bins, on the DESI SV3 + M2P<sub>1D,α</sub> measurement. Each line is associated to a systematic considered in this article. The left panels show the absolute uncertainties, and the right panels, their relative values with respect to statistical uncertainties showed in Fig. 16.</figcaption></figure>

described by a uniform distribution between 0 and 1. The standard deviation of the distribution, equal to 0.30, quantifies the spread among the possible values, leading to a systematic uncertainty equal to 30 per cent of the correction.

Opposite to eBOSS (Chabanier et al. 2019), we chose to not account for the incompleteness of the BAL catalog in our analysis, as it is one of the weaker contaminants. The study of BAL catalogue completeness will be performed on further studies.
<figure><img src="image_19.png" /><figcaption>**Figure 17.** Value of the systematic uncertainties $\sigma_{\text{syst}}$ in Å, for different redshift bins, on the DESI SV3 + M2P<sub>1D,α</sub> measurement. Each line is associated to a systematic considered in this article. The left panels show the absolute uncertainties, and the right panels, their relative values with respect to statistical uncertainties showed in Fig. 16.</figcaption></figure>

described by a uniform distribution between 0 and 1. The standard deviation of the distribution, equal to 0.30, quantifies the spread among the possible values, leading to a systematic uncertainty equal to 30 per cent of the correction.

Opposite to eBOSS (Chabanier et al. 2019), we chose to not account for the incompleteness of the BAL catalog in our analysis, as it is one of the weaker contaminants. The study of BAL catalogue completeness will be performed on further studies.
The general trends in Fig. 17 are similar to those of the eBOSS measurement (Chabanier et al. 2019). Given our limited statistics, most of the systematic errors are smaller than the statistical uncertainties for all redshift bins and all scales. However, this is expected to change for future DESI measurements, which will offer unprecedented statistics, thus reducing the statistical errors.

There is room for improvement for the major source of systematic uncertainties presented above. The noise modelling can be improved by understanding and correcting the source of unaccounted noise. Regarding the resolution modelling, the mathematical model and its verification with the relatively new CCD mocks presented in Section 4.1 can be improved with additional tests and larger data sets. Decreasing statistical error on the side-band power spectrum will directly reduce the associated systematic error. Concerning the spectrum pixel masking, especially for the atmospheric lines, thinner masks can be applied considering the improvement of atmospheric emission line subtraction in DESI (see Guy et al. 2022). Furthermore, a more complex but analytical mathematical correction could be derived for this regular masking. For the DLA completeness, a more advanced study as the one performed in Chabanier et al. (2019) is needed to reduce the associated systematic uncertainties. Additionally, as shown in Rogers et al. (2017), the impact and correction of objects with lower column density (sub-DLAs, Lyman limit systems, etc.) should be accounted for in a future, more developed study. Finally, the systematics are defined from a simplified assumption of uncertainty propagation, along with conservative choices. We plan to improve the modeling of each systematic and obtain more reliable uncertainties by modeling them directly into large samples of synthetic data.

## 7 DESI MEASUREMENT

We apply the methodology and corrections described in previous sections to the combination of SV3, and M2 data sets, noted SV3 + M2. We choose to remove the SV1 data set from this measurement due to considerations on the noise power spectrum shown in Section 4.2 and measurements performed in the Appendix A.

The $P_{1D,α}$ measurement is done using the pipeline and parameters presented in Section 3. Considering all the corrections defined in the previous section, the final $P_{1D,α}$ estimator is defined by

$$ P_{1D,α}(k) = A_{line}(z, k) \cdot A_{dla}(z, k) \cdot A_{cont}(z, k) \cdot A_{res}(k) \cdot \left( \langle [P_{raw}(k) - P_{pipeline}(k) - \alpha] \cdot \mathbf{R}^{-2}(k) \rangle - P_{SB1,m}(k) \rangle \right). \quad (25) $$

Fig. 18 presents the normalized $P_{1D,α}$ measurement such that $\Delta_{1D,α}^2 = kP_{1D,α}/\pi$. This observable is shown for 9 redshift bins ranging from 2.2 to 3.8, for wavenumbers $0.145 \le k \le 2 \text{ Å}^{-1}$, and using a total of 73 839 subforests extracted from 26 330 quasar spectra. The represented error bars are the statistical and systematic uncertainties added in quadrature. The details of subforest properties for each redshift bin are given in Table 4.

### 7.1 Comparison with other measurements

We perform a comparison with past measurements. The comparison of the DESI $P_{1D,α}$ with the last moderate-resolution measurement (eBOSS survey Chabanier et al. 2019) is shown in Fig. 19 (left). We also compare our measurement with the last high-resolution measurement obtained using the combination of KODIAQ, SQUAD, and XQ-100 surveys (Karaçaylı et al. 2022) in Fig. 19 (right).

Both moderate- and high-resolution measurements were expressed in velocity units (s km$^{-1}$). The conversion between this unit system

and the one used for DESI (Å$^{-1}$) is defined by

$$ k [\text{s} \cdot \text{km}^{-1}] = k [\text{Å}^{-1}] \times \lambda_{\alpha}(1+z)/c. \quad (26) $$

This conversion is performed at the stage of the FFT estimation by converting the terms in equation (14) in velocity units before performing the ensemble average to compute $P_{1D,α}$. We found that doing this conversion on the averaged $P_{1D,α}$ significantly shift $P_{1D,α}$ because of the redshift term in the equation (26). This is caused by the non-uniform redshift distribution in each z-bins, pointed out by Fig. 2. Similarly, we converted all corrections terms appearing in equation (25) in velocity units. For the comparison with eBOSS measurement, we compute $P_{1D,α}$ in velocity units with the same binning. For the high-resolution measurement, we rebin our DESI measurement to the same wavenumber binning, and account this rebinning in the calculation of error bars, to obtain a fair comparison.

The comparison with the eBOSS measurement in Fig. 19 (left) yields a ~15 per cent difference at small wavenumber ($k < 0.01 \text{ km s}^{-1}$). To investigate this discrepancy, we performed a detailed investigation by varying most of the parameters of our analysis, which are susceptible to impact small scales. The discussion concerning those tests is detailed in Appendix D. On the other scales ($0.01 \text{ km s}^{-1} < k < 0.02 \text{ km s}^{-1}$), our measurement agrees with eBOSS considering the error bars.

The major improvement lays at large wavenumber ($k > 0.015 \text{ km s}^{-1}$) where the improved DESI resolution and noise modelling allows us reaching much smaller scales than eBOSS, especially for high redshifts. We are able to conservatively reach the wavenumber $k_{max} = 2.0 \text{ Å}^{-1}$ for all redshift bins. In comparison, the eBOSS measurement (Chabanier et al. 2019) achieved maximal wavenumber of $k_{max} = 1.54 \text{ Å}^{-1}$ for $z = 2.2$ and $k_{max} = 1.03 \text{ Å}^{-1}$ for $z = 3.8$. At large wavenumbers, the eBOSS measurement is highly contaminated by noise and resolution. We consider that our measurement is more suitable to probe those scales, considering the improvement in resolution and noise estimation.

We compare our measurement with the high-resolution measurements in Fig. 19 (right). This measurement (Karaçaylı et al. 2022) is performed with a statistically smaller sample of quasars but with a very high spectral resolution (5000 $\le R \le 60 000$) and SNR. The high-resolution $P_{1D,α}$ measurements thus allow reaching very small scales with large error bars. Our agreement shows a 20 per cent disagreement at small scales ($k > 0.02 \text{ km s}^{-1}$), mostly for the lower redshifts measured. It indicates that there is still room for improvement in our measurement's noise and resolution modeling. We note that considering error bars, both measurements agree on intermediate scales.

## 8 CONCLUSION AND PROSPECTS

We performed the first measurement of the 1D power spectrum ($P_{1D,α}$) with DESI data. The main objective of this paper is to carefully characterize the different contaminants of $P_{1D,α}$ with regard to DESI instrument. In particular, we modelled the noise and spectral resolution of DESI. In comparison to the previous eBOSS measurement (Chabanier et al. 2019), we improved the analysis of side-band power spectrum and atmospheric emission lines.

We used adapted synthetic data to correct the impact of spectrum pixel masking, continuum fitting and spectral resolution modelling. We performed a complete review of the systematic uncertainties linked to the $P_{1D,α}$ pipeline and compared the DESI measurement with previous moderate- and high-resolution measurements. We find a relatively good agreement, except for a slight difference at
The general trends in Fig. 17 are similar to those of the eBOSS measurement (Chabanier et al. 2019). Given our limited statistics, most of the systematic errors are smaller than the statistical uncertainties for all redshift bins and all scales. However, this is expected to change for future DESI measurements, which will offer unprecedented statistics, thus reducing the statistical errors.

There is room for improvement for the major source of systematic uncertainties presented above. The noise modelling can be improved by understanding and correcting the source of unaccounted noise. Regarding the resolution modelling, the mathematical model and its verification with the relatively new CCD mocks presented in Section 4.1 can be improved with additional tests and larger data sets. Decreasing statistical error on the side-band power spectrum will directly reduce the associated systematic error. Concerning the spectrum pixel masking, especially for the atmospheric lines, thinner masks can be applied considering the improvement of atmospheric emission line subtraction in DESI (see Guy et al. 2022). Furthermore, a more complex but analytical mathematical correction could be derived for this regular masking. For the DLA completeness, a more advanced study as the one performed in Chabanier et al. (2019) is needed to reduce the associated systematic uncertainties. Additionally, as shown in Rogers et al. (2017), the impact and correction of objects with lower column density (sub-DLAs, Lyman limit systems, etc.) should be accounted for in a future, more developed study. Finally, the systematics are defined from a simplified assumption of uncertainty propagation, along with conservative choices. We plan to improve the modeling of each systematic and obtain more reliable uncertainties by modeling them directly into large samples of synthetic data.

## 7 DESI MEASUREMENT

We apply the methodology and corrections described in previous sections to the combination of SV3, and M2 data sets, noted SV3 + M2. We choose to remove the SV1 data set from this measurement due to considerations on the noise power spectrum shown in Section 4.2 and measurements performed in the Appendix A.

The $P_{1D,α}$ measurement is done using the pipeline and parameters presented in Section 3. Considering all the corrections defined in the previous section, the final $P_{1D,α}$ estimator is defined by

$$ P_{1D,α}(k) = A_{line}(z, k) \cdot A_{dla}(z, k) \cdot A_{cont}(z, k) \cdot A_{res}(k) \cdot \left( \langle \left[ P_{raw}(k) - P_{pipeline}(k) - \alpha \right] \cdot \mathbf{R}^{-2}(k) \rangle - P_{SB1,m}(k) \right). \quad (25) $$

Fig. 18 presents the normalized $P_{1D,α}$ measurement such that $\Delta_{1D,α}^2 = kP_{1D,α}/\pi$. This observable is shown for 9 redshift bins ranging from 2.2 to 3.8, for wavenumbers $0.145 \le k \le 2 \text{ Å}^{-1}$, and using a total of 73 839 subforests extracted from 26 330 quasar spectra. The represented error bars are the statistical and systematic uncertainties added in quadrature. The details of subforest properties for each redshift bin are given in Table 4.

### 7.1 Comparison with other measurements

We perform a comparison with past measurements. The comparison of the DESI $P_{1D,α}$ with the last moderate-resolution measurement (eBOSS survey Chabanier et al. 2019) is shown in Fig. 19 (left). We also compare our measurement with the last high-resolution measurement obtained using the combination of KODIAQ, SQUAD, and XQ-100 surveys (Karaçaylı et al. 2022) in Fig. 19 (right).

Both moderate- and high-resolution measurements were expressed in velocity units (s km$^{-1}$). The conversion between this unit system

and the one used for DESI (Å$^{-1}$) is defined by

$$ k [\text{s} \cdot \text{km}^{-1}] = k [\text{Å}^{-1}] \times \lambda_{\alpha}(1+z)/c. \quad (26) $$

This conversion is performed at the stage of the FFT estimation by converting the terms in equation (14) in velocity units before performing the ensemble average to compute $P_{1D,α}$. We found that doing this conversion on the averaged $P_{1D,α}$ significantly shift $P_{1D,α}$ because of the redshift term in the equation (26). This is caused by the non-uniform redshift distribution in each z-bins, pointed out by Fig. 2. Similarly, we converted all corrections terms appearing in equation (25) in velocity units. For the comparison with eBOSS measurement, we compute $P_{1D,α}$ in velocity units with the same binning. For the high-resolution measurement, we rebin our DESI measurement to the same wavenumber binning, and account this rebinning in the calculation of error bars, to obtain a fair comparison.

The comparison with the eBOSS measurement in Fig. 19 (left) yields a ~15 per cent difference at small wavenumber ($k < 0.01 \text{ km s}^{-1}$). To investigate this discrepancy, we performed a detailed investigation by varying most of the parameters of our analysis, which are susceptible to impact small scales. The discussion concerning those tests is detailed in Appendix D. On the other scales ($0.01 \text{ km s}^{-1} < k < 0.02 \text{ km s}^{-1}$), our measurement agrees with eBOSS considering the error bars.

The major improvement lays at large wavenumber ($k > 0.015 \text{ km s}^{-1}$) where the improved DESI resolution and noise modelling allows us reaching much smaller scales than eBOSS, especially for high redshifts. We are able to conservatively reach the wavenumber $k_{max} = 2.0 \text{ Å}^{-1}$ for all redshift bins. In comparison, the eBOSS measurement (Chabanier et al. 2019) achieved maximal wavenumber of $k_{max} = 1.54 \text{ Å}^{-1}$ for $z = 2.2$ and $k_{max} = 1.03 \text{ Å}^{-1}$ for $z = 3.8$. At large wavenumbers, the eBOSS measurement is highly contaminated by noise and resolution. We consider that our measurement is more suitable to probe those scales, considering the improvement in resolution and noise estimation.

We compare our measurement with the high-resolution measurements in Fig. 19 (right). This measurement (Karaçaylı et al. 2022) is performed with a statistically smaller sample of quasars but with a very high spectral resolution ($5000 \le R \le 60 000$) and SNR. The high-resolution $P_{1D,α}$ measurements thus allow reaching very small scales with large error bars. Our agreement shows a 20 per cent disagreement at small scales ($k > 0.02 \text{ km s}^{-1}$), mostly for the lower redshifts measured. It indicates that there is still room for improvement in our measurement’s noise and resolution modeling. We note that considering error bars, both measurements agree on intermediate scales.

## 8 CONCLUSION AND PROSPECTS

We performed the first measurement of the 1D power spectrum ($P_{1D,α}$) with DESI data. The main objective of this paper is to carefully characterize the different contaminants of $P_{1D,α}$ with regard to DESI instrument. In particular, we modelled the noise and spectral resolution of DESI. In comparison to the previous eBOSS measurement (Chabanier et al. 2019), we improved the analysis of side-band power spectrum and atmospheric emission lines.

We used adapted synthetic data to correct the impact of spectrum pixel masking, continuum fitting and spectral resolution modelling. We performed a complete review of the systematic uncertainties linked to the $P_{1D,α}$ pipeline and compared the DESI measurement with previous moderate- and high-resolution measurements. We find a relatively good agreement, except for a slight difference at
<figure><img src="image_19.png" /><figcaption>**Figure 18.** Normalized 1D Ly $\alpha$ forest power spectrum ($\Delta_{1D,\alpha}(k)$) using the SV3 + M2 data set, for redshift bins from $z = 2.2$ to $z = 3.8$. All the corrections given in equation (25) are applied to perform this measurement. As an illustration, wavenumbers in velocity space for different redshifts are represented at the top of the figure. Error bars are systematic and statistical uncertainties added in quadrature.</figcaption></figure>

<table><caption>**Table 4.** Number of subforest, average redshift, and signal-to-noise ratio for each redshift bins in the final data set sample used in this measurement.</caption><thead><tr><th></th><th>z bin</th><th>2.2</th><th>2.4</th><th>2.6</th><th>2.8</th><th>3.0</th><th>3.2</th><th>3.4</th><th>3.6</th><th>3.8</th></tr></thead><tbody><tr><th>#</th><td>25 846</td><td>16 371</td><td>11 501</td><td>8205</td><td>5343</td><td>3267</td><td>1800</td><td>1006</td><td>500</td></tr><tr><th>⟨z⟩</th><td>2.2</td><td>2.4</td><td>2.6</td><td>2.8</td><td>2.99</td><td>3.19</td><td>3.39</td><td>3.59</td><td>3.79</td></tr><tr><th>SNR</th><td>2.79</td><td>3.0</td><td>3.08</td><td>3.13</td><td>3.28</td><td>3.27</td><td>3.14</td><td>3.07</td><td>3.21</td></tr></tbody></table>

large scales with the eBOSS measurement, partially due to the different residual correction we apply. Our estimation of $P_{1D,\alpha}$ is also compared with the QMLE method in a companion paper Karaçaylı et al. (2023a). Measurements with FFT and QMLE methods agree at 1 per cent level precision up to half the Nyquist frequency.

The DESI spectral resolution is approximately two times better than SDSS. Consequently, our $P_{1D,\alpha}$ measurement is of high scientific interest to probe the small scales of the intergalactic medium. However, the data sets we exploited remains inferior to eBOSS in terms of statistics. If we apply the same SNR cut as in Chabanier et al. (2019), our sample contains 17 333 subforests compared to 94 558 for eBOSS. However, we expect future DESI data sets to provide high increase of statistics (up to 1 million Ly $\alpha$ forest, thus almost 3 million subforests). This unprecedented data

set will allow obtaining a subpercent precision measurement. The resulting $P_{1D,\alpha}$ measurement will provide stringent constraints on the sum of neutrinos masses, warm dark matter models, and on the parameters of the intergalactic medium (DESI 2016a; Valluri et al. 2022).

We plan to improve our analysis to keep the level of systematic error close to the statistical one. First, applying stricter constraints on a larger Ly $\alpha$ forest sample will be beneficial to reduce systematic uncertainties. Furthermore, we also plan to improve the treatment of contaminants presented in Section 6. In particular, we plan to test extensively the resolution and noise estimations on pixel-level simulations of the CCD camera.

High-resolution hydrodynamical simulations, associated with Gaussian processes emulator such as in Pedersen et al. (2020b)
<figure><img src="image_20.png" /><figcaption>**Figure 18.** Normalized 1D Ly $\alpha$ forest power spectrum ($\Delta_{1D,\alpha}(k)$) using the SV3 + M2 data set, for redshift bins from $z = 2.2$ to $z = 3.8$. All the corrections given in equation (25) are applied to perform this measurement. As an illustration, wavenumbers in velocity space for different redshifts are represented at the top of the figure. Error bars are systematic and statistical uncertainties added in quadrature.</figcaption></figure>

<table><caption>**Table 4.** Number of subforest, average redshift, and signal-to-noise ratio for each redshift bins in the final data set sample used in this measurement.</caption><thead><tr><th></th><th>z bin</th><th>2.2</th><th>2.4</th><th>2.6</th><th>2.8</th><th>3.0</th><th>3.2</th><th>3.4</th><th>3.6</th><th>3.8</th></tr></thead><tbody><tr><th>#</th><td>25 846</td><td>16 371</td><td>11 501</td><td>8205</td><td>5343</td><td>3267</td><td>1800</td><td>1006</td><td>500</td></tr><tr><th>⟨z⟩</th><td>2.2</td><td>2.4</td><td>2.6</td><td>2.8</td><td>2.99</td><td>3.19</td><td>3.39</td><td>3.59</td><td>3.79</td></tr><tr><th>SNR</th><td>2.79</td><td>3.0</td><td>3.08</td><td>3.13</td><td>3.28</td><td>3.27</td><td>3.14</td><td>3.07</td><td>3.21</td></tr></tbody></table>

large scales with the eBOSS measurement, partially due to the different residual correction we apply. Our estimation of $P_{1D,\alpha}$ is also compared with the QMLE method in a companion paper Karaçaylı et al. (2023a). Measurements with FFT and QMLE methods agree at 1 per cent level precision up to half the Nyquist frequency.

The DESI spectral resolution is approximately two times better than SDSS. Consequently, our $P_{1D,\alpha}$ measurement is of high scientific interest to probe the small scales of the intergalactic medium. However, the data sets we exploited remains inferior to eBOSS in terms of statistics. If we apply the same SNR cut as in Chabanier et al. (2019), our sample contains 17 333 subforests compared to 94 558 for eBOSS. However, we expect future DESI data sets to provide high increase of statistics (up to 1 million Ly $\alpha$ forest, thus almost 3 million subforests). This unprecedented data

set will allow obtaining a subpercent precision measurement. The resulting $P_{1D,\alpha}$ measurement will provide stringent constraints on the sum of neutrinos masses, warm dark matter models, and on the parameters of the intergalactic medium (DESI 2016a; Valluri et al. 2022).

We plan to improve our analysis to keep the level of systematic error close to the statistical one. First, applying stricter constraints on a larger Ly $\alpha$ forest sample will be beneficial to reduce systematic uncertainties. Furthermore, we also plan to improve the treatment of contaminants presented in Section 6. In particular, we plan to test extensively the resolution and noise estimations on pixel-level simulations of the CCD camera.

High-resolution hydrodynamical simulations, associated with Gaussian processes emulator such as in Pedersen et al. (2020b)