Hong-Xin Zhang
email: hzhang18@ustc.edu.cn

Rory Smith

Se-Heon Oh

Sanjaya Paudel

Pierre-Alain Duc

Alessandro Boselli

Patrick Côté

Laura Ferrarese

Yu Gao

Deidre A Hunter

Thomas H Puzia

Eric W Peng

Yu Rong

Jihye Shin

Yinghe Zhao

Typeset using L A T E X twocolumn style in AASTeX63

Keywords: galaxies, evolution -galaxies, dwarf -galaxies, starburst -galaxies, interactionsgalaxies, ISM -galaxies, kinematics and dynamics -galaxies, individual(VCC 848)

A clear link between a dwarf-dwarf merger event and enhanced star formation (SF) in the recent past was recently identified in the gas-dominated merger remnant VCC 848, offering by far the clearest view of a gas-rich late-stage dwarf-dwarf merger. We present a joint analysis of JVLA H i emission-line mapping, optical imaging and numerical simulations of VCC 848, in order to examine the impact of the merger on the stellar and gaseous distributions. VCC 848 has less than 30% of its H i gas concentrated within the central high-surface-brightness star-forming region, while the remaining H i is entrained in outlying tidal features. Particularly, a well-defined tidal arm reaches N (H i) comparable to the galaxy center but lacks SF. The molecular gas mass inferred from the current SF rate (SFR) dominates over the atomic gas mass in the central ∼ 1.5 kpc. VCC 848 is consistent with being a main-sequence starforming galaxy for its current stellar mass and SFR. The H ii region luminosity distribution largely agrees with that of normal dwarf irregulars with similar luminosities, except that the brightest H ii region is extraordinarily luminous. Our N -body/hydrodynamical simulations imply that VCC 848 is a merger between a gas-dominated primary progenitor and a gas-bearing star-dominated secondary. The progenitors had their first passage on a near-radial non-coplanar orbit more than 1 Gyr ago. The merger did not build up a core as compact as typical compact dwarfs with centralized starburst, which may be partly ascribed to the star-dominated nature of the secondary, and in a general sense, a negative stellar feedback following intense starbursts triggered at early stages of the merger.

INTRODUCTION

In the standard Λ cold dark matter (ΛCDM) paradigm, dark matter halos, together with visible galaxies sitting close to their centers, assemble hierarchically under the influence of gravity. Galaxy mergers, especially those involving gas-rich galaxies of comparable masses (i.e. "wet" major mergers with primary-tosecondary mass ratios 4), can dramatically change galaxy morphologies (e.g. Toomre & Toomre 1972; Barnes & Hernquist 1991;Mihos & Hernquist 1994), enhance star formation activities (e.g. Sanders & Mirabel 1996;Mihos & Hernquist 1996;Duc et al. 1997;Zhang et al. 2010;Luo et al. 2014;Cibinel et al. 2019) and trigger active galactic nuclei (e.g. Ellison et al. 2011;Treister et al. 2012;Weston et al. 2017) in short timescales. Observational studies in the past decade have established that the galaxy merger rates increase steadily with redshift at least up to z ∼ 2-3 (e.g. Lotz et al. 2011;Lopez-Sanjuan et al. 2015;Mundy et al. 2017;Duncan et al. 2019;Ventou et al. 2019).

Although galaxy merging events in the local universe are not as common as in the early universe, nearby galaxy mergers provide unique laboratories for detailed investigations of the influence of merging process on galaxy evolution. A vast majority of studies of galaxy mergers so far have focused on relatively massive galaxies, while mergers between dwarf galaxies (M < 10 9 M ) received little attention until very recently, which is partly due to the limited depth of most existing surveys. Being close to the bottom of the hierarchical structure formation, dwarf galaxy mergers are expected to be more common than their massive counterparts. Recent cosmological zoom-in simulations suggest that about 10%-20% of Local Group satellite galaxies with M > 10 6 M have experienced a major merger event since z = 1 (Deason et al. 2014).

The first systematic study of gas-rich dwarf-dwarf interacting pairs was presented by Stierwalt et al. (2015), who found that most of their isolated dwarf pairs are as gas-rich as unpaired dwarfs of similar stellar masses, and that the current star formation rate (SFR) of dwarf pairs is enhanced by a factor of ∼ 2 on average at relatively small projected separations, similar to what has been found for massive galaxy pairs (e.g. Ellison et al. 2013;Silva et al. 2018;Pearson et al. 2019). Pearson et al. (2016) further showed that dwarfs involved in interacting pairs tend to have more extended atomic gas distribution than their unpaired analogues. By searching for low surface brightness merger signatures, Paudel et al. (2018a) compiled a catalog of 177 relatively low mass merger candidates (M < 10 10 M ), most of which turn out to be star-forming galaxies. Kado-Fong et al. (2020) found that 15%-20% of nearby dwarf galaxies with extreme starburst activities show signs of tidal debris. In addition to the above relatively systematic studies, a few more case studies of gas-rich dwarf pairs have also been carried out (Annibali et al. 2016;Privon et al. 2017;Paudel & Sengupta 2017b;Paudel et al. 2018b;Makarova et al. 2018;Johnston et al. 2019).

Previous studies of gas-rich dwarf pairs are biased toward systems at relatively early merging stages (e.g. Stierwalt et al. 2015). It remains unclear how dwarfdwarf merging events impact the overall star formation activities and re-shape the gaseous and stellar distributions of the merger remnants. Particularly, it has long been conjectured that many blue compact dwarf galaxies (BCD), which are characterized by having bluer colors, more intense star formation activities and unusually higher central surface brightness compared to ordinary dwarf irregular galaxies (e.g. Gil de Paz et al. 2003), might be formed through gas-rich dwarf-dwarf merging (e.g., van Zee et al. 1998;Noeske et al. 2001;Ostlin et al. 2001;Bekki 2008;Lelli et al. 2012Lelli et al. , 2014;;Watts & Bekki 2016). Other plausible formation mechanisms for the compact stellar distributions of BCDs include inspiraling of giant star formation clumps driven by dynamical friction (Elmegreen, Zhang & Hunter 2012) and central starburst sustained by pristine gas accretion from the local environs or even the cosmic web (e.g. Johnson et al. 2012;Lopez-Sanchez et al. 2012;Ashley et al. 2014;Verbeke et al. 2014). BCD galaxies are the closest local analogues to UV-luminous galaxies such as Lyman break galaxies and Lyα emitters detected at high redshift (e.g., Gawiser et al. 2007;Finkelstein et al. 2011;Shibuya et al. 2019). It is generally difficult to differentiate the above mentioned formation mechanisms for BCDs and other starburst dwarf galaxies, partly due to the often morphologically irregular appearance of starforming dwarf galaxies.

In a recent work (Zhang et al. 2020, hereafter Paper I), we report a discovery of remarkably extended stellar shells around the BCD galaxy VCC 848 located in the outskirts of the Virgo cluster, based on deep optical imaging data from the Next Generation Virgo Cluster Survey (NGVS; Ferrarese et al. 2012). This discovery confirms that VCC 848 is the remnant of a merger be-tween dwarf galaxies. VCC 848 is perhaps the clearest known example of a star-forming dwarf galaxy resulting from a gas-rich merger. Paper I shows that VCC 848 is likely formed by merging between two dwarf galaxies with comparable masses within a factor of a few, and the merging event has significantly enhanced the star formation activities in the past ∼ 1 Gyr. In the rest of the current paper, we present interferometric observations of H i gas in VCC 848, study the recent star formation activities and discuss the merging process by invoking numerical simulations. Some relevant properties of VCC 848, which are either from the literature or derived in this work, are summarized in Table 1. Throughout this paper, we adopt a distance of 16.5 Mpc for VCC 848 (Mei et al. 2007;Blakeslee et al. 2009) and use the Schlafly & Finkbeiner (2011) Galactic extinction map to correct the photometry.

OBSERVATIONS AND DATA REDUCTION

Optical and Far-infrared Images

Broadband u-, g-, i-, and z-band images of VCC 848 were obtained by NGVS with the MegaCam instrument on the Canada-France-Hawaii Telescope. With a dedicated data acquisition strategy and processing pipeline, the NGVS reaches a 2σ surface brightness limit of µ g 29 mag arcsec -2 . The processed NGVS images have a pixel scale of 0. 186. The full width at half maximums (FWHM) of the point spread function (PSF) are 0. 78, 0. 71, 0. 53, and 0. 62, respectively, for the u, g, i, and z passbands of the VCC 848 field. An in-depth analysis of these broadband images of VCC 848 has been presented in Paper I.

Narrow-band Hα imaging of VCC 848 was obtained with the Las Campanas 100-inch du Pont telescope by Gil de Paz et al. (2003), as part of a large imaging observational campaign for nearby BCDs. The total exposure time for VCC 848 was 45 minutes and the seeing FWHM was ∼ 1. 4 (pixel scale = 0. 26). Gil de Paz et al. (2003) did continuum subtraction of their narrow-band images by using the R broadband images. We correct the continuum-subtracted fluxes for the contamination of [N ii] lines by adopting a [N ii]λ6584/Hα line ratio of 0.069 determined by Vilchez et al. (2003).

Archival far-UV (FUV) image from the Galaxy Evolution Explorer (GALEX) All-sky Imaging survey is available for VCC 848. With an exposure time of 112s, the FUV emission is detected only in the central high surface brightness region of VCC 848, and the integrated FUV magnitude is 16.9 ± 0.05 (Voyer et al. 2014). We do not use the shallow FUV image in the work, except that we calculate an integrated SFR by combining the FUV (corrected for the Galactic extinction) and infrared (see below) photometry (Table 1).

Herschel Photodetector Array Camera and Spectrometer (PACS) 100 µm, 160 µm and Spectral and Photometric Imaging Receiver (SPIRE) 250 µm, 350 µm, 500 µm observations of VCC 848 were presented by Grossi et al. (2015), as part of the Herschel Virgo Cluster Survey (Auld et al. 2013). Among these far-infrared (FIR) images, the PACS 100 µm has the highest spatial resolution (FWHM = 9. 4), and more importantly has been shown to be the most reliable monochromatic estimator of total IR (TIR) luminosity among all of the Spitzer or Herschel wavebands, with a scatter of 0.05 dex (Galametz et al. 2013). We will therefore use the 100 µm image to evaluate the contribution of obscured star formation to the total star formation budget at local scales of VCC 848 .

JVLA H i Emission Line Mapping

Observations

VCC 848 was observed in B, C, and D array configurations of the Jansky Very Large Array (JVLA; Perley et al. 2011) radio interferometer. The observations were performed in the 1-2 GHz L-band (primary beam 30 ) through multiple scheduling blocks (SBs) between 2016 May and 2017 Apr (Project ID: 16A-074; PI: Zhang; Table 2). The JVLA WIDAR correlator was configured with two 1-GHz baseband IF pairs (A0/C0 and B0/D0) and the 8-bit sampler. The correlator integration time is 3 seconds. Each baseband was divided into 22 sub-bands with dual polarization products (RR and LL) each. Only one 4 MHz wide sub-band of the A0/C0 baseband was used for H i emission line observations presented in this paper. The 4 MHz sub-band was tuned to a central frequency corresponding to the systemic velocity of VCC 848. By using a recirculation factor of 4, the 4-MHz wide sub-band consists of 512 spectral channels, with a channel width of 7.81kHz (∼1.7 km s -1 for a rest frequency of 1420.4 MHz).

Each SB of our observations began and ended with a 10-to 15-minute scan on the primary flux calibrator 3C286. The primary flux calibrator was used to calibrate the flux scale and bandpass. In between the observations of the primary flux calibrator, the phase calibrator J1254+1141 was observed for 5-6 minutes before and after each scan (10-25 minutes) on our target galaxy. VCC 848 was observed for a total on-source time of 4.0 hrs in the B array, 2.0 hrs in the C array, and 0.9 hrs in the D array (Table 2).

Calibration

The multi-configuration data sets from the spectral line sub-band were calibrated and imaged using the Mei et al. 2007;Blakeslee et al. 2009; (3) Vilchez et al. 2003;Lee et al. 2003;Zhao et al. 2013;(4) Before combining the calibrated multi-configuration MSs using the CASA task concat, the (relative) weight of visibilities in each calibrated MSs was first determined based on line-free channels (statwt) and then the task cvel2 was used to transform channel coordinates to a common heliocentric velocity reference frame with a 1.7 km s -1 channel width. Spectral continuum was subtracted from the combined MS (uvcontsub) with a first-order polynomial fit to visibilities in the line-free channels.

The continuum subtracted MS was deconvolved using tclean, with a cell size of 1. 5 × 1. 5, a robust parameter of 0.7, and a circular clean region of 4 in radius around the target. A robust parameter of 0.7 leads to a restoring beam size of 7. 1 × 6. 6 (PA = 56.78 • ) and rms noise of 0.64 mJy beam -1 . The multi-scale clean algorithm was used for recovering emission structures on different spatial scales. After several tests, a set of scale sizes of 6 , 12 , 24 and 36 and a "smallscalebias" parameter of 0.7 were adopted for cleaning, down to a global stopping threshold of 2 times the rms noise. The residual maps were checked to make sure all the flux has been retrieved.

The H i data cube was blanked in order to separate features of genuine emission from pure noise, following similar procedures adopted in previous studies (e.g. Walter et al. 2008;Hunter et al. 2012). In particular, the original data cube was first smoothed to a spatial resolution of 15 ×15 , then pixels that belong to spatial features appearing in at least 3 consecutive channels and above 2× the rms noise of the smoothed cube were masked as area of real emission (imager.drawmask). The blanking mask thus generated was applied to the original data cube (immath). The blanked data cube at the original cm -2 beam -1 ch -1 ) 1.5

resolution will be used for the following analysis in this paper.

Comparison with ALFALFA Observation

The spatially integrated JVLA H i line profile of VCC 848 is plotted together with that from the single-dish Arecibo Legacy Fast ALfa Survey (ALFALFA; Haynes et al. 2011) observation in Figure 1. The JVLA line profile has been smoothed to the 11 km s -1 spectral resolution of ALFALFA data. The rms noise σ rms of the ALFALFA spectrum is about 2.6 mJy. The two profiles are in good agreement, except that the JVLA line profile is slightly (< 2σ rms ) lower than the ALFALFA profile on the high-velocity side. It is worth noting that the JVLA H i detection covers virtually the same velocity range as the ALFALFA detection. The total H i flux density measured on the JVLA line profile is 6.30 ± 0.06 Jy km s -1 , accounting for ∼ 94% of that (6.68 ± 0.09 Jy km s -1 ) from the ALFALFA observation.

H i GAS DISTRIBUTIONS

In contrast to the stellar component, the gaseous component of galaxies is collisional and can dissipate its orbital energy via shocks. This may lead to diverged dynamical evolution of gas and stars during a violent merging process. The H i 21-cm emission has proven to provide indispensable clues to the evolutionary status of merging systems (e.g. Hibbard et al. 2001).

Integrated H i Line Profile

The ALFALFA H i velocity profile of VCC 848 has a ∼ 160 km s -1 velocity width 1 measured at 20% of the peak intensity W 20 (Haynes et al. 2011). This velocity width is larger than that (∼ 130 km s -1 ) expected for its baryonic mass (approximated as M +1.36M HI = 7.4×10 M ), based on the baryonic Tully-Fisher relation (e.g., McGaugh 2012). The H i velocity profile is asymmetric with respect to the systemic velocity (1532 km s -1 ). As we will show later, the approaching half (V < 1532 km s -1 ) of the H i distribution follows a more regular rotation pattern than the receding half (V > 1532 km s -1 ). So we mirror the H i profile of the approaching half to the receding side, and compare the mirrored profile with the observed profile at the receding side in Figure 1. It can be seen that the receding side shows a deficit of H i at V 1570 km s -1 whereas an excess of H i at V 1570 km s -1 , indicating significant tidal forces that disturbed the gas velocity distribution. W 20 measured on the mirrored line profile is 140 km s -1 , which is close to that predicted by the Tully-Fisher relation.

H i Maps and Comparison with the Stellar Light Distribution

Figure 2 shows the column density map, intensityweighted velocity field and velocity dispersion maps of VCC 848. To help visualise the spatial association between stellar light and H i gas, the g-band surface brightness contours are overlaid on the H i maps in the upper row, while the H i distributions are contoured on the g-band image in the lower row. The g band image has been adaptively smoothed to achieve a minimum signalto-noise ratio of 5 pixel -1 (Paper I).

The H i gas spatial distribution is asymmetric with respect to the photometric center of the stellar main body, and it appears to be confined within the inner edges of the two outermost stellar shells. Moreover, the H i distribution is bounded on the southeast by a near-circular outer edge of the stellar tidal arm that wraps nearly 180 • around the galaxy center (Figure 2 in Paper I). 2 Clearcut edges of the H i distribution in the southeastern half are manifested by the closely crowded H i intensity contour lines. H i gas distribution in the northwestern half of VCC 848 extends slightly further north than the detected stellar light distribution and has typical N (H i) 10 20 cm -2 (0.8 M pc -2 ). High column density H i gas (N (H i) 10 21 cm -2 ), corresponding to the yellowish contours in the lower left panel of Figure 2), is mainly concentrated in two distinct spatial areas. One is broadly associated with the central (twisted) stellar disk of VCC 848 and accounts for 15% of the total H i flux of the system, while the other one is associated with the above mentioned stellar tidal arm and accounts for 13% of the total H i flux. We point out that the general association between the H i and stellar features, together with the gas-dominated nature of VCC 848, suggest that ram pressure stripping, if any, does not have an appreciable influence on the H i morphology.

The overall H i velocity gradient is approximately along the direction of the major axis of the stellar main body, but there exist remarkable local deviations from this overall trend. In particular, while the velocity field at the approaching side (V los < 1532 km s -1 ) is largely consistent with a solid-body rotation (see below), velocity contours of the receding side are closely crowded (i.e. steep velocity gradient) and twisted towards the southwest edge of the stellar main body, giving rise to the apparent high line-of-sight velocity dispersion (up to ∼ 25 km s -1 ) of the H i gas there. The velocity contours at the receding half is consistent with being affected by tidal forces largely along the plane of sky. The above mentioned high column density H i associated with the tidal arm is mainly at the receding side, where the arm appears to be launched from.

By performing multi-gaussian line profile fitting to the H i data cube of VCC 848 (see Section 4), we find that, except for sporadic locations with high column densities 2 The spatial association of atomic gas with the stream-like feature reported in Paper I suggests that this feature is a tidal arm emanating from the primary progenitor rather than a stellar shell or stream stripped from the secondary progenitor, because the H i gas of the secondary, if any, should have been completely detached from the secondary's stellar component and assimilated by the primary after the first two to three periapsis passages, as will be demonstrated by computer modeling in Section 6.

ROTATION CURVES AND MASS PROFILES

Derivation of Rotation Curves

We use the 2D Bayesian Automated Tilted-ring fitter (2dbat) recently developed by Oh et al. (2018) to extract rotation curves from the H i velocity field. Instead of using the first moment map which could be affected by the disturbed motions from the merging process, we extract a bulk velocity field for the underlying circular rotation of the galaxy in an iterative manner. In particular, we perform profile decomposition analysis of the individual velocity profiles of the H i data cubes as described in Oh et al. (2019) and extract the decomposed components whose central velocities are close to a reference velocity field (e.g., single Gaussian velocity field). We call this a bulk velocity field, and classify the other ones deviating from the bulk motions as disturbed, noncircular motions. We then improve the reference velocity field for profile decomposition with a new model velocity field from the 2dbat analysis of the first reference velocity field and repeat the profile decomposition described above. In this way, we extract the bulk velocity field of VCC 848 for the rotation curve analysis. Next, we create a mask for the bulk velocity field to isolate a more or less elliptical region encompassing the main body of VCC 848 (left panel of Figure 3). This quasi-elliptical region has a semi-major axis length of R maj 30 and accounts for ∼ 30% of the total H i flux of the system. The approaching half of this region largely follows the characteristic spider-shaped velocity field expected for an inclined rotating disk.

We fit the velocity field of the masked quasi-elliptical region with 2dbat by fixing the center position to the one obtained from our stellar isophotal analysis (Table 1; Paper I). The systemic velocity V sys , kinematic PA and inclination angle are fitted as free parameters but are kept constant with radius. 2dbat iteratively searches for the tilted-ring parameters that give the best fit to the observed velocity field. The rotation velocities (filled circles in the right panel of Figure 3) are extracted with a 7. 5 ring width that is chosen to be slightly larger than the beam size. The best-fit kinematic parameters are indicated in the caption of Figure 3. Lastly, we apply asymmetric drift correction to the extracted rotation velocities, following the method described in Bureau & Carignan (2002) (see also Oh et al. 2011). As is shown in Figure 3, the asymmetric drift correction is < 2 km s -1 across the probed radial range.

Rotation Curves

Rotation velocities extracted from the approaching side increase nearly linearly with radius at R maj 20 , beyond which the rotation curve appears to flatten out. At the receding side, where the velocity field is substantially disturbed by the ongoing merging, the derived circular velocities are about half that of the approaching side at R maj 25 , beyond which the radial trend steepens and reaches velocities comparable to that of the approaching side at R maj ∼ 30 . Projected radial velocities derived based on the extracted rotation velocities are over-plotted on the H i position-velocity diagrams in Figure 4, where the largely regular rotation pattern in the inner 30 is clearly illustrated.

In order to quantify the radial trend of rotation velocities extracted from the approaching side, we parameterize the rotation curve with the same functional form as that used in Barrera-Ballesteros et al. (2018) (see their Equation 4). The functional form is defined by three parameters: the maximum velocity (V max ), the transition radius (R turn ) where the rotation curve transitions from a solid-body inner part to a flat outer part and the sharpness (α) of the transition. We obtain a best-fit V max = 41.5±0.5 km s -1 , R turn = 19.6 ±0.9 (∼ 1.6 kpc) and α = 130.1±1.0. The large α value indicates a sharp transition at R turn . The corresponding best-fit curve is over-plotted in the right panel of Figure 3.

Based on the above parametrized fitting, we obtain a slope of 26.5±3.9 km s -1 kpc -1 (V max /R turn ) for the solid-body rising part of the rotation curve, where the uncertainty is determined based on repeat fitting to randomly disturbed rotation velocities for 1000 times. This slope of rotation curve falls within the range of ordinary dwarf irregular galaxies or BCDs with off-centered starburst that have inclination-corrected central brightness (µ g,corr = 22.6 mag arcsec -2 for a kinematic inclination = 45.2 • ) comparable to VCC 848 (Lelli et al. 2014). Because the rotation velocity gradient slope is proportional to the square root of local matter density for a solid-body rotation curve, we can infer that the merging event of VCC 848 has not built up an exceptionally compact core.

Dynamical and Baryonic Mass Profiles

We use the rotation curve of the approaching half to determine the cumulative dynamical mass profiles (= V rot (R) 2 R/G), which is shown in Figure 5. We also overplot mass profiles of the baryonic components, including stars, atomic and molecular gas in Figure 5. All of the profiles are extracted using the same geometric parameters as determined from the 2dbat kinematic fitting. The stellar mass is estimated by using the (gi)-massto-light relation calibrated for Local Group dwarf galaxies (Zhang et al. 2017). Details about the procedure of stellar mass estimation are given in Paper I. We adopt a logarithmic stellar mass uncertainty of 0.2 dex which is applicable to color-based mass estimation for low-dustextinction cases (Zhang et al. 2017). The atomic gas mass is equal to the H i gas mass multiplied by a factor of 1.36 to account for helium. The molecular gas mass is indirectly inferred from the SFR estimated with Hα luminosities (Catalan-Torrecilla et al. 2015, Equation 10), by adopting a constant molecular gas consumption timescale of (5.25 ± 2.5) × 10 -10 yr -1 found for nearby disk galaxies at (sub-)kpc scales (Bigiel et al. 2008;Leroy et al. 2008). 3 As is illustrated in Figure 5, the stellar and atomic gas components together account for 20% of the dynamical mass enclosed within the central 30 , and the fraction increases to ∼ 24% if including molecular gas in the baryon budget. Stars constitute more than half (54%) of the baryonic mass, which is in contrast to the overwhelmingly gas-rich nature of the system as a whole.

3 Another method of estimating the molecular gas mass is based on a metallicity-dependent gas-to-dust ratio (G/D). According to the G/D-metallicity relation from Rémy-Ruyer et al. ( 2014), VCC 848 is expected to have a log(G/D) of 3.00 ± 0.37 for 12+ log(O/H) = 8.03 (Table 1). In addition, Grossi et al. (2015) determined a total dust mass of 10 5.01 M (the value has been adjusted to reflect our adopted distance) within the central ∼ 30 of VCC 848. Therefore, the total gas mass (atomic plus molecular) within the central 30 is estimated to be 10 8.01±0.37 M and hence the most likely molecular gas mass be ∼ 10 7.5 M by subtracting the atomic gas contribution (Figure 5). This molecular gas mass is in reasonable agreement with that (10 7.6 M ) inferred from the SFR-based method shown in Figure 5. We choose to use the SFR-based method in this paper. This is because, besides their low spatial resolution, the Herschel images of VCC 848 generally have signal-to-noise ratios that are too low (except for the brightest star-forming sites) to allow for constructing a robust radial profile of dust mass densities.

This might be partly attributed to a highly efficient star formation toward the stellar main body and a less severe tidal stripping of the stellar component than the gas during the merging. Lastly, it is also notable that the molecular gas appears to dominate over the atomic gas in mass within the central ∼16 -20 (1.3-1.6 kpc; accounting for the uncertainties of molecular gas mass estimate) in radius. Similarly high fraction of molecular gas was also indirectly inferred for some nearby BCD or starburst dwarf galaxies (Hunter et al. 2019).

STAR FORMATION ACTIVITIES

Specific SFR and Star Formation Efficiencies

The global logarithmic specific SFR (SFR per unit stellar mass, sSFR) of VCC 848 is -9.87 (see Table 1 for the total SFR and stellar mass), which is 0.19 dex higher than that expected for normal star-forming galaxies of similar stellar masses (Shin et al. 2019, Equation 19). However, given a 0.27 dex scatter of the Shin et al. (2019) relation, VCC 848 is consistent with being a main sequence star-forming galaxy for its stellar mass. The (gi) and Hα maps of the central 1.5×1.5 arcmin are shown in Figure 6. The Herschel 100 µm continuum emission, which traces the dust obscured star formation, is contoured in red color on the Hα map. A vast (Figure 3). The four pseudo-slits (centered at the galaxy center) are indicated as long magenta lines overlaid on the H i column density map (left panel). Every 20 along the pseudo-slit direction is marked with a 7 -long line perpendicular to the slit direction. Also overlaid on the column density map are the intensity-weighted H i velocity contours, with 15 km s -1 intervals.

In the right panels, the dot-connected red curves are the projected radial velocities based on rotation curves determined for the central 30 of VCC 848 (Figure 3).

majority ( 90%) of current star formation traced by

Hα and 100 µm is confined to the central high surface brightness regions enclosed by the N (H i) = 10 21 cm -2 contour. Although the tidal arm in the southeast reaches similar N (H i) to the central region, the SFR there accounts for merely 6% of the total. The most intense star-forming site of the system is located ∼ 10 offset toward the northeast of the photometric center, and is associated with the brightest star cluster (Paper I) and H ii region (Section 5.2) in VCC 848.

It is remarkable that the youngest stellar populations in the central ∼ 10 (in radius) are confined to a straight narrow strip that appears to trace the symmetry axis of the twisted disk. The age spread of star clusters detected along this narrow strip suggest that it has existed for at least ∼ 1 Gyr (Paper I). We note that this linear stellar structure is not likely a bar, because the associated H i velocity field does not show the characteristic S-shaped isovelocity contours expected for a bar structure (e.g., Athanassoula 1992).

To further illustrate the relative strength of star formation activities across the system, we show spatial distributions of sSFR and star formation efficiency (SFE) with respect to the atomic plus molecular gas (see Section 4 for the inference of molecular gas mass) in Figure 7, where the involved input images, including the g, i, Hα and H i, have been smoothed to match the beam size of H i map and rebinned to 7 ×7 (560×560 pc) pixel sizes. The SFR is estimated from the Hα luminosities. Note that we have ignored the contribution of dust obscured star formation in Figure 7 for a high resolution view of the overall distribution of star formation. This practice is justified by the fact that the TIR inferred from the 100 µm (Galametz et al. 2013) accounts for at most ∼ 10% of the local SFR at spatial scales comparable to the resolution of 100 µm, according to the hybrid Hα+FIR recipe for SFR estimation presented in Catalan-Torrecilla et al. (2015). This is probably not surprising given the low metallicities of VCC 848.

As is shown in Figure 7, the sSFR and SFE reach their logarithmic maximum of -9.1 yr -1 and -9.4 yr -1 respectively at the most intense (and also bluest) star-forming site of the system. The H i tidal arm has typical sSFR and SFE that are more than one order of magnitude lower than the central high surface brightness regions. In Figure 8, the SFR surface densities of individual pixels shown in Figure 7 are plotted as a function of the H i mass surface densities. Above a detection limit of 10 -4.6 M yr -1 kpc -2 for Σ SFR(Hα) , the individual regions fall into two groups which are separated by a ∼ 0.3 dex gap in log Σ SFR(Hα) from -3.6 to -3.3. The regions in the higher-Σ SFR(Hα) group belong to the central high surface brightness regions, whereas those in the lower-Σ SFR(Hα) group and with logΣ HI 0.88 be- Figure 5. Accumulated radial profiles of the dynamical mass as well as the stellar, atomic and molecular gas masses within the central 30 . The same geometric parameters derived from the 2dbat fitting of the H i velocity field (Figure 3) are used for extracting all of the profiles. The atomic gas mass is calculated as the H i gas mass multiplied by a factor of 1.36 to account for the contribution of helium. Estimation of the molecular gas mass is described in Section 4.3

The gray-shaded region represents the uncertainties of the estimation of molecular gas mass.

long to the H i concentration along the tidal arm (Figure 2). The lack of a direct connection between star formation and H i gas and the apparent relevance of stellar surface densities are vividly demonstrated in VCC 848 (see, e.g., Ekta et al. 2008 for similar examples).

Hα Luminosity Distribution of H ii Regions

Detection of H ii Regions

H ii regions trace massive stars ( 10 M ) formed in the recent (∼ 10-20 Myr) past, and their spatial and luminosity distributions reflect the mode and strength of star-forming activities. We use the Sextractor software (Bertin & Arnouts 1996) to detect H ii regions based on the Hα image of VCC 848. With a detection threshold of 2.5 times the rms noise, we detect 23 H ii regions in total. To determine the detection completeness, we repeatedly add artificial point sources (one at a time) that are convolved with a 1.8 FWHM Gaussian (Section 5.2.2) and span a uniform range of flux densities at 1000 randomly selected positions within the central main body of VCC 848 (delineated by the dotted ellipse overplotted on the Hα image in Figure 9), and then run Sextractor for detection with the same parameters as for the original image. This test suggests that our H ii region detection reaches ≥ 90% recovery rate (our completeness limit) at log L(Hα) [erg s -1 ] 37.7. Among the 23 H ii regions detected in VCC 848, 10 have L(Hα) above the 90% limit. The following analysis in this section will focus on the 10 brightest H ii regions.

Aperture Photometry of H ii Regions

An inspection of the radial light distributions of the 10 brightest H ii regions suggests that they have FWHMs ranging from 1.7 to 1.9 and nearly all of them are reasonably described by single Gaussian profiles 4 . We perform photometry for each region with a 2.1 (8 pixels) diameter circular aperture and subtract a local median background using an 8-pixel wide annulus that is sufficiently far away from the center of the H ii region in question. The aperture size is chosen to be larger than the FWHM but at the same time to avoid significant contamination from neighboring regions. Pixels belonging to neighboring H ii regions are excluded when determining local background and noise. For the two H ii regions close to the brightest H ii region (left panel of Figure 9), we manually select background regions that are at the same distance from the brightest H ii region but in different directions. The total Hα flux of each region, as expected for a Gaussian profile with FWHM of 1.8 , is obtained by applying a multiplicative aperture 4 The similar radial profiles might not be very surprising given that our H ii regions are only marginally resolved in the Hα image and that the intrinsic light profiles of H ii regions are found to have nearly invariant gradients at log L(Hα) [erg s -1 ] 38.6 (Rozas et al. 1998). The same H i column density contours shown in Figure 2 are overlaid in both panels here. The Herschel 100µm is overlaid as red contours (2σ, 3σ and 4σ above the sky background) on the Hα flux density map. Note that the dust IR emission contributes at most ∼ 10% of the total SFR budget (Hα+TIR) across the system at spatial scales comparable to the resolution (9. 4) of the 100 µm image.

correction factor of 1.77 to the background-subtracted aperture photometry. The derived Hα flux densities and uncertainties of the 10 H ii regions are given in Table 4. We note that the marginally-resolved nature of our detections means that more sophisticated procedures for H ii region photometry such as HIIphot (Thilker et al. 2000) cannot play their advantages.

H ii Region Distributions

The H ii region distributions are shown in Figure 9, where the left panel marks the spatial distribution while the right panel presents the cumulative Hα luminosity distribution above the 90% completeness limit. No binning in luminosity is used for constructing the cumulative distribution. The variation of detection completeness with luminosities are also overplotted in the right panel. To provide some perspective on the H ii regions detected in VCC 848, the 90% completeness limit of L(Hα) is ∼ 5 times that of the Orion nebula (Kennicutt 1984), and for a Case B recombination, this corresponds to a Lyman continuum photon emission rate of 10 49.6 s -1 , which is equivalent to ∼ 13 O9V stars (Sternberg et al. 2003). The brightest H ii region is associated with the brightest star cluster in VCC 848 (Paper I) and has a Hα luminosity (10 38.94 erg s -1 ) slightly higher than the brightest H ii region N66 (10 38.78 erg s -1 , Kennicutt 1984) in the Small Magellanic Cloud (SMC). We note that the Hα luminosity of the brightest H ii region may be excited by the associated (and the brightest) star cluster which has an age of ∼ 6 -8 Myr and mass of 10 5 M (Paper I).

We perform a least-squares fitting to the completeness-corrected cumulative Hα luminosity distribution of VCC 848 by adopting a power-law form for the (probability) luminosity distribution (dN (L)/dL ∝ L -α ) that has been widely used in the literature. We find a best-fit power-law index α = 2.6 +0.7 -0.5 , where the uncertainties correspond to the range of α with reduced χ 2 ≤ χ 2 min +1. As is illustrated in Figure 9, the luminosity distribution is well described by a single power-law form, except for the brightest H ii region which is ∼ 0.2 dex more luminous than that expected for the power-law fitting.

Considering the uncertainties, the power-law index determined above for VCC 848 falls in the broad range of values found for nearby spiral galaxies (∼ 2.0 ± 0.5; e.g. ), but appears to be at the upper bound for dIrr galaxies (∼ 1.5 ± 0.5; e.g. Strobel et al. 1991;Miller & Hodge 1994;Youngblood & Hunter 1999;van Zee 2000). However, it has been shown that measurements of α are in many cases dependent on the L(Hα) range used in the power-law fitting, in the sense that α measured at higher luminosities tends to be similar to or larger than that at lower luminosities (e.g., Youngblood & Hunter 1999). The H ii region detection limit for VCC 848 is relatively high compared to most previous studies of dIrr galaxies in the literature, so the large α measured for VCC 848 may not be surprising. We shall return to this point when comparing VCC 848 with ordinary dIrr galaxies in the next section.

Comparison with Isolated dIrr Galaxies

To probe the impact of galaxy merging, we compare the Hα luminosity distribution of VCC 848 to that of ordinary late-type dwarf galaxies with similar luminosities. Youngblood & Hunter (1999) and van Zee (2000) studied H ii region luminosity distributions for the largest samples of relatively isolated gas-rich dIrr galaxies so far. In particular, the samples of Youngblood & Hunter (1999) and van Zee ( 2000) include 4 and 6 dIrr galaxies, respectively, that have M B within ± 0.5 mag of VCC 848 and at the same time have H ii regions with maximum log L(Hα) ≥ 37.7. We note that all of the 4 galaxies from Youngblood & Hunter (1999) happen to be fainter than VCC 848 while all of the 6 from van Zee (2000) are brighter than VCC 848. These 10 galaxies constitute our comparison sample. Cumulative Hα luminosity distributions of the 10 galaxies are shown in Figure 9. In addition, we construct a composite luminosity distribution by co-adding the H ii region catalogs of the 10 galaxies at log L(Hα) ≥ 37.7, and show the resultant cumulative luminosity distribution in Figure 9. For the sake of comparison, the co-added distribution is normalized such that the total number of H ii regions at log L(Hα) ≥ 37.7 is equal to that of VCC 848.

The cumulative Hα luminosity distribution of VCC 848 is in reasonable agreement with the co-added distribution of the comparison sample. This general agreement is remarkable and may imply that the galaxy merging does not significantly affect the birth mass distribution of star clusters, at least at the current stage of the merger and for the most massive ones. We note that, once the full catalogs of H ii regions down to their re- 6 and7, and are color-coded according to their stellar mass surface densities. Note that only data points that fall within the 5×10 20 cm -2 contours (Figure 2) are plotted. The diagonal dotted lines represent lines of constant SFE(H i) and thus constant H i gas depletion (by star formation after transitioning to the molecular phase) times τHI (10 8 , 10 9 and 10 10 yr). Above the detection limit of SFR(Hα), the data points fall into two groups which are separated by a 0.3 dex gap in log Σ SFR(Hα) from -3.6 to -3.3. Data points of the group with higher log Σ SFR(Hα) belong to the central high surface brightness star-forming regions, whereas those of the lower-Σ SFR(Hα) group with log ΣHI 0.88 are from the southeast H i tidal arm (Figure 2). spective detection limit are used for the power-law fitting, the galaxies in the comparison sample have a median α = 1.5 ± 0.5, which is in fair agreement with the typical values found for dIrr galaxies and suggests that the Hα luminosity distributions tend to flatten toward fainter luminosities.

Besides exploring the overall shape of the Hα luminosity distribution above the detection limit, we can also compare the fractional luminosity contribution of the detected H ii regions in VCC 848 and the comparison sample. In particular, at log L(Hα) ≥ 37.7, the H ii regions collectively contribute 47% of the total Hα emission of VCC 848, which is higher than any of the galaxies in the comparison sample (∼ 15% ± 10%). The extraordinary high fractional contribution in VCC 848 is largely driven by the single brightest H ii region that contributes 22% of the total Hα emission. As has been shown in Paper I, the brightest star cluster, which is associated with the brightest H ii region in VCC 848, is also exceptionally luminous for the current SFR of VCC 848.

NUMERICAL SIMULATIONS

We use hybrid N -body/hydrodynamical simulations to gain insight into the merging process of VCC 848. In particular, we make use of a Treecode-smoothed particle hydrodynamics (SPH) algorithm that largely follows the techniques described in Hernquist & Katz (1989), where the Treecode allows rapid calculation of gravitational accelerations and the SPH code allows for modeling the collisional gas component. The relevant Treecode-SPH parameters used in our simulations, such as the gravitational softening length, time-step scheme, Treecode opening angle and viscosity prescription, are the same as that adopted in Smith, Davies & Nelson (2010). This Treecode-SPH code has been extensively used for modeling galaxy interactions and environmental effects acting on dwarf galaxies (e.g., Smith et al. 2012aSmith et al. ,b, 2013aSmith et al. ,b, 2015)).

We build multi-component models of a pair of otherwise isolated dwarf galaxies using the publicly available code DICE [START_REF] Perret | DICE: Disk Initial Conditions Environment v4[END_REF]. In particular, each galaxy consists of a spherical dark matter halo with a NFW density profile (Navarro et al. 1996), an exponential stellar disk and, if there is gas, an exponential gas disk. The dark matter halo component consists of 400,000 particles, while the stellar and gas components consist of 50,000 particles each. The gas component is treated isothermally, with a constant sound speed of 7.5 km s -1 (roughly equivalent to a temperature of 10 4 K). This has been shown to be a valid approximation of a rigorous treatment of warm interstellar atomic gas that has a highly efficient (inefficient) cooling efficiency at T > (≤) 10 4 K (e.g., Barnes 2002). The gas-to-stellar disk scalelength ratio is fixed to 1.5, which can be compared to the typical ratio of ∼ 1.3 for nearby star-forming dwarf galaxies (e.g., Hunter et al. 2012). The scale height is 0.3 times the scale length for both the stellar and gas disks. Lastly, the dark matter halos have a concentration parameter c = 10 and their masses (M h ) are linked to M following the relation derived from abundance matching technique (e.g., Guo et al. 2010).

Given the above simulation ingredients, there is a large parameter space that involves the progenitor mass ratios, gas fractions, disk sizes, orbit geometries, etc. It is beyond the scope of the current work to conduct a complete search for all possible combinations of these parameters to match every detail of the observations. As elaborated below, our choice of the parameters is largely guided by the observed stellar distribution of VCC 848. The primary goal of our simulation is to infer the merger In the left panel, the green and red circular apertures mark the H ii regions brighter and fainter than the 90% completeness limit (> 10 37.7 erg s -1 ), respectively. The big dotted ellipse encloses the region used for completeness estimation. In the right panel, the cumulative Hα luminosity distribution of VCC 848 above the 90% completeness limit is shown as red filled circles while the completeness-corrected distribution as red open circles. The red solid curve corresponds to a power-law fit to the completeness-corrected Hα luminosity distribution of VCC 848, and the red shaded region represents the uncertainties of the power-law fit. The green dotted curve delineates the completeness limit (marked in the right y-axis) as a function of Hα luminosities, and the gray shaded region represents the regime below the 90% completeness limit. Cumulative Hα luminosity distributions for 10 relatively isolated dIrr galaxies whose MB are within ± 0.5 mag of VCC 848 are plotted as gray curves. The black dashed curve corresponds to the normalized co-added distribution of the 10 dIrr galaxies.

stage by qualitatively reproducing the most prominent stellar features. While our simulations are not meant to match the observed H i gas distribution of VCC 848, we will try to discuss the response of the gas component and its impact on the merging process that may generally apply to gas-dominated mergers of dwarf galaxies. 6.1. Parameter Setup 6.1.1. Stellar Mass and Size of Progenitor Galaxies An upper limit of the primary-to-secondary progenitor stellar mass ratio ( 5) of VCC 848 has been estimated in Paper I. This upper limit was obtained by ignoring the contribution of the secondary progenitor to the central exponentially declining part of the stellar light distribution of VCC 848, and the mass estimate for the secondary only takes into account of the stellar light excess (mainly) caused by stellar shells towards large radii. Since the outlying stellar shells of the merger remnant consist of the least bounded stars stripped from the outskirts of the secondary progenitor (e.g., Quinn 1984), it is plausible that at least a similar fraction of stars that initially inhabited smaller radii of the secondary has been deposited into the inner part of the remnant. Therefore, the true primary-to-secondary stellar mass ratio may be something close to 2 (i.e., 5-1:1+1), which is adopted in our simulations, with M ,primary = 1.4×10 8 M and M ,secondary = 0.7×10 8 M . The corresponding dark matter halo mass ratio is 1.3, with M h,primary = 5.3×10 10 M and M h,secondary = 4.2×10 10 M . We note that the stellar to dark matter halo mass ratio decreases steeply with decreasing stellar (or halo) mass towards the low mass end and the dark matter halo accounts for ≥ 99% of the total mass budget (and hence the gravitational force) of each galaxy in our simulations. This means that the halo mass ratio and merging timescale are not very sensitive to the exact choice of stellar mass ratio. For example, an increase of the stellar mass ratio by a factor of 2 leads to merely a ∼ 20% increase in the halo mass ratio. Lastly, the exponential stellar disk scalelength of the simulated primary progenitor is set to 1.0 kpc while the scalelength of the secondary is set to be 0.8 kpc, in general accord with the luminosityscalelength relation followed by nearby dwarf irregular galaxies (e.g., Hunter et al. 2006).

Orbit Geometry

The setup for the orbit geometry of the merging pair is guided by the observed stellar light distribution of VCC 848 (Paper I). Firstly, the remarkably extended outlying stellar shells that are largely aligned along the east-west direction implies a nearly radial encounter that is largely perpendicular to the line of sight direction. Secondly, the stellar main body of VCC 848 is twisted by ∼ 10

• from the center to larger radii, which implies a moderately off-center collision. Thirdly, the major axis of the stellar main body is highly inclined with respect to the direction of shell alignment, which implies a non-coplanar collision. Given these observational hints, we adopt the following initial orbit geometry and disk orientation in our simulations. The two progenitor galaxies are initially 50 kpc apart, and they approach each other with a relative velocity of 50 km s -1 along the collision direction and a sideways velocity of 10 km s -1 . Regarding the relative orientation of the galaxy disks, we ran several tests and achieve the best match to the observations when the disk plane of the primary (secondary) progenitor is inclined at 60

• (90 • ) with respect to the orbital plane and the disk planes of both progenitors are parallel with the sideways direction. An illustration of the 3D view of the star particle distribution at a pre-encounter moment is shown in Figures 10. We note that the relative approaching velocity is chosen to be smaller than the primary's escape velocity (∼ 77 km s -1 at a virial radius of ∼ 75 kpc). A larger approaching velocity delays the instant when the coalescence occurs and shells start to form, but makes little difference to the number of passages after coalescence nor to the final appearance. In addition, after testing with several different sideways velocities, we find that, while a sideways velocity is needed to make off-center collision and asymmetrical distribution of shells along the collision direction, the choice of its exact value is not crucial, due to a rapid orbital decay of the secondary caused by dynamical friction.

Results of the Simulations

We invoke two representative simulations to discuss the progenitor properties and the merging process of VCC 848. Both simulations adopt a 2:1 gas-to-stellar mass ratio for the primary progenitor, but they differ in the gas content of the secondary progenitor. Specifically, the secondary progenitor is gas-free in the first simulation (hereafter S1) but has a 1:1 gas-to-stellar mass ratio in the second simulation (hereafter S2). As we will show, it is not possible to constrain the exact gas richness of the secondary, but we can use a 1:1 ratio in the S2 simulation to demonstrate whether the secondary is expected to be star or gas dominated. Figures 11 and13 show the projected spatial distributions of star/gas particles at selected snapshots of the S1 and S2 simulations, respec-tively. In the following subsections, we will first invoke the S1 simulation to discuss the merging stage of VCC 848 and then the S2 simulation to discuss the impact of a gas-bearing secondary progenitor on the merger.

Reproducing the Stellar Shells of VCC 848

The most prominent tidal features in VCC 848 are the three outlying interleaved stellar shells, among which the innermost one has significantly higher stellar surface densities, sharper outer edges and smaller opening angle than the other two. Based on the S1 simulation, the projected star particle density distributions around 1.4 Gyr since the first close encounter of the two simulated galaxies match the above-described tidal features of VCC 848 best (Figure 11). We however note that the two outermost shells of VCC 848 have somewhat less smooth outer edges than do the ones in the S1 simulation. This may be ascribed to the secondary progenitor of VCC 848 being gas-bearing rather than gas-free as assumed in the S1 simulation (see the next two sections for more discussion). By the best-match moment, the progressively disassembled secondary progenitor has crossed through the primary for more than 7 times. The secondary starts to be dismantled after making the 3rd passage through the primary, and since then, one new stellar shell is produced every time the secondary crosses though the primary. While the outlying extended stellar shells are largely aligned along the east-west direction, as is observed in VCC 848, the direction of alignment of the shells produced near or after the best-match moment appears to be tilted toward the major-axis direction of the primary. This may partly explain the observed isophotal twist of the stellar main body of VCC 848.

Response of the Gas Component

The gas and stellar components of the primary progenitor respond to the gravitational disturbances in a similar manner, except that the gas component develops much narrower tidal features than the stellar component, in general agreement with previous simulations of star-dominated galaxy mergers (e.g., Barnes & Hernquist 1996). By the moment that the simulation matches the stellar light distribution of VCC 848, the gas component has largely settled into the warped central disk, with few gas particles associated with the extra-planar stellar tidal arms emanating from the primary, which should be attributed to orbital energy dissipation (through shock heating and radiative cooling) of the gas entrained in (self-)intersecting tidal arms. The fraction of gas particles with relatively high local (volume) number densities, as is illustrated in Figure 12, increases steadily during the course of merging. About 7% of the 50,000 gas particles end up in the central 0.1 kpc in radius of the merger remnant. The enhancement of dense gas fraction may eventually results in an enhanced star formation activities in reality. In addition, the curve of temporal increase of dense gas fraction flattens with time, especially after the 3rd passage ( 0.6 Gyr). The flattening corresponds to a gradual weakening of the gravitational distortion of the merger remnant, due to the progressive disintegration of the secondary progenitor. No gas or stellar bar forms during the merging. which may be attributed to the gas dominant nature of the primary's disk (see also Athanassoula et al. 2013).

There are noticeable differences between the resultant gas distribution in our simulation and the observed H i distribution in VCC 848. In particular, nearly all of the gas particles end up being confined to a central disk in the simulation, whereas in VCC 848, besides the high surface density H i gas confined to the disk plane, a significant fraction of H i is associated with the extraplanar stellar tidal arms as well (Section 3). We have run simulations with a initial gas disk scale length ∼ 2 times our fiducial choice, but did not succeed in producing extra-planar gas arms as significant as observed in VCC 848. It appears that the observed H i distribution corresponds to an earlier merging stage (i.e., more violent tidal disruption) than that indicated by the stellar distributions. Future simulations with a comprehensive parameter-space exploration may hold the promise to achieve a better match to observations. It also remains to be seen whether starburst-driven galactic winds can significantly affect the dynamical evolution of the gas component. With that said, we note that neither of the two previous simulations of gas-rich dwarf galaxy mergers by Bekki (2008) and Starkenburg et al. (2016) which include star formation and feedback produced prominent gas tidal tails or arms. In spite of the discrepancy, our discussion about the simulated gas component, including the tendency of settling towards the central disk, enhancement of gas densities and the impact on the stellar component (see below), should be largely valid in a qualitative sense.

Gas Content of the Secondary Progenitor

In the S2 simulation illustrated in Figure 13, the merging proceeds ∼ 80 Myr faster, and the contrast of newly formed shells relative to preceding generation of shells is significantly lower than in the S1 simulation. The faster merging in the S2 simulation reflects a faster orbital decay due to an additional mass contribution of the gas component to the secondary progenitor. During the first three passages, the secondary progenitor experiences a nearly complete gas removal and the removed gas quickly settles into the primary's gas disk. The gas removal from the secondary, which results from a dissipative collision with the primary's gas component, leaves behind a stellar system that is out of dynamical equilibrium and experiences a rapid expansion and dissolution. The contrast and decay rate of stellar shells are mainly determined by the phase space volume of the accreted galaxy (Hernquist & Quinn 1988), in the sense that a larger phase space volume leads to lower contrast and quicker decay. This explains why the shell system in the S2 simulation has much lower contrast and faster decay than that in the S1 simulation. We note that the dynamical response of the secondary's stellar system to rapid gas removal is analogous to the dispersal of star clusters via "infant mortality" (e.g., Hills 1980).

The near absence of shell contrast in the S2 simulation does not match the observations, which implies that the disk of the secondary progenitor of VCC 848 was dominated by stars rather than gas. Nevertheless, the fact that the two outermost shells of VCC 848 have less smooth outer edges than that in the S1 simulation (Section 6.2.1) probably indicates that the shells in VCC 848 decay faster than in the S1 simulation and hence the secondary progenitor may not be completely devoid of (atomic) gas. It is however not possible to infer the exact gas fraction of the secondary progenitor based on the observed gas distribution, because the gas component originally belonging to the secondary, if any, quickly loses its identity and is phase-mixed with the primary progenitor's gas component before stellar shells start to be produced.

In Paper I, we have shown that the average cluster formation rate of VCC 848 increased by a factor of ∼ 7-10 in the recent ∼ 1 Gyr, which implies that the merger event has triggered star formation that accounts for ∼ 40% of the total stellar mass (2.1×10 8 M ) by assuming a constant cluster mass function. Given the current gas mass of 4.2×10 8 M , the triggered star formation has consumed ∼ 17% of the total pre-merger gas, which is 3 times the gas inflow or high-density gas in the S1 simulation. In addition, an indirect constraint on the secondary progenitor's gas content may be placed by its location on the color-magnitude diagram. By adopting a (gi) color of 0.6 -0.7 mag as observed for the outlying stellar shells (Paper I), which is likely an upper (redder) limit for the secondary progenitor as a whole, we find that the secondary is 0.1 mag bluer than the blue edge of the (gi) distribution of Virgo early-type dwarf galaxies (Janz & Lisker 2009;Roediger et al. 2017) for a plausible g band magnitude -15 (Paper I). Above all, although being dominated by stars, the secondary progenitor of VCC 848 might contain a non-negligible amount of gas. As illustrated in Figure 12, gas from the secondary progenitor can dramatically affect the gas distribution of the merger remnant. Particularly, the S2 simulation ends up with an order of magnitude higher fraction of high density gas than does the S1 simulation, with ∼ 50% of the gas particles flowing into the central 0.1 kpc in radius. This dramatic difference must be ascribed to the hydrodynamic interaction (e.g., ram-pressure sweeping) between the two progenitor gas disks. Blumenthal & Barnes (2018) pointed out that inward gas migration induced by the hydrodynamic interaction results in more effective coupling of the stellar and gaseous components, which favors more efficient gas inflow to the center driven by gravitational torques exerted by the distorted stellar disk. The case for VCC 848 is expected to be in between the S1 and S2 simulations.

SUMMARY AND DISCUSSION

Summary of VCC 848

Understanding the impact of dwarf-dwarf merging on star formation and morphological evolution is of funda-mental importance in the hierarchical structure formation paradigm. In Paper I, we reported the first study of VCC 848 that reveals the tell-tale signature of galaxy merger: extended stellar shells, and a significant enhancement of star cluster formation rate of the system during the past ∼ 1 Gyr. VCC 848 is probably by far the clearest example of star-forming galaxies formed as remnants of mergers between gas-bearing dwarf galaxies (M < 10 9 M ). In the current paper, we have studied the H i gas distribution and star formation activities of VCC 848. We also invoke numerical simulations to probe the evolution of the stellar and gaseous components in response to gas-dominated merging between dwarf galaxies. Based on these analysis, we can draw the following summary for VCC 848.

VCC 848 is the remnant of a merger between two gasbearing dwarf galaxies with a primary-to-secondary stellar mass ratio of ∼ 2 -5 and dark matter halo mass ratio of 2. The primary progenitor is expected to be gas-dominated whereas the secondary progenitor be gasbearing but star-dominated. The two progenitor galaxies experienced their first periapsis passage more than 1 Gyr ago, and since then a burst of star formation, as traced by the formation of super star clusters, has been triggered near the center of the system several hundred Myr ago (Paper I). SFR in the last ∼ 10 -100 Myr, as traced by the Hα emission and FUV continuum, has dropped to a level that is consistent with being a main-sequence star-forming galaxy. The luminous H ii regions in VCC 848 follow a Hα luminosity distribution similar to that of ordinary star-forming dwarfs, except that the brightest H ii region (and the associated brightest star cluster as well) in VCC 848 has extraordinarily high luminosity compared to its counterparts in ordinary dwarfs.

Although VCC 848 is globally dominated by neutral atomic gas, less than 30% of the atomic gas is associated with the stellar main body which contains nearly all of the total SFR. Molecular gas mass inferred from the SFR seems to dominate over the atomic gas mass in the central ∼ 1.5 kpc. It is conceivable that a large amount of atomic gas has been transformed to the molecular phase that fuels the active star formation activities in the recent past. By contrast, a significant fraction of the atomic gas entrained in tidal arms reaches column densities as high ( 10 21 cm -2 or 10 M pc -2 ) as that in the main body but with negligible associated star formation.

The observed stellar distribution of VCC 848 implies a nearly radial and non-coplanar encounter between two dwarf galaxies whose disks are highly inclined with respect to the orbital plane. Early simulations of spiralspiral mergers (e.g., Mihos & Hernquist 1996) that non-coplanar encounters do not result in gas inflow and starburst activities as strong as coplanar encounters do. More recent simulations by Blumenthal & Barnes (2018) find an intricate dependence of gas inflow on both encounter parameters and progenitor properties. In particular, they find that a higher gas content or larger gas-to-stellar disk size ratio generally leads to less effective gas inflow. One may therefore be tempted to deduce that dwarf-dwarf mergers are not as effective as spiral-spiral mergers in driving gas inflow and centralized starburst in general. Nevertheless, our simulations tailored to VCC 848 suggests that, while the merger of a gas-dominated primary with a gas-free secondary results in weak gas inflow (∼ 5%), the merger with a gas-bearing secondary can drive a substantial gas inflow (up to ∼ 50%) to the center of the remnant, even in the seemingly disfavored encounter geometry. The substantial gas inflow seen in our simulations with a gas-bearing secondary is in agreement with SPH-based simulations by Bekki (2008) and Watts & Bekki (2016). Our dynamical analysis of the H i velocity field suggests that the merger event has not built up a central region as compact as typical compact dwarfs with centralized starburst (Section 4.2). This may be explained either by a relatively low gas content in the secondary progenitor (thus ineffective gas inflow) according to our simulations, or by a substantial gas outflow driven by stellar feedback following centralized starburst in early stages of the merger. Stellar feedback is expected to be more effective in driving galactic outflows in dwarf galaxies than in massive galaxies. Kim et al. (2009) ran the first adaptive mesh refinement (AMR) simulation of a gas-rich dwarf-dwarf merger by incorporating a high-resolution treatment of the multiphase interstellar medium and shock-induced star formation, and they found that more than 90% of the gas is expelled by stellar feedback. In VCC 848, the centralized starburst triggered at early stages of the merger may have driven substantial galactic outflows that evacuate the gas that once flowed into the central region.

Reflection on the outcome of dwarf-dwarf mergers

A general conclusion from early simulations (e.g., Mihos & Hernquist 1996) is that mergers of galaxies lack-ing compact bulges (as is generally true for normal starforming dwarf galaxies) develop substantial gas inflows and induce intense starburst activities at much earlier merging stages than those with bulges. This is because bulges act to suppress the development of strong nonaxisymmetric structures (e.g., bars) that helps driving gas inflows until the final coalescence of mergers. Our result that the most intense starburst activities in VCC 848 happened near the center at earlier stages of the merger appear to be in line with these simulations. It is noteworthy that a recent study by Privon et al. (2017) found widespread starburst activities in a well-separated pair of dwarfs that appear to be at an early stage of merging. In addition, the lack of a compact core in VCC 848 is also in line with the general expectation that stellar feedback may drive substantial galactic outflows that act to suppress the formation of stellar bulges in lowmass galaxy mergers (e.g., Ubler et al. 2014).

If these findings are verified to be generally true for dwarf-dwarf mergers, it would negate the possibility that mergers between (normal) low surface brightness starforming dwarfs can result in compact dwarf remnants. The observed compact dwarfs with centralized starburst activities, if being involved in mergers, may be formed from progenitors that already had bulge-like compact cores prior to the merger event. Mechanisms such as gas accretion from the cosmic web (e.g., Sánchez Almeida et al. 2014) and disk instabilities (e.g., Elmegreen, Zhang & Hunter 2012) may have contributed to formation of pre-existing compact cores. Pre-existing compact cores may help sustain strong starburst activities when their host galaxies are involved in galaxy mergers.

It has long been recognized that gas plays a disproportionately large role in interacting and merging galaxies. The usually gas-dominated nature of star-forming dwarf galaxies means that dwarf-dwarf mergers may not be simply the scale-down counterparts of (star-dominated) spiral-spiral mergers in the local universe. Although star-forming dwarf galaxies in the local universe can reach similarly high gas richness as the observed highredshift massive star-forming galaxies, gas in the disks of local dwarfs is primarily in a warm diffuse atomic phase whereas it is mainly in a cold clumpy molecular phase for high-redshift massive galaxies (e.g., Tacconi et al. 2010). As has been pointed out by Fensch et al. (2017), local gas-rich dwarf mergers may induce more significant enhancement of gas inflow and star formation than mergers of high-redshift massive clumpy galaxies in a relative sense. A complete understanding of the outcome of dwarf-dwarf mergers requires detailed analysis of representative samples of dwarf systems at both early (e.g., Stierwalt et al. 2015) and late stages of merging. On the simulation side, a more realistic treatment of the multiphase interstellar medium, star formation and feedback process, will be indispensable in unraveling the physical processes governing the evolution of dwarf-dwarf mergers in general.

a

  ergs s -1 b 1 Total IR luminosity . . . . . . . . . . . . . . . . . . . . 4.0 × 10 7 L b,c 5 SFRFUV+IR . . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.025 M yr -1 d This work SFRHα+IR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.023 M yr -1 d This work SFRIR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.002 M yr -1 d This work Central surface brightness and color measured within the central 1 arcsec in semi-major axis b Value has been adjusted to our adopted distance c Total IR luminosity derived by combining the Herschel 100 µm, 160 µm and 250 µm photometry of Grossi et al. (2015), following the recipe of Galametz et al. (2013) d SFR derived from the observed FUV, Hα or IR luminosity by adopting equations 8 and 10 of Catalan-Torrecilla et al. (2015) References: (1) Gil de Paz et al. 2003; (2)

Figure 2 .

 2 Figure2. H i gas column density map (left), intensity-weighted velocity field (middle) and velocity dispersion map (right) of VCC 848. The g-band image is contoured onto the H i distributions in the top panels, while the H i distributions are contoured onto the g-band image in the bottom panels. The g-band contour levels run from 19.8 to 25.3 mag arcsec -2 with a 0.5 mag interval, and from 25.3 to 28.3 mag arcsec -2 with a 1.0 mag interval. Only pixels with velocity-integrated intensities ≥ 3 times the rms noise are displayed in the velocity field and ≥ 5 times the rms noise in the velocity dispersion map. The H i column density contours are at levels of (2, 5, 10, 15) × 10 20 cm -2 . The velocity field contours are drawn with 10 km s -1 intervals, and the systemic velocity of 1532 km s -1 is drawn as the yellowish contour lines. The red plus symbol in each panel marks the photometric center of the system. The beam size (7. 1 × 6. 6, PA = 56.78 • ) is indicated in the bottom-left of the top panels. near the central disk region, nearly all the other spatial locations have H i line profiles that each are adequately fitted by a single Gaussian component. This suggests that the intensity-weighted velocity field shown in Figure 2 is largely a reasonable representation of the large scale motion (along the line of sight) of the H i gas. The H i channel maps (binned by a factor of 3 along the velocity dimension) are presented in Figure A1 of the Appendix.

Figure 3 .

 3 Figure 3. Results of the H i velocity field modeling with the 2dbat algorithm. The left panel shows the model velocity field constructed using the best-fit tilted-ring parameters, where Vsys = 1532 km s -1 , PA = 21.8 • and inclination angle = 45.2 • . The same velocity contours shown in Figure 2 for the observed velocity field are overlaid on the model velocity field. The right panel shows the rotation velocities determined for the approaching (blue) and receding (red) sides separately. Note that the receding side (to the southwest) is severely distorted by the ongoing merging event. The approaching side is largely consistent with a solid-body rotation that can be reasonably modeled. The filled symbols represent the extracted rotation velocities, and the open circular symbols represent the asymmetric drift corrected rotation velocities. The blue curve represents a parametric fitting to the rotation curve of the approaching side by adopting the same functional form used in Barrera-Ballesteros et al. (2018). See Section 4.2 for details.

Figure 4 .

 4 Figure 4. H i position-velocity diagrams (right panels) extracted along four 7 -wide pseudo-slits oriented at PA = 38.3 • (PV1), 83.3 • (PV2), 128.3 • (PV3) and 173.3 • (PV4), respectively, where 38.3 • is the PA of the best-fit kinematic major axis of VCC 848 (Figure3). The four pseudo-slits (centered at the galaxy center) are indicated as long magenta lines overlaid on the H i column density map (left panel). Every 20 along the pseudo-slit direction is marked with a 7 -long line perpendicular to the slit direction. Also overlaid on the column density map are the intensity-weighted H i velocity contours, with 15 km s -1 intervals. In the right panels, the dot-connected red curves are the projected radial velocities based on rotation curves determined for the central 30 of VCC 848 (Figure3).

  Figure 6. gi color map (left) and Hα flux density map (right) of the central 1.5×1.5 arcmin of VCC 848. The color bar for the Hα map is labeled with the corresponding SFR surface density, by adopting the recipe from Catalan-Torrecilla et al. (2015).The same H i column density contours shown in Figure2are overlaid in both panels here. The Herschel 100µm is overlaid as red contours (2σ, 3σ and 4σ above the sky background) on the Hα flux density map. Note that the dust IR emission contributes at most ∼ 10% of the total SFR budget (Hα+TIR) across the system at spatial scales comparable to the resolution (9. 4) of the 100 µm image.

Figure 7 .

 7 Figure 7. Maps of specific SFR (sSFR, left) and star formation efficiency (SFE) with respect to atomic plus molecular gas (right) of the central 1.5×1.5 arcmin of VCC 848. As in Figure 6, the SFR is estimated based on the observed Hα flux intensities.The Hα and optical images have been smoothed to match the spatial resolution of H i data and re-binned to a 7 ×7 pixel size before deriving the maps shown here. The same H i column density contours shown in Figure2are overlaid in both panels. See Section 5.1 for details.Kennicutt et al. 1989;Helmboldt et al. 2005;Liu et al. 2013), but appears to be at the upper bound for dIrr galaxies (∼ 1.5 ± 0.5; e.g.Strobel et al. 1991;Miller & Hodge 1994;Youngblood & Hunter 1999;van Zee 2000). However, it has been shown that measurements of α are in many cases dependent on the L(Hα) range used in the power-law fitting, in the sense that α measured at higher luminosities tends to be similar to or larger than that at lower luminosities (e.g.,Youngblood & Hunter 1999). The H ii region detection limit for VCC 848 is relatively high compared to most previous studies of dIrr galaxies in the literature, so the large α measured for VCC 848 may not be surprising. We shall return to this point when comparing VCC 848 with ordinary dIrr galaxies in the next section.

  Figure 8.Σ SFR(Hα) is plotted as a function of ΣHI at 560 pc resolution. The data points correspond to the individual pixels plotted in Figures6 and 7, and are color-coded according to their stellar mass surface densities. Note that only data points that fall within the 5×10 20 cm -2 contours (Figure2) are plotted. The diagonal dotted lines represent lines of constant SFE(H i) and thus constant H i gas depletion (by star formation after transitioning to the molecular phase) times τHI (10 8 , 10 9 and 10 10 yr). Above the detection limit of SFR(Hα), the data points fall into two groups which are separated by a 0.3 dex gap in log Σ SFR(Hα) from -3.6 to -3.3. Data points of the group with higher log Σ SFR(Hα) belong to the central high surface brightness star-forming regions, whereas those of the lower-Σ SFR(Hα) group with log ΣHI 0.88 are from the southeast H i tidal arm (Figure2).

Figure 9 .

 9 Figure 9. Central 1×1 arcmin of the Hα image of VCC 848 (left) and cumulative Hα luminosity distribution of H ii regions (right).In the left panel, the green and red circular apertures mark the H ii regions brighter and fainter than the 90% completeness limit (> 10 37.7 erg s -1 ), respectively. The big dotted ellipse encloses the region used for completeness estimation. In the right panel, the cumulative Hα luminosity distribution of VCC 848 above the 90% completeness limit is shown as red filled circles while the completeness-corrected distribution as red open circles. The red solid curve corresponds to a power-law fit to the completeness-corrected Hα luminosity distribution of VCC 848, and the red shaded region represents the uncertainties of the power-law fit. The green dotted curve delineates the completeness limit (marked in the right y-axis) as a function of Hα luminosities, and the gray shaded region represents the regime below the 90% completeness limit. Cumulative Hα luminosity distributions for 10 relatively isolated dIrr galaxies whose MB are within ± 0.5 mag of VCC 848 are plotted as gray curves. The black dashed curve corresponds to the normalized co-added distribution of the 10 dIrr galaxies.

Figure 10 .

 10 Figure10.Illustration of the viewing direction (60 • of azimuth angle in the orbital plane, 75 • of elevation angle above the orbital plane) that leads to the best match between our simulation results and the stellar light distribution of VCC 848. The disk planes of the two progenitors (at an angle of 60 • to each other) are parallel with the sideways axis. This 3D view shows projected distributions of star particles of the primary (blue) and secondary progenitors (red) at a snapshot 0.05 Gyr before the first close passage. For clarity, only 5000 randomly selected star particles of each progenitor are plotted.

Figure 11 .

 11 Figure11. Projected surface number density distribution of star particles (grayscale in a logarithmic stretch) at four selected snapshots of the S1 simulation, where the primary progenitor has a 2:1 gas-to-stellar mass ratio while the secondary is gas-free. The adaptively smoothed g-band image of VCC 848 is shown in the top left panel for comparison purposes. For the rest of panels that show the simulation results, the time elapsed since the first close passage is indicated for each column. The surface number densities are determined at a 0.2×0.2 kpc spatial resolution. The top panels show the primary+secondary combined stellar maps, the middle panels show the stellar maps of the primary, and the bottom panels show the stellar maps of the secondary. The gas particle surface number density distributions are overplotted in the top and middle panels as green contours with logarithmic intervals. The viewing (or projection) direction is as illustrated in Figure10.

Figure 12 .

 12 Figure12. Temporal variation of the fraction of gas particles with relatively high local number densities ngas in the S1 (solid curves) and S2 (dashed curves) simulations, respectively. ngas is calculated as the average within cubic cells with 50 pc side lengths. Two different ngas thresholds (in the regime of cold neutral medium phase) are used to illustrate the temporal variation.

Figure 13 .

 13 Figure13. Same as Figure11, but for selected snapshots of the S2 simulation, where the primary and secondary progenitors have gas-to-stellar mass ratios of 2:1 and 1:1 respectively. Note that the first close passage in the S2 simulation occurs ∼ 80 Myr earlier than the S1 simulation.

Table 1 .

 1 Properties of VCC 848 Property Value Reference Other name . . . . . . . . . . . . . . . . . . . . . . . . . . . i) color . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.46 mag This work 12+log(O/H) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20.8 mag arcsec -2 a This work This work Stellar mass . . . . . . . . . . . . . . . . . . . . . . . . . . . Vmax . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41.5 km s -1 This work Radius of Vmax . . . . . . . . . . . . . . . . . . . . . . . .

	Morphological classification	A1223+06 iI,M BCD	• • • 1
	Distance . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .	16.5 Mpc	2
	Heliocentric radial velocity. . . . . . . . . . . . .	1532 km s -1	This work
	Angular distance from M87 . . . . . . . . . . . . Angular distance from M49 . . . . . . . . . . . .	6.7 2.4	•	• • •
	-16.05 mag b -16.46 mag (g 8.03 Absolute g magnitude . . . . . . . . . . . . . . . . .	1 This work 3
	i-band scale length of stellar main body	0.67 kpc	This work
	µ0,g

•

• • • Absolute B magnitude . . . . . . . . . . . . . . . . . µ0,g -µ0,i . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 0.3 mag a 2.1 × 10 8 M This work H i mass . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 × 10 8 M b 4 1.6 kpc This work Hα luminosity . . . . . . . . . . . . . . . . . . . . . . . . . 3.9 × 10 39

Table 2 .

 2 The JVLA H i Observations of VCC 848

	Haynes et al. 2011;

(3) problematic data points were located (plotms) and flagged (flagdata); (4) derive the ITRF antenna position corrections (gencal); (5) set the model visibility amplitude and phase of the flux calibrator (setjy), and

Table 3 .

 3 H i mapping parameters of VCC 848 Multiscale parameter. . . . . . . . . . . . . . . . . . . . . . . . . 6 , 12 , 24 , 36 Pixel size . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

	Parameters	Value
	Briggs robustness parameter . . . . . . . . . . . . . . . . .	0.7
		1. 5
	Synthesized beam . . . . . . . . . . . . . . . . . . . . . . . . . . . . Beam position angle . . . . . . . . . . . . . . . . . . . . . . . . .	7. 1 × 6. 6 56.8 •
	rms noise (mJy beam -1 ch -1 ) . . . . . . . . . . . . . . .	0.64
	N (H i) rms noise (×10

  19

Table 4 .

 4 Hα luminosities of the brightest H ii re-

	gions in VCC 848
	ID	R.A.	Decl.	logLHα	σ(logLHα)
		(degree) (degree) (erg s -1 ) (erg s -1 )
	(1)	(2)	(3)	(4)	(5)
	1	186.46985 5.81168	38.94	0.01
	2	186.46889 5.81087	38.36	0.02
	3	186.46735 5.80558	38.16	0.02
	4	186.46820 5.80808	38.06	0.03
	5	186.46973 5.81301	37.99	0.04
	6	186.46762 5.80660	37.96	0.04
	7	186.46974 5.81087	37.92	0.02
	8	186.46677 5.80580	37.87	0.03
	9	186.46882 5.81001	37.82	0.03
	10 186.46925 5.81512	37.74	0.05

It is not straightforward to apply an inclination correction to this observed velocity width, as we will show that the H i gas in VCC 848 appears not to have a simple planar spatial distribution.

ACKNOWLEDGMENTS

We thank the anonymous referee for her/his constructive comments that improved the manuscript. We acknowledge support from the National Key R&D Program of China (2017YFA0402702), the NSFC grant (Nos. 11973039 and 11421303), and the CAS Pioneer Hundred Talents Program. SP acknowledges support from the New Researcher Program (Shinjin grant No. 2019R1C1C1009600) through the National Research Foundation of Korea. SHO acknowledges a support from the National Research Foundation of Korea (NRF) grant funded by the Korea government (Ministry of Science and ICT: MSIT) (No. NRF-2020R1A2C1008706). YG's research is supported by National Key Basic Research and Development Program of China (grant No. 2017YFA0402700), National Natural Science Foundation of China (grant Nos. 11861131007, 11420101002), and Chinese Academy of Sciences Key Research Program of Frontier Sciences (grant No. QYZDJSSW-SLH008). Software: casa (McMullin et al. 2007), 2dbat (Oh et al. 2018), SExtractor (Bertin & Arnouts 1996), DICE (Perret 2016)

APPENDIX

A. H i CHANNEL MAPS OF VCC 848

The H i channel maps of VCC 848 are contoured on the g band image in Figure A1.