Shock Excitation in Narrow-line Regions Powered by AGN Outflows

Outflows in an active galactic nucleus (AGN) are considered to play a key role in the evolution of the host galaxy through transfer of a large amount of energy. A narrow-line region (NLR) in the AGN is composed of ionized gas extending from parsec to kiloparsec scales. It has been suggested that shocks are required to ionize the NLR gas. If AGN outflows generate these shocks, they will sweep through the NLR, and the outflow energy will be transferred into a galaxy-scale region. In order to study the contribution of the AGN outflow to the NLR-scale shock, we measure the [Fe ii]λ12570/[P ii]λ11886 line ratio, which is a good tracer of shocks, using near-infrared spectroscopic observations with the Warm INfrared Echelle spectrograph to Realize Extreme Dispersion and sensitivity (WINERED) mounted on the New Technology Telescope. In 13 Seyfert galaxies we observed, the [Fe ii] and [P ii] lines were detected in 12 and 6 targets, respectively. The [Fe ii]/[P ii] ratios in 4 targets were found to be higher than 10, which implies the existence of shocks. We also found that the shock is likely to exist where an ionized outflow, i.e., a blue wing in [S iii]λ9533, is present. Our result implies that the ionized outflow present in an NLR-scale region sweeps through the interstellar medium and generates a shock.


INTRODUCTION
The influence of the Active Galactic Nuclei (AGN) on the host galaxy evolution has been actively studied in recent decades (Harrison et al. 2018 for review).
Corresponding author: Misaki Mizumoto mizumoto-m@fukuoka-edu.ac.jpOne example is the correlation between the mass of the central Super-Massive Black Hole (SMBH) and velocity dispersion of the galactic spheroids (i.e., bulges and elliptical galaxies; Ferrarese & Merritt 2000;Gebhardt et al. 2000).The fact that there is a positive correlation between them implies their co-evolution, that is, SMBHs control the evolution of their host galaxies and vice versa.Theoretical studies predict that AGN outflows, which originate from the radiation pressure due to accretion on to the black hole, offer "a plausible physical origin" (King & Pounds 2015) for the co-evolution.Impact of AGN on the host galaxy has also been implemented by many cosmological hydrodynamical simulations, such as the EAGLE project (Booth & Schaye 2009;Crain et al. 2015;Schaye et al. 2015), the SIMBA project (Davé et al. 2019;Thomas et al. 2019;Appleby et al. 2020), and the IllustrisTNG project (Weinberger et al. 2017;Pillepich et al. 2018;Weinberger et al. 2018;Nelson et al. 2019).
In this work, we focus on physical processes in and around a Narrow Line Region (NLR).The NLR extends from pc-scales to kpc-scales and emits forbidden lines with a line width of several hundred km s −1 (Peterson 1997).The main ionization mechanism in the NLR is photo-ionization from the AGN radiation (Ho et al. 1993).However, it is known that the photo-ionization cannot explain all the observed line features (e.g., Simpson et al. 1996) and instead or in addition, shocks are likely to be responsible for some features (e.g., Knop et al. 1996;Wilson & Raymond 1999;Mouri et al. 2000;Rodríguez-Ardila et al. 2004).Here, we came up with a scenario that a shock triggered by the AGN outflow contributes to the NLR ionization; the fast velocity of the AGN outflow will generate the shock by interaction with the interstellar medium (ISM) (e.g., Richings & Faucher-Giguère 2018a,b).The shock may make a critical contribution to the NLR ionization, transferring energy from AGN to the NLR-scale cloud (i.e., the spheroid-scale cloud).In these manners, the AGN outflow can carry its energy to the spheroid.Some studies have partially supported this scenario; for example, Joh et al. (2021) found that the gas density and velocity dispersion of the NLR region are larger when the AGN is more active, which implies that NLR gas clouds are brought from the AGN outflows.
To investigate the existence of shocks in the NLR, Terao et al. (2016, T16) studied the flux ratio of [Fe II]λ12570 and [P II]λ11886 emission lines observed in the near-infrared (NIR) J band.Oliva et al. (2001) proposed that the line flux ratio of [Fe II] to [P II] is a good tracer of shocks.The two lines share similar critical densities and ionization potentials, and are expected to be excited in similar physical environments.By contrast, dust depletions are very different (Hobbs et al. 1993); iron is abundant in dust whereas phosphorus is almost absent.When a shock breaks dust, [Fe II] appears much stronger than [P II].Indeed, it is known that [Fe II] is much stronger in supernova remnants than [P II], where shocks are present, whereas both have about the same flux in the Orion bar, where there are no known shocks (Oliva et al. 2001).T16 investigated the line ratios of 44 Seyfert galaxies (including some data taken from literature) and found that three of them (NGC 2782, NGC 5005, and Mrk 463) have [Fe II]/[P II] ratios larger than 10, whereas more than half of them have a smaller one (∼ 2).They also found that the line ratios are not correlated with the radio loudness or starburst intensity, suggesting that the possible mechanism to trigger the shock is the AGN outflow.
The goal of this paper is to study the validity of our scenario that the AGN pc-scale outflow generate the shock in the NLR.We describe our sample selection, observation setup, and data reduction in Section 2. The obtained line profiles are explained in Section 3. We discuss how the shocks are generated in NLRs in Section 4 and summarize our result in Section 5.

Target selection and observations
The observations were carried out with the WINERED 1 (Warm INfrared Echelle spectrograph to Realize Extreme Dispersion and sensitivity) spectrograph (Ikeda et al. 2016(Ikeda et al. , 2022) ) mounted on the New Technology Telescope (NTT) at the La Silla observatory as a visitor instrument.WINERED has a bandwidth of 9100-13500 Å.The slit length is 16.34 arcsec and the pixel scale is 0.27 arcsec/pixel.The slit width was set to 1.08 arcsec (=4.0 pixel; WIDE mode) with R = 18000 (∆v = 17 km s −1 ).The observations were carried out for five nights between 2017 November 30 and 2018 March 3, with seeing of 0.7-1.5 arcsec.The X-ray UFO catalog papers (Tombesi et al. 2010(Tombesi et al. , 2014;;Gofford et al. 2013) listed 36 objects.Within them, 11 objects can match our selection thresholds; they are located in the southern sky and have Jmag < 15.5 and z < 0.07, for which the spectral lines of interest ([Fe II]λ12570 and [P II]λ11886) fall in the observable wavelength range.In addition, two Seyfert galaxies, 1H 0707-495 and IRAS 13224−3809, have been reported to host strong UFOs (Hagino et al. 2016;Parker et al. 2017).Thus we also observed them; the total number of targets is 132 .We also note that all of them are radioquiet (log R < 1, where R is the ratio of the flux density at 5 GHz to 4400 Å), that is, we can exclude potential energy inputs from radio jets and to study the influence of multi-scale outflows.Tables 1 and 2 show the target list and the observation log, respectively.
When the host galaxy emission was anticipated to potentially occupy most of the slit length of 16.34 arcsec, we used an Object-Sky-Object mode ("OSO" mode), that is, the blank sky was observed between the object observations.Any bright sources were not observed immediately before the targeted observations to avoid persistence (see section 5.3 in Ikeda et al. 2022).Where necessary, we waited until the persistence gets sufficiently low before we start observing the targets.
When a target did not occupy the slit length, we adopt the "ABBA" mode.In this mode we position the object at a location offset from the center of the slit (designated as A-position), and then, for the next exposure, position it at a location shifted in the opposite direction (designated as B-position).By subtracting the spectra obtained from these two positions (that is, A-B or B-A), the contribution of the sky background is removed and we get an image with two spectra, one positive and one negative.In order to discuss NLR in AGN with kpcscale spatial extent, only those with the throw distance (distance between the centroids of A and B exposures) more than 2 kpc are treated in this mode.Data reduction was performed using IRAF3 and PyRAF4 .The detailed explanation of the data reduction (e.g., wavelength calibration, dispersion calibration, flat fielding) in WINERED are described in Section 8 in Ikeda et al. (2022).Hereafter data reduction unique to this paper are explained.

Data reduction
The source region was selected by the following steps.First, a spatial profile along the slit length was plotted and the peak was searched (red vertical lines in Figure 1).Next, the spatial profile around the peak was fitted with a single Gaussian to calculate a Full Width at Half maximum (FWHM).Last, the source region was selected from −1×FWHM to +1×FWHM (red-shaded regions in Figures 1).The sky spectrum was set to have the same width as the source one.We used Skycorr (Noll et al. 2014) to subtract the sky spectrum under the situation that the infrared background can change substantially.The sky-subtracted spectra were normalized to unity for the continuum, using a continuum command in PyRAF.In the continuum command, the spectra were normalized by a cubic Legendre polynomial, without emission/absorption-line-dominant bands.As a next step, telluric absorption lines were corrected by a telluric command, using A0-type standard stars with similar airmass observed before or after a suit of observations (Table 2).In the telluric command, we determined the optimal value of the scale and wavelength shift of the telluric standard spectrum after removing intrinsic lines (see Sameshima et al. 2018 for detailed method of removing the intrinsic lines) and divided it from the target spectrum.
It is imperative to prevent the extended signal of a galaxy from contaminating the sky background, which is then subtracted from each of the pair of A and B positions.The targets presented in this study possess low redshift; therefore, it is imperative to prevent the extended signal of a galaxy from contaminating the sky background, which is then subtracted from each of the pair of A and B positions.Column (5) in Table 3 provides insights into the size of the spectral extraction region, a product of Gaussian fitting as illustrated in Figure 1.The spatial extension of the emission region along the slit length falls within the range of 9-13 pixels for the "ABBA" mode.Although it is worth noting that this spatial extent is determined using the NIR continuum and may not be entirely congruent with the NLR extent, which should ideally be assessed via narrow lines, we maintain confidence that over-subtraction was effectively mitigated in the A-B sky subtraction process.This confidence is supported by the fact that most of the emission resides within this defined region, and the spatial extent (9-13 pixels for A-B subtraction) remains significantly smaller than the throw distance of 45 pixels by a factor exceeding three.Figure 2 visually exemplifies the spectral extraction region and the sky region, reaffirming that extended emission in the Aposition does not contaminate the B-position, and vice versa.
Since we are interested in the line intensity ratio and the shape of the line profile, we did not observe the flux standard star for measuring the absolute value of flux.Then, all the exposure were stacked and the continuum of the total spectrum was normalized again.As the last step, the barycentric correction was performed.We set the number of independent chains (walkers) exploring the parameter space (nwalker) to be 512, and the number of iterations performed by each walker in the MCMC process (niter) to be 3000.The fitting procedure is as follows: (1) One positive Gaussian was introduced.The line width (standard deviation; σ v ) was set to be larger than the velocity resolution of 17 km s −1 .Because the line center may be slightly shifted due to gas motion, the central velocity was set to be free within ±200 km/s from the rest frame.We summarize the line intensity in Table 4 as an equivalent width (EW).Its error (1σ) was calculated from the Bayesian posterior probability (See Figure 11 in Appendix B for corner plots of the MCMC fitting).The adjacent [S IX] emission line was simulta-neously modeled for the fitting of [Fe II] in IC 4329A, because the fitting model was distorted without it.(2) If the line intensity did not exceed three times its error, the line was considered undetectable and the three times of the error was adopted as the upper limit.(3) If some residuals were seen after introducing the first Gaussian and they were overlapped with the line of interest, additional Gaussians were introduced ([Fe II] in NGC 4507 and [P II] in NGC 6240).To determine whether the additional Gaussian was necessary, an F test (Bevington 1969;Eadie et al. 1971) was performed based on the chi squared and degrees of freedom before and after the Gaussian was added.The additional one was introduced only when it was judged to be necessary at a 99.7% confidence level.This procedure was not performed when the residuals are apart enough to be unaffected by fitting of the line of interest (e.g., [Fe II] in NGC 3783).
The EW ratios of [Fe II] to [P II], as well as the obtained EW and σ v of each line, are listed in Table 4.We detected the [Fe II] line in 12 out of 13 AGNs and [P II] in 6 of them.Neither of the lines was detected in NGC 1365.
While T16 calculated the [Fe II]/[P II] ratios using the line fluxes, we use the EWs.Our ratios must be F λ, [Fe II] /F λ, [P II] times larger than those of T16, where F λ, k refers to the continuum flux density at the k line.The NIR power-law slope in the wavelength unit (α λ , where F λ ∝ λ α λ ) in the quasar composite spectrum is known to be α λ = −0.19(Glikman et al. 2006).If the slope in Seyfert galaxies is similar to those in quasars, the factor is F λ, [Fe II] /F λ, [P II] = (1.257/1.188)−0.19 = 0.99.Therefore the difference of the line ratio definitions by us and T16 hardly affects the results.
We found that four sources (Ark 120, ESO 323-G77, MCG −5-23-16, and NGC 4507) have the [Fe II]/[P II] ratios larger than 10; notably, for NGC 4507, [P II] was detected and the ratio exceeded 10, which is the first report.They are likely to have shocks in the NLR.This result doubles the number of the known AGN samples with the line ratio larger than 10, where the most major previous work was done by T16.

The [S III] line
We use the [S III]λ9533 line as a tracer of the ionized outflow, instead of the conventional optical [O III]λ5007 line.The two lines share the same transition ( 3 P 2 -1 D 2 ) and belong to the same family in the periodic table.It is predicted that they are emitted from the similar region in the AGN (see e.g., Liu et al. 2015).We also note that sulfur is an almost perfect volatile element since it is not locked in dust and thus is suitable for tracing the pure outflow component.Figure 4 shows the obtained   which is shown only when the line is modeled by a single Gaussian.( 4) and ( 5), ( 7) and ( 8), ( 9) and ( 11) are the same as ( 2) and ( 3), but for the [P II] line, the central peak of [S III], and the blue wing of [S III], respectively.(10) Velocity shift for the line center of [S III] blue wing from the central peak.[S III]λ9533 line profiles.The fitting procedure is almost the same as the one in Figure 3, but the line center is allowed to be ±1000 km s −1 because the gas motion traced by the [S III] may be larger.The first Gaussian to be introduced, whose line velocity is almost zero in the rest frame, is referred to as the "main" one, and the additional ones are as "sub".If the center velocity of the "sub" Gaussian was negative (within a range of −1000 km s −1 to 0 km s −1 ) at a 99.7% confidence level, it was interpreted as the "blue wing" and we categorized that the object has the "ionized outflow".Consequently, in 6 out of the 13 AGNs (Ark 120, IRAS 13224−3809, MCG −5-23-16, NGC 1068, NGC 4507, NGC 7582), the "blue wing" is detected.Table 4 lists the fitting results, including the velocity shifts, σ v , and EWs.

Absorption features around the [P II] line
Some absorption features were apparent in some AGNs around the [P II] line.In particular, those at v = −1400 km s −1 and +200 km s −1 (corresponding to 11830 Å and 11895 Å, respectively), were observed in multiple sources: the 11830 Å line in MCG −5-23-16, MCG −6-30-15, NGC 1068, NGC 1365, NGC 4507, and NGC 7582, and the 11895 Å line in ESO 323-G77, NGC 1365, and NGC 7582.These lines cannot be the remaining counts of the telluric absorption because some of them are seen in the wavelength range in which any potential remaining effect is small.They were broader than the telluric lines and instrumental resolution.Therefore they are not artificial.They were located at a similar wavelength in the rest frame of each AGN, implying that they are associated with AGNs and not the Galactic absorption.Interestingly, they were present mostly in Seyfert 2 galaxies; the exceptional cases are MCG −6-30-15 (which has complex absorbers close to the central SMBH; Lee et al. 2001) and NGC 1365 (whose broad Hα line is extremely faint and classified as Seyfert 1.9 Trippe et al. 2010.This implies that most of the broad line region may be covered).Consequently, we suggest that these absorptions are due to molecular gas in the torus or some obscuring material (e.g., polar dust; Hönig et al. 2013) in the vicinity of AGNs.The line identification and search for other features will be presented in our forthcoming paper.

Spatial distribution of the [Fe II]/[P II] ratio in NGC 1068
NGC 1068 is the brightest and most nearby targets and thus suitable for studying position dependence on the [Fe II]/[P II] ratio; it was reported that, whereas the [Fe II]/[P II] ratio is small (∼ 1) near the central core, it increases up to ∼ 7 at a distance of 2-4 arcsec from the center (Hashimoto et al. 2011;Riffel et al. 2014).In order to confirm such spatial dependence in our observation, we created spectra divided in the slit-length direction and calculated the EWs of the two emission lines.The continuum peak position was set as the origin, and the spectra were extracted from every 5 pixels (= 1.35 arcsec; almost equivalent to the seeing).The position angle was 220 deg and the slit was oriented roughly from southeast to northwest.Eleven spectra were derived, covering 7.5 arcsec on each side in total.Whereas the [Fe II] emission lines were detected in all of them, the [P II] emission lines were detected only in the three central ones.The EW values or their upper limits were calculated in the same way as in section 3.1.The width of the undetected [P II] lines used for the upper limit was the same as that for the spectrum at the central bin.The obtained spatial dependence of the [Fe II]/[P II] ratio is shown in Figure 5.As in the previous studies (Hashimoto et al. 2011;Riffel et al. 2014), we confirmed that the ratio is smaller at the center and higher at a projected distance of ∼ 3 arcsec.This means that photo-ionization is dominant in the vicinity of the core in NGC 1068, while the shock ionization appears farther from the core.an EW of ∼ 1.8 Å (n.b. this value is not very accurate because it is read from the spectrum figure).On the other hand, our study can detect the emission line with EW = 0.49 Å (NGC 4507), which is 3.6 times fainter than the detection limit of T16.For example, if the detection limit of the [P II] line were 1.8 Å in NGC 4507 like T16, then the [P II] line could not be detected and the ratio value had only a lower limit of >6.7.In this case, we could not conclude whether the shock exists.Consequently, our observation provides a more robust conclusion on the existence of shock than the previous study.
Figure 7 compares the [Fe II]/[P II] ratios with the radio loudness.The result of T16 is also superposed.T16 found that no significant difference exists between the average value of the [Fe II]/[P II] ratios in the radioloud AGNs (R > 10; the average value is 3.67 ± 1.09) and that in the radio-quiet AGNs (R < 10; the average value is 2.87 ± 1.38) for the [P II]-detected objects, and stated that the radio jet is unlikely to be the main contributor to the shock in the NLR.We performed a Kolmogorov-Smirnov (KS) test to quantitatively compare the distributions of the [Fe II]/[P II] ratios between this study and T16.Only the [P II]-detected sources were used.The obtained p-value was 0.562, so that We cannot discuss any significant differences in this sample.The notations are the same as in Figure 6.The results of T16 are superposed in gray.We find no correlation between the existence of shocks and radio loudness.

UFO contribution
To study the UFO contribution to the shock, we compare the UFO mass loss rate and kinetic power with the [Fe II]/[P II] ratios in Figure 8.Using the four objects for which both UFO parameters and the [Fe II]/[P II] ratios were determined (NGC 7582, NGC 1068, IRAS 13224-3809, NGC 6240), we examined their correlations.The Pearson correlation coefficient values and p-values were -0.465 and 0.535 for log mass loss rate, and -0.368 and 0.632 for log UFO power, respectively.Thus, no correlations were detected between them.
While the UFO is considered to sweep through the ISM and produce a kpc-scale outflow (e.g., Tombesi et al. 2015;Marasco et al. 2020), the efficiency of the kinetic power transfer from the UFO to the kpc-scale outflow has a wide variety, depending on the flow type; notably, it can be in the energy conserving or momentum conserving modes (see, e.g., Mizumoto et al. 2019).If these two types exist in the sample, no correlation is observed even if UFOs have some (indirect) contribution to kpc-scale outflows and shocks.There is a possibility that if only AGNs with energy-conserving flows were selected, a correlation between UFOs and NLRscale shocks might be identified.Our sample size is, however, too small to conduct this type of study.Determination of the flow type from sub-pc to kpc will be required for a direct assessment of the effect of UFOs to shocks.

Correlation with the [S III] outflow
Figure 9 shows the obtained relation between the [S III] outflow velocity, that is, the velocity shift in the center of the blue wing from the rest frame if the source has an outflow, and the line ratio.We use the off-centered value for NGC 1068.There is a hint that targets associated with [S III] outflows have a larger [Fe II]/[P II] ratio, that is, they are more likely to produce shocks.
As discussed in section 3.4, the [Fe II]/[P II] ratio of NGC 1068 is higher (∼7) at an angular distance of ∼ 3 arcsec, whereas small (< 1) within 1 arcsec from the center.In addition, integral field unit (IFU) observations revealed an annulus of r ∼ 2 arcsec with high [Fe II]/[P II] ratio (∼5-8), whereas an inner circular region of r < 1 arcsec has a small ratio (≲ 2).This kind of off-centered, high [Fe II]/[P II] regions have been also found in other AGNs (e.g., Storchi-Bergmann et al. 2009;Schönell et al. 2019;Riffel et al. 2020) 5).The other notations are the same as in Figure 6.ure 9).Consequently, all the AGNs with the NLR-scale outflows are found to be high [Fe II]/[P II] ratio of ≳ 7 in our sample, and thus the ionized outflows are considered to play a dominant role to the shock ionization.
In addition, we notice that the [Fe II]/[P II] peak region in NGC 1068 reported by Riffel et al. (2014) and the [O III] spatial distribution shown by Das et al. (2006) share a similar size scale.The [O III] emission region is extended toward the northeast direction with a radius of 0-2 arcsec (Das et al. 2006).This means that the ionized outflow is blown from the AGN core in this direction.On the other hand, the annulus with high [Fe II]/[P II] ratio (r ∼ 2-3 arcsec) has asymmetry and the shock is more prominent in the northeast part (Riffel et al. 2014).These geometry matches our scenario that the AGN ionized outflow triggers the shock ionization, shown in Figure 10.Based on this picture, the kinetic energy or the momentum flux of the ionized outflow triggered by AGN are transferred to the gas in a several-hundred-pc-scale region via the shock.In other words, the AGN feedback (i.e., AGN affects to the host galaxy via the outflow) occurs in the NLR via the ionized outflow and the shock.
The exception is ESO 323-G77, which is not associated with the [S III] outflows but has a ratio larger than ten.The [S III] emission line of ESO 323-G77 shows an asymmetric profile, which is present only in this source among our sample (Figure 4).If the asymmetry originates in the outward motion of the NLR cloud in a dusty region (Peterson 1997), such a bulk motion may be responsible for some shock ionization.Another possi- ble scenario is that a potential dusty wind, the presence of which was suggested based on the detection of dust emission extending axially in this source (Leftley et al. 2018), generates the shocks.In any case, some outflowing gas may exist in ESO 323-G77.

Future work
In our observations, some AGNs exhibit [P II] lines too weak to be detected, necessitating the calculation of upper limits on the [Fe II]/[P II] ratios only.In these cases, firm conclusions regarding the presence of shocks for these targets are elusive.Therefore, it is neccesary to employ multiple diagnostic methods to confirm the existence of shocks.For instance, Figure 13 et al. 2009;Riffel et al. 2014;Fazeli et al. 2019;Schönell et al. 2019;Riffel et al. 2020Riffel et al. , 2021)).For instance, Schönell et al. (2019) presented [Fe II]/[P II] maps for six nearby AGNs, highlighting that NGC 3227 and NGC 5899 exhibit peak values offset from the nuclear position, aligning with the high velocity dispersion of [Fe II].How-ever, the sample size remains limited, underscoring the significance of our study in expanding the roster of suitable candidates for such observations.
The schematic diagram depicted in Figure 10 closely aligns with theoretical shock models (King 2010;Zubovas & King 2012;Faucher-Giguère et al. 2012;Faucher-Giguère & Quataert 2012;King & Pounds 2015;Richings & Faucher-Giguère 2018a,b;Zubovas & King 2019).In these models, the inner radiative-driven wind interacts with the ISM, giving rise to a shock front.This shock front subsequently sweeps up the ISM, generating a region of shocked ISM.The nature of this shock hinges on whether cooling, typically radiation, significantly dissipates energy from the hot shocked gas on a timescale shorter than its flow time.If cooling exerts a strong influence in this context, most of the pre-shock kinetic energy is lost via radiation, characterizing a momentumdriven flow.Intense cooling results in substantial compression of the shocked gas, yielding a geometrically narrow post-shock region.Conversely, in scenarios where cooling is negligible, the post-shock gas retains all the mechanical energy imparted by the shock and undergoes adiabatic expansion into the ISM.In contrast to the momentum-driven (isothermal) case, the post-shock gas in energy-driven flows is characterized by a geometrically elongated structure.Given that energy-driven flows are significantly more intense than momentumdriven flows, the assessment of radiative outflow's impact on the surrounding environment (in this case, the NLR) hinges on the type of gas flow described above.To unravel these complexities, it is essential to measure momentum fluxes before and after the shock, necessitating an accurate assessment of spatial scales.Therefore, spatially-resolved observations, as outlined in the preceding paragraph, assume paramount importance.

CONCLUSION
We performed NIR slit spectroscopy of 13 Seyfert galaxies for which UFOs had been observed in X-rays, and measured their [Fe II]/[P II] ratio, with the aim of determining whether AGN outflows can generate shocks in the NLR.The [Fe II] and [P II] lines were detected in 12 and 6 sources, respectively, among them.From the derived ratios, we found that shocks are likely to be present in 4 objects, Ark 120,and NGC 4507.No direct relation between the presence/absence of the shocks and UFOs is found.On the other hand, we found that when the [Fe II]/[P II] ratio is higher than ∼ 7, the [S III] emission line shows signs of the ionized outflow.We favor a scenario that the ionized outflows trigger the NLR-scale shock.The scenario had been evaluated so far only for NGC 1068, using the IFU observation, before this work.This work revealed that the scenario is likely to be applicable for a wider range of sample.
wing (v = −481 ± 4 km s −1 and σ v = 583 ± 16 km s −1 ), and the [S IX] line.The [S III] line also has the blue wing at −100 ± 20 km s −1 with σ v = 420 ± 120 km s −1 .Considering that the [P II] line has no hint of blue wing, some shock due to the outflow is expected in this target.We also note that a strong absorption is present in the [P II] band (−1400 km s −1 ).
A.12. NGC 6240 NGC 6240 is a merging luminous infrared galaxy and has two cores.In this study, the brighter southern core was observed in the slit.
(3) Airmass (minimum-maximum) (4) Observation mode.A/B shows the A-postion and B-position, and O/S shows the object and sky observation, respectively.(5) (Exposure time in sec for each snapshot)×(the number of snapshots).(6) Name of the telluric standard star.(7) Airmass of the tellluric standard star.(8) Seeing of the telluric standard star.

Figure 1 .
Figure 1.Spatial profiles of the emission along the slit.The upper panel shows the NGC 1068 case, in which the blank sky was observed between the object observations for the sky subtraction.The profile shows the one in the object observation.The source region at the center is automatically selected with Gaussian fitting.The lower panel shows the MCG-5-23-16 case, in which the image of A-position was subtracted from the image of B-position and the onedimensional spatial profile around the slit was calculated.The positive peak at a positive distance (right-hand side in the panel) is the signal at the A-position, whereas the negative one in the left hand is for B. The spectral extraction region (red) is automatically selected, as in the upper panel.The peak at the center is the remaining counts from a bright telluric standard star observed shortly before this target observation.

Figure 2 .
Figure 2. The slit viewer image of MCG-5-23-26 (Aposition).Please note that the left and right sides are reversed from Figure 1.The white, yellow-shaded, and orangeshaded rectangles show the WINERED slit, the spectral extraction region (A-position), and the sky region (B-position), respectively.We note that the yellow-shaded rectangle corresponds to the red-shaded region in Figure 1.The signal on the left is a serendipitous point source (2MASS J09473861-3057050).
The [Fe II] and [P II] lines The obtained line profiles of [Fe II] and [P II] are shown in Figure 3.The line fitting was performed by the Markov chain Monte Carlo methods (MCMC) using the emcee module (Foreman-Mackey et al. 2013) in Python.

Figure 3 .
Figure 3. Velocity profiles of the [Fe II] and [P II] emission lines.Positive velocities mean redshift, and vice versa.The red line in each panel shows the model (continuum + Gaussian) when the parameter with the maximum likelihood in the sampling is adopted.The orange shade shows the 1σ posterior spread from the sampling median for each wavelength bin.The middle and bottom panels show the residuals and the telluric absorption, respectively.The blue bar at the upper left shows the typical error on the continuum.The gray shows the spectrum in the adjacent Echelle order.The adjacent [S IX] line is labeled.In the line detection case, its significance is written at the upper right.(a) 1H 0707−495, (b) Ark 120, (c) ESO 323-G77, (d) IC 4329A, (e) IRAS 13224−3809, (f) MCG −5-23-16, (g) MCG −6-30-15, and (h) NGC 1068.

Figure 3 .
Figure 3. Continued.In NGC 4507, the adjacent [S IX] line is introduced as the gray solid line.(i) NGC 1365, (j) NGC 3783, (k) NGC 4507, (l) NGC 6240, and (m) NGC 7582.When an additional Gaussian is necessary to introduce, each Gaussian component is shown in the blue dotted line.

Figure 4 .
Figure 4. Velocity profiles of the [S III] emission lines.The gray solid line shows the spectrum in the adjacent Echelle order.Gaussians whose central velocity is v = 0 within 3σ are shown in the green dashed lines, whereas those on the blue side ("blue wing") are in the blue dotted lines.The exceptional case is MCG −6-30-15, which shows an excess on the red side (yellow dot-dashed line).The other notations are the same as those in Figure 3.

Figure 5 .
Figure 5.The [Fe II]/[P II] ratio in NGC 1068 as a function of projected distance from the core, i.e., the continuum peak.The positive direction was determined to be northwest.
Figure 6 shows the relation between the EWs of the two emission lines, [Fe II] and [P II].All the sources for which [P II] is not detected show only a weak [Fe II] line.This suggests that non-detection of [P II] is simply because the narrow lines are weak regardless of the species, rather than due to the existence of shock.The detection of weak [P II] emission lines is critical for accurate measurement of the [Fe II]/[P II] ratio.For example, the weakest [P II] line detected in T16 (NGC 5506) has

Figure 6 .
Figure 6.Relation between the EWs of [Fe II] and [P II].Red points show the targets from which both the [Fe II] and [P II] lines are detected, whereas blue points do those with no detection of [P II].Dashed lines show [Fe II]/[P II]=2, 10, and 20.The yellow-shaded area shows the shock region.

Figure 7 .
Figure 7. Relation between the [Fe II]/[P II] ratio and radio loudness.The notations are the same as in Figure6.The results of T16 are superposed in gray.We find no correlation between the existence of shocks and radio loudness.

Figure 8 .
Figure 8. Relation between the [Fe II]/[P II] ratio and UFO mass-loss rate (in upper panel) and kinetic power (in lower panel).The notations are the same as in Figure 6.No correlation between the existence of shock and the UFO parameter.
. If we use the [Fe II]/[P II] ratio in the off-centered region (see Figure 5), the ratio is similar to the ones in the other objects with the [S III] outflows (see the open red circle in Fig-

Figure 9 .
Figure 9. Relation between the [Fe II]/[P II] ratio and the [S III] outflow velocity.The outflow velocity of the targets with no [S III] outflow is set to 0. The open red circle is for NGC 1068 but off-centered (where the projected distance is ∼ −3 arcsec; see Figure5).The other notations are the same as in Figure6.

Figure 10 .
Figure 10.Schematic picture of the ionized outflow, NLR, and shock region.In the NGC 1068 case, the ionized outflow is traced in the [O III] and [S III] lines.The outer radius of the ionized outflow matches the radius of the shock-dominant annulus.
in Calabrò et al. (2023) introduces various diagnostic diagrams for shock identification.Among these, the [S III] 9530 Å/Paγ vs [Fe II] 1.257 µm/Paβ diagram can be constructed based on our WINERED observations.Our subsequent publications will provide comprehensive spectral data, including full spectra, line profiles, and equivalent widths (EWs), encompassing Paβ and Paγ.Spatially-resolved observations of [Fe II], [P II], and [S III] play a pivotal role in elucidating the ionization mechanism within the NLR and facilitating direct observations of energy transfer from the AGN to kpc-scale regions.Some prior studies have generated [Fe II] and [P II] maps within the NLR context (Storchi-Bergmann Both the [Fe II] and [P II] lines are strong and broad.In particular, the [P II] line profile required two Gaussian component with v = 400 ± 80/ − 440 ± 30 km s −1 and σ v = 850 ± 50/370 ± 40 km s −1 , respectively.The [S III] line has a similar line broadening to the [Fe II] line.A.13. NGC 7582 The absorption features are observed at around the [P II] line at −1400 km s −1 , +200 km s −1 , and +1000 km s −1 .The [S III] blue wing is detected at −120 km s −1 .B. CORNER PLOTS FOR THE [FE II] AND [P II] LINE FITTING We fit the [Fe II] and [P II] lines with Gaussians.Figure 11 displays the corner plots in the MCMC fitting shown in Figure 3.The Gaussian amplitudes, central velocity, and velocity dispersion are shown.For [Fe II] in NGC 4507 and [P II] in NGC 6240, the results for the Gaussian with the largest EW are only shown.

Figure 11 .
Figure 11.Corner plots for the [Fe II] (upper panel) and [P II] (lower panel) line fitting in Figure 3. See the caption of Figure 3 for the target name.

Table 1 .
Target list

Table 3 .
Pixel scale, spectral extraction region, and throw distance of each target

Table 4 .
Results of the NTT/WINERED observations peak EW [S III] peak σv [S III]wing EW [S III]wing velocity [S III]wing σv