Panchromatic HST/WFC3 Imaging Studies of Young, Rapidly Evolving Planetary Nebulae. II. NGC 7027

The iconic planetary nebula (PN) NGC 7027 is bright, nearby (D ~ 1 kpc), highly ionized, intricately structured, and well observed. This nebula is hence an ideal case study for understanding PN shaping and evolution processes. Accordingly, we have conducted a comprehensive imaging survey of NGC 7027 comprised of twelve HST Wide Field Camera 3 images in narrow-band and continuum filters spanning the wavelength range 0.243--1.67 microns. The resulting panchromatic image suite reveals the spatial distributions of emission lines covering low-ionization species such as singly ionized Fe, N, and Si, through H recombination lines, to more highly ionized O and Ne. These images, combined with available X-ray and radio data, provide the most extensive view of the structure of NGC 7027 obtained to date. Among other findings, we have traced the ionization structure and dust extinction within the nebula in sub-arcsecond detail; uncovered multipolar structures actively driven by collimated winds that protrude through and beyond the PN's bright inner core; compared the ionization patterns in the WFC3 images to X-ray and radio images of its interior hot gas and to its molecular outflows; pinpointed the loci of thin, shocked interfaces deep inside the nebula; and more precisely characterized the central star. We use these results to describe the recent history of this young and rapidly evolving PN in terms of a series of shaping events. This evolutionary sequence involves both thermal and ram pressures, and is far more complex than predicted by extant models of UV photoionization or winds from a single central progenitor star, thereby highlighting the likely influence of an unseen binary companion.

A planetary nebula (PN) is generated when an intermediate mass star ( 0.8-8 M ) evolves off the asymptotic giant branch (AGB), having ejected its dusty AGB envelope at mass loss rates exceeding ∼ 10 −5 M yr −1 (Höfner & Olofsson 2018). Within a few tens or hundreds of years after AGB envelope ejection, depending on the star's initial mass (e.g., Miller Bertolami 2016a), the increasing wind speed of the rapidly transitioning, mass-losing star -from ∼10 km s −1 on the AGB to ∼1000 km s −1 , post-AGB -sweeps the AGB material into a dense shell, even as the UV from the newly exposed, hot (T ef f ∼ 200 kK) AGB core begins to ionize the AGB ejecta, forming the PN. The "classical" result of the foregoing post-AGB evolutionary processes is a spherically symmetric, thickwalled PN with semi-evacuated interior that is bright in atomic emission lines across the UV through IR wavelength range (Kwok et al. 1978;Schönberner et al. 2018Schönberner et al. , 2005Toalá & Arthur 2014).
However, a significant fraction of PNe display profoundly axisymmetric and/or point-symmetric structures (see, e.g., Parker et al. 2006, and references therein). Furthermore, it has been known for decades that many of these highly non-spherical objects harbor large masses of residual molecular gas and dust that are directly descended from AGB envelope material (e.g. Kastner et al. 1996;Bublitz et al. 2019). The present consensus is that the non-spherical shapes and large residual molecular masses of these PNe are the direct consequence of the presence and influence of binary companions during AGB evolution (e.g., De Marco & Izzard 2017, and references therein). Binary companions can influence mass-launching and collimation mechanisms, which then shape PNe by way of ram pressure in collimated flows or jets (e.g., Balick et al. 2019). The details of these binary-induced PN shaping processes remain elusive, however; fast mass loss driven by a binary companion's accretion disk as well as common envelope evolution are leading hypotheses to explain the origin of collimated outflows (e.g., García-Segura et al. 2018;Zou et al. 2020Zou et al. , 2022. Hubble Space Telescope (HST) emission-line imaging studies of the most recently formed, and hence most rapidly evolving, PNe provide particularly effective means to understand PN ionization and shaping processes (e.g., Sahai & Trauger 1998). In HST's Cycle 27, NGC 7027 was one of the first two PNe to be targeted for a comprehensive, contemporaneous set of emission-line images, from near-UV through optical to near-IR, with the Wide Field Camera 3 (WFC3) (Kastner et al. 2020). We selected NGC 7027 for such a study because it is among the youngest and most rapidly evolving PN within ∼1 kpc of the Sun (D = 0.89 kpc; Masson 1989). Indeed, NGC 7027 represents an "industry standard" for emission-line searches and surveys, thanks to the sheer number of bright lines it displays across the electromagnetic spectrum, from radio to X-ray (e.g., Zhang et al. 2005Zhang et al. , 2008Montez & Kastner 2018;Neufeld et al. 2020).
Previous HST imaging using its Wide Field Planetary Camera 2 (WFPC2) and Near-Infrared Camera and Multi-object Spectrometer (NICMOS) have established the remarkably complex structure of the nebula (e.g., Latter et al. 2000;Schönberner et al. 2018;Guerrero et al. 2020, and references therein). Its distinguishing features include a bright, elliptical inner region, a dusty equatorial belt, and multiple, point-symmetric collimated outflows, with the latter confirmed kinematically via ground-based near-infrared and radio interferometric imaging spectroscopy (e.g., Cox et al. 2002;Nakashima et al. 2010). Concentric ring-like or spiral structures are observed to surround the core region of NGC 7027 in HST/WFCP2 Hα and V-band images (Guerrero et al. 2020, and references therein). Radio-regime proper motion analyses of the inner elliptical shell indicate a kinematic age of ∼600 yrs (Masson 1989;Zijlstra et al. 2008), while nebular expansion analyses that are based on multi-epoch HST imaging of the extended ring system suggest a dynamical age in the range 1080-1620 yrs (Schönberner et al. 2018;Guerrero et al. 2020). Consistent with such a recent departure from the AGB, the nebula harbors an exceedingly hot (∼200 kK) proto-white dwarf (pWD) at its center (Latter et al. 2000).
First results from our Cycle 27 HST/WFC3 imaging programs targeting NGC 7027 and a second, well-studied PN -the similarly young (expansion age ∼2000 yr), nearby (D ∼ 1.0 kpc) bipolar PN NGC 6302 -were presented in Kastner et al. (2020). In Kastner et al. (2022), we presented a detailed overview of the HST/WFC3 imaging survey of NGC 6302. Here, we present the full suite of Cycle 27 HST/WFC3 images of NGC 7027, and we highlight compelling results obtained from the image suite. In § 2, we describe the observations and image processing that was performed to obtain the final suite of a dozen (mostly narrow-band) WFC3 images, which are presented and described in § 3. In § 4, we present key line ratio images that reveal the extinction and ionization structure of the nebula. In § 5, the focus is shifted onto the central star (CS) where its extinction, luminosity, and present-day mass are evaluated. In § 6, we develop and discuss extinction, ionization, and other lines of evidence for the presence and structure of nebula-shaping shocks. The implications of the results for the ongoing structural evolution of NGC 7027 are described in § 7. A summary of the results of our HST/WFC3 imaging study of NGC 7027 is presented in § 8.

OBSERVATIONS
Images of NGC 7027 presented in this paper were obtained with both the UVIS and NIR units of HST's WFC3 during Cycle 27 in 2019 September. The UVIS (CCD sensor) channel provides a field of view of 2.7 x 2.7 and pixel scale of 0.04 , and the NIR channel uses a HgCdTe array with field of view of 2.27 x 2.27 and pixel scale of 0.13 . Over the range of wavelengths imaged here, the point-spread function (PSF) varies from ∼ 0.025 (near-UV) to ∼ 0.16 (NIR). Table 1 summarizes the WFC3 filters used to obtain images of NGC 7027 as well as the targeted emission line, the date the image was obtained, and the total exposure time of imaging through each filter. The UVIS images were obtained in 2-point GAP-LINE dither mode (DITHER-LINE for quad filter FQ243N), and the IR images were obtained in 2-point DITHERBLOB dither mode. UVIS and IR images were obtained at a position angle of -16.9 • and -17.2 • measured east of north, respectively. Processing of images is described in the comprehensive HST/WFC3 imaging targeting NGC 6302 ) and is briefly described here. The standard pipeline calibration and processing, using CALWF3 v3.5.0, included bias correction, dark current subtraction, flat field and shutter shading corrections, cosmic-ray rejection, and (for UVIS images) CTE correction. Additional postpipeline processing was performed using the DrizzlePac package software 1 , particularly the tweakreg and astrodrizzle modules. Fine image registration was accomplished by defining the World Coordinate System of the F160W filter image via the tweakreg module, using Gaia Data Release 2 catalog stars as positional references. The other WFC3 images were then aligned (using tweakreg) with the F160W image serving as the reference. This positional calibration also included geometric distortion corrections and additional cosmic-ray corrections. The dithered exposures of each long exposure image were then merged using the astrodrizzle module. This process was implemented to correct for any misalignment that remained after the initial pass through HST's image reduction pipeline. We estimate that, following this image registration procedure, the WFC3 images of NGC 7027 are aligned to a relative accuracy of ≤0.05 . Note-a) Filter pivot wavelength and effective bandwidth. b) Potential contamination from O iii λ3312, λ3341, and λ3444 (Zhang et al. 2005). c) Potential contamination from He ii λ4859 and λ6560, respectively (Zhang et al. 2005). d) Short exposure images used for pixel grafting of saturated pixels in the long exposures.
Evidence for pixel saturation is present in the long-exposure F502N and F656N images at the regions of highest surface brightness within the nebula, specifically ∼10 to the northwest of the central star. Pixels that exceeded the saturation threshold of 65,000 counts (see WFC3 Instrument Handbook; Dressel 2021) were replaced with count rates from short-exposure image pixels. This "pixel grafting" procedure was confined to displacements of 0 to −10 (RA) and 0 to 10 (dec) with respect to the position of NGC 7027's central star.
3. THE HST/WFC3 IMAGE SUITE Fig. 1 presents the complete suite of HST/WFC3 images obtained for NGC 7027, spanning a range from near-UV (∼245 nm) to near-IR (∼1665 nm). These images highlight the overall dramatic variation in the morphology of NGC 7027 over this wavelength range. Much of these differences can be attributed to the varying effects of extinction, but also reflect the intensity of the specific emission lines and continuum isolated by each of the filter passbands ( Table 1). The structures observed in these images will be further discussed in the following sections; here, we briefly summarize the main features of the nebular morphology that is highlighted by the HST/WFC3 image suite.
In all images, the overall SE-NW extension of the nebula, which represents its main fast outflow axis (e.g., Cox et al. 2002), is apparent. At the shortest (near-UV) wavelengths, the nebula displays the most profound asymmetry along this (SE-NW) direction. Specifically, in both the F243N and F343N images, the regions immediately to the NW of the central star appear brighter than the regions SE of the central star; this asymmetry is a consequence of the smaller dust extinction toward the (forward-facing) NW regions (Montez & Kastner 2018). The prominent extinction "hole" located ∼3 NW of the central star discussed in the Chandra X-ray imaging study of Montez & Kastner (2018) is especially apparent in both of these near-UV images. The nebula's equatorial dust "belt," as well as a system of dusty filamentary structures located ∼5 directly north of the central star, appear as prominent features of the F343N, F487N, F502N images and, to a lesser extent, the F656N and F673N images.
The near-IR images, covering the range from 1.15 µm (F110W) to 1.664 µm (F167N), are less affected by extinction due to intranebular dust, and (as a result) are dominated by the inner elliptical shell and immediately surrounding ("cloverleaf") structures within the nebula that were previously imaged by HST/NICMOS (Latter et al. 2000). The two images with the highest surface brightness, F502N ([O iii]) and F656N (Hα), best reveal the extended system of ring-like dust structures surrounding the nebula that was previously imaged by HST using WFPC2 (Schönberner et al. 2018;Guerrero et al. 2020).
In all of the WFC3 images (to greater and lesser degrees), jet-like protrusions at position angles (PAs, measured counterclockwise from N) of ∼120 • (ESE) and ∼300 • (WNW) are seen superimposed on the aforementioned nebular structures, with the WNW protrusion particularly prominent. Lesser protrusions appear at PAs of ∼10 • (nearly due N) and ∼200 • (SSW). In addition, a complex of dusty ejecta is seen oriented nearly due S (most notably in the F343N, F487N, F502N, and F656N images); this south-directed ejecta complex has a faint counterpart at PA ∼340 • (NNW) that is most apparent in the F502N and F656N images.

Nebular extinction maps from H recombination lines
The three H recombination line images obtained as part of the WFC3 image suite -i.e., Hβ, Hα, and Paβ, as obtained through WFC3 filters F487N and F656N and from the F128N−F130N difference image, respectively -provide a means to map dust extinction across the nebula at the ∼0.1 resolution of HST/WFC3. To do so, we followed the methodology described in Kastner et al. (2022); i.e., we constructed Hα/Hβ and Paβ/Hβ line ratio images, and compared these per-pixel line ratios with the theoretically predicted ratios, so as to obtain the spatial distribution of the extinction parameter c(Hβ) for a given choice of reddening law. To construct c(Hβ) images from the Hα/Hβ and Paβ/Hβ line ratio images, we adopt intrinsic relative intensities for Hα/Hβ and Paβ/Hβ of 2.85 and 0.162, respectively, which correspond to a "canonical" values of nebular electron temperature (10 4 K) and electron density (10 4 cm −3 ), with very little sensitivity to either parameter (e.g., Osterbrock 1989). We then adopt a uniform value of R V = 3.1 to generate the c(Hβ) images, with the caveat that R depends on grain size, and grain size would not realistically be uniform across the entire PN.
The resulting c(Hβ) images are shown in the bottom panels of Fig. 2. Ideally, if the R λ value used is correct, the c(Hβ) images obtained from the two line ratios should be identical. Indeed, the c(Hβ) images, shown in the bottom panels of Fig.4, are nearly identical in morphology. The key differences in values occur in the nebula's waist and the eastern region outlining the north lobe, with the c(Hβ) map derived using Paβ/Hβ having systematically larger values than the c(Hβ) map obtained from Hα/Hβ. Specifically, values across the c(Hβ) map obtained from Paβ/Hβ range from  ∼ 0.9 (lobes) to ∼ 2.1 (nebula waist), whereas the values in the c(Hβ) map obtained from Hα/Hβ lie in the range ∼ 0.9 − 1.7. Other, more subtle differences between the extinction maps could be due to the relatively large variation in the PSFs between Paβ and Hβ, which would affect the resolution of structures across the Paβ/Hβ ratio image moreso than in the Hα/Hβ ratio image, and any dust scattering of Hα emission that may be present.
Overall, the results for c(Hβ) are consistent with previous studies that mapped extinction across NGC 7027 (Kastner et al. 2002;Montez & Kastner 2018, and references therein). In particular, previous extinction maps also feature extinction maxima along the waist of the nebula, as well as a prominent gap or "hole" in extinction located ∼5-10 to the WNW of the central star, a region that is also bright in soft X-rays (Montez & Kastner 2018). The apparent lack of obscuring nebular material at this location in NGC 7027 is likely due to the disruptive effects of a fast, collimated outflow, as is discussed in detail in § 6.1.
However, in contrast to earlier extinction studies, the ∼0.1 resolution of the line ratio images and (hence) c(Hβ) maps obtained from the HST/WFC3 images ( Figure 2) reveal the dust structures within NGC 7027 in unprecedented detail. North of the CSPN lies dusty lobe structure consisting of a series of ripple-like filaments that terminate at the lobe's northernmost extension, ∼10 from the CSPN. In the southern regions of the inner nebula, opposite the CSPN from the northern filament structures, the c(Hβ) map reveals dust rims and arcs that appear to outline bow-shock-like structures. These features, along with the large knot along the waist of the nebula, appear more prominent in the c(Hβ) image derived from Paβ/Hβ (Figure 2, lower right panel). All of the dust structures revealed in the c(Hβ) maps in Figure 2 appear to have direct counterparts in near-IR H 2 emission (see, e.g., Fig. 3 in Latter et al. 2000), indicating that these structures constitute dense regions of dust and molecular gas that have survived the processes of photoionization and shock-driven ionization that are ongoing in the inner nebula (see next).
In Appendix A, we describe how these c(Hβ) images are used to perform extinction corrections for several of the HST/WFC3 images.  (Fig. 4) does not require extinction correction because of the small wavelength separation between the two emission lines; this ratio image displays structure essentially identical to that seen in the extinction corrected [O iii]/Hα ratio image (as displayed in Fig. 13    ratio images are displayed alongside an image of the ratio of [Fe ii] 1.64 µm forbidden line emission with respect to neighboring 1.67 µm continuum emission (we use this F164N/F167N ratio image in preference to a F164N−F167N difference image, due to its superior signal-to-noise ratio). The 1.67 µm continuum emission is likely due to free-free emission and hence (like Hα) traces ionized H. As discussed in Kastner et al. (2022), both the 0.676 µm [S ii] and (especially) 1.64 µm [Fe ii] emission lines are potential tracers of shock-driven excitation. Fig. 4   of O + and, hence, to [O iii] emission, which would be consistent with studies indicating that [O iii] emission can trace shock zones (in addition to photoionization) in PNe (Guerrero et al. 2013).

High and low excitation lines: ionization gradients
The [Ne v]/[Ne iv] ratio image in Fig. 3 should also provide a map of the reach of highly ionizing radiation, with the caveat that the wide bandwidth of the F343N [Ne v] filter encompasses several permitted O iii lines that are relatively bright in NGC 7027 (see Table 1). However, it is apparent that the signature of intranebular extinction remains strong in this ratio image, despite both component images being subject to the extinction correction procedure described in the Appendix. Specifically, the dark features around the nebula's waist and toward its northern reaches in the [Ne v]/[Ne iv] image correspond to regions where extinction maps reveal large amounts of dust and dust scattering (see Fig. 2). Intriguingly, the overall morphology of this image resembles that of the broad-band Chandra X-ray image of NGC 7027 (Montez & Kastner 2018, see below). This suggests that the regions of enhanced [Ne v] emission (relative to [Ne iv]) might trace a zone of heat conduction from the shocked, X-ray-emitting (∼ 10 6 K) plasma to the (∼ 10 4 K) photoionized gas. Alternatively, the [Ne v]/[Ne iv] image could be affected by scattering, especially in the outer regions Note-a) Estimated using stellar evolution models in Miller Bertolami (2016b). b) Blackbody effective temperature estimate derived in Schönberner et al. (2005). c) Blackbody effective temperature estimate derived in Latter et al. (2000).
of the point-symmetric structures, given that scattering is more efficient at smaller wavelengths, and the scattered UV photons can most easily escape through the low-extinction "holes" formed by the blowouts associated with regions of shocks (see below).

EXTINCTION TOWARD THE CENTRAL STAR
The emission line ratio-based approach to measuring intranebular extinction described in § 4.1 cannot be applied to determine the extinction along the line of sight directly toward the CSPN of NGC 7027, as the star is a bright continuum source. Such an extinction estimate is essential in order to estimate the CSPN's present effective temperature T ef f and luminosity L, and thereby ascertain the present mass and evolutionary state of the CSPN via comparison with theoretical models (e.g., Miller Bertolami 2016a). Hence, we estimated the extinction toward the CSPN via comparison of the star's spectral energy distribution with stellar atmosphere models calculated for hot central stars of PNe. To this end, measurements of the CSPN flux, corrected for local nebular background, were made for 9 HST/WFC3 images spanning a wavelength range from 0.3435 µm to 1.6677 µm. Fig. 5 displays these observed fluxes.
The resulting spectral energy distribution (SED) was then compared with a Tübingen NLTE model atmosphere 2 (Werner et al. 2012) for a model with T ef f = 200 kK (based on the CS effective temperature inferred by Latter et al. 2000). A pure H model with log(g) = 7.0 was selected, for simplicity. Including near-IR and optical points from Latter et al. (2000), the extinction is estimated to be A V = 2.63 using the Cardelli et al. (1989) reddening law with R λ =3.1, A V = 2.54 using the Whittet (1992)  To correct the CS fluxes, we adopted the mean of these extinction estimates, A V = 2.69. In Fig. 5, uncorrected and extinction corrected flux points are overlaid with a 198 kK blackbody, as suggested by Latter et al. (2000), for comparison. Since the extinction correction results for the CS are consistent with the 198 kK blackbody approximation, the blackbody model was used to approximate the total flux of the CS for purposes of luminosity estimation. The total bolometric luminosity obtained for the CS of NGC 7027 is ∼ 6.2 × 10 3 L , which agrees well with previous results, all falling within the range ∼ 5.5 × 10 3 − 2 × 10 4 L (Beintema et al. 1996;Masson 1989). Table 2 summarizes the estimated age and present-day mass of NGC 7027's CS obtained from the preceding results for luminosity and effective temperature as well as the CS parameters derived by these estimates all suggest that the CS is descended from a progenitor of mass ≥3 M (e.g., Miller Bertolami 2016a).

EVIDENCE FOR SHOCKS IN NGC 7027
The presence of strong, wind-collision-generated shocks in NGC 7027 was revealed by the Chandra X-ray Observatory detection of extended, luminous X-ray emission from the nebula extending over the energy range ∼0.3-3.0 keV (Kastner et al. 2001;Montez & Kastner 2018). Montez & Kastner (2018) estimate that the velocities of the outflows that generate the X-ray-emitting regions of NGC 7027 may be as high as v wind ∼ 500 km s −1 , on the basis of the inferred characteristic temperature of the X-ray emitting plasma (∼ 3 × 10 6 K). In Fig. 6 we present an overlay of the hard (1.0−3.0 keV) Chandra X-ray image from Montez & Kastner (2018) on the HST/WFC3 [Ne v] and Hα images. It is apparent from this Figure that the hard X-ray emission closely traces the edge of the inner, highly ionized ([Ne v]-bright) region of the nebula, with the exception of the equatorial region to the ENE and WSW of the central star. These regions correspond to directions of the largest dust extinction, as seen in Fig. 2, indicating that the X-ray emission in these regions is likely being heavily absorbed by cooler, surrounding, CNO-rich gas (see discussion in Montez & Kastner 2018). Furthermore, the X-ray surface brightness protrusions to the NW and SE of the central star closely correspond to regions of high-velocity gas detected in near-IR H i (Brγ) emission (Cox et al. 2002).
In this section, we present additional evidence for shocks in NGC 7027, through line ratio diagnostics obtained from the [O iii]/Hα and [S ii]/Hα ratio images along with a comparison of nebula cross sections in these line ratios and those of [Fe ii] and X-ray surface brightness.

Line ratio diagnostics
Excitation diagrams obtained from diagnostic line ratios, long used to classify AGN (e.g., Kewley et al. 2006), can, in principle, also be used to distinguish between photoionized regions and fast lowionization structures (LISs), or shocked material, in PNe (Danehkar et al. 2018;Montoro-Molina et al. 2022). Here, we follow the methodology described in Danehkar et al. (2018) Danehkar et al. (2018) defines this nebular photon-shock dividing line to be parallel to a classification line that differentiates between Seyfert and low-ionization nuclear emission-line region (LINER) galaxies. The definition was modified to take into account diagnostic diagrams for axisymmetric simulations of fast LISs and observations of photoionized gas in (Raga et al. 2008).
The excitation diagram was then used to classify each pixel and then to color code an H-alpha intensity image of the nebula. The result is shown in the middle panel of Fig. 7, where blue represents photoionized regions and red represents fast LISs. The results indicate that the elliptical interior region of NGC 7027 is dominated by photoionization, while the most extensive shocked regions, according to this analysis, correspond to the same directions as the X-ray protrusions (Montez & Kastner 2018) and the fast flows seen in Brγ emission Cox et al. (2002). Most PN studies have, by necessity, employed spatially integrated line ratios to determine where nebulae lie on a given excitation diagnostic diagram such as that in the left panel of Fig. 7 (e.g., Frew & Parker 2010;Raga et al. 2008 Fig. 7. Note that these mean ratios are within the "photoionized" zone of the diagram, as a consequence of the relatively small surface area of the nebula that displays LIS-like ratios. Thus, if only spatially integrated line ratios were available, one might conclude there were no shocks present in NGC 7027. We caution, however, that this young and rapidly evolving PN may represent an unusual case study in the application of spatially resolved line ratio diagnostics to distinguish between shocked and photoionized regions (see discussion in Montoro-Molina et al. 2022).  (Fig. 7) -serves as clear evidence that the shocks due to fast outflows impinging on the nebular material are strongest in this particular direction.
To further elucidate this structure, we present in Fig. 9 a cross-section of the nebula in X-rays, [S ii], [Fe ii] as well as [O iii] emission along this same direction. This cross-section plot confirms the nested structure of these shock tracers, and reveals the detailed structure of the particularly strong shocks that are propagating along this direction within the nebula. As discussed below ( §7), the clarity of these shocks compared to the nebula's other outflows indicates that the jets emanating along the SW-NE direction could dominate the future shaping process of NGC 7027, perhaps eventually forming a bilobed, pinch-waisted PN that appears as a classical "butterfly" nebula, like NGC 6302 .
The nested spatial distributions of the [S ii] and [Fe ii] emission lines also should aid in constraining the densities within these regions of NGC 7027. In a study of the symbiotic system BI Crucis, whose excitation is largely collisional, Contini et al. (2009, and references therein) point out that [Fe ii] emission should dominate over Fe ii permitted lines if the gas density of a nebula falls between 10 2 −10 4 cm −3 , and both [Fe ii] and Fe ii lines should be excited if the gas density lies between 10 6 −10 8 cm −3 . Similarly, the presence of [S ii] emission constrains the gas density to be <10 4 cm −3 (Acker & Jaschek 1995). Since NGC 7027's optical spectrum displays fairly weak permitted Fe ii and S ii lines compared to their forbidden emission counterparts (Zhang et al. 2005), we infer a nebular density of <10 4 cm −3 in the NW and SE lobe structures. However, these constraints also suggest that the gas from which [Fe ii] emission arises should have a higher density than regions where [S ii] emission is observed, since regions of [S ii] emission are collisionally de-excited (hence quenched) at a density <10 4 cm −3 , but regions of [Fe ii] emission can reach densities of ∼10 5 cm −3 without having visible Fe ii lines. The expectation that [Fe ii] should trace denser gas than [S ii] is indeed consistent with the nested emission structure apparent in Fig. 9.

THE SHAPING OF NGC 7027
Previous sections of this paper have developed the deepest-ever look into the structure and current state of NGC 7027 and its central star. In this section, we use these data to attempt to peek back into the source and evolution of mass injection and nebular shaping.  Bublitz et al. 2022); a schematic diagram of NGC 7027 outlining the primary structural features is presented in the right panel, overlaid on a color composite constructed from the WFC3 image suite by STScI staff 3 . In the schematic diagram, each of the three primary previously identified outflows (Cox et al. 2002;Nakashima et al. 2010;Bublitz et al. 2022) is traced by arrows, with color-coded red and blue extremes to signify their relative redshifts and blueshifts, respectively. The elliptical shell seen in the WFC3 NIR images and radio continuum image is traced by a dashed ellipse. The outer ring system, which is most clearly visible Figure 10. Left panel, three color composite of NGC 7027 generated from Paβ/Hβ (red), extinction corrected Hα (green), and 1.3 mm radio continuum (blue; from Bublitz et al. 2022). Right panel, schematic diagram outlining the structural components discussed in § 7.1, overlaid on a color composite constructed from the WFC3 image suite by STScI staff. Three pairs of jets, labeled as Outflows 1, 2, and 3 (see Cox et al. 2002), are traced with arrows that are color-coded to indicate the relative redshift/blueshift for each outflow component. The nebula's outer ring system, central elliptical shell (long-dashed ellipse), and equatorial belt (short-dashed arc) are also indicated.

A series of shaping events
in the Hα and [O iii] images ( Fig. 1), as well as the equatorial belt of NGC 7027 are also indicated in the schematic.
We propose that this complex set of structures is most likely the result of a series of PN shaping events resulting from binary system interactions involving NGC 7027's AGB star progenitor and an unseen compact companion. Binary interactions while the primary star was still on the AGB likely resulted in the spirals seen surrounding NGC 7027 (in scattered light) in the [O iii] and Hα images (at least 9 have been identified; Corradi et al. 2004). Such ring or spiral patterns in AGB ejecta are potential signatures of the presence of a wide (separation ∼1-10 au) binary system (e.g., Mastrodemos & Morris 1999;Corradi et al. 2004;Chen et al. 2020, and references therein) and, indeed, have been detected in the ejected envelopes of AGB stars that are suspected binary systems (e.g., R Scl; Maercker et al. 2012, and references therein). The estimated post-AGB age of NGC 7027's central star, < 1 kyr ( § 5), is significantly less than the estimated age of the rings, ∼ 2 kyr (Su 2004;Guerrero et al. 2020), supporting the idea that these outer rings formed before the CS left the AGB. The nebula's prominent equatorial belt might then trace the binary's orbital plane, where the bulk of the AGB star's ejected envelope mass would have accumulated (e.g., Mastrodemos & Morris 1999).
Of the three main outflows noted in Fig. 10, Outflow 1 (the focus of the discussion in § 6.2) displays the most pronounced evidence for strong shocks, at the present epoch. The fact that Outflows 2 and 3 appear less prominent or absent in the shock tracers [S ii] and [Fe ii] relative to Outflow 1 may indicate that Outflows 2 and 3 are older than (were ejected prior to) Outflow 1, despite their similar projected extents. Indeed, Outflows 2 and 3 could have occurred during NGC 7027's post-AGB but pre-planetary nebula (PPN) phase (Sahai & Trauger 1998;Sahai 2000;Akashi & Soker 2021), perhaps during repeated close binary encounters that culminated in the unveiling of the central star (see below). These two pairs of jets have triangular shapes that may be bow shocks generated by narrow, collimated jets. Furthermore, Outflow 3 displays a striking misalignment relative to Outflows 1 and 2, wherein Outflow 3 appears to be oriented close to NGC 7027's equatorial plane (Bublitz et al. 2022, and references therein). Some PPNe display a similar combination of recent, collimated outflows -some with similarly surprising misalignments -impinging on older, more slowly expanding ring/spiral structures (AFGL 2688 perhaps being the best-studied example of this phenomenon; e.g., Cox et al. 2000;Balick et al. 2012).
The central elliptical shell of NGC 7027 (Fig. 10) likely consists of the ejected material from the former AGB star that has been swept up by a fast, quasi-spherical post-AGB wind. The age of the CS since leaving the AGB has been estimated consistently to be ∼100-1000 yrs ( § 5). However, since this post-AGB CS age pertains to a single star rather than a central star subject to binary system interactions that would, presumably, truncate AGB evolution, this estimate is most likely an upper limit. Given that the youngest outflow, Outflow 1, may have been ejected as recently as ∼100 yrs ago (Cox et al. 2002), it is possible that this outflow and the unveiling of the CS occurred essentially simultaneously. A thin shell of [Fe ii] emission is seen surrounding the elliptical shell ( Fig. 10), connecting the shell to the shocked regions of Outflow 1 (as seen in [S ii] and [Fe ii]), perhaps providing evidence for this simultaneity.

Did NGC 7027 undergo an ILOT?
The foregoing sequence of events in the shaping history of NGC 7027 appears consistent with a model of formation for bipolar PNe formation via Intermediate-Luminosity Optical Transient (ILOT) ejection events (Soker & Kashi 2012). This (ILOT) model offers an explanation for the formation of bi-lobed PNe that display evidence for multiple, paired mass ejections. Soker & Kashi (2012) set out the following criteria to identify PNe that may have been formed by ILOT events: the PN (1) contains structures or components that exhibit a linear relationship between velocity and distance; (2) is bipolar in structure, indicating the presence of a binary system with orbital separation of ∼ 1 AU; (3) has lobe expansion velocities exceeding ∼ 100 km s −1 ; and (4) has a total fast outflow kinetic energy in the range ∼ 10 46 -10 49 erg. Soker & Kashi (2012) further point to NGC 6302 as an exemplar of a bipolar PN that was generated by an ILOT event some ∼2000 yr ago (see also Kastner et al. 2022).
Indeed, NGC 7027 would appear to satisfy all four ILOT criteria. With regard to criteria 1-3: its cloverleaf-shaped H 2 emitting region exhibits a velocity field proportional to radius (Cox et al. 2002); the nebula displays (nascent) bipolar structure, in the form of its dusty equatorial belt and collimated outflows (especially the youngest, Outflow 1; Fig. 10); and the presence of [Fe ii] and X-ray emission, as well as high-velocity Brγ emission (Cox et al. 2002), along Outflow 1 indicates outflow speeds in excess of ∼ 100 km s −1 , at least along this (SE-NW) direction (see § 6.2). NGC 7027 also appears to meet the fourth (energy) criterion: adopting the mass of 0.1 M estimated by Santander-García et al. (2012) for their "high-velocity molecular blobs", and an assumed velocity of 100 km s −1 , we arrive at a total energy of ∼ 2 × 10 46 erg for its collimated outflows. With a post-AGB age of just ∼100-1000 yr ( §5) and at least one collimated outflow (Outflow 1) that is of similar dynamical age (Cox et al. 2002), NGC 7027 may therefore be the youngest known example of an ILOT-generated bipolar PN.

SUMMARY
The young, rapidly evolving, nearby NGC 7027 is among the best objects for studying PN shaping and excitation mechanisms. We have obtained and analyzed a comprehensive suite of narrow-band HST/WFC3 images for this nebula, spanning a range from near-UV (∼245 nm) to near-IR (∼1665 nm). Our analysis and main results can be summarized as follows.
• Hα/Hβ and Paβ/Hβ emission line ratio images are constructed and a per-pixel comparison is made with the theoretically predicted ratios, so as to obtain the spatial distribution of the extinction parameter c(Hβ). Overall, the results for c(Hβ) are consistent with previous studies that mapped extinction across NGC 7027, featuring extinction maxima along the waist of the nebula, as well as a prominent gap in extinction located ∼ 5 − 10 to the WNW of the central star; however, the extinction maps reveal the dusty substructures within NGC 7027 in unprecedented detail. The c(Hβ) images are used to perform extinction corrections for several of the HST/WFC3 images (see Appendix A).
• Aperture photometry performed on the WFC3 images and near-IR and optical measurements from Latter et al. (2000) are used to construct the spectral energy distribution of the CSPN from the near-UV to near-IR. Analysis of these data yields an estimate of extinction of A V = 2.69 toward the central star. The resulting CSPN luminosity is 6.2 × 10 3 erg s −1 for an assumed effective temperature of 198 kK (Latter et al. 2000). These parameters indicate that the central star is descended from a progenitor of ≥3 M , has a present-day mass of ∼0.6-0.7 M , and left the AGB within the last ∼1000 years. . This structure appears best explained as a juxtaposition of photoionized and shock-ionized zones within the nebula. An excitation diagram is used to assess the spatial regimes of photoionized vs. shocked nebular material ( § 6.1). This analysis reveals that the outflow structures oriented SW-NE -Outflow 1 in Fig. 10 -display the strongest shocks. This conclusion is bolstered by the consistent orientations of the extended X-ray, [S ii], and [Fe ii] emitting regions along this same direction. The nested spatial distributions of these tracers (and [Oiii]) reveal the detailed shock structure (and constrain the gas densities) along Outflow 1.
• A timeline for the shaping of NGC 7027 is proposed, in which a close binary system has produced the nebula's ring system and equatorial belt, and subsequent binary interactions, possibly in a series of short-lived eruptions (the ILOT scenario), then generate the multipolar jets observed. This model awaits confirmation through a study of the proper motions of the complex NGC 7027 jet system, via these WFC3 and archival WFPC2 HST images, so as to establish more precise dynamical ages.