Improved Lyman Alpha Tomography using Optimized Reconstruction with Constraints on Absorption (ORCA)
DDraft version February 25, 2021
Typeset using L A TEX twocolumn style in AASTeX63
Improved Lyman Alpha Tomography using Optimized Reconstruction with Constraints onAbsorption (ORCA)
Zihao Li,
1, 2
Benjamin Horowitz, and Zheng Cai Department of Astronomy, Tsinghua University, Beijing 100084, People’s Republic of China Department of Aerospace Engineering, Sichuan University, Chengdu 610207, People’s Republic of China Department of Astrophysical Sciences, Princeton University, Princeton, NJ, USA
Submitted to APJABSTRACTIn this work we propose an improved approach to reconstruct the three dimensional intergalacticmedium from observed Lyman- α forest absorption features. We present our new method, the Opti-mized Reconstruction with Constraints on Absorption (ORCA), which outperforms the current baselineWiener Filter (WF) when tested on mock Lyman Alpha forest data generated from hydrodynamicalsimulations. We find that both reconstructed flux errors and cosmic web classification improve sub-stantially with ORCA, equivalent to 30-40% additional sight-lines with the standard WF. We use thismethod to identify and classify extremal objects, i.e. voids and (proto)-clusters, and find improvedreconstruction across all summary statistics explored. We apply ORCA to existing Lyman Alpha forestdata from the COSMOS Lyman Alpha Mapping and Tomography Observations (CLAMATO) Surveyand compare to the WF reconstruction. Keywords: cosmology: observations — galaxies: high-redshift — intergalactic medium — quasars:absorption lines — techniques: spectroscopic - methods: numerical INTRODUCTIONThe galactic superclusters, sheets, filaments and voidsare a significant part of our observable Universe, whichare also referred to as the large scale structure of the uni-verse. Cosmography is the science of mapping and de-scribing those features using a variety of probes. Whilethe local universe can be mapped out with galaxy sur-veys like SDSS(Eisenstein et al. 2011), all sky galaxysurveys are restricted to universe at low redshift sincethe surface brightness scales with redshift as ∝ (1+ z ) − ,and we even no longer resolve a galaxy disk at longer dis-tances, making it increasingly difficult to map the largescale structure at higher redshift.An alternative probe of large scale structure at z > α forest absorption in spectra of backgroundquasars and galaxies, which is complementary to lowredshift galaxy surveys. The Ly α forest traces the neu-tral hydrogen density (Gunn & Peterson 1965) which Corresponding author: Zihao Lileezihao [email protected] can correspond to the underlying dark matter den-sity through fluctuating Gunn-Peterson approximation(FGPA) (Croft et al. 1998). While the flux along a sin-gle line-of-sight towards a background source (quasar orgalaxy) only provides one dimensional information, 3-D map can be reconstructed using inversion methodsgiven a set of lines of sight (LOSs) toward a group ofbackground objects (Caucci et al. 2008). A long-usedmethod for this purpose is Wiener Filter (Pichon et al.2001), which has been validated to recover the dark mat-ter field by Caucci et al. (2008) and the observationalrequirements for implementing such method is discussedby Lee et al. (2014a). Since this early work, the WienerFilter has been widely used in many surveys includingCLAMATO (Lee et al. 2018), eBOSS (Ravoux et al.2020), and LATIS (Newman et al. 2020). Going for-ward, there is significant interest in expanding Ly α to-mography as a probe of the IGM over more sky areaand cross-correlate the properties of the IGM with over-lapping galaxy samples, such as in the already plannedSubaru Prime Focus Spectrograph (PFS) (Takada et al. a r X i v : . [ a s t r o - ph . GA ] F e b Li, Horowitz, and Cai a posteriori ini-tial density field which gives rise to the observed field.A unique advantage of TARDIS is that, owing to gravi-tational evolution, it yields both information on the un-derlying dark matter density field and on velocity al-lowing us to deconvolve redshift space and real spacequantities and provide a more accurate reconstruction.However, such a method depends on N-body simulationand FGPA model which inherently rely on cosmologi-cal and astrophysical assumptions. These methods re-constructed field all within known frames and this mayomit some unknowing physics (e.g., the underlying trueuniverse may not be well modeled by the current cos-mological model we use).One important feature of large scale structure arevoids, which occupies the majority of volume of cosmicweb. Voids are regions lacking large galaxy populationsand their density are below mean cosmic density (VanDe Weygaert & Platen 2011). However, understand-ing and cataloging these regions has been found to beimportant for a number of cosmological analysis includ-ing constraining neutrino properties (Kreisch et al. 2019;Zhang et al. 2020; Liu et al. 2020), dark energy (Lee &Park 2009; Bos et al. 2012), and modified gravity (Pericoet al. 2019; Contarini et al. 2020). Cosmic voids catalogshave been constructed from spectroscopic surveys (Maoet al. 2017) and photometric surveys (Fang et al. 2019).However, identifying cosmic voids at high redshift withstatistical significance is difficult using galaxies them-selves due to the low number density and a high samplevariance.The counterpart of voids are galactic (proto)-clusters,which are important for studying star formation historyand the origin of large scale structure. Observationsof low redshift galaxies in cluster environments haveshown older stellar population and lower star formationrates than those located in the field (Wake et al. 2005;Skibba et al. 2009), indicating that these cluster envi- ronments underwent a cycle of significant star formationand quenching at high redshift ( z > .
0) (Tran et al.2010), so-called “Cosmic Dawn.” Deep galaxy redshiftsurveys provide a promising way to study these (proto)-clusters, for example in the COSMOS field (Chiang et al.2014). However, similar to void discovery, identifying(proto)-clusters at high redshift is difficult using galaxiesresulting in finding only the most massive protoclustersin the deepest fields (such as COSMOS). Observationswith Ly α forest provide a useful alternative for iden-tification of proto-cluster environments, either throughtomographic reconstruction (e.g.Stark et al. (2015b)) orthrough a group of spectra analysis within a protoclusterscale (Cai et al. 2016; Cai et al. 2017).The properties of cosmic voids and proto-clusters arethe result of non-linear structure formation from theGaussian early universe to the time of observations.However, the Wiener Filter only provides an unbiasedminimum-variance estimate in the limit that the under-lying field is Gaussian (Tegmark 1997). While this isa valid approximation for power-spectra analysis of theCMB (Bond et al. 1998) or galaxy surveys (Vogeley &Szalay 1996; Tegmark et al. 1998), when studying thestructure of non-linear features one expects there to beadditional information that is not captured or is other-wise smoothed by the filtering. With next generationfacilities coming online later this decade, Ly α tomogra-phy will be possible across large regions of the sky atcomparatively spatial resolution ( ∼ M = 0 . , Ω Λ = 0 .
69 and H = 70 km s − Mpc − . RECONSTRUCTION METHODS2.1.
Wiener Filter
The Wiener Filter is a special case of the maximumlikelihood method, which analytically reconstructs thefield with maximum likelihood given observed data withknown priors (e.g., the underlying distribution of the
RCA: Optimized Reconstruction with Constraints on Absorption M = C MD · ( C DD + N ) − · d, (1)where C MD and C DD are the map-data and data-datacovariance matrix, N is the noise covariance matrix,which is diagonal assuming noise is uncorrelated, and d is the input data. The shape of covariance matrix can becalculated in different approaches with cosmological pri-ors, and in a widely accepted ad hoc approach, a Gaus-sian random linear filed prior is assumed (Pichon et al.2001). As is implemented in Caucci et al. (2008); Leeet al. (2014a), covariance of two points r and r eitherin map and data is Gaussian, so that C MD and C DD can have the same definition: C MD = C DD = C ( r , r )and C ( r , r ) = σ F exp (cid:34) − (∆ r (cid:107) ) L (cid:107) (cid:35) exp (cid:20) − (∆ r ⊥ ) L ⊥ (cid:21) (2)where ∆ r (cid:107) and ∆ r ⊥ are the distance between r and r along, and transverse, to the line-of-sight, respec-tively, and σ F is the the priori expected variance of 3-DLy α Forest flux fluctuations in a volume of order L ⊥ L (cid:107) ,while L (cid:107) and L ⊥ are correlation lengths along and per-pendicular to the LOSs. L ⊥ is often chose to be onthe order of LOSs mean separation (cid:104) d LOS (cid:105) to avoid fic-titious structures while L (cid:107) depends on the specific sce-nario (e.g., Lee et al. (2014a) set L (cid:107) to FWHM of the as-sumed instrumental resolution while Caucci et al. (2008)set it to be of the order of the Jeans length in order toavoid information loss for small scales along the LOSs).The choice of all parameters we use in this paper are dis-cussed in Section 3.2. We use the Wiener Filter codes dachshund developed by Stark et al. (2015b).In practice, d is a column vector containing observedflux contrast from all lines of sight, δ F = F/ (cid:104) F (cid:105) − C DD + N ) and C MD are two large matrices con-taining correlation information. Although the noise ma-trix N is usually assumed to be diagonal, the covariancematrix C DD is complicated as the signals are correlatedto each other, which makes it tremendously computa-tional intensive to inverse ( C DD + N ) as our surveysbecome larger. There are attempts to reduce computa-tion like dividing the box into small overlapping chunks(Lee et al. 2014a), but there is still significant costs as-sociated with the matrix to inverse operation.Note that Wiener filtering only matches the maximuma posteriori estimator in the case of a Gaussian randomfield whose properties are described solely by the signal http://github.com/caseywstark/dachshund covariance. Even without taking into account the non-linear evolution of the matter density of the universe,this will not describe the flux field at z ∼ . ORCA
As standard Wiener Filter requires the inversion ofa large matrix, we can approach the reconstruction asan optimization problem using an L-BFGS optimizerto avoid intensive computation, which is also used inHorowitz et al. (2019). While Horowitz et al. (2019) re-construct the initial density field using Gunn PetersonApproximation, we are directly reconstructing the fluxfield without such an assumption, allowing generaliza-tion to different cosmological models without residualbiases.An estimate of the underlying flux field, s , can befound by minimizing the loss function L , L = k ( S m ( s ) − s ) + ( R ( s ) − d ) T N − ( R ( s ) − d )+ k (cid:88) clip ( s, , + ∞ ) + k (cid:88) clip ( s, , α ) (3)where S m is the smoothing operator which Gaussiansmooths with kernel size m the output flux field s , d isthe input (obersved) flux on each skewers, R is a skewer-selector function which maps the 3-D field to the ob-served skewers, N is the noise covariance matrix and clip is the function clipping values outside the interval to theinterval edges. These clipping functions have the effectof penalizing extremal and/or non-physical values, and α is an empirical constant between 0 and 1 dependingon the mean flux and smoothing scale (which dependson spectra resolution and sightline spacing).It is useful to note that the first two terms of thisoptimization procedure alone will reconstruct the stan-dard Wiener Filter result for the case of L (cid:107) = L ⊥ = m ;both terms have the same effect of penalizing the result-ing likelihood for small-scale variations below that scale.However, it yields significant improvement over existingWF implementations in computational cost and mem-ory by virtue of performing a local smoothing operationrather than a full non-sparse matrix inversion (for addi-tional general discussion, see Horowitz et al. (2019)).We have adjusted our optimization to perform amultiscale, annealed optimization as implemented inHorowitz et al. (2020) within a L-BFGS framework. Ateach step in the optimization, we smooth the result-ing flux field in steps, progressing from 1.25 h − Mpcsmoothing to 0.5 h − Mpc smoothing, with the step of
Li, Horowitz, and Cai h − Mpc. The smoothing range we use was foundthrough empirical testing to give the best reconstruc-tion (minimal difference between the real map and re-constructed map). It should be noted that if the finalsmoothing scale is too small or zero, there will be manyfake small structures in the reconstructed map, because,in such case, the first term of L is close or equivalent tozero and then the second term of L simply tries to makepixels on skewers equal to input data regardless of noise,and 0.5 h − Mpc smoothing in the last step is sufficientto avoid such problems. This anneal scheme helps us toapproach the global minimum of L , without getting theoptimizer stuck in local minima. With sufficient tests,we get evident smaller mean squared error (MSE) in fluxcompared to that in optimization with a single smooth-ing scale 0.5 h − Mpc even though both of them reachnumerical convergence. MOCK SURVEYS3.1.
Mock datasets
We use a hydrodynamical simulation with Nyx code(Almgren et al. 2013) for our mock survey in this paper,which has a 100 h − Mpc box size with particle resolu-tion 4096 . The simulation uses a flat ΛCDM cosmologywith Ω m = 0 . , Ω b = 0 . , h = 0 . , n s = 0 . σ = 0 .
8. We downsample the simulation to particleresolution 200 as the true field in our mock survey.Following Krolewski et al. (2018), We randomly se-lect skewers with a mean sightline separation (cid:104) d LOS (cid:105) =2 . h − Mpc, comparable to CLAMATO (cid:104) d LOS (cid:105) =2 . h − Mpc, and further, within the predicted range ofthe upcoming Subaru Prime Focus Spectrograph (Sub-aru/PFS) high redshift tomography program. We alsoemulate the predicted observations of Thirty Meter Tele-scope (TMT) assuming (cid:104) d LOS (cid:105) = 1 h − Mpc and samenoise properties. We hereafter use N-PFS (PFS-likemock survey using Nyx simulation) and N-TMT (TMT-like mock survey using Nyx simulation) respectively forthe two different mock surveys in this paper.We apply the procedure provided by Horowitz et al.(2019) and Horowitz et al. (2020) to add pixel noise toeach skewer and model the continuum-fitting error. Thenoise level on each skewer is determined by drawing aS/N ratio from a distribution between a minimum andmaximum S/N. From Stark et al. (2015b), the distribu-tion follows a power-law: dn skewer /d ( S/N ) =
S/N − α ,based on LBG luminosity function (Reddy et al. 2008)and observed distribution in Lee et al. (2014b). We use α = 2.7 which provides the best approximation of S/Ndistribution of CLAMATO and PFS.To account for continuum misclassification error, theflux values within each skewer is offset such that the final observed flux is F obs = F δ c , (4)with δ c being a value drawn from a Gaussian distributionwith mean 0 and width σ = 0 . S/N + 0 . . (5)where the constants are fitted from CLAMATO data.The S/N ratio for N-PFS ranges from 1.4 to 10and that for N-TMT ranges form 2.8 to 10 followingHorowitz et al. (2019).3.2. Flux Reconstruction
We apply ORCA and Wiener Filter to the mock skew-ers and get the reconstructed flux fields. We use σ F =0 . , L ⊥ = (cid:104) d LOS (cid:105) = 2 . h − Mpc , and L (cid:107) = 2h − Mpcin Equation(2) for Wiener Filter reconstruction and k = 5 , k = 0 . , k = 0 .
025 in Equation(3) for theORCA reconstruction. We use the same ORCA pa-rameters for both mock survey and CLAMATO surveydiscussed in Section 4. While Lee et al. (2018) uses σ F = 0 .
05 for WF in CLAMATO data, we adjust σ F to0.082 in our mock survey to match the PDF of CLAM-ATO δ rec F . We apply Gaussian-smoothing to the outputfield with a Gaussian kernel of 2 h − Mpc in the followinganalysis except for the void finding discussed in Section4.1.Figure 1 shows the probability density distribution ofsmoothed flux in true map, Wiener filtered map andORCA reconstructed map. The smoothed true field isquite non-Gaussian, indicating that there is more cos-mological information beyond the two-point correlationfunction (Krolewski et al. 2018). We find that the dis-tribution of ORCA is more consistency with the truedistribution and that of the Wiener filter, which has alarger deviation. It is also notable that part of the fluxvalues from Wiener filter are beyond one which is non-physical and this is a ubiquitous pattern for Wiener filtersince it has no additional constraints. ORCA improvesthe flux PDF on both high and low flux tail, providingmore realistic flux values in the reconstructed field.We show scatter plot of the reconstructed flux againstthe true flux in Figure 2. We find that the relationsbetween Flux recon and Flux true are biased, which is alsofound by Lee et al. (2014a) and Ozbek et al. (2016). Wenote that they plot with flux contrast, but we use fluxinstead to better illustrate points with flux value beyondone. ORCA has an obvious improvement in reconstruct-ing flux at 0.8 to 1.After applying linear regression to the scatteringpoints, we find better slope and Pearson coefficient for
RCA: Optimized Reconstruction with Constraints on Absorption Flux P D F N-PFS TrueN-PFS ORCAN-PFS WF Figure 1.
Frequency distributions of flux in true map,Wiener filtered map and ORCA reconstructed map. Thesmall plots inside the figure show low counts regions in thelog scale.
Flux true F l ux r ec ORCA k = 0.677r = 0.711y=x
Flux true WF k = 0.653r = 0.680y=x Figure 2.
Scatter plot of the reconstructed flux against thetrue flux in Nyx simulation. The red dotted line representsthe Flux recon = Flux true relation while the blue dash line isthe best linear fit of the points with k slope and r Pearsoncoefficient of the fit. These metrics show that ORCA pro-vides a less biased, smaller variance estimate, compared tothe WF solution.
ORCA, which are 0.677 and 0.711 compared to 0.653and 0.680 for Wiener Filter respectively. While Ozbeket al. (2016) linearly corrected the flux according toslope of the regression (subtract flux value by fitting y-intercept and divide them by slope), we do not correctthe bias in the following analysis, because our cosmicweb classification procedure (see Section 3.3, Eq. 6) isthe second order in nature and our results do not changewith a first order correction. In addition, we find thatthe linear correction will cause additional flux values ex-ceeding one as the correction rotates the points in Figure2. It also causes the flux distribution deviating from thetrue distribution in Figure 1, which meets our expecta-tion since ORCA has already optimized the field to getminimal difference from the true field and such a simplecorrection would break the optimized state. x( h Mpc) y ( h M p c ) True NodesFliamentsSheetVoid x( h Mpc) y ( h M p c ) ORCA x( h Mpc)WF
Figure 3.
Cosmic Web Structures from a single slice. Void,sheet, filaments and nodes are marked in dark blue, lightblue, green and yellow respectively. All maps are smoothedwith a 2 h − Mpc Gaussian kernel. N od e s F il a m e n t s S h ee t V o i d Predicted
NodesFilamentsSheetVoid T r u e N od e s F il a m e n t s S h ee t V o i d Predicted0.267 0.533 0.190 0.0100.053 0.487 0.416 0.0440.007 0.184 0.624 0.1850.000 0.049 0.420 0.531WF
Figure 4.
Confusion matrix for cosmic web classificationfor the CLAMATO mock survey. Numbers on the diagonalthe fraction of each true cosmic structure correctly identifiedby the reconstruction, while those off the diagonal reflect thefaction of each structure misidentified.
Cosmic Web Classification
We choose the Tidal Shear Tensor (T-web) method(Forero-Romero et al. 2009) to classify the cosmic struc-tures to test how well our method reconstructs the field.The T-web method uses the deformation tensor, theHessian of the underlying gravitational potential, D ij = ∂ Φ ∂x i ∂x j (6) Li, Horowitz, and Cai cos e N-PFS WFN-PFS ORCAN-TMT WFN-TMT ORCA 0.0 0.2 0.4 0.6 0.8 cos e cos e Figure 5.
Histogram of the cosine of the angle between the reconstructed eigenvectors and true eigenvectors for each eigenvector.Black dashed line would correspond to random orientations of reconstructed eigenvectors, while the furthest right bin correspondsto complete agreement of eigenvectors (i.e. cos θ ∼ h − MpcGaussian kernel.
Table 1.
Cosmic Web RecoveryMock Data Method Pearson Coefficients Volume overlap (%) λ λ λ Node Filament Sheet VoidN-PFS WF 0.473 0.545 0.659 0.268 0.487 0.624 0.531ORCA 0.522 0.613 0.714 0.318 0.557 0.628 0.582N-TMT WF 0.856 0.902 0.948 0.649 0.781 0.827 0.804ORCA 0.882 0.925 0.962 0.720 0.811 0.851 0.843 which is numerically practical to compute in Fourierspace, ˜ D ij = k i k j k δ k (7)where δ k is the density field and we obtain D ij by inverseFourier-transforming ˜ D ij .The eigenvectors of the deformation tensor relate tothe principle curvature axes of the density field at eachpoint, corresponding in the Zel’dovich approximationwith the principle inflow/outflow directions. The corre-sponding eigenvalues determine if the net flow is inwardor outward. Points with three eigenvalues above thresh-old value λ th are classified as nodes, two values above λ th are filaments, one value above λ th are sheets, andzero values above λ th are voids. While the density fieldis related to the flux field with high flux indicating lowdensity and low flux indicating high density, we use fluxfield in the deformation tensor computation as describein Lee & White (2016). Thus, the relation between cos-mic structures and eigenvalues above λ th is reversed (i.e.three eigenvalues above λ th are voids). Following Horowitz et al. (2020), we define our thresh-old value λ th for each reconstruction such that the voidsoccupy 22% of the total volume. In the true flux field, wefind that [22.0, 49.9, 25.4, 2.7]% of the volume is occu-pied by voids, sheets, filaments, and nodes, respectively.In ORCA and WF maps of N-PFS mock survey, thosenumbers are [22.0, 49.3, 25.8, 2.9]% and [22.0, 51.5, 24.1,2.4]%, respectively. The volume fraction of all the fourstructures in ORCA map are closer to that in true mapcompared to WF. A random slice of the classificationis shown in figure 3, where we can visually find a no-table improvement of ORCA for a better recovery ofvoids, e.g, ORCA recovers the big void at the lower-leftof the slice while WF breaks it down into small voids.The accuracy of the cosmic web classification is furtherquantitatively measured by confusion matrix in Figure4. The volume overlap of ORCA map with true fieldfor node, filament, sheet, and void, is [31.8, 55.7, 62.8,58.2]%, compared with [26.8, 48.7, 62.4, 53.1]% of WFmap. ORCA outperforms WF in cosmic web classifica-tion through node, filament and void, and it’s compara- RCA: Optimized Reconstruction with Constraints on Absorption ≈ e , ˆ e , and ˆe of the pseudo-deformation tensor D ij . Figure 5 illustrates the distribution of cos θ of theangle between the reconstructed eigenvectors and trueeigenvectors. In Figure 3, we could find that ORCA hasan improved performance in recovering all three eigen-vectors. We further increase the density of sightlines sothat (cid:104) d LOS (cid:105) = 2 h − Mpc and then we find WF has sim-ilar quality in cosmic web alignments. It indicates thatreconstructing by ORCA has equivalent effects of in-creasing sightlines. And for N-TMT survey, ORCA stillnotably outperforms WF. We also quantify the agree-ment between reconstructed field and true field in termsof Pearson coefficients of reconstructed and true eigen-values in Table 1. It shows that we get a stronger cor-relation between three reconstructed and true eigenval-ues [ λ , λ , λ ] either for N-PFS or N-TMT surveys withORCA, which ranges from [0.522, 0.613, 0.714] to [0.882,0.925, 0.962] in contrast with that ranging from [0.473,0.545, 0.659] to [0.856, 0.902, 0.948] with WF. APPLICATION TO CLAMATO DATAWe apply our technique to the first data release ofthe COSMOS Lyman Alpha Mapping And TomographicObservations (CLAMATO) survey (Lee et al. 2018).This data includes reduced Ly α forest signatures from240 galaxies and quasars with redshifts ranging 2 . 00, allowing reconstruction at 2 . < z < . ∼ . 157 square degree. Standard Wiener filter re-constructions of this data have been used to detect alarge number of cosmic voids and proto-cluster regions.4.1. Void and cluster finding We adopted the void finding procedure presented inStark et al. (2015a) and compared void catalogs in map The data is available here: https://doi.org/10.5281/zenodo.1292459. We use map 2017 v4.bin and map 2017 v4 sm2.0.binas Wiener Filter map and pixel data v4.bin as input to ORCA. reconstructed by Wiener Filter and that by ORCA. Inthe flux field, we begin by finding all points with δ F above a threshold (SO threshold). Spheres then growcentered on all those points until the average δ F insidethe sphere reaching a second threshold (SO average).Due to the pixel noise and continuum error, the PDFof WF δ rec F is broadened, especially in high flux regionwhere we find voids, and we use a different thresholdsfor true flux field δ F and WF δ rec F . For both true andWF field in mock survey, we use those thresholds pro-vided by Krolewski et al. (2018). The process of choos-ing thresholds can be summarized as follows (also seeTable 2):1. Find voids in the true redshift-space density field ofthe Nyx simulation. • SO thresh = 0 . ρ and SO avg = 0 . ρ 2. Find voids in the true flux field of the Nyx simula-tion. • Choose SO thresh and SO avg to best match voidsfound in 1. • SO thresh = 0 . . • Choose SO thresh and SO avg to best matchvoids found in mock survey with voids found in1. • SO thresh = 0 . . 175 for bothmock survey and CLAMATO.(b) ORCA • Use the same thresholds in 2. • SO thresh = 0 . . 152 for bothmock survey and CLAMATO.In the true field, the voids are firstly identified inredshift-space density field using SO threshold= 0 . ρ and SO average= 0 . ρ . The thresholds for underlyingflux field δ F and WF reconstructed field δ rec F are chosento be [0.192,0.152] and [0.220,0.175] respectively, match-ing the void fraction in redshift-space density field. Nev-ertheless, we did not use the same threshold for ORCA δ rec F as that of WF, because ORCA has constraints onabsorption and we find ORCA working well around thethresholds of the true field. Thus, we simply use thesame thresholds for ORCA δ rec F .Qualitatively, we find a good agreement between thetwo CLAMATO maps . While Krolewski et al. (2018) For Wiener Filter map, we use the void catalog here: https://doi.org/10.5281/zenodo.1295839. Li, Horowitz, and Cai O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F Void Center O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F O RC A W F O RC A Line-of-sight ( h Mpc) W F Redshift 0.40.20.00.20.4 F D ec ( h M p c ) Figure 6. Comparison of voids found in CLAMATO map reconstructed by Wiener Filter and ORCA. Each strips are the stackof 4 slices with 2 h − Mpc in thickness through the RA direction. The black circles represent the voids intersecting the slices andthe + marks all voids centered at one of four slices of a strip. The strips connected by staples on the right are reconstructionswith ORCA and Wiener Filter at the same RA marked by ORCA and WF on the left. RCA: Optimized Reconstruction with Constraints on Absorption h Mpc)10 n ( r )( h M p c ) CLAMATO WFCLAMATO ORCAN-PFS TrueN-PFS WFN-PFS ORCA Figure 7. Void radius function for our mock catalog (NyX)and the CLAMATO data. Due to the geometry of theCLAMATO volume, there are few large voids ( r > . h − Mpc) identified. identified 355 r > h − Mpc voids, including 48 higher-quality r ≥ r > h − Mpc voids and55 r ≥ h − Mpc voids. We find that there is 70.52% ofvolume of voids found in the Wiener Filtered map thatare also found in ORCA reconstructed map. This highoverlap fraction indicates a good agreement of the twomethods. We plot the comparison of voids found in twomaps in Figure 6. The figure shows the stack of everyfour slices with 2 h − Mpc in the RA or DEC direction.We also compare the void radius function in CLAMATOto that in the N-PFS mock survey of the two methods inFigure 7. Due to edge effects, voids are more likely to befound near survey boundaries, so we exclude all the voidswith distance from center to the boundary smaller thanthe void radius. Following Krolewski et al. (2018), wecompute the void radius function with weights to eachvoid by the effective volume over which it could havebeen observed (e.g., for the geometry of CLAMATO,the effective volume is (30 − r )(24 − r )438 Mpc h − with r the void radius). We find a good agreement ofvoid radius function in CLAMATO and mock survey,either with ORCA or WF.To test the improvement of ORCA in void recovery,we define the void volume overlap completeness as thefraction of voids found in the true flux map from ourmock catalog that are also found in reconstructed map,and the volume overlap purity as the fraction of voidsfound in reconstructed map that are also found in thetrue flux map, i.e. Purity ≡ V true ∩ V rec V rec , (8)Completeness ≡ V true ∩ V rec V true . (9) r void ( h Mpc) V o l u m e ov e r l a p fr ac ti on ORCA Volume overlap completenessWF Volume overlap completenessORCA volume overlap purityWF volume overlap purity Figure 8. Purity (solid line) and completeness (dashed line)of volume overlap fraction of N-PFS mock survey, as definedin Eq, 8 and Eq. 9. Each bin is 1 h − Mpc. In the N-PFS mock survey, we find that ORCA recon-structed map has a higher volume overlap completenesscomputed using all voids, which is 34.50% compared to30.72% in the WF reconstructed map. In Figure 8, weplot the completeness and purity of the volume overlapfraction compared between voids in mock survey recon-structed by ORCA and WF, and the true flux field in theNyx simulation as a function of void radius. While forsmall voids ( r < r > r ≥ h − Mpc, nor-malizing each void to its void radius and stacking inunits of the void radius r/r void . We could see a goodagreement between void profiles in mock survey andCLAMATO data. The void profiles for r > r void inmock surveys trace well the true underlying flux field δ F . However, the large deviation between void profilesin WF δ rec F and δ F means that void profiles in the recon-struction can hardly trace void profiles in matter, whilethe void profiles in ORCA δ rec F match better with theprofile of Nyx and it provides us a possibility to tracethe matter inside the void.0 Li, Horowitz, and Cai r / r void F CLAMATO WFCLAMATO ORCAN-PFS TrueN-PFS ORCAN-PFS WF r / r cluster + Figure 9. Solid lines represent the mean of void (left) and cluster (right) profiles in Wiener reconstructed map from CLAMATOdata (red), ORCA reconstructed map from CLAMATO data (blue), Nyx map (black), ORCA reconstructed map from Nyxmock observations (magenta) and Wiener reconstructed map from Nyx mock observations (cyan) and errorbars show standarddeviation σ of the stacks. Table 2. Volume fraction for different thresholds of void and cluster insimulated and CLAMATO catalogs.Type Data Field SO thresh SO avg Vol. frac.Void N-PFS δ F δ rec F δ rec F δ rec F δ rec F δ F -0.304 -0.271 2.71%ORCA δ rec F -0.304 -0.271 2.11%WF δ rec F -0.304 -0.271 1.72%CLAMATO ORCA δ rec F -0.304 -0.271 1.96%WF δ rec F -0.304 -0.271 1.20% Note —Comparison of volume fraction of voids and clusters found in Nyx mock survey (100 h − Mpc box) and CLAMATOsurvey with two methods. All the maps for finding clusters are smoothed with a 2 h − Mpc Gaussian kernel, while that forfinding voids are without additional smoothing. We also apply SO method to find clusters using a sim-ilar procedure. We smooth both true and reconstructedmaps with a 2 h − Mpc Gaussian kernel following (Starket al. 2015b). The volume fraction of nodes in true nyxmap classified by T-web method (see Section 3.3) is 2.7%and we choose thresholds for SO cluster finder to makevolume fraction of clusters in the nyx map match thisfraction. We use the same thresholds for both true andreconstructed maps. It can be seen at right panel of Fig-ure 9 that ORCA tends to recover overdensities better,while Wiener Filter underestimates the density insidethe clusters. CONCLUSION In this work we have introduced a new tomographicflux reconstruction technique to use on Ly α Forest ob-servations. Testing our approach on mock catalogs fromhydrodynamical simulations, we have shown improvedcosmic web reconstruction vs. standard Wiener filteringapproaches. This improvement can be seen both in clas-sification accuracy as well as reconstruction of profilesand number statistics of voids and (proto)clusters. Inaddition to testing on mock catalogs, we have also ap-plied our technique to data from the CLAMATO Sur-vey and have found good agreement with void profilesfrom simulations. In our simulations, we have also foundthat our method can reconstruct void profiles more con- RCA: Optimized Reconstruction with Constraints on Absorption α Tomography is expected to play a majorrole in upcoming spectroscopic surveys, such as SubaruPrime Focus Spectrograph (PFS) (Takada et al. 2014),it is important to gain the maximum information fromthe limited time available. We have found that ORCAprovides large scale structure constructions comparableto WF with 30-40% more sightlines, depending on met-ric used. With the increasing of survey volume andsightline density in upcoming surveys the computationalcosts of the WF reconstructions becomes more apparent.We find ORCA reconstructs 10-100 times faster thanWiener Filter ( dachshund ) when ORCA runs on oneNVIDIA Tesla V100S GPU 32 GB and dachshund runson one thread of Intel(R) Xeon(R) Gold 6226R CPU @2.90GHz.Going forward, this technique provides a complemen-tary tool with forward model density reconstruct tech-niques (e.g. Horowitz et al. (2020)), which rely on strongassumptions about IGM physics, and provides a usefultool to cross-correlate with galaxy properties. A regime that is of particular interest for galaxy evolution studiesis the center of proto-cluster regions, where one expectssignificant deviations from FGPA. We hope to furtherexplore ORCA and similar extensions to WF in thisregime in future works.ACKNOWLEDGMENTSZL and ZC are supported by the National Key R&DProgram of China (grant No. 2018YFA0404503). BH issupported by the AI Accelerator program of the SchmidtFutures Foundation. We thank Khee-Gan Lee and AlexKrolewski for helpful discussions. We thank BenjaminZhang for help on technical issues relating to Dachshundimplementation.This research used resources of the National EnergyResearch Scientific Computing Center, a DOE Officeof Science User Facility supported by the Office of Sci-ence of the U.S. Department of Energy under ContractNo. DEC02-05CH11231. We also acknowledge the com-puting resources from Department of Astronomy at Ts-inghua University.APPENDIX A. ERROR ANALYSIS ON ORCA RECONSTRUCTIONIn order to propagate uncertainties to whatever analysis the tomographic map will be used for we need to understandthe error properties of our map. Since we have a complete likelihood, it is easy to test the relative likelihood of agiven flux realization. The errors will be correlated, as the signal covariance term forces the reconstructed map to besmooth. If we are interested in a pixel by pixel flux error we can calculate this via a response formalism; i.e. by varyingeach pixel we can study the change in the resulting likelihood. In order to know how the loss function responds, weseparate the L in Equation (2.2) into 4 terms ( L , L , L , L ): L = k ( S m ( s ) − s ) L = ( R ( s ) − d ) T N − ( R ( s ) − d ) L = k (cid:88) clip ( s, , + ∞ ) L = k (cid:88) clip ( s, , α ) (A1)We add a small increment to flux at every pixel and see how value of each term change to penalize the optimization.Figure 10 shows the flux field and the change of loss function L and its four components at a slice perpendicular to theLOSs, added 0.1 to flux at each pixel. We can see in (b) and (c) that the pixels on skewers are the most importantones with the biggest impact on the first two terms compared to other pixels off the skewers, which is expected instandard Wiener Filter as it uses information around skewers. We see the effect of ORCA in (e) and (f ) where we addtwo constraints to the optimization. In (e) we see that L works when flux values exceeding one and the yellow andlight blue regions can correspond to high flux region in (a) . With such a constraint, we could avoid the non-physicalvalues in our final optimized map. In (f ) we see that L works at low flux region and it penalizes the optimizationwhen we lose some low flux values, which helps us better recover overdensity. It can also prove that the optimizationhas reached the minimum of L where any increment to flux will increase the total loss function, shown in (d) whereall ∆ L are positive.2 Li, Horowitz, and Cai (a) Flux (b) (c) (e) (f) (d) Figure 10. (a) shows the original flux field from the ORCA reconstruction. (d) shows the change of total loss function ∆ L inresponse to an increment to the flux, while (b) , (c) , (e) and (f ) show the change of each separated term. Black circles indicateposition of skewers. REFERENCES Almgren, A. S., Bell, J. B., Lijewski, M. J., Luki´c, Z., &Andel, E. V. 2013, The Astrophysical Journal, 765, 39,doi: 10.1088/0004-637x/765/1/39Bond, J. R., Jaffe, A. H., & Knox, L. 1998, PhRvD, 57,2117, doi: 10.1103/PhysRevD.57.2117Bos, E. P., van de Weygaert, R., Dolag, K., & Pettorino, V.2012, Monthly Notices of the Royal AstronomicalSociety, 426, 440Cai, Z., Fan, X., Peirani, S., et al. 2016, ApJ, 833, 135,doi: 10.3847/1538-4357/833/2/135Cai, Z., Fan, X., Bian, F., et al. 2017, The AstrophysicalJournal, 839, 131, doi: 10.3847/1538-4357/aa6a1aCaucci, S., Colombi, S., Pichon, C., et al. 2008, MNRAS,386, 211, doi: 10.1111/j.1365-2966.2008.13016.xChiang, Y.-K., Overzier, R., & Gebhardt, K. 2014, TheAstrophysical Journal, 782, L3,doi: 10.1088/2041-8205/782/1/l3Contarini, S., Marulli, F., Moscardini, L., et al. 2020, arXive-prints, arXiv:2009.03309.https://arxiv.org/abs/2009.03309Croft, R. A. C., Weinberg, D. H., Katz, N., & Hernquist, L.1998, ApJ, 495, 44, doi: 10.1086/305289 Eisenstein, D. J., Weinberg, D. H., Agol, E., et al. 2011,The Astronomical Journal, 142, 72,doi: 10.1088/0004-6256/142/3/72Fang, Y., Hamaus, N., Jain, B., et al. 2019, MNRAS, 490,3573, doi: 10.1093/mnras/stz2805Feng, Y., Chu, M.-Y., Seljak, U., & McDonald, P. 2019,FastPM: Scaling N-body Particle Mesh solver.http://ascl.net/1905.010Forero-Romero, J, E., Hoffman, Y., Gottl¨ober, S. J.,Klypin, A., & Yepes, G. 2009, MNRAS, 396, 127,doi: 10.1111/j.1365-2966.2009.14885.xGunn, J. E., & Peterson, B. A. 1965, ApJ, 142, 1633,doi: 10.1086/148444Horowitz, B., Lee, K.-G., White, M., Krolewski, A., & Ata,M. 2019, The Astrophysical Journal, 887, 61,doi: 10.3847/1538-4357/ab4d4cHorowitz, B., Seljak, U., & Aslanyan, G. 2019, JCAP, 2019,035, doi: 10.1088/1475-7516/2019/10/035 RCA: Optimized Reconstruction with Constraints on Absorption13