Early Planet Formation in Embedded Disks (eDisk). I. Overview of the Program and First Results

We present an overview of the Large Program, “ Early Planet Formation in Embedded Disks ( eDisk ) , ” conducted with the Atacama Large Millimeter / submillimeter Array ( ALMA ) . The ubiquitous detections of substructures, particularly rings and gaps, in protoplanetary disks around T Tauri stars raise the possibility that at least some planet formation may have already started during the embedded stages of star formation. In order to address exactly how and when planet formation is initiated, the program focuses on searching for substructures in disks around 12 Class 0 and 7 Class I protostars in nearby ( < 200 pc ) star-forming regions through 1.3 mm continuum observations at a resolution of ∼ 7 au ( 0 04 ) . The initial results show that the continuum emission, mostly arising from dust disks around the sample protostars, has relatively few distinctive substructures, such as rings and spirals, in marked contrast to Class II disks. The dramatic difference may suggest that substructures quickly develop in disks when the systems evolve from protostars to Class II sources, or alternatively that high optical depth of the continuum emission could obscure internal structures. Kinematic information obtained through CO isotopologue lines and other lines reveals the presence of Keplerian disks around protostars, providing us with crucial physical parameters, in particular, the dynamical mass of the central protostars. We describe the background of the eDisk program, the sample selection and their ALMA observations, and the data reduction, and we also highlight representative ﬁ rst-look results.


INTRODUCTION
Protoplanetary disks are ubiquitously found around young stellar objects (YSOs), particularly around Class II YSOs (e.g., Mann et al. 2014;Ansdell et al. 2016;Pascucci et al. 2016;Cieza et al. 2019;Cazzoletti et al. 2019;van Terwisga et al. 2019;Grant et al. 2021).Recent direct imaging of protoplanets in the disk around the T Tauri star PDS 70 (Keppler et al. 2018;Müller et al. 2018;Haffert et al. 2019;Isella et al. 2019;Benisty et al. 2021;Casassus & Cárcamo 2022) lends strong support to the notion that these disks are the sites of planet formation.The detection of rings and gaps in the HL Tau disk (ALMA Partnership et al. 2015;Yen et al. 2016Yen et al. , 2019) ) revolutionized the field of planet formation and made the use of disk substructures as signposts for young planets a reality (e.g., Dong et al. 2015;Zhang et al. 2018).Although it is still debated whether such substructures are in fact always created by planets (e.g., Zhang et al. 2015;Flock et al. 2015), they are important in their own right.For example, regions of density enhancement, such as rings (either complete or partial, e.g., van der Marel et al. 2013) and spirals (e.g., Pérez et al. 2016), have a local pressure maximum that is conducive to the trapping of relatively large grains, which, in turn, facilitates the formation of planetesimals and ultimately planets (Youdin & Goodman 2005;Johansen et al. 2007;Johansen et al. 2009).
Given their fundamental importance, there have been a number of ALMA programs aiming at detecting and characterizing disk substructures (e.g., Akiyama et al. 2016;Andrews et al. 2016;Long et al. 2018).The Cycle 4 Large Program, Disk Substructures at High Angular Resolution Project (DSHARP; see Andrews et al. 2018), surveyed 20 nearby Class II disks in dust emission at 1.25 mm at an angular resolution of 0. 04 and showed that substructures are ubiquitous, mostly in the form of concentric rings and gaps.Similar high-resolution observations of Class II disks were also performed for 7 sources in Ophiuchus at resolutions of 3-5 au as a part of the Ophiuchus DIsc Survey Employing ALMA (ODISEA) project (Cieza et al. 2021), finding ring/gap structures in all of them.These surveys strongly suggest that the important disk substructures are well developed in the Class II stage of star formation.This raises the obvious question: when do substructures start to develop in disks?
If planets are indeed responsible for producing the rings and gaps in Class II disks, planet formation could have started earlier, i.e., during the Class I or even Class 0 stages.Here, we consider that the dust growth forming millimeter-or centimeter-sized dust grains at least should occur prior to the Class II stage.The general lack of sufficient solids in Class II disks to form giant planets also points to an earlier epoch of planet formation (Manara et al. 2018;Tychoniec et al. 2018;Ward-Duong et al. 2018).Theoretical studies have suggested that Class 0/I disks are more plausible sites for planetesimal formation than Class II disks because of their large amount of gas and dust (Tsukamoto et al. 2017).A systematic survey focused on protostellar disks can directly test this possibility by searching for structural signatures of early planet formation.For example, if such early planets are formed through core accretion (e.g., Helled & Bodenheimer 2014), dust traps that facilitate the formation of planetesimals and planetary cores should be seen.Alternatively, if planets formed through gravitational instabilities (e.g., Boss 1997), prominent spiral arms are expected (e.g., Tomida et al. 2017) .
In order to investigate the possibility of early planet formation, it is obviously necessary for us to investigate disks around Class 0/I protostars.However, observations of Class 0/I disks are more difficult than those of Class II disks because they are still deeply embedded in envelopes.In particular, the inner parts of envelopes within a few 1000 au often show flattened structures and velocity gradients and it is not straightforward to distinguish embedded disks from these structures.Careful analysis, including their kinematics, is required to identify rotationally-supported embedded disks.
In this paper, we present the ALMA Large Program "Early Planet Formation in Embedded Disks (eDisk)", in which 19 Class 0/I embedded protostars in nearby star-forming regions are systematically surveyed at resolutions of ∼ 5 au (0. 04) for the first time.The main scientific objective of the eDisk program is to investigate whether or not substructures exist in disks around embedded, Class 0/I, protostars that might be indicative of early planet formation.This, of course, also involves the important task of securely identifying the nature of the embedded disks, e.g., by searching for the kinematical signatures of rotationally supported (Keplerian) disks through molecular line observations.In turn, this will provide fundamental (global) parameters of the young stellar objects, such as the disk radii and dynamical stellar masses based on the Keplerian motions.The first homogeneous samples of embedded disks obtained in the eDisk program can be directly compared with the corresponding samples of Class II disks, enabling us to shed light on the disk evolution from the Class 0/I to Class II phases.
The paper is organized as follows: an overview of previous studies of embedded disks around protostars is provided in Sect. 2. Section 3 describes the selection of the targeted sample of sources and the overall observing strategy, while the data reduction procedures are presented in Section 4. In Sect.5, the first results of the program are presented with emphasis on observations of the 1.3 mm continuum of the sample of sources as well as an example of the C 18 O 2-1 emission toward one, R CrA IRS7B.These parts of the continuum and line observations are key to achieving the main objectives mentioned above.In Sect.6, we briefly discuss implications of the first results from the survey, and brief descriptions of the data release from the eDisk program are provided in Sect.7 Finally, a summary of the paper is provided in Sect.8.

PREVIOUS STUDIES OF EMBEDDED DISKS AROUND PROTOSTARS
Even though disks around embedded protostars are expected theoretically to be a natural consequence of the collapse of a rotating dense molecular core (e.g., Terebey et al. 1984) their firm observational detections are challenging as such disks are still deeply embedded in the envelopes surrounding them.Envelopes often show flattened structures with typical sizes of a few thousand au, which is roughly 10-50 times larger than typical disks, and also have rotating as well as infalling motions (e.g., Ohashi et al. 1997a;Momose et al. 1998).High resolution and sensitivity are therefore required to carefully disentangle emission from disks and envelopes.
The confusion between the envelope and disk emission also makes it difficult to unambiguously address the disk sizes from the continuum emission by them-selves and thus whether significant disks may develop early in protostellar evolution.However, some attempts suggest that very large (> 100 au) disks are rare around deeply embedded protostars (Maury et al. 2019).

Kinematic identification of embedded disks
One approach to distinguish embedded disks and their envelopes is to utilize their kinematic difference; embedded disks are expected to have Keplerian rotation similarly to Class II disks, while envelopes are expected to have rotation conserving their specific angular momenta.As a result, rotational velocities of embedded disks are proportional to r −0.5 , while those of envelopes are proportional to r −1 , where r is the cylindrical radius from the rotational axis (e.g., Ohashi et al. 1997b;Simon et al. 2000;Belloche 2013) This approach was taken to distinguish the embedded disk and its surrounding envelope around the protostar HH111 (Lee 2011) and also those around 6 other protostars (Yen et al. 2013) in C 18 O 2-1 with the Submillimeter Array (SMA) at angular resolution of ∼1 -4 .In some cases (L1527 IRS and L1448-mm), only rotation proportional to r −1 was found because of insufficient angular resolution.However, L1527 IRS was also observed with the Combined Array for Millimeter-wave Astronomy (CARMA) at an angular resolution of ∼ 1 in 13 CO 2-1, displaying Keplerian rotation instead of that proportional to r −1 (Tobin et al. 2012).
Further observations following up these pioneering works were enabled by ALMA with its high angular resolution and sensitivity.L1527 IRS was observed again as a benchmark at sub-arcsecond resolutions in C 18 O 2-1, clearly demonstrating that its infalling envelope transforms into a Keperian disk at a radius of ∼70 au (Ohashi et al. 2014;Aso et al. 2017).Similar findings have been found in other low-mass Class 0/I protostars (e.g., Yen et al. 2014Yen et al. , 2017;;Aso et al. 2015;Tobin et al. 2016;Harsono et al. 2018;Sai et al. 2020;Maret et al. 2020).C 18 O has often been adopted to investigate kinematics of envelopes and disks, as well as other more rare isotopes, such as C 17 O or C 34 S (e.g., Artur de la Villarmois et al. 2019).Important products from detected Keplerian motions are dynamical masses of protostars and mass accretion rates.Notably, direct measurements of masses for optically invisible protostars are usually difficult because their stellar photospheres are not observable.These physical parameters could shed light on the formation and evolution of protostars (e.g., Yen et al. 2017;Artur de la Villarmois et al. 2019) based on more systematic observations of embedded disks.

Chemical identification of embedded disks
The presence of disks around deeply embedded protostars may also affect their chemical structures and consequently, the molecular line emission observed toward them.For example, such embedded disks have been suggested to be associated with distinct features such as slow accretion shocks that in turn may reveal the presence of the disks (e.g., Sakai et al. 2014Sakai et al. , 2016;;Oya et al. 2016;Okoda et al. 2018).Also, grain growth to millimeter-or centimeter-sized grains in disks may cause the continuum emission to become optically thick affecting the line emission signatures observed from the embedded disks (e.g., Harsono et al. 2018).
We should, however, note that consistent chemical indicators across larger samples of disks are still lacking.For example, it has been suggested that SO is enhanced at the disk-envelope surface due to slow accretion shocks in L1527 IRS (Sakai et al. 2014), while in L483 and Elias 29 SO appears to trace the Keplerian rotating disks themselves (Oya et al. 2017;Oya et al. 2019).Oya et al. (2017) also suggested that the emission from CS and complex organic molecules are arising from the disk around L483, while Jacobsen et al. (2019) demonstrated that those species show velocity structures that can be explained better by infalling and rotating material rather than a Keplerian disk.These observations indicate that the molecules tracing the disks may vary among objects, depending on their physical properties (e.g.stellar luminosities) and chemical history in prestellar stage (e.g., Oya et al. 2019).

Substructures of embedded disks
Although embedded disks have been successfully observed around a number of protostars in the past decade, only a few have shown detectable substructures suggestive of ongoing planet formation to date because higher angular resolutions of 0. 1 or less are required to look for substructures in embedded disks.Sheehan & Eisner (2017) have found a central cavity and ring toward the Class I protostar WL 17, in 3 mm continuum at an angular resolution of 0. 06.Sheehan & Eisner (2018) have also found concentric rings/gaps within the disk of the Class I protostar GY 91, also known as ISO-Oph 54 (see Cieza et al. 2021) in the 0.9 mm continuum emission at an angular resolution of 0. 13.More recently, Sheehan et al. (2020) presented seven protostars including two Class 0 sources from the Orion Very Large Array/ALMA Nascent Disk and Multiplicity (VANDAM) survey that exhibit ring/gap substructures in their disks, although detailed modeling showed an offset in the inner disk for a few sources, which may suggest that the substructures could be due to close binary formation.In addition, Segura-Cox et al. (2020) have found four annular struc-tures in the Class I protostar Oph IRS 63 in 1.3 mm continuum at an angular resolution of ∼0.04 although their contrast is relatively low.
In addition to ring structures, spiral structures have also been observed in embedded disks.Lee et al. (2020) have found a pair of symmetric spiral structures in the disk around the Class I protostar HH111 VLA1 at a resolution of 0. 04 (16 au).The disk is relatively massive (0.33-0.5 M ), with Toomre's Q parameter near unity in the outer part of the disk, where the spiral structures are detected, supporting the idea that mass accretion from the envelope is driving the outer disk to be gravitationally unstable.Obviously, it is definitely crucial to observe more embedded sources at angular resolutions higher than 0. 1 to search for substructures, particularly around Class 0 protostars, which show less substructures in their disks to date.

DESCRIPTION OF THE EDISK PROGRAM
The foundation of eDisk is the ALMA large program (2019.1.00261.L: PI N. Ohashi) with supplemental data obtained through ALMA DDT observations (2019.A.00034.S: PI J. Tobin).We also used some ALMA archival data.In this section we present the sample ( §3.1) and the details of the observations ( §3.2).

Sample
For the eDisk program we focus on sources located in nearby (d < 200 pc) star forming regions that are relatively bright (L bol > 0.1 L ).The first criterion is required to spatially resolve observed disks at a given angular resolution (∼0.04; see below).The second criterion is required to detect disks around sample protostars within a reasonable observation time.We checked how many such protostars there are based on available catalogs of young stellar objects (YSOs): the largest catalogs of YSOs in nearby star-forming regions were provided through large survey observations at mid-and far-infrared wavelengths using the Spitzer Space Telescope.Rebull et al. (2010) counted the number of Class 0/I protostars in Taurus based on the near-to mid-IR slopes, α obtained in the Spitzer observations, finding that there are 21 Class 0/I protostars with L bol > 0.1 L 1 , excluding potential Taurus members newly identified.Dunham et al. (2015) also counted the number of Class 0/I protostars based on Spitzer observations and found 37 Class 0/I protostars in the regions of Chamaeleon, 1 Since T bol and L bol were not provided for these sources in Rebull et al. (2010), we calculated them based on 2MASS, Spitzer IRAC and MIPS, Herschel data, except for one source, for which IR photometry data from Kryukova et al. (2012) was used instead of Spitzer data.
Corona Australis, Lupus, Musca, and Ophiuchus.Note that two of the 37 sources show unreasonably high T bol (> 3000 K) for embedded protostars.Based on these previous works, there would be ∼60 protostars within a distance of 200 pc and with luminosities more than 0.1 L .Note that protostars in isolated star-forming regions are not included here.Within a reasonable observation time, it is impossible for us to observe all these protostars from a data volume point of view as well as the restrictions of the maximum time allocated to large program in a given LST range.In the Cycle 7 semester 645 hours were available in total for Large Programs, of which we aimed for a required observation time for the eDisk program of ∼100 hours.According to our previous ALMA observations of embedded disks around protostars and also our simulations of observations, we set a target sensitivity for molecular line observations at ∼2.7 mJy beam −1 with a velocity resolution of ∼0.17 km s −1 .This sensitivity was able to be achieved with an on-source integration time of ∼2.4 hours, requiring ∼6 hours observation time, including overhead.Note that this observation time can provide us with a sensitivity of ∼13 µJy beam −1 for the continuum emission, which is sufficiently high to detect continuum emission around protostars.With these estimations of the required observation time, we found that ∼17 sources could be observed within ∼100 hours of observation time.
When 17 sources were selected we tried to avoid any bias as much as possible; we selected 17 "representative" protostars, which have been relatively well studied before, are located in differnet star-forming regions (i.e., Chamaeleon, Corona Australis, Lupus, Ophiuchus, Taurus, and also isolated Bok globules), and present a wide distribution of the bolometric luminosities (L bol ) and temperatures (T bol ).Selecting sources spreading out over various star-forming regions can also help to carry out our observations more efficiently within a limited period of the observations.In addition to these 17 sources, two sources, B335 and TMC-1A, for which necessary data has been already obtained with ALMA before, were added to the sample.The 19 sample protostellar systems selected for the eDisk program are listed in Table 1.
All of the selected sources have complete spectral energy distributions (SEDs) measured from the nearinfrared to submillimeter that are compiled in Appendix A. Based on those SEDs we re-derived their bolometric temperatures (T bol ) and luminosities (L bol ) also utilising the most recent distances from Gaia measurements listed in Table 1.Note that none of the 19 targets were detected with Gaia directly due to their embedded natures, but distances to nearby sources, presumably related to the same star-forming complexes, were typically estimated based on Gaia observations and provided the distance estimates quoted in literature.The final sample comprises 12 Class 0 and 7 Class I protostars and represents various evolutionary stages of protostars with a roughly uniform distribution of L bol and T bol (Figure 1a).For comparison purposes, Figure 1b shows the same T bol -L bol diagram, but for the 56 protostars with d <200 pc and L bol > 0.1 L discussed above, except for two sources with unreasonably high T bol , indicating that the final 19 protostars reasonably represent the nearby and relatively bright protostars.

Observations
Observations of the eDisk program using the ALMA 12 m array in band 6 were conducted from April 2021 to July 2022.We aimed to achieve high angular resolutions of ∼0.04, corresponding to ∼6 au at a distance of 140 pc.Such high angular resolution is necessary for us to investigate substructures of embedded disks and make direct comparisons with previous maps of Class II disks obtained at similar angular resolutions, such as those obtained by DSHARP or ODISEA.
The correlator was configured to observe both 1.3 mm continuum emission as well as CO (J = 2 − 1) isotopologues and other lines providing additional information about the physical and chemical structures of the gas on the scales of the circumstellar disks surrounding the embedded protostars.Baseband 1 of the correlator was configured in the lower sideband with 4 spectral windows (spw) to include our primary molecular lines, C 18 O (2-1), 13 CO (2-1), and SO (6 5 −5 4 ), at a velocity resolution of ∼0.17 km s −1 .The 12 CO (2-1) emission was included in baseband 4 in the opposite sideband to study protostellar outflows.The bandwidth of baseband 4 was set to be 937.5 MHz to secure the continuum bandwidth, and as a result, the velocity resolution of the 12 CO data is 4 times worse than that of the C 18 O data.The remaining two basebands are located in the upper and lower sidebands, respectively, with the maximum bandwidth of 1875 MHz.While these two basebands were adopted for the 1.3 mm continuum emission, the central frequency of baseband 2 was adjusted to incorporate as many important molecular lines as possible, such as CH 3 OH, SiO, DCN, H 2 CO, and c-C 3 H 2 .Even though the velocity resolution of these observations, ∼1.3 km s −1 , is not ideal to study the kinematics from these molecular transitions, their detection can still provide useful constraints on circumstellar structures based on, e.g. the integrated intensities.All the major emission lines observed in our observations are listed in Table 2.In order to achieve an angular resolution of ∼0.04, the C43-8 antenna configuration was adopted for our observations, providing long baselines up to ∼12.6 km.In addition to the C43-8 antenna configuration, the C43-5 antenna configuration was also adopted to cover shorter baselines with a minimum of ∼15 m.This additional antenna configuration covering shorter baselines is crucial for us to investigate relatively extended structures up to ∼ 2 − 3 in the envelopes surrounding our sampled protostars.We should note that a large fraction of the emission from the envelopes around the sampled protostars is still filtered out in our observations.However, this missing flux has less impact on the main objectives of the eDisk program focusing on dust and gas on the scales of the disks in the vicinity of the protostars.The Field of View (FoV) of the 12 m array at the observed frequency is ∼ 26 .Logs of the observations are provided in Appendix B (see Table 5).
For two sources in our sample, TMC-1A and B335, comparable observations were obtained with ALMA prior to our eDisk program and we collected those data from the ALMA archive.In addition, due to scheduling constraints only short baseline data were obtained for Oph IRS 63 as a part of eDisk.Hence we relied on ALMA archival observations of Oph IRS 63 for long baseline data.The information on the ALMA archival data we use in the eDisk program is summarized in Appendix B (see Table 6).After delivery, the data were restored to a calibrated state running scriptForPI.py in the same CASA version that the calibration pipeline ran with.For four sources (BHR71 IRS1/IRS2, IRAS 04166+2706, and IRAS 04169+2702), there were Execution Blocks (EBs) that were given a QA0 status of 'semi-pass' (see Table 5 in Appendix B), which in the case of these datasets was because of excessive phase decorrelation during the observations.For those we downloaded the raw data and ran the ALMA calibration pipeline manually to inspect the datasets and determine if they could be salvaged using self-calibration.The QA0 semi-pass data that are used in the final images are denoted in the Observing Logs included in Appendix B.
Outside of standard calibration, we needed to perform additional processing to improve the S/N of the continuum data through self-calibration, to remove the continuum emission from the line data, and to combine data taken in two configurations in order to achieve the balance between flux recovery and angular resolution needed to achieve the science goals of the project.Nearly all of the eDisk targets were sufficiently bright that images created using the standard ALMA calibration were dynamic range limited, meaning that the im-  Note- † For these sources, we used only ALMA archival data without new observations.Some of the archival data we used have slightly different phase centers from those shown in this table (see Table 6).Column 1: Target name.ages were not reaching the thermal noise limit, and were instead limited by phase and/or amplitude variations during the science target observations.To both combine the two configurations, self-calibrate, continuum subtract, and image the data, we constructed a set of semi-interactive scripts for each target.Our methodology was based on the DSHARP data reduction procedure (Andrews et al. 2018) with substantial modification.One particular choice we made was to operate on each measurement set individually for calculation and application of phase and amplitude gain solutions and not concatenate our data to a single measurement set.
This avoids issues that can be encountered when operating on concatenated measurement sets within various CASA tasks; for example, interpolation of gain solutions between distinct observations that might be far apart in time, spectral window numbers being offset, and possible issues with combined metadata to name a few.Many of the issues with concatenated measurement sets can be worked around, but we avoid the known issues and the possible unknown issues by operating on individual measurement sets.The scripts used for each source are provided at http://github.com/jjtobin/edisk.The reduction procedure starts with the construction of a spectrally-averaged continuum dataset.We identified line-free regions of the individual spectral windows to include in the continuum dataset in two ways: first, we made use of the cont.datfile produced by the ALMA imaging pipeline, and edits to the continuum ranges in that file were made if necessary, based on inspection of the hif findcont weblog.Second, we manually flagged the spectral lines that are expected to be present in all datasets regardless of whether they were already included in the cont.datfile.We then translated these frequency ranges defined in the kinematic local standard of rest frame (LSRK) to the channel ranges for each EB and flagged those channels.Then the separate spws from each EB were spectrally averaged individually.The four high-resolution spws were averaged over 480 channels and the three low-resolution spws were averaged 60 channels, resulting in ∼30 MHz channels, approximating the spectral resolution of the low-spectral resolution (time domain mode; TDM) of the ALMA correlator.
We next checked the relative flux density scaling of each individual EB using their azimuthally-averaged visibility amplitude data.To do this, we first imaged each EB individually and shifted the brightest source in the image to the phase center using the CASA task fixvis.We then computed the azimuthally-averaged visibility amplitudes for each EB, plotted them together, and computed the average scale factors between each EB.The scale factors were estimated between the minimum overlapping uv -distance out to 800 kλ; this outer uvdistance was chosen to avoid a possibility that decorrelation at longer baselines could significantly impact the scale factors.The offsets were typically at the 10-15% level between EBs and we scaled each EB to the 'best' EB which was judged by having the least amount of scatter in its data.Thus, our goal was a relative flux normalization of all the EBs, given that we were not improving on the overall flux density scale accuracy.We also note that there is possibility that some scale factors include intrinsic source variability in addition to calibration differences.While the ALMA Technical Handbook indicates that the accuracy of the flux density scale is expected to be ∼5%, the scaling that we have had to apply to datasets indicates that the flux density scale can have uncertainties up to the rescaling factors that we apply to the data. 2 .The scale factor applied to each EB is provided in Appendix C.
There were some cases where this flux normalization could not be performed at this stage due to decorrelation, which caused the visibility amplitude profile of one or more EBs to have a different behavior as a function of uv -distance than others.In those cases, we did not determine scale factors at this stage, but we completed the self-calibration steps that followed and determined the scale factors after self-calibration.We then re-ran the whole procedure from the start, applying these scale factors a priori ; we referred to this as the 'two-pass' method.The 'two-pass' method was applied to Oph IRS43, IRAS 04166+2706, IRAS 04169+2702, and BHR71 IRS1 and IRS2.The relative scaling of each EB was then applied to the continuum datasets that did not have the brightest source shifted to the phase center because the unshifted data were what we ultimately used in the self-calibration process.
Following flux normalization of each EB we began the self-calibration process.For the combined processing of the compact and extended configuration data, we first self-calibrated the compact configuration data and then we self-calibrated the pre-self-calibrated compact configuration data with the extended configuration data.Self-calibration was performed using the data that did not have the phase center shifted.It was not neces-2 There is the Jorsater & van Moorsel effect in clean (see Section 4.2) that may also add uncertainties to the measured flux densities in the cleaned images sary to align our data prior to self-calibration because the time difference between observations was less than 1 year and all the long baseline data, where proper motion would have the greatest effect, were all taken within a 3 month time period.We create our model for selfcalibration using the data from all the measurement sets; however, it is only the calculation and application of the phase and amplitude gain solutions that happen on a per-measurement set basis.
To prepare for self-calibration, we first calculated a ranked list of reference antennas for each EB, equally weighing the amount of flagging on a particular antenna and its geometric position within the array.Following this step, we made an initial image using all the continuum datasets to determine the peak S/N of the image from standard calibration; the creation of all images during the self-calibration process made use of the tclean auto-multithresh algorithm (Kepley et al. 2020) to automatically create clean masks for deconvolution.Following the generation of this initial image, we analyzed the scan lengths within an observation and determined calibration solution intervals that most evenly divide a scan by the number of integrations per scan.This approach will avoid the uneven division of scans which can result in some solutions being calculated from significantly less data.The list of solution intervals to attempt always includes a solution that spans the length of an EB, referred to as 'inf EB' where gaincal is run with the parameter combine='scan,spw'.The 'inf EB' solution interval corrects systematic effects, in particular, antenna position offsets that translate to systematic phase offsets between the calibrator and science target.Then for shorter solution intervals, we are correcting the atmospheric phase variations.We include an 'inf' solution, where the 'scan' setting is removed from the gaincal parameters and solutions break at scan boundaries, then for shorter solutions, we divide the median scan length by 2 and all subsequent solution intervals divide the previous interval by 2 until the final solution interval of 'int' (one integration) is reached.
Then for each solution that is to be attempted, we also determine the clean depths in terms N× σ to use for each solution interval in order to progressively clean deeper during each round of self-calibration.We begin with the assumption that we would trust features with intensities > (Peak S/N)/15, which means that for a Peak S/N = 150, we would trust features > 10σ.Note that these numbers are empirical, but small changes do not appreciably affect the outcome of self-calibration.So we divide the S/N of the initial image by 15 to obtain our starting clean threshold, and decreased the thresh-olds for subsequent solution intervals logarithmically until the clean depth is equivalent to 3σ.
With a set of solution intervals to attempt and their accompanying clean thresholds, the self-calibration process begins.We collect all the model creation with tclean, calibration solution generation with gaincal, application of the calibration solutions with applycal, and checking the outcome of a self-calibration iteration with tclean into a single Python function to simplify to process.We always use the combine='spw' parameter in gaincal to use all the available bandwidth for selfcalibration and we use the parameter interp='linearPD' which will scale the phase solutions with frequency to compensate for the fact that the optimal phase solutions will be different between the low and high frequency end of the bandwidth.
At the start of each self-calibration iteration, the selfcalibrated data from the previous iteration (corrected data column) is split out into a new measurement set such that each self-calibration solution interval is incremental upon the previous iteration.We compare the image S/N before self-calibration and after self-calibration of each iteration is applied to ensure that the S/N is increasing.We do this by running tclean again but using the self-calibration model that was created at the start of a self-calibration iteration as the startmodel for tclean.This approach enables us to effectively determine the improvement in S/N just due to the gain corrections derived from that model.When the image S/N increases within an iteration, we consider that iteration successful, if the S/N decreases, that iteration has failed and we consider the last successful iteration as the final selfcalibration iteration.The above steps consider phaseonly self-calibration and we do attempt amplitude selfcalibration (solving for amplitude and phase simultaneously) following the phase-only self-calibration using the same criteria for success or failure of amplitude selfcalibration.Following the completion of self-calibration for the successful solution intervals, the self-calibrated data for each EB are saved to new measurement sets and are considered the final continuum data.
Once self-calibration is completed on the continuum data, we process the line data.First, we apply the flux normalization to each full measurement set using the scaling factors derived from the continuum.Then we apply the self-calibration solutions to the full measurement sets for each EB.However, the S/N improvement for the line data are typically not significant because the line data are not dynamic range limited like the continuum.For each EB, we apply all the gain tables for each successful self-calibration iteration, which is necessary because each table is incremental, so applying the full correction requires all the tables.Once self-calibration solutions are applied, we subtract continuum emission from the line data using the uvcontsub task.To subtract the continuum data, we select the line-free regions of the data that will have their continuum data fitted and subtracted.To select the frequency ranges to fit, we select the complement of the frequency ranges that were flagged to produce the continuum dataset, that is the line-free regions of the bandwidth.Following uv continuum subtraction, we have the final spectral line data set and it is ready for imaging.
We first create the continuum images using tclean and the auto-multithresh algorithm to create the clean masks.Then we create images using Briggs weighting, for the robust parameters of -2, -1, -0.5, 0, 0.5, 1, and 2 in order to compare the features present with different weighting.We also created images using uv -tapering, to reduce the influence of the long baselines and bring out larger-scale structure, and we created tapered images with the tapers starting at 1000, 2000, and 3000 kλ and used robust parameters of 1 and 2. The typical synthesized beam size for the robust=0.5images was ∼0.05.The images for most sources use 6000 pixels with a pixel size of 0. 003; however, some sources require image sizes of ∼14000 in order to include sources toward the edge of the primary beam.
Then we create the spectral line images for each line listed in Table 2.We typically created images using Briggs weighting with robust values of 0.5 and 2, and uv -tapering starting at 2000 kλ.The resulting images have synthesized beams of ∼0. 1 and ∼0.15 for the two robust values and applied tapering, respectively.The tapering was essential for the spectral line images because the surface brightness sensitivity was otherwise too low at the native resolution of the data to detect the spectral lines with high S/N.The spectral line images typically used 4000 pixels and a pixel size of 0. 01.The spectral line images were all created with velocity-defined channels around the spectral lines of interest using their associated rest frequencies (Table 2).
We then save the non-primary beam corrected images (flat noise), the primary beam corrected images, the final clean mask, and the primary beam response for each individual image.As such, there are ∼52 images created per target when all the combinations of robust parameters, lines, and tapering are considered.The images presented in this overview paper and subsequent papers from the eDisk program are selected based on the balance between resolution and sensitivity, and also the structures discussed.

Jorsater & van Moorsel effect
Recently it has been pointed out that data imaged with the clean algorithm may suffer systematic errors of the flux scale due to the so-called Jorsater & van Moorsel (JvM) effect (Jorsater & van Moorsel (1995); see, in particular, Czekala et al. (2021)).The reason for this is that the dirty beams in the CLEAN processes deviate from the used Gaussian restoring beams.As discussed in Czekala et al. (2021), the JvM effect is more significant when data taken from different antenna configurations are combined like in the case of the eDisk program.The impact due to the JvM effect to the flux scales of our eDisk maps differs from source to source: maps with more extended diffuse emission are more impacted by the JvM effect.However, it has also been pointed out that the correction method to fix the flux scale issue, implemented by Czekala et al. (2021), may artificially manipulate the noise levels of maps after the correction (Casassus & Cárcamo 2022).As the main discussions of the first set of papers from eDisk are based primarily on the morphology of maps and velocity structure of molecular emission rather than the flux scale of maps, we did not apply any correction to our maps in this overview paper and the related first set of papers.

FIRST RESULTS
In this section, we present the first results of our eDisk program, an overview of the observations of the 1.3 mm continuum for the sample of sources as well as an example of analysis of the C 18 O 2-1 emission toward R CrA IRS7B.The structures of embedded disks can be traced well in the 1.3 mm dust continuum emission, while C 18 O is the key molecular line to study the kinematics of the gas on small scales, e.g.searches for the presence of Keplerian disks and estimates of the dynamical masses of the central protostars.In this sense, the presentation here provides a glance at the type of scientific results that can be obtained from the eDisk data.More details of the first results about each targeted source can be found in a series of subsequent papers (see the references in Table 3).

Continuum emission
We present representative maps of the continuum images obtained in the eDisk program on the same spatial scale in Figure 2. The same continuum images but enlarged are also presented in Figure 3. Robust parameters used to create these continuum images and their angular resolutions are shown in Table 3.The rms noise level of these continuum maps is ∼25 µJy beam −1 on average.Continuum emission is clearly detected at all source positions, and four out of the 19 sources also show continuum emission arising from companions; three of them, i.e., Ced110 IRS4, R CrA IRS7B, and R CrA IRAS 32, are newly found to be binaries in our observations.These newly found binaries are named based on the IAU convention; Ced110 IRS4A and Ced110 IRS4B, R CrA IRS7B-a and R CrA IRS7B-b, and R CrA IRAS 32A and R CrA IRAS 32B.
All the continuum emission images are spatially resolved, exhibiting elongated structures.These elongated structures strongly suggest that the dust continuum traces disks around the targeted sources.Note that for some sources envelope contamination may still play a role for both dust and gas emission: to assess that more detailed radiative transfer modeling will be needed on a source-by-source basis.The extent of the continuum emission, measured through 2-D Gaussian fittings and listed in Table 3, varies significantly from source to source; L1489 IRS has the largest emission with an angular extension of ∼ 3. 9 × 1. 3, while the companion of Oph IRS43 (VLA 2) is the most compact, marginally resolved (2σ) with an angular extent ∼ 0. 019 × 0. 016.We note that some of the continuum emissions were not well fitted with 2-D Gaussian, suggesting that the extents of the emissions provided in Table 3 for those sources may not reflect their actual extents.In such a sense, the extents of the continuum emissions provided here should be considered as a quick reference of their extents measured by a simple procedure.More accurate measurements of the extents of the dust emissions through more sophisticated methods will be presented in future papers.A comparison between the extents of the continuum emission from the different sources and their bolometric temperatures is shown in Figure 4.In agreement with the systematic studies of protostars in Orion with the VLA (Tobin et al. 2020), no clear correlation between size and evolutionary stage is seenalthough the two largest disks in our sample are both associated with Class I protostars (L1489 IRS and IRAS 04302+2247).Most of the disks in the eDisk sample are smaller than about 100 au but we caution that these are estimates of the sizes of the dust disks and that the gas disks may be larger as it is often inferred for more evolved disks (see, e.g., Miotello et al. 2022).
From a comparison between the major and minor axes of the continuum emission, the inclination of the system can be also estimated on the assumption that the continuum emission is a vertically thin disk.We find that 7 out of 19 sources, excluding companions, have inclination angles of more than 70 degrees (see Table 3), roughly consistent with the expectation from a random distribution.We caution that the estimate of the inclination is a lower limit to the true inclination because the disks have most likely some vertical thickness.Still, these are different from most Class II disks that have been ob- served at high angular resolution to date which tend to be more face-on (Andrews et al. 2018;Cieza et al. 2021).3. Images are presented in an order of T bol listed in Table 1.The beam size and the scale of 20 au are shown at the bottom left and right corners in each panel, respectively.
Figure 5 shows a comparison between the inclination angle and the bolometric temperature, demonstrating no correlation between them.This suggests that the disk inclination has less impact on our classification of the 19 protostellar systems based on the bolometric temperature even though it was also suggested that higher inclination angles of disks could misclassify protostars as Class 0 sources (Tamura et al. 1996;Nakazato et al. 2003;Tobin et al. 2008;Fischer et al. 2017).We also note that some of the continuum structures show relatively high brightness temperatures of more than 100 K.Such high brightness temperature may suggest that the detected continuum emission is at least partially optically thick.Basic physical parameters of the continuum emission are summarized in Table 3.We note that statistical errors are not given in this table because the ac-  tual uncertainties of these measurements are more likely determined by systematic errors, such as the flux calibration error and the gain calibration error.
When we carefully inspect the continuum maps, it is found that some of them show ring-like structures that can be visually identified by naked eye.The most obvious source is L1489 IRS.Oph IRS63 also weakly shows ring-like structures.More details of the continuum maps obtained in the eDisk program will be presented in subsequent papers.Ring-like structures in continuum emission around these two sources have been suggested by previous studies as well (Ohashi et al. 2022;Segura-Cox et al. 2020) although our map of L1489 IRS seems to show ring-like structures more clearly than in the previous work (Ohashi et al. 2022).The continuum emission around IRAS 04169+2702 seems to show a bean-like structure with brighter emission on its southwest side, which could also be a ring-like structure.We should stress that a more detailed and systematic analysis of continuum emission is definitely required to identify their substructures, including those that cannot be visually identified by the naked eye, as will be discussed in a forthcoming paper.Nevertheless, the first analysis on the observations suggests that continuum emission around our Class 0/I sample does not show very sharp ring or spiral structures that are often seen in Class II continuum disks (e.g., DSHARP, ODISEA).This striking difference will be discussed in more detail in Sec.6.1.On the other hand, several disk continuum images show brightness asymmetries, particularly along their minor axes.IRAS 04302+2247, GSS 30 IRS3, IRAS 16253-2429, IRAS 16544-1604, R CrA IRS7B-a, R CrA IRAS 32A, and R CrA IRAS 32B show clear asymmetry that can be identified by eye, while asymmetry can be also identified with more careful inspection of continuum intensity distributions for some additional sources, such as L1527 IRS.For these sources showing brightness asymmetries, the geometric centers of the continuum emission derived from single 2-D Gaussian fittings listed in Table 3 are shifted from their actual continuum peak positions.The possible origin of the brightness asymmetry is discussed in Section 6.1.

C 18 O 2-1
The C 18 O 2-1 emission is also detected toward all 19 sources in our sample.Here, we show the C 18 O emission detected in one of our sources, R CrA IRS7B, as an example of what can be seen toward the eDisk targets.Lindberg et al. (2014) have carried out ALMA observations of IRS7B in Band 7 at an angular resolution of ∼ 0. 4 and detected continuum emission and several molecular lines, including C 17 O 3-2.Both continuum and C 17 O emission have extended structures elon- Note- † For these sources, we used only ALMA archival data without new observations.Column 1: target names.C 18 O channel maps with a velocity resolution of 0.334 km s −1 (two channels combined) in R CrA IRS7B.Solid contours are drawn from 3σ to 15σ every 3σ for positive emissions, while dashed contours are drawn from -3σ to -15σ every -3σ and less every -5σ.1σ is 1.2 mJy beam −1 .The center of each channel map is the position of R CrA IRS7B-a shown in Table 3.The black and red crosses in each panel indicate the positions of IRS7B-a and IRS7B-b, respectively.The beam size and a scale bar of 50 au are shown at the bottom-left and the bottom-right corners in the bottom-left panel.
gated in the southeast to northwest direction.The C 17 O emission also shows a clear velocity gradient along the elongation, which was interpreted as Keplerian rotation (Lindberg et al. 2014).Although a systemic velocity of 6.2 km s −1 was adopted for the analysis of the velocity structures of C 17 O in Lindberg et al. (2014), we adopt ∼6.0 km s −1 as the systemic velocity of IRS7B in this paper (see more details below).C 18 O 2-1 was detected at more than 3σ in the LSR velocity range between −3.15 km s −1 and 14.22 km s −1 at an angular resolution of ∼ 0. 13 × 0. 11 as shown in Figure 6.Note that the velocity resolution of the C 18 O data shown here is ∼0.33 km s −1 , which is double the original velocity resolution to achieve a higher sensitivity.Since the systemic velocity of this system is ∼6 km s −1 (see Sec. 6.2), LSR velocities smaller and larger than 6 km s −1 correspond to blue-and red-shifted, respectively.At higher blue-shifted velocities (V LSR = −3.15 to 2.86 km s −1 ) and red-shifted velocities (V LSR = 8.20 to 14.22 km s −1 ), C 18 O is mainly detected on the southeast side and the northwest side of IRS7B-a within ∼0. 5, respectively.At lower blue-shifted (V LSR = 3.19 to 4.20 km s −1 ) and red-shifted velocities (V LSR = 7.20 to 7.87 km s −1 ), more extended C 18 O emission is detected within the field of 2 ×2 , demonstrating that our observations can also trace extended emission as well.At velocities close to the systemic velocity, absorption of C 18 O is detected.In particular, very deep absorption is seen at the position of IRS7B-a at V LSR between 5.20 and 5.87 km s −1 .The extent of the deep absorption is almost the same as that of the continuum emission associated with IRS7B-a, suggesting that the deep absorption is due to cold, foreground C 18 O emission against the bright continuum emission.Compact and weak C 18 O emission is detected at the position of IRS7B-b at several channels of red-shifted velocities, such as V LSR = 9.21 and 11.21 km s −1 , and also absorption at V LSR = 5.20 and 5.53 km s −1 .The absorption at the position of IRS7B-b may be also due to cold foreground C 18 O against the continuum emission associated with IRS7B-b.
To investigate the overall distribution and velocity structures of C 18 O around IRS7B-a, moment 0 and 1 maps of IRS7B are shown in Fig. 7b.For comparison, the continuum map of IRS7B is also shown in Fig. 7a.The moment 0 map shows two emission components with semi-elliptic shapes, one located to the southeastern side of IRS7B-a and the other located on the northwestern side of IRS7B-a, forming a disk-like structure divided into two along its minor axis as a whole.The C 18 O disk-like structure is basically similar to that seen in the continuum associated with IRS7B-a, although it extends slightly larger than the continuum emission, suggesting that C 18 O probably traces gas associated with the protoplanetary disk around IRS7B-a, but further investigation of the kinematics of C 18 O is required for us to have a firm conclusion (see below).Although both continuum and C 18 O show similar disk-like structures, their intensity distributions are quite different, e.g., C 18 O is stronger in outer parts of the disk-like structure and it gets weaker toward the position of IRS7B-a, where the continuum emission shows a peak.This discrepancy would be most likely due to absorption of the warm continuum emission by colder C 18 O gas.
The moment 1 map presented in Figure 7b shows a velocity structure of the C 18 O emission having blueshifted velocities with respect to the systemic velocity of ∼6 km s −1 on the southeastern side and red-shifted velocities on the northwestern side.This clear velocity gradient along the major axis of the disk-like structure is naturally interpreted as rotation.The nature of the velocity structures can be investigated further with a position-velocity (PV) diagram cutting along the major axis, as discussed in Sec.6.2.The diagram shows that the velocity of C 18 O increases toward the position of IRS7B-a, a clear sign of differential rotation.Further discussions on the nature of the rotation are presented in Section 6.2.

Substructures and asymmetry of continuum emission
When we compare our continuum maps from Fig. 2 with those toward Class II sources, from, e.g., the DSHARP and ODISEA programs, the absence of clear ring or spiralstructures identifiable by the naked eye with few exceptions is remarkable.The noticeable difference between protostars and Class II sources could suggest that substructures are not well developed in the embedded phase yet, and they could be quickly formed when Class 0/I protostars evolve into Class II sources.Interestingly, the two sources, L1489 IRS and Oph IRS 63 showing visually identified ring-like structures and another source, IRAS 04169+2702, possibly showing a ring-like structure in the eDisk sample are all Class I protostars.In addition, two protostars previously found to have ring-like structures in their continuum maps are also Class I protostars (WL 17 and GY 91).Note that substructures found in dust continuum around GY 91 are less sharpe than those seen around Class II sources (see Sheehan & Eisner 2018;Cieza et al. 2021) while WL 17 does show a sharp ringlike structure (Sheehan & Eisner 2017).Although seven sources from the Orion VANDAM survey exhibit rela- shown in contour and moment 1 map shown in color.The moment 0 map was created by integrating emission having intensity more than 3σ (1 σ = 1.2 mJy beam −1 ) in the velocity range of -3.15 km s −1 to 14.22 km s −1 to avoid the influence of the deep absorption features.The moment 1 map was created over the same velocity range as the moment 0 map with the emission having an intensity more than 5σ.The contours of the moment 0 map are drawn every 2σ from 3σ (1σ = 2.9 mJy beam −1 km s −1 ).The black and red crosses indicate the positions of IRS7B-a and IRS7B-b, respectively.The beam size is indicated in the bottom-left corner with the black ellipse.
tively clear ring structures, these ring structures could be due to close binary formation (Sheehan et al. 2020).
We should note, however, that a reason behind the lack of obvious sharp substructures in continuum emission around the Class 0/I protostars could also be that the continuum emission at 225 GHz may be optically thick.In fact, some of the continuum emission observed in the eDisk program may be optically thick as indicated by their higher brightness temperatures (see Table 3).To address this possibility, it will be necessary to observe continuum emission at lower frequencies, where dust continuum emission is expected to be optically thinner.One of the sources in the eDisk sample, L1527 IRS, has been observed in 7 mm continuum with VLA, finding no clear gaps (Sheehan et al. 2022).
We should also note that some of the disk structures seem to be more inclined (see Table 3), particularly as compared with Class II disks observed at high angular resolution to date, which could make substructures more difficult to be observed in disks around our sample.In addition, some of the continuum emissions are also more compact, and our angular resolution may not be high enough to detect substructures.More detailed analyses and discussions of substructures in the continuum emission, including the visibility analysis, will be presented in forthcoming papers.
In contrast to substructures, several of the target sources show brightness asymmetries, particularly along their minor axes (see Figure 3).A possible explanation for these brightness asymmetries along the minor axes might come from a geometrical effect, namely, the farside of a highly inclined disk with a finite vertical thickness is brighter than the near-side.Such a near-far side brightness asymmetry is an indication that the dust has yet to completely settle onto the disk midplane.In more edge-on systems such as HH 212, the non-settled dust produces a dark lane near the mid-plane sandwiched between two brighter features (Lee et al. 2017;Lin et al. 2021).In contrast to embedded disks, Class II disks rarely show any asymmetry along their minor axes (e.g., Andrews et al. 2018;Long et al. 2019).The difference could suggest that dust is not quite settled in the embedded phase but quickly settles as embedded disks evolve into Class II disks.

Identifying Keplerian rotation
In addition to investigating substructures in continuum emission, another important purpose of the eDisk program is to identify Keplerian rotation around our target protostars.Here we discuss the case of R CrA IRS7B-a as an example showing how we identify Keplerian rotation.
One of the methods to identify Keplerian rotation, which is also adopted by the eDisk program, is fitting to position-velocity (PV) diagrams.Our fitting process includes two independent methods.One method uses the ridge in 1D intensity profiles along the positional axis (the velocity axis) at a given velocity (position).The other method uses the outer edge in the 1D intensity profiles.We call the former the ridge method, while we call the latter the edge method.In each method, PV diagrams are fitted with a power-law function.The function is a single-or double-power function, which is defined with a power-law index, p in .The systemic velocity, v sys , is also included as a free parameter in the fitting process.More details of the fitting method and process are explained in Appendix D.
As was shown in Section 5.2, the C 18 O emission in R CrA IRS7B-a shows a velocity gradient along its disk-like elongated structure centered at the position of R CrA IRS7B-a, suggestive of rotation.The C 18 O PV diagram cutting along the major axis of the dust emission around R CrA IRS7B-a presented in Figure 8 shows that the rotation can be described as differential rotation, and it is fitted in both the edge and ridge methods described above.The fitting results are presented in Table 4.For demonstration purposes, only fittings with a single-power law function are presented here.Both edge and ridge fittings with a single powerlaw show that p in is close to 0.5, suggesting that the rotation can be explained as Keplerian rotation.The dynamical mass of the central protostar was estimated to be ∼ 3.2 M with the edge method, while it was estimated to be ∼ 2.1 M with the ridge method.As explained in Appendix D, the actual dynamical mass of the central star is considered to be between ∼ 2.1 M and ∼ 3.2 M .Discussions on further details of the Keplerian disk around R CrA IRS7B-a are beyond the scope of this overview paper, and they will be presented in the subsequent paper.

DATA RELEASE
As one of ALMA large programs, the eDisk program will release a suite of data products that have legacy values for the community.The released data sets, which will be available at https://almascience.org/alma-  Contours are drawn from 3σ to 15σ in every 2σ for positive emission, while dashed contours are drawn from -3σ to -15σ in every -3σ and less every -5σ.1σ is 1.2 mJy beam −1 .Cicles and triangles show data points obtained through the ridge and edge methods, respectively, while the solid and dashed curves represent fitting results obtained through the ridge and edge methods, respectively.Different colors show data points obtained along different direction: blue and red show data points obtained along the position axis (xcut), while cyan and magenta represent data points obtained along the velocity axis (vcut).
data/lp/edisk, include (1) CASA scripts and associated python modules used to calibrate and image the data; (2) calibrated visibility data from each execution block, including those from the ALMA archive; (3) continuum images created from combined visibility data with a robust value of 0.5; (4) line cube images created from combined visibility data with a robust value of 0.5; (5) fits file of the primary beam for each image mentioned above.These released data sets enable the community to reproduce all the images presented in this paper and forth coming papers.All the continuum emission is spatially resolved at angular resolutions of ∼ 0. 04, revealing elongated disklike structures.The detected continuum emission is considered to trace embedded disks around the sampled Class 0/I protostars although they could also partially arise from inner parts of the envelopes surrounding the protostars.More importantly, the continuum emission does not show any sharp ring-like or spiral-like structures with high contrast, which are often seen in Class II disks.This could suggest that substructures are formed in disks quickly when Class 0/I protostars evolve to Class II sources or otherwise that high optical depth of dust emission may obscure such features.Note that some of the continuum emission show more compact or inclined disk-like structures, which could prevent us from finding substructures there.Asymmetries in the brightness distributions along minor axes seen by eye for some sources may also reflect the dust optical thickness and orientations of the disks.
C 18 O 2-1 emission detected in R CrA IRS7B-a also shows an elongated disk-like structure similar to that seen in its dust continuum emission.There is a clear velocity gradient along the elongation.The analysis of the velocity structures in the PV-diagram cutting along the elongation suggests that the velocity gradient seen in the C 18 O elongation can be explained as Keplerian rotation with a central dynamical mass of 2.1-3.2M .
The first results here already demonstrate the great potential of the data from the eDisk program.A set of accompanying papers go into more detail about the first results for individual sources setting the stage for future more synergistic studies.These analyses will provide new, important constraints on these pivotal stages in the formation and early evolution of protostars and their disks.ridge points in a PV diagram and fitting the obtained points with a power-law function.First, the edge and ridge radii are obtained in the 1D intensity profile along the positional axis at the given velocity (xcut).The edge radius is defined as the outermost position with a given threshold level (3σ in this paper) of emission.The uncertainty of the edge radius is calculated from the emission gradient and the observational noise level.The ridge radius is defined as the intensity-weighted mean position in the 1D profile or the center derived from the Gaussian fitting.This fitting also calculates an uncertainty of the ridge radius from the observational noise level.Similarly, the edge and ridge velocities are obtained in the 1D profile along the velocity axis at a given position (vcut).Then, in each of edge and ridge methods, the xcut points at velocities higher than a "middle" velocity and the vcut points at radii larger than a "middle" radius are combined as the final pairs of (r, v).The middle velocity and radius are determined from the closest pair of points between the xcut and vcut, as explained in Aso & Machida (2020) in more detail.Note that 1D intensity profiles along only one direction (position or velocity) would be enough to obtain the edge and ridge points if a PV-diagram is mainly elongated along the velocity or positional axis, respectively.Second, the obtained pairs of (r, v) are fitted with a power-law function, in the edge and ridge methods separately.The function is a single-or double-power function as follows: where v b is the velocity at the radius of r b , p in is the power-law index, dp is the difference of the index between the inner and outer parts, and v sys is the systemic velocity.dp is limited to be ≥ 0 to produce a steeper or the same index outside r b , such as the inner Keplerian rotation versus the outer rotation where the angular momentum is conserved.For the single-power fitting, v b is fixed at a middle value and dp = 0. Equation D2 provides a function of radius, defined here as V fit (r).When dp ≥ 0, the inverse function of velocity R fit (v) can be defined.Using these functions, χ 2 to be minimized is defined as where σ is the uncertainty of each radius/velocity obtained above.The power-law fitting provides the pair of (r b , v b ).From this pair, the central stellar mass can be calculated as M * = v 2 b r b /G/ sin 2 i, where G is the gravitational constant and i is an inclination angle.The inclination angle is often estimated from the aspect ratio of the associated continuum emission or independent modeling of the associated outflow, which is beyond the scope of this fitting with PV diagrams.If the power-law index p in is significantly different from the Keplerian index 0.5, the estimated M * does not have the meaning of the central stellar mass.

Figure 1 .
Figure 1.(a) T bol -L bol diagram of the eDisk sample sources.(b) T bol -L bol diagram of the protostars with d < 200 pc and L bol > 0.1 L .Two dashed vertical lines in each diagram indicate the boundaries between Class 0 and Class I (70 K; see Evans et al. 2009) and Class I and Class II (650 K; see Evans et al. 2009).
Detailed procedureThe eDisk data were all calibrated using the ALMA calibration pipeline packaged with the Common Astronomy Software Application (CASA)(McMullin et al. 2007) version 6.2.1 and pipeline version 2021.2.0.128.

Column 2 :
Alternative common name also used in previous works.Column 3: Right ascension of the observation phase center in the International Celestial Reference System (ICRS).Column 4. Declination of the observation phase center in ICRS.Column 5. Evolutionary stage classification based on T bol .Column 6. Distance to the target based on recent estimations from Gaia measurements.Column 7: Bolometric temperature re-derived from the newly compiled SED.See Appendix A. Column 8: Bolometric luminosity re-derived from the newly compiled SED.See Appendix A. Column 9: Reference for the distance estimations.

Baseband 3 (
233.3125-235.1875GHz)   Mostly for the 1.3 mm continuum emission Note- † These two lines are blended.‡ In the DDT observations, spw 1 and spw 3 in the baseband 1 were swapped, and the frequency range of the spw 3 in the DDT observations was set at 220.369383-220.427977GHz.

Figure 2 .
Figure 2. Gallery of the continuum maps obtained in the eDisk program.All maps are shown in the same linear scale (large panels: 640 au × 640 au, medium panels: 320 au × 320 au, small panels: 160 au × 160 au).All the maps are shown in individual asinh-stretch intensity scales, except for IRAS 04169+2702 that is shown in a linear intensity scale.All the map centers are source positions listed in Table 3, except for Ced110 IRS4, R CrA IRAS 32, and Oph IRS43, where the map centers are shifted from the source positions by (0. 75, 0. 0), (−0.5, 0. 5), and (0. 0, −0.3), respectively to show the binary.The beam size and the scale of 20 au are shown at the bottom left and right corners in each panel, respectively.

Figure 3 .
Figure 3. Enlarged continuum images.Images of IRAS 04302+2247 and L1489 IRS are not enlarged from those in Figure 2. The intensity scales are the same as those in Figure 2 except for the 4 companions.All the map centers are source positions listed in Table3.Images are presented in an order of T bol listed in Table1.The beam size and the scale of 20 au are shown at the bottom left and right corners in each panel, respectively.

Figure 4 .
Figure 4. Dust disk size in diameter is plotted as a function of bolometric temperature of their central protostellar system.Sources located in different regions are shown in different symbols.Open symbols show the 4 binary systems (Ced110 IRS4A/B, Oph IRS43 VLA1/VLA2, R CrA IRS7B-a/b, and R CrA IRAS 32A/B).The vertical dashed line shows the boundary between Class 0 and I sources.

Figure 5 .
Figure 5. Dust disk inclination angle is plotted as a function of bolometric temperature of their systems.Four binary companions (Ced110 IRS4B, Oph IRS43 VLA2, R CrA IRS7B-b, and R CrA IRAS 32B) are excluded.The vertical dashed line shows the boundary between Class 0 and I sources.
Column 2: Right ascension of the continuum emission measured with single 2-D Gaussian fittings.For L1489 IRS 2-D Gaussian fittings with twocomponent were performed because single 2-D Gaussian fittings did not provide a solution.Column 3: Declination of the continuum emission measured with single 2-D Gaussian fittings.For L1489 IRS 2-D Gaussian fittings with two-component were performed because single 2-D Gaussian fittings did not provide a solution.Column 4: Robust parameter used to create the continuum image shown in Figures 2 and 3. Column 5: synthesized beam FWHM and position angle of the representative map shown in Figure 2. Column 6: Beam deconvolved size and position angle of the continuum emission measured with single 2-D Gaussian fittings.For L1489 IRS 2-D Gaussian fittings with two-component were performed because single 2-D Gaussian fittings did not provide a solution.Column 7: Peak intensity of the continuum emission and its brightness temperature estimated with the full Plank function.Column 8: Integrated flux density of the continuum emission measured with single 2-D Gaussian fittings.For L1489 IRS 2-D Gaussian fittings with two-component were performed because single 2-D Gaussian fittings did not provide a solution.Column 9: Inclination angle of the continuum emission estimated from the ratio of the major and minor axes of the continuum emission listed in the column 5. Column 10: References from the eDisk program.1. Yamato et al. (2023); 2. Phuong et al. (2023); 3. Han et al. (2023); 4. Lin et al. (2023); 5. van 't Hoff et al. (2023); 6. Sai et al. (2023); 7. Gavino et al. (2023); 8. Thieme et al. (2023); 9. Santamaría-Miranda et al. (2023); 10.Narayanan et al. (2023); 11.Aso et al. (2023); 12. Flores et al. (2023); 13.Kido et al. (2023); 14.Sharma et al. (2023); 15.Ohashi et al. (2023); 16.Encalada et al. (2023); 17.This work

Figure 7 .
Figure 7. (a)  Continuum emission map in IRS7B.IRS7B-a is located at the center of the map, while IRS7B-b is located to the northwest side of IRS7B-a.The beam size is shown in the bottom-left corner as the white ellipse.(b) C 18 O moment 0 map shown in contour and moment 1 map shown in color.The moment 0 map was created by integrating emission having intensity more than 3σ (1 σ = 1.2 mJy beam −1 ) in the velocity range of -3.15 km s −1 to 14.22 km s −1 to avoid the influence of the deep absorption features.The moment 1 map was created over the same velocity range as the moment 0 map with the emission having an intensity more than 5σ.The contours of the moment 0 map are drawn every 2σ from 3σ (1σ = 2.9 mJy beam −1 km s −1 ).The black and red crosses indicate the positions of IRS7B-a and IRS7B-b, respectively.The beam size is indicated in the bottom-left corner with the black ellipse.

Figure 8 .
Figure 8. C 18 O position-velocity diagram of R CrA IRS7Ba with a velocity resolution of 0.334 km s −1 .Contours are drawn from 3σ to 15σ in every 2σ for positive emission, while dashed contours are drawn from -3σ to -15σ in every -3σ and less every -5σ.1σ is 1.2 mJy beam −1 .Cicles and triangles show data points obtained through the ridge and edge methods, respectively, while the solid and dashed curves represent fitting results obtained through the ridge and edge methods, respectively.Different colors show data points obtained along different direction: blue and red show data points obtained along the position axis (xcut), while cyan and magenta represent data points obtained along the velocity axis (vcut).
provided an introduction to the ALMA large program, eDisk, including the background of the program, the sample selection of the 19 Class 0/I protostars, the ALMA observations at Band 6, and the data reduction processes.The first look results of the 1.3 mm continuum emission detected around all the 19 targets for the eDisk program and C 18 O 2-1 emission detected around one of the 19 targets, R CrA IRS7B, are presented to demonstrate what kind of results we would expect to see in the eDisk program.

Figure 9 .
Figure 9. Spectral energy distributions from the near-infrared to millimeter ranges of the eDisk sources.Photometric data and spectroscopic data are shown in filled squares and empty circles, respectively.Flux data taken from the IRSA database are shown in green and Spitzer IRS spectra are shown in blue.Herschel PACS/SPIRE data are shown in red.Deconvolved GSS30 IRS3 data taken from Je et al. (2015) are shown in orange.B335 flux data shown in purple are taken from Evans et al. (2023) and Herschel PACS/SPIRE spectral data shown in red are taken fromGreen et al. (2016).Photometry data for BHR71 IRS1 and IRS2 data taken fromTobin et al. (2019) are shown in magenta.The black line represents the spline fit of the SED values used to calculate each source's L bol and T bol .

Table 2 .
Main molecular lines covered in the spectral setups for the eDisk program

Table 3 .
Physical parameters of continuum emission

Table 4 .
Results of the best fitting to the C 18 O PV diagram

Table 7 .
Scale factors