© M. M. Sardaneta et al. 2022

1. Introduction

The environment plays a major role in shaping galaxy evolution (e.g., Peng et al. 2010). Rich clusters are dominated by quiescent systems (morphology-segregation effect; Dressler 1980), and the few spirals that inhabit these environments are generally deficient in all the components of their interstellar medium (ISM; HI – Cayatte et al. 1990; Solanes et al. 2001; Gavazzi et al. 2005; Chung et al. 2009; H₂ – Fumagalli et al. 2009; Boselli et al. 2014a; dust – Cortese et al. 2010, 2012b; Longobardi et al. 2020). Their star formation is also reduced with respect to that of similar objects in the field (Gavazzi et al. 1998, 2002a,b, 2006, 2010, 2013; Lewis et al. 2002; Gómez et al. 2003; Peng et al. 2010; Boselli et al. 2014b, 2016b; Cortese et al. 2021).

Different physical mechanisms have been proposed in the literature to explain these results, as extensively reviewed in Boselli & Gavazzi (2006, 2014). These include gravitational interactions with other cluster members (Merritt 1983) or with the gravitational potential well of the cluster (Byrd & Valtonen 1990) – including their combined effect, generally called galaxy harassment (Moore et al. 1998) – or the hydrodynamic interaction of the ISM of galaxies with the hot and dense intracluster gas (thermal evaporation, Cowie & Songaila 1977; starvation, Larson et al. 1980; ram pressure stripping, Gunn et al. 1972). Several observations seem to indicate that while gravitational perturbations are dominant in relatively small systems such as groups thanks to the low velocity dispersion of galaxies within the dark matter halo, hydrodynamic interactions are dominant in massive clusters, where the high velocity dispersion and the density of the gas make them particularly efficient (e.g., Vollmer et al. 2001; Boselli et al. 2008a,b, 2014a; Gavazzi et al. 2013). It is nevertheless clear that all these different perturbing mechanisms might be acting jointly on galaxies in rich environments, with different effects on their evolution.

Despite this evolutionary picture is becoming clearer, there are still several aspects of the galaxy transformation in rich environments that are poorly known. Among these is the fate of the ISM removed during a ram pressure stripping event. The stripped ISM is often observed in spectacular tails of atomic (e.g., Chung et al. 2007), ionised (e.g., Gavazzi et al. 2001; Yagi et al. 2010; Boselli et al. 2016a), or hot gas (e.g., Sun et al. 2007). Indeed, it has been observed that while in some objects the stripped gas collapses within the tail to form giant molecular clouds (Jáchym et al. 2014, 2017, 2019; Moretti et al. 2020) and stars (Fossati et al. 2016; Poggianti et al. 2019), in other systems it remains diffuse (e.g., Boissier et al. 2012). In these objects the cold gas stripped from the disc does not form stars because it changes phase, becoming first ionised then hot (e.g., Boselli et al. 2016a). The detailed analysis of a few nearby systems with high quality multi-frequency spectroscopic and imaging data seems to suggest that the gas collapses into giant molecular clouds only when its velocity dispersion is sufficiently low (e.g., Fossati et al. 2016; Boselli et al. 2018b). Otherwise, the cold gas changes phase for heat conduction, turbulence, magneto-hydrodynamic waves, or shocks. It thus becomes first ionised and then hot gas visible in the X-ray spectral domain (Gavazzi et al. 2001; Sun et al. 2007; Fossati et al. 2016; Boselli et al. 2016a).

The kinematics of the stripped gas is thus a crucial parameter for understanding the process of star formation within the tails. Observations and simulations of ram pressure stripped galaxies consistently indicate that a ram pressure stripping process is able to drastically perturb the velocity field of a spiral galaxy, displacing its kinematic centre with respect to the photometric one and affecting its rotation curve (RC; Vollmer et al. 2001, 2004, 2006, 2008; Kronberger et al. 2008; Merluzzi et al. 2013; Fumagalli et al. 2014; Consolandi et al. 2017; Sheen et al. 2017; Bellhouse et al. 2017, 2019; Boselli et al. 2022). Despite the advent of extraordinary new instrumentation in the optical (Multi Unit Spectroscopic Explorer (MUSE) at the Very Large Telescope (VLT)) and radio (Atacama Large Millimeter/submillimeter Array – ALMA) domain, multi-frequency data at high angular and spectral resolution of nearby ram pressure stripped galaxies, where the proximity allows us to study the star formation process down to the scale of individual HII regions, are only sporadically available. To understand the effects of the stripping process on the kinematics of the gas, observations and simulations of galaxies spanning a wide range in stellar mass and with different impact parameters are necessary. A particularly interesting case are edge-on galaxies, which allow us to observe as the gas is stripped out from the disc plane in the z direction. Among these, there are only two systems undergoing a ram pressure stripping event that have been studied in detail: UGC 6697 in the cluster A1367 (Consolandi et al. 2017) and JO204 in the farther cluster, A957 (z ∼ 0.052; Gullieuszik et al. 2017). The kinematic properties of the stripped gas are significantly different in the two galaxies given their impact parameters, UGC 6697 suffering an edge-on stripping event and JO204 a face-on one.

The Virgo Environmental Survey Tracing Ionised Gas Emission (VESTIGE; Boselli et al. 2018c), a blind Hα narrow-band imaging survey of Virgo, allows us to identify a large number of galaxies undergoing a ram pressure stripping event in the closest cluster of galaxies. The presence of ionised gas tails without any stellar counterpart is the clearest demonstration that a system is undergoing a hydrodynamic interaction with the surrounding environment. Indeed, gravitational perturbations act indifferently on the stellar and gaseous components, also producing tidal structures visible in deep broad-band imaging data. Several galaxies undergoing a ram pressure stripping event have been identified so far (Boselli et al. 2016a, 2018a, 2022; Fossati et al. 2018; Longobardi et al. 2020; Junais et al. 2021). Among them, NGC 4330 is an edge-on galaxy with evident tails of atomic (Chung et al. 2007, 2009) and ionised (Abramson et al. 2011; Fossati et al. 2018) gas, which indicates that the galaxy is entering the cluster from the south and impacting the intracluster medium (ICM) almost face on. It is thus an ideal candidate for studying in great detail the effects of the perturbation on the kinematics of the perturbed gas on the disc of the galaxy and within the tail.

In this paper we use two different sets of high spectral resolution (up to R ≃ 10000) Fabry-Perot (FP) data to study the kinematics of the ionised gas in this representative object. We combine these new sets of observations with tuned hydrodynamic simulations developed to reproduce the structural and spectrophotometric properties of this galaxy (Vollmer et al. 2021). Unfortunately, most of the gas in the tail has a very low surface brightness, and thus it is out of reach for our dedicated observations, which detect only a few high-surface-brightness star-forming regions. The proximity of the galaxy (16.5 Mpc; Gavazzi et al. 1999; Mei et al. 2007), however, and the adopted spectral resolution allow us to study in detail the kinematics of the gas of these compact regions down to scales of ∼100 pc in angular resolution and velocity dispersion σ ∼ 13 km s⁻¹. Such resolution is not reached for the more distant UGC 6697 and JO204 despite their extraordinary MUSE data.

The paper is structured as follows: in Sect. 2 we describe the MegaCam Canada-French-Hawaii Telescope (CFHT) narrow-band imaging and the FP high-resolution spectroscopic observations along with the multi-frequency data used in the analysis. In Sects. 3 and 4 we analyse the imaging data and construct the RC of the galaxy, while we study the two-dimensional kinematics of the gas in Sect. 5 and the specific kinematics of the gas out of the plane in Sect. 6. The results of the analysis are compared to the hydrodynamic simulations of Vollmer et al. (2021) in Sect. 7. The discussion and conclusions are given in Sect. 8.

2. Observations and data reduction

2.1. VESTIGE narrow-band imaging

Deep Hα narrow-band imaging of NGC 4330 has been gathered with MegaCam at the CFHT by the VESTIGE survey (Boselli et al. 2018c), and the first imaging results have been presented in Fossati et al. (2018). Briefly, the galaxy has been observed using the narrow-band filter MP9603 (λ_c = 6591 Å; Δλ = 106 Å), which includes the Hα line and the two [NII] lines at λ 6548 Å and λ 6584 Ź. The galaxy has been observed during the main VESTIGE survey but also during a pilot run. The data used in this work are the combination of these two datasets, and are thus deeper than those analysed in Fossati et al. (2018), which were limited to the pilot observations. To secure the subtraction of the stellar continuum, the galaxy has been also observed in the broad r-band filter with shorter exposures to avoid saturation within the nucleus. The different exposures were gathered using a large dither pattern (15 arcmin in RA and 20 arcmin in Dec) to minimise the effects of scattered light and reflections in the flat-fielding, thus optimising the detection of low-surface-brightness features in the ionised gas component. A slightly different pattern was used for the pilot observations. The total integration time is 21060 sec in the narrow-band filter and 1506 sec in the broad-band r filter². These deep observations have a typical sensitivity of f(Hα)≃2 × 10⁻¹⁷ erg s⁻¹ cm⁻² (5σ) for point sources and Σ(Hα) ≃10⁻¹⁸ erg s⁻¹ cm⁻² arcsec⁻² (1σ after smoothing the data to ∼3 arcsec resolution) for extended sources, thus approximately a factor of 1.7 deeper than the mean sensitivity of VESTIGE. NGC 4330 has been observed during excellent seeing conditions (FWHM = 0.73 arcsec).

The data have been reduced using Elixir-LSB (Ferrarese et al. 2012), a data reduction pipeline especially designed to detect the diffuse emission of extended low-surface-brightness structures formed during the interaction of galaxies with their surrounding environment (e.g., Boselli et al. 2016b). The astrometric and photometric calibrations have been done using the Sloan Digital Sky Survey (SDSS) and the Panoramic Survey Telescope and Rapid Respons System (Pan-STARRS) surveys by means of the MegaPipe pipeline (Gwyn 2008). The typical photometric uncertainty in the data is of 0.02-0.03 mag in both bands. The emission of the ionised gas has been derived after subtracting the stellar continuum as described in Boselli et al. (2019). This is done by combining the r-band image with the g-band frame obtained during the Next Generation Virgo Cluster Survey (NGVS; Ferrarese et al. 2012) to take the colour of the stellar continuum into account (Spector et al. 2012).

2.2. Fabry-Perot spectroscopy

High-resolution Integral Field Unit (IFU) spectroscopy was obtained using two different FP interferometers, PUMA attached at the 2.1 m telescope in San Pedro Mártir, Baja California, Mexico, and GHASP, mounted on the 1.93 m telescope at the Observatoire de Haute Provence (OHP), France. A detailed description of the instruments, of the observations, and of the data reduction is given in Appendix A. Briefly, PUMA is an FP spectrograph with a field of view of ∼10 × 10 arcmin² coupled with a 512 × 512 charge-coupled device (CCD) camera, with pixel size of 1.27 × 1.27 arcsec² and has a spectral resolution R ∼ 6500. GHASP has a field of view of 5.8 × 5.8 arcmin², is coupled with a 512 × 512 Image Photon Counting System (IPCS) with pixel size of 0.68 × 0.68 arcsec², and has a spectral resolution R ∼ 10 000.

2.3. Multi-frequency data

The FP data were combined to those gathered during other surveys of the Virgo cluster for the following analysis, as extensively described in Fossati et al. (2018). Namely, the distribution of the old stellar population has been derived using the very deep i-band image obtained during the NGVS survey of the cluster (Ferrarese et al. 2012). HI data are also available from the Imaging of Virgo Spirals in Atomic Gas (VIVA), gathered at the Very Large Array (VLA) in C-short array configuration (Chung et al. 2009). These data have been obtained with a velocity sampling of ∼10 km s⁻¹ and a spatial resolution of ∼15 arcsec. The same set of data has been combined with other VLA data obtained in D configuration to increase the sensitivity to low column density, extended features by Abramson et al. (2011) with a final velocity resolution of 10.4 km s⁻¹ and an angular resolution of full width at half maximum (FWHM) of 26.4 × 24.0 arcsec². ¹²CO(2 − 1) observations obtained with the Submillimeter Array (SMA), with a synthesised beam of 6.35 × 4.47 arcsec² and a spectral resolution of 5 km s⁻¹, have been presented in Lee et al. (2017).

3. Imaging analysis

3.1. Stellar distribution

The NGVS i-band image (see Fig. 1) traces the bulk of the stellar population and can thus be used to derive the main morphological and structural parameters of the galaxy. It can also be used as reference to identify the different ionised gas structures removed from the disc plane during the interaction. This image shows an edge-on disc with a symmetric and unperturbed stellar distribution down to the surface brightness limit of the NGVS survey (μ_i = 27.4 mag arcsec⁻²), proving that the galaxy is not undergoing any major gravitational perturbation. The galaxy has a small bulge without any evident bar (e.g., Abramson et al. 2011).

Fig. 1. Photometric properties of NGC 4330. Left panels: NGVS i-band image at low contrast (upper panel) to see the distribution of the absorbing dust within the disc, and at high contrast in logarithmic scale (lower panel) to see the boxy shape of the disc at its edges. The red and blue ellipses indicate the i-band mean surface brightness μ(i) = 23.5 and the effective surface brightness μe = 21.5 mag arcsec−2, respectively. The green bullet and the red cross indicate the position of the photometric and kinematical centres, respectively. Right panels, from top to bottom: radial variation in the i-band surface brightness (μ), position angle (PA), ellipticity (ε), morphological inclination (iph), and the fourth degree cosine coefficient of the Fourier series (B4) parameters. The red and blue vertical lines correspond to the ellipses drawn in the upper left panel and indicate the radius R23.5(i)∼145.0 arcsec (∼11.6 kpc) and the effective radius Re(i)∼66.5 arcsec (∼5.3 kpc), respectively.

We derived the photometric properties of the stellar disc through the elliptical isophote fitting method using the IRAF task ELLIPSE (Jedrzejewski 1987). In this task, ellipses are drawn to match the isophotes. The surface brightness along the ellipse is expanded in Fourier series:

$$\begin{array}{c}\hfill \mu (\theta )={\mu }_0+{\mathrm{\Sigma }}_{n=1,\infty }(A_{n}sinn\theta +B_{n}cosn\theta ),\end{array}$$(1)

where μ₀ is the surface brightness averaged over the ellipse, A_n and B_n the higher order Fourier coefficients. If an isophote is a perfect ellipse, all the A_n and B_n coefficient would be exactly zero. The term A₁, B₁, A₂, and B₂ indicate errors in the fitting procedure. The terms A₃ and B₃ give the ‘egg-shaped’ isophote but the most interesting is B₄. If B₄ > 0, the galaxy has discy-shape isophotes, while if B₄ < 0 boxy-shaped isophotes. The ELLIPSE task provides a number of parameters that describe the surface brightness of the galaxy as a function of the semi-major axis, including the ellipticity (ε), the position angle (PA), and it reconstructs the growth curve that gives the total apparent magnitude inside each isophote (e.g., Barbosa et al. 2015; Korsaga et al. 2019). All these parameters were allowed to vary as a function of the radial distance from the galaxy centre.

The isophotal radius derived using the NGVS i-band image is R_23.5(i)∼145.0 arcsec (∼11.6 kpc). This value is slightly larger than the one derived by Cortese et al. (2012a) using the SDSS-i-band image (R_23.5(i)∼140.3 arcsec, ∼11.2 kpc). The difference with this previous measurement might be due to the presence of dust and to the adopted sigma clipping selection in the ELLIPSE task. The effective radius is R_e(i)∼66.5 arcsec (∼5 kpc) and the effective surface brightness μ_e(i)∼21.5 mag arcsec⁻², consistently with those given in Cortese et al. (2012a). The photometric centre is located at RA(J2000) = 12^h23^m17.26^s, Dec = +11° 22′07″.9, with an uncertainty of ∼0.9 arcsec (∼76 pc). The sudden jump of the parameters observed in Fig. 1 at a radius of ∼5 arcsec (∼0.4 kpc) is due to the dusty region around the centre.

Consistently with, for example, Barbosa et al. (2015), we estimate the morphological inclination of the galaxy using the relation

$$\begin{array}{c}\hfill cos{i}_{\mathrm{ph}}=\frac{{(1-\epsilon )}^2-q_0^2}{1-q_0^2},\end{array}$$(2)

where q₀ = 0.1 is the intrinsic flattening of a Scd galaxy (Haynes & Giovanelli 1984). Figure 1 indicates that the morphological inclination of the galaxy converges to i_ph = 89.5° and the photometric position angle to PA_ph = 59 ± 1° for R(i) ≥ 75 arcsec (∼6 kpc) (μ_i ≥ 22.75 mag arcsec⁻²). The B₄ parameter, which is negative, indicates that the galaxy disc has a boxy shape outside of the effective radius, R_e shown in Fig. 1. This particular shape becomes more evident in the saturated i-band image and is significantly stronger in the south-west direction. Tuned simulations suggest that this particular shape could be produced after a ram pressure stripping event. The gas removed out from the disc plane in the z direction can alter the shape of the gravitational potential well, modifying the stellar orbits, thus producing a thicker disc (Farouki & Shapiro 1980; Clarke et al. 2017; Safarzadeh & Scannapieco 2017; Steyrleithner et al. 2020). This phenomenon might explain the peculiar morphology of NGC 4330, with a boxy shape particularly evident in the south-west direction, where the stripping is more efficient. Unfortunately, we lack of stellar kinematical data necessary to further probe this scenario.

3.2. Warm gas distribution

The ionised gas distribution of NGC 4330 has been extensively described in Fossati et al. (2018). Here we summarise the main characteristics, referring the reader to that work for a more detailed description. The ionised gas emission mainly comes from the inner regions of the stellar disc, where most of the HII regions are situated (see Fig. 2). The Hα disc is truncated with respect to the stellar disc, indicating that the activity of star formation has been quenched at its edges (Fossati et al. 2018). The north-eastern edge of the star-forming disc is bent to the south in a ‘hook-shaped’ structure (e.g., Abramson et al. 2011), while a few compact and bright HII regions are present just south to the stellar disc at the south-western edge. There is also a very extended, diffuse and low-surface-brightness emission (Σ(Hα)∼ 3 × 10⁻¹⁸ erg s⁻¹ cm⁻² arcsec⁻²) towards the south of the disc, produced by the gas stripped during the ram pressure episode, indicating that the galaxy is moving northwards on the plane of the sky within the cluster. Associated with the diffuse emission, the VESTIGE Hα image shows a few weak compact regions with a surface brightness over Σ(Hα) ∼ 10⁻¹⁷ erg s⁻¹ cm⁻² arcsec⁻², which corresponds to the GHASP detection limit where star formation is probably taking place. Not all the features located out of the plane of the galaxy and above this threshold limit are clearly detected in the monochromatic FP image of NGC 4330, neither in the PUMA data nor in the deeper GHASP data. The bent region at the north-east edge of the disc (labelled A in Fig. 2f) is clearly detected as well as a few regions on the south-west extension of the disc (labelled K, L and M). A few ionised regions have also been detected at the limit of the stellar disc (labelled E, G and F).

Fig. 2. Images of NGC 4330. Panel a: full resolution VESTIGE continuum-subtracted Hα narrow-band image of NGC 4330 in linear scale. Panel b: NGVS i-band image. Panels c,d: PUMA and GHASP FP Hα monochromatic images, respectively. Panel e: NGVS i-band image. The dashed black ellipse and the inner ellipse indicate the surface brightness, μ(i) = 23.5, and the effective surface brightness, μe(i) = 21.5 mag arcsec−2, respectively, and the green bullet indicates the photometric centre (see Fig. 1). Panel f: NGVS i-band image with marks on the zones of interest labelled by letters (A to M) and two parallel dashed lines at r0 ∼ 30 arcsec (∼2.4 kpc), which divide the galactic disc into three main regions: north-eastern, central, and south-western (see Sect. 4). The cyan, red, and blue contours in panels a to d show the VESTIGE surface brightness Σ(Hα) = 0.3, 1.0, and 12.5 × 10−17 erg s−1 cm−2 arcsec−2, respectively. The red contours from panels e and f show the VESTIGE surface brightness Σ(Hα) = 10−17 erg s−1 cm−2 arcsec−2. The blue contours in panel e indicate the surface brightness limit of PUMA, Σ(Hα)∼3.5 × 10−17 erg s−1 cm−2 arcsec−2, while those in panel f show the surface brightness limit of GHASP, Σ(Hα)∼10−17 erg s−1 cm−2 arcsec−2.

4. Rotation curves

Rotation curves of galaxies are generally derived using the intensity-peak method (IPM). However, this method hardly works for galaxies with high inclination because the velocity profiles in a given position integrate the whole flux along the line of sight (LoS), thus at different radii on the galactic plane (Sancisi & Allen 1979; Swaters et al. 1997; García-Ruiz et al. 2002; Rosado et al. 2013). To overcome this problem, two different methods are generally used: the envelope-tracing method (ETM; e.g., Sofue et al. 1997, 1999a) and the iteration method (ITM; e.g., Swaters et al. 1997; Takamiya & Sofue 2000; Heald et al. 2006, 2007). We applied these three methods to NGC 4330 following the methodology described in detail in Appendix B. In particular, we invite readers to refer to Appendix B.1 for explanations of position-velocity diagrams (PVDs), which are extensively used in this section on RCs as well as in Sect. 6.

4.1. Comparison of the three methods

Figure 3 shows the RCs derived using three different methods: (1) the IPM, (2) the ETM, and (3) the ITM. The RCs obtained with the IPM using both GHASP and PUMA datasets show that NGC 4330 is rotating as a solid-body inside a radius of ∼30 arcsec. The output of the fit of Eq. (B.4) (ETM) provides similar terminal velocities reached around ∼30 arcsec (∼2.4 kpc) as well as similar transition radii between the rising and the flat part of the RCs compared to the IPM method. After this transition, the RCs computed with the ETM method are flat up to the photometric radius R_23.5 ∼ 140 arcsec (∼11.6 kpc), which is represented by the red line in panels a and c of Fig. 3. The approaching and receding sides of the RCs measured with the IPM method start to differ at a radius of ∼30 arcsec (∼2.5 kpc) and do not have a single point of matching where the two sides overlap until the photometric radius. The rotation velocity obtained for each dataset with the ETM differs by 8 km s⁻¹. This number is within the uncertainties due to the differences in resolution between the two instruments and could be partly due to a different threshold for the PVD extraction.

Fig. 3. PUMA (panels a and b) and GHASP (panels c and d) RCs (panels a and c) and PVDs (panels b and d) of NGC 4330. Panels a,c: RCs derived using three different methods: the IPM, the ETM, and the ITM (see Appendix B.3 for a further explanation on the methodology). In the case of the IPM (black squares), an azimuthal sector of ±5° in the sky plane has been taken around the major axis to compute the rotation velocities; filled and empty black squares indicate the receding and the approaching sides, respectively. The upper-left insert shows the azimuthal sector and the location of the regions pointed on the RC. For the ETM, the rotation velocities have been derived from the PVD along the major axis (z = 0), shown in panels b and d. The terminal velocities are represented by empty blue circles for the receding side and filled red circles for the approaching side; the pink line corresponds to the best fit to these points using η = 0.3, and the shaded area shows the rotation velocity amplitude for η = [0.2, 0.5] (see relation (B.3)). The ITM based on a cylindrical model (described in Appendix B.4) is represented by the black line. Panels b,d: PVDs showing the LoS velocities computed with the ETM from which the rotation velocities are computed in the upper panels. On top of the PVD, the filled red and open blue circles indicate the extremal LoS velocities and the pink curve the arctan model fitting those velocities. The red and blue vertical dotted lines in panels a and c show the photometric and effective radius, and the grey dotted lines in the lowest panels show the position of r0 (see Sect. 4).

Both datasets show that the RCs derived using the ETM are steeper than those derived using the IPM in the inner regions. This difference is due to the fact that here the contribution of regions along the same LoS but located at different radial distances from the galaxy centre is large, while this is not the case at the edges of the disc where the ETM measures the rotation. Close to the galaxy centre, the width of the PVD is supposed to narrow because along the minor axis the rotation velocity vector lies in the sky plane, and thus its projected component along the LoS is null. Since we do not observe a large decrease in the PVD width in the centre, we conclude that it is mainly due to radial velocity dispersion.

For the GHASP dataset, the RCs derived from the ETM and the ITM have a similar inflexion shape along the transition radius (r_t) between the rising and the flat part of the RC (see Fig. 3). Nevertheless, the fit shows that the parameter r_t differs by ∼5 arcsec (∼0.4 kpc) using the two methods. For the PUMA dataset, r_t occurs at a radius ∼9 arcsec (∼0.7 kpc) larger with the ETM than with the ITM, which makes the RC steeper with this last method (see Fig. 3). Despite the difference between r_t computed using the ETM and ITM, an acceptable discrepancy between the maximal rotational velocities of ∼8 km s⁻¹ is observed at the photometric radius in the RCs of each dataset.

The three methods consistently show that the rotation of two-thirds of the ionised gas disc of NGC 4330 (i.e. inside a diameter of ∼60 arcsec, ∼4.8 kpc) can be modelled as a solid rotating body. The same methods have also revealed that the transition between the rising and flat part of the RCs occurs in between r₀ ∼ 30 arcsec (∼2.4 kpc) and slightly before the photometric radius.

The IPM allows us to analyse the rotation on both the receding and approaching sides of the galaxy using the pixels within an angular sector of ±5° from the kinematic major axis. The RCs derived using this method show that at the photometric radius the rotation starts to be dominated by streaming motions. This occurs mainly at the south-west of the disc, where the gas has non-circular motions and where the ionised disc is distorted, and hence the projection parameters may be no longer valid in these outer regions. These peculiar motions are not modelled with the other methods. Nevertheless, the ETM and ITM withdraw the effects of LoS to the rotational velocity on the inner disc and allow us to obtain a reliable measure of the maximal rotational velocity of the galaxy, we fix it to the average value of the four velocities given in Table B.1 for the two datasets: $V_{\mathrm{rot}}^{\max }$ = 158 ± 12 km s⁻¹, the uncertainty being estimated from their standard deviation.

4.2. The dynamical mass

The RC of NGC 4330 is not flat because the gas that is used as a mass-tracer is subject to ram pressure stripping outside the transition radius. Indeed, taking into account its size, morphological type, and maximum rotation velocity, NGC 4330 is an intermediate-mass galaxy that is supposed to have a flat RC (e.g., Sofue et al. 1999b), as previously described in this section by different models. Despite the fact that NGC 4330 is located at ∼0.4 virial radius from M87 (e.g., Chung et al. 2007), its mass distribution is probably not affected by the cluster environment (Amram et al. 1992) and can be considered as unperturbed. Following Lequeux (1983), the actual dynamical mass of a spiral galaxy that contains a stellar disc, a spherical stellar bulge, and a spherical dark matter halo ranges in the mass interval given by $M(R)=\alpha R{V}_{\mathrm{rot}}^2(R)/G$, where G is the universal gravity constant, V_rot the flat rotation velocity, and α a parameter spanning between the lowest value of α = 0.6 for a galaxy dominated by a flat disc component and the higher value of α = 1.0 if the spherical halo components dominate. Inside the radius R_23.5(i−band)∼11.6 kpc, using $V_{\mathrm{rot}}^{max}=158\pm 12$ km s⁻¹ and using α = 0.8 ± 0.2 to be conservative, the dynamical mass of NGC 4330 within R_23.5(i−band) is M_dyn(R_{23.5, i − band}) = 4.4 ± 0.6 × 10¹⁰ M_⊙. For the same maximum rotation velocity, the mass within the HI radius R_HI = 11.3 kpc is M_dyn(R_HI) = 4.9 ± 1.2 × 10¹⁰ M_⊙. These estimates can be compared to the stellar mass of the galaxy, M_star = 5 × 10⁹ M_⊙ (see Table 1).

4.3. Comparison with HI and CO kinematical data

In Fig. C.4 we show the HI velocity fields from Chung et al. (2009) and the RCs we derive from the PVDs. As expected, due to beam smearing effects, the RCs obtained using the minimum and maximum HI surface brightnesses strongly differ. These data are used to compute the ETM RC, which matches fairly well with the Hα RC, except maybe within the first 2 kpc where the HI velocities have probably been overestimated (see Fig. 4). The comparison between the Hα, HI, and CO velocity fields are shown in Figs. C.5 and C.6 respectively. The HI velocity field is the most extended with a diameter of ∼282 arcsec (∼22.6 kpc) from the north-eastern top to the south-western lowest corner. The Hα velocity field has been measured on a diameter of ∼205 arcsec (∼16.4 kpc), and exceeds by more than a factor of two the one derived using CO data (∼94 arcsec, ∼7.5 kpc). Despite these differences, which might be partly due to the different sensitivity in the data combined with their different angular resolution (15 arcsec for HI and 6 arcsec for CO), the three velocity fields are comparable indicating a similar kinematics of the three components. We applied the ETM to the ¹²CO (2 − 1) data of Lee et al. (2017) and to the HI data of Chung et al. (2009, see Sect. 2.3). From the PVD extracted with a pseudo-slit along the kinematic major axis we computed the terminal velocity using η = 0.3, where η is a fraction of maximum intensity to set the envelope intensity (see Eq. (B.3)), then we fitted the data to the arc-tangent model (Eq. (B.4)) using a least-squares minimisation. The HI, CO, and Hα ETM RCs are plotted in Fig. 4, which also includes the Hα RCs computed with the ITM and the IPM.

Fig. 4. RCs of NGC 4330 derived using the 12CO(2 − 1) at λ = 1.3 mm data from Lee et al. (2017) (green), the HI (λ = 21 cm) data from Chung et al. (2009) (purple), and the FP Hα data from PUMA (black) and GHASP (blue). The RCs from CO and HI data were computed using the ETM. The solid lines correspond to the best fit to the points of the terminal velocity using η = 0.3, and the shaded areas show the rotation velocity amplitude for η = [0.2, 0.5] (see relation (B.3)); the dashed lines show the RC computed from the PVD maximum intensity. For comparison, the Hα RCs shown in Fig. 3 are also reported on this plot. The green, purple, and red vertical lines indicate the CO extension (∼47 arcsec, ∼3.8 kpc), the HI radius (∼141 arcsec, ∼11.3 kpc), and the outermost photometric radius (R23.5 ∼ 145 arcsec, ∼11.6 kpc), respectively. The top-right insert indicates the angular resolution of each dataset.

For the HI RC, the ETM fit gave a coefficient of determination R² = 0.97, with r₀ = 18 ± 2 arcsecond r_t = 6 ± 2 arcsec, providing a RC with the steepest inner slope, represented by the purple line in Fig. 4. This cuspy inner slope is due to the low HI surface brightness, which is reached thanks to the low-resolution VLA HI data (∼25 arcsec; see Sect. 2.3). This low resolution induces strong beam smearing effects and an overcorrection at the lowest HI levels. The ETM provides also a way to estimate of the lowest limit of the RC using only the maximum intensity I_max (in setting η = I_min = 0; see Eq. (B.3)). This curve is represented by the purple dashed line in Fig. 4. The interval between the two purple lines delimits all the possible HI RCs that could be computed and gives inner slopes compatible with those obtained at other wavelengths. The maximum HI rotation velocity is ∼10 km s⁻¹ lower than those obtained with Hα data.

Since the CO is detected only in the inner regions, we extrapolated its ETM RC, using the parameters of the arctan function fitted to the observational data (see Appendix B.3), to compare it with those computed at other wavelengths using the same method. The CO RC is represented by the green dashed line in Fig. 4. The fit of the CO RC, which was obtained with a coefficient of determination R² = 0.76, gave r₀ = 19 ± 10 arcsec, and the transition between the rising and flat part of the RC at r_t = 33 ± 26 arcsec. These values are in agreement with those determined for the Hα RCs using the same method. At large radii, the CO and the HI RCs match and have values ∼10 km s⁻¹ lower than those derived using the Hα data. We recall, however, that due to its limited extension, we do not have any observational constrain on the CO RC outside ∼50 arcsec (∼0.4 kpc).

5. Kinematical analysis from two-dimensional maps

5.1. Velocity fields

Figure 5 shows the velocity field of the ionised gas of NGC 4330 derived from the PUMA and GHASP FP datasets. The velocity fields of both datasets are consistent, they display a rigid body shape inside of a diameter of ∼60 arcsec (∼4.8 kpc) where the isovelocity contours are parallel to the minor axis of the galaxy. We traced on these maps two parallel dashed lines at r₀ ± 30 arcsec (±2.4 kpc) enclosing this region. Outside this diameter the isovelocity contours become almost parallel to the major axis in the NW hook (region A) and on its counterpart, between regions I and J. However, these maps differ in the southern tail (i.e. in regions K and L; see point 5 hereafter). These last features are the signature of non-circular motions that indicate strong perturbations in the kinematics of the regions affected by the stripping process. The two-dimensional velocity field obtained with high-resolution FP interferometry allows us also to measure the kinematic of the gas in very low-surface-brightness regions. These regions, labelled from A to M, were previously identified in the monochromatic image and they are now identified on the velocity fields of both datasets in Fig. 5. Our data allowed us to measure the LoS velocities of the stripped gas in the southern regions, in region F located on the compression front on the north-eastern side of the disc and in region I on the disc plane. The LoS velocities of all these regions are given in Table 2 and will be further discussed in Sect. 5.4.

Fig. 5. Velocity fields of the ionised gas of NGC 4330 (left PUMA, right GHASP). The general features of the main disc regions (north-east, central, and south-west), identified by two parallel dashed lines, and the regions with particular features, identified with letters A to H, are discussed in the text. The black bullet and the cross indicate the position of the photometric and kinematical centre, respectively. The black arrow shows the wind direction as derived from the models of Vollmer et al. (2012, 2021).

We made seven observations. First, the detection and estimation of the LoS velocity of the gas over the galaxy disc, and more specifically inside of regions A, E, F, G, H, I, and J (westernmost part), are robust since derived from independent sets of data with consistent results. Second, the gas inside region A has a radial velocity of ∼1650 km s⁻¹, which is ∼25 km s⁻¹ lower than that of the disc at the same galactocentric distance.

Third, if we compare the radial velocities in the south-eastern regions of the galaxy with respect to the radial velocities in the plane of the disc, either with respect to their position perpendicular to the plane of the disc or with respect to the direction of the wind as deduced from the direction of the tails and predicted by simulations (e.g., Vollmer et al. 2021), we notice that the gas inside regions B, C, D, and H has radial velocities significantly different (∼±50 km s⁻¹) from those measured in the disc from where the gas has been stripped. The radial velocities of the gas regions B, D, and C is lower than that in the disc, while the radial velocity measured in region H is higher by the same amount.

Fourth, despite the fact that they are above the FP detection limit observed in other regions, the detection of ionised gas in regions B, C, D, and M and in the western tail of region J could have been affected by night skyline contamination. Indeed, the night-line subtraction was very difficult because the OH-radical has a strong emission line at ∼6597 Å (e.g., Osterbrock et al. 1996), roughly at the same wavelength that the redshifted Hα emission. On the other hand, the LoS velocity of the gas in region H, which has been probably originated in the approaching side of the disc, is lower and different from the LoS velocity of the gas measured in regions B, C, and D, in the receding side. This difference suggests that these estimates are not largely affected by the residual night skyline emission.

Fifth, ionised gas is detected in regions K and L with similar LoS velocities. A velocity difference of ∼50 km s⁻¹ is observed between both datasets, which could also be due to parasitic night skylines. Thus, despite this difference in the two datasets, these velocities are consistent with streaming motions caused by the stripping process.

Sixth, the gas detected in regions J, K, and H has a velocity similar to that measured at the edge of the rotating disc (1450 km s⁻¹), suggesting a flat RC.

Finally, Fig. 5 also shows two dashed lines at ∼30 arcsec (∼2.4 kpc) from the galactic centre. These dashed lines enclose the disc section where the ionised gas shows a rigid-body rotation as indicated by the velocity contours parallel to the minor axis. To the south-west, region I marks the limit where the velocity contours start to become parallel to the major axis. Its counterpart to the north-east is the beginning of the hook shape where the isovelocities become not parallel to the minor axis. Region I might represent the outer layer of the vertical disc stratification, at the interface with the surrounding medium, where ram pressure begins its action. Since the stripping pressure is not perfectly vertical on the disc, the I region suffers a radial component of the stripping force.

In order to check the highly uncertain velocity of region D gathered with the FP data, we obtained a long-slit spectrum on February 17, 2021, using the low-resolution (R = 700) MISTRAL spectrograph installed on the 193 cm OHP (Adami et al. 2018). We took four independent exposures of 15 min each, for a total of 60 min. The exposures were gathered using a slit of width 1.9 arcsec oriented north-south to include at the same time region D and the disc of the galaxy. The low-resolution spectrum of region D is shown in Fig. C.2. Because of the night skyline contamination, the Hα line is marginally detected with a signal-to-noise ratio (S/N) ∼2 and a mean LoS velocity of 1458 ± 92 km s⁻¹, only marginally consistent with the FP measurements of 1580 ± 30 km s⁻¹. Figure C.2 also shows the corresponding MISTRAL spectrum obtained on the galaxy disc north of region D, where the average velocity is 1629 ± 46 km s⁻¹, in very good agreement with the FP data. As for the other extraplanar regions detected in the FP data, the velocity in the stripped gas is lower than that measured in their corresponding regions to the north on the galaxy disc.

The gas stripping process removes the gas from the disc and transfers and dilutes it into the ICM. The only velocity component that we can measure is the one along the LoS. For edge-on galaxies it is straightforward to interpret the LoS velocity as rotation velocity in the galaxy disc plane. The intergalactic wind might change the game, in adding a velocity component related to the direction of motion of the galaxy within the cluster. Nevertheless, in the specific case of NGC 4330, the wind is almost perpendicular to the galaxy disc and has a large component in the sky plane, as suggested by the numerical simulations (see Sect. 7). For these reasons the observed velocities of the gas expelled from the disc are not expected to dramatically change. Thus, the gas continues to rotate outside of the disc and shows the same LoS velocities than those it had before being stripped. Since regions B, C, and D (H) have LoS velocities lower (higher) than their corresponding region in the disc, this suggests that here the gas partly lost its angular momentum during the ram pressure stripping event.

5.2. Velocity dispersion

Figure 6 shows the velocity dispersion of the galaxy as measured from the PUMA and GHASP datasets. Following Boselli et al. (2022) we checked that the velocity profiles have a single peak, (i.e. that the barycentric method applied to measure the velocity dispersions does not introduce any strong systematic bias in the results).

Fig. 6. Velocity dispersion fields of NGC 4330 (left PUMA, right GHASP). The general features of the main disc regions (north-east, central, and south-west), identified by two parallel dashed lines, and the regions with particular features, identified with letters A to H, are discussed in the text. The black bullet indicates the position of the photometric centre. The black arrow shows the wind direction as derived from the models of Vollmer et al. (2012, 2021).

At first order, the velocity dispersion of the gas is fairly constant over the disc of the galaxy, with σ ∼ 35 ± 10 km s⁻¹ and σ ∼ 21 ± 7 km s⁻¹ for the PUMA and the GHASP data, respectively. The PUMA velocity dispersion is higher than the one measured in the GHASP data because of its reduced spectral resolution. It is indeed difficult to measure velocities far below the instrumental resolution (σ_LSF = 19.5 and 13.1 km s⁻¹, respectively) as it has been recently discussed in Boselli et al. (2022) in comparing VLT/MUSE data with FP data obtained with GHASP.

At higher spatial frequencies, fluctuations in the velocity dispersion are observed on scales comparable to or larger than the ones set by the seeing. These fluctuations are smoother in the GHASP data than in the PUMA ones. This is due to the difference in the method used to compute the map (Gaussian smoothing versus Voronoi tessellation, respectively).

Panels c and e of Fig. C.3 show intensity-velocity dispersion diagnostic diagrams and panels a and d the histograms for the intensities and of the velocity dispersions. These diagnostic diagrams are generally used to characterise the dynamics of the ISM in extragalactic giant HII regions (Munoz-Tunon et al. 1996) and also to investigate the dynamic of HII galaxies (Moiseev & Lozinskaya 2012; Carvalho & Plana 2018). Moiseev & Lozinskaya (2012) show a sketch where different features and regions in the plot can be interpreted in terms of (case 1) HII regions (constant velocity dispersion region with high Hα intensity), (case 2) shell (inclined bands in the plot), (case 3) low density turbulent ISM (triangular shape with a base spanning from low to high velocity dispersion and a triangle height ranging from low to moderate Hα intensity).

Shells (case 2) are not observed in the diagnostic diagrams shown in panels c and e, probably because of the lack of spatial resolution for this relatively distant galaxy. The triangular shape distribution of the diagnostic diagram (case 3) is asymmetric, with a deviation towards the lowest velocity dispersions at low intensity and towards larger velocity dispersions at intermediate intensities, indicating the supersonic behaviour of the ISM in NGC 4330. In panel e we selected four areas, represented by ellipses, of different pixel density (from high to low mean pixel density: red, green, dark, and light blue), and identified their locations in the galaxy in panel b. The red pixels with subsonic velocities are mainly located on the outskirts of the galaxy while the dark blue ones, with supersonic velocities, are a little bit on the inner side (case 3). The supersonic green pixels are located along the middle plane of the galaxy while the light blue with the highest intensities (case 1) match with the nucleus of the galaxy and the brightest HII regions. Their velocity dispersion ranges between 20 km s⁻¹ and 26 km s⁻¹.

5.3. The residual velocity fields

Two-dimensional residual velocity fields are generally used to optimise the free parameter determination of RCs that are computed from a two-dimensional velocity fields. In general, because the RCs are based on axisymmetric models, a pattern observed on the residual velocity field would be due to incorrect parameter determination. As a consequence, when the free parameters of the RC are optimised, the residual velocity field is pattern-free and the dispersion in the two-dimensional residual velocity is minimum (see Warner et al. 1973; Fuentes-Carrera et al. 2004; Erroz-Ferrer et al. 2012).

The PUMA and GHASP residual velocity maps (V_res = V_LoS − V_mod) shown in Fig. 7 were obtained after subtracting the iteration model velocity field (V_mod) from the observed one (V_LoS). Once we identified the regions presenting a large scatter with the mean velocity dispersion of the residual velocity field, which are not supposed to be in rotation in the plane of the disc, they were masked from the observed velocity field and the axisymmetric model was recomputed. This is because the axisymmetric model is not supposed to be observationally constrained by non-axisymmetric motions.

Fig. 7. Residual velocity fields obtained after subtracting the velocity field computed from the ITM model data cube (left PUMA, right GHASP). The model, based on a solid-body rotation, is shown in the inset. The purple parallel lines show the isovelocity contours of the edge-on rotating disc model. The bullet and the cross indicate the position of the photometric and kinematical centre, respectively. The arrow indicates the wind direction as derived from the models of Vollmer et al. (2012, 2021).

The residual velocity maps given in Fig. 7 do not show any systematic pattern in the galactic disc and have a velocity dispersion of σ_res ∼ 19 km s⁻¹ for the GHASP dataset and σ_res∼ 24 km s⁻¹ for the PUMA dataset. The residual velocity is ∼0 km s⁻¹, suggesting that the adopted cylindrical model (described in Appendix B.4) reproduces the axisymmetric rotation of the disc modelled by the RC and its associated parameters well.

An accurate examination of Fig. 7 indicates how the kinematical properties of the various extraplanar ionised gas features created after a ram pressure stripping event differ from those properties predicted by a purely axisymmetric cylindrical model. The most discrepant regions in terms of relative velocity of the gas are regions G, H, and L for the PUMA data and G, H, and from J to M on the south-western side of the galaxy, for the GHASP data. Here the mean gas velocity is ∼50 km s⁻¹ higher than that observed in the model. A less pronounced difference in the same direction is observed in region E (∼20 km s⁻¹), which is located in projection close to the kinematic centre of the disc. On the contrary, the gas inside the regions D (PUMA) and B, C, D (GHASP) has a mean velocity ∼40 km s⁻¹ higher than that of the cylindrical model. Additional analysis about these individual regions is given in Sect. 5.4.

5.4. Peculiar motions

Table 2 gives the main kinematical properties of the gas in regions A to M. The GHASP data indicate that the gas inside most of these particular regions has a velocity dispersion near to the mean value measured over the galactic disc, and this is probably due to the Gaussian smoothing used to create the map of the velocity dispersion. The PUMA data rather indicate that the gas inside all the same regions has a velocity dispersion higher than the one measured within the galaxy disc, with the higher velocity dispersions measured in regions E, G, and K. On the other hand, the residual velocity field derived from both sets of data shows that the gas in regions A, E, F, and I follows the rigid-body rotation with a low residual velocity. This result is expected for the gas in regions E, F, and I because they are still associated with the galactic disc but not for the region A, which has been perturbed by the stripping process.

In region A, the velocity contours from both datasets indicate that the gas is rotating with a slight bending with respect to the minor axis. The PVD obtained with the PUMA dataset shows that the gas at the north-east, where the Hα monochromatic map has a hook feature, has higher radial velocities than those expected from the model. For this reason, we could expect that here the gas has a higher velocity dispersion due to non-circular motions, corresponding to high values on the residual velocity map. However, the RCs computed with IPM, which show the gas rotating as a rigid body in region A, do not confirm this result. A decrease in the velocity rotation with respect to the rigid-body rotation is instead observed in the RC and in the residual velocity maps in the south-western side of the galaxy.

Of particular interest is region I, where the GHASP data indicate that the velocity dispersion is ∼5 km s⁻¹ higher than the mean value estimated over the disc of the galaxy. The gas inside region I still follows the rigid-body rotation as indicated by a low residual velocity in both maps. Region I is located at the boundary between the rising and the flat part of the RCs, where the velocity contours in the velocity fields starts to become parallel to the major axis. The PVD shows that here the ionised gas emission start to decrease significantly, and that outside this region the gas starts to rotate slower than expected by the model. Its counterpart on the north-east side of the disc also shows a drop in surface brightness in the PVDs, but outside this region the rotation velocity of the gas seems to increase with respect to the model. In this symmetric region at the north-east of the disc, the gas has no high velocity dispersion and no high residual velocity values. Therefore, the increase in the line width of the gas observed in region I might be due to its acceleration, which is able to remove it from the galaxy disc in the south-west of the galaxy during the stripping event.

Figure 5 shows that the gas inside regions B, C, and D has a lower LoS velocity with respect to those measured in their corresponding positions on the disc plane while region H has a higher one. This means that these four regions rotate slower than their counterparts on the galactic disc. This can also be seen on the residual velocity field (Fig. 7) where the velocities are given in the rest frame of the rotating disc and in which negative (positive) velocities on the approaching (receding) side indicate a lower amplitude rotation outside the plane of the disc. The gas in regions B, C, D, and H does not have a high velocity dispersion, meaning that the wind does not bring a substantial quantity of kinetic energy to those regions. This probably indicates that the gas of these regions had lost angular momentum when it was detached from the galactic disc.

At the southern edge of the galactic disc, the velocity contours are bent with respect to the minor axis. This is particularly evident for the gas associated with region G, where the velocity dispersion and residual velocity increase smoothly with galactocentric distance, reaching their largest value in the detached cloud at the south-western edge of the region. The same cloud at the edge of region G has also a higher LoS velocity than its counterpart on the disc, which means a lower rotation velocity because this cloud is located on the receding side of the disc. Regions B, C, D, and H show a similar behaviour, they display lower LoS velocities, thus lower rotation velocities because they are approaching to the observer with respect to the galaxy centre. Then, we can infer that the gas inside regions J to M, with a low velocity dispersion but a high residual velocity, has evolved in a similar way as the gas in regions B, C, D, and H (i.e. having lost angular moment once stripped from the disc). Finally, the gas in region L has a higher velocity dispersion that the one observed around regions K and M. This does not mean that it gets kinetics energy from the wind because the line broadening is probably due to the presence of a star-forming region (see Fig. 2).

6. Kinematics of the gas out the plane

In this section we focus on the PVDs at different locations in the galaxy disc, parallel and perpendicular to the major axis of the galaxy, in order to study the kinematics of the gas out of the plane of the stellar disc. Figures 8 and 9, in which we have also located regions A to M, show the PVDs extracted from the PUMA data cube measured along the major and the minor axis, respectively, at different distances from disc plane (e.g., Epinat et al. 2008a; Rosado et al. 2013).

Fig. 8. PVDs parallel to the kinematic major axis extracted from the PUMA dataset. The distance of the pseudo-slit from the disc plane along the z axis (positive to the north) is given in the upper-left corner of each panel. The positions of all the pseudo-slits are shown on the Hα monochromatic map in the right panel. From the top at +3.7 to −5.5 arcsec (from ∼0.3 to ∼0.4 kpc), the intensity of the contour level is at 1.5 σ, and from −7.4 to −11 arcsec (from ∼0.6 to ∼0.9 kpc) the intensity level is at 1.2σ. In the left panels, coloured lines overlaid on the PVD represent the best fit of the terminal velocity (Relation (B.3) with η = 0.3) according to the arc-tangent model (Eq. (B.4)). There is a correspondence in colour between the terminal velocities and the pseudo-slits plotted on the monochromatic map in the right panel. In the middle panels, the PVDs derived from the data cube model, which was built in the context of the ITM to compute the RC (see Appendix B.4), are overlaid on the observed PVD in grey. All the PVDs extracted from the data cube model are at 8% of the maximum intensity of the model. The two parallel dashed lines are at r0 = ±30 arcsec (±2.4 kpc) and indicate the three main disc regions: north-east, central, and south-west. The different features of interest are indicated with letters A to H.

Fig. 9. PVDs parallel to the kinematic minor axis obtained from the PUMA dataset. The locations of the pseudo-slits are in the upper-left corner of each panel, and they are traced on the Hα monochromatic map in the left panel with colour correspondence. The distance increases from the south to the north of the disc, with zero at the kinematic centre. The black dotted ellipse on the Hα monochromatic map indicates the i-band μe(i) = 23 mag arcsec2 isophote. The intensity of the contour levels is 1.3, 1.4, and 1.6 σ from −45.5 to +63 arcsec (∼−3.6 to 5 kpc) and 1, 1.1, and 1.25σ from −130.9 to −52 arcsec (∼−10.5 to −4.2 kpc). The PVDs extracted from the data cube model (black contours) are at 1, 10, and 50% of the maximum intensity of the model.

6.1. PVDs along the major axis

Figure 8 shows the PVDs along the major axis at different distances from the disc plane. The pseudo-slits are placed along the z axis from +3.7 arcsec (∼0.3 kpc) in the northern direction (upwind direction) of the kinematic major axis, up to −11 arcsec (∼0.9 kpc) in the southern direction (tail direction), and they are separated from each other by ∼1.8 arcsec (∼0.1 kpc; see Appendix B.1). Figure 8 shows that the rotating gas is less extended along the disc upwind (+1.8 arcsec), where it reaches 50 arcsec (∼4 kpc) in both the north-east and south-west directions, than in the tail (70 arcsec to the north-east and 60 arcsec to the south-west). This result is expected in a ram pressure stripping scenario, where the gas is removed outside-in (e.g., Quilis et al. 2000; Vollmer et al. 2001; Boselli et al. 2006, 2022; Roediger & Hensler 2005; Tonnesen & Bryan 2009).

To identify the main kinematic features resulting from the ram pressure stripping event, we compare the predicted values by theoretical fittings with the observed data. We applied the ETM to each PVD measured along the major axis at different distances from the disc plane. In the left panels of Fig. 8 the best fit to the terminal velocity of each PVD is overlaid on its corresponding plot. In addition, in the context of the ITM, we built a first order approximation data cube model to compute the RC of NGC 4330 (see in Sect. 4). Here we compare the PVDs from the observational data with the PVDs at the same position obtained from that data cube model in the central panels of Fig. 8. We notice that the cylindrical model follows the emission intensity of the galaxy in velocity and matches with the kinematical behaviour of the major axis PVD.

Some low-surface-brightness plumes of Hα emission were detected on top of the dust lane inside region F. These plumes, which can be seen in the upper two panels of Fig. 8, show an increment in velocity by ∼30 km s⁻¹ breaking the symmetry of the PVDs. If the Hα emission, after being compressed on the front side of the galaxy during the stripping event, is pulled as a whole to the south-west of the stellar disc, these plumes represent gas dragged that is still resisting to change, giving us a clue of the initial state of the disc.

The PVD along the major axis presents a kinematic symmetry with respect to the kinematic centre without any plume or knot with important variations in velocity with respect to the model. The kinematic major axis matches the bent dust lane in the central region, within a radius of ∼30 arcsec (∼2.4 kpc), enclosed by the vertical dashed lines. As is suggested by the bending of the velocity contours (see Fig. 5), the radial velocities tend to be constant ∼150 km s⁻¹ to the south-west, while to the north-east the PVDs show velocities that change by ∼ ± 25 km s⁻¹ perpendicular to the disc. To the north-east, the gas in region A has a radial velocity of ∼200 km s⁻¹. As indicated by the residual velocity map (see Sect. 4), although the gas inside region E has a peculiar shape, it rather follows the rotation of the galaxy. On the contrary, the gas in region G has a lower velocity than in the disc.

6.2. The lag gradient

Kinematical studies of edge-on galaxies have revealed the presence of a vertical velocity gradient in the RC moving away from the disc plane. This gradient has been observed at optical (e.g., Heald et al. 2006; Rosado et al. 2013; Bizyaev et al. 2017; Levy et al. 2019) and radio (e.g., Miller & Veilleux 2003; Fraternali & Binney 2006; Zschaechner & Rand 2015) wavelengths, including in the Milky Way (e.g., Kalberla et al. 2014). This phenomenon has been explained by a combination of vertical movements due to energy input from galactic fountains, radial movements of the gas due to the pressure gradient of the halo and a declining rotation velocity due to the conservation of the angular momentum (Fraternali & Binney 2006; Heald et al. 2007). This vertical velocity gradient can be modified whenever a galaxy such as NGC 4330 is suffering an external perturbation. In NGC 4330, the amplitude of the RC does not decrease sufficiently with the vertical distance from the disc plane to explain the observed lags. In this galaxy, which is suffering an almost face-on ram pressure stripping event, galactic fountains were probably suppressed, in particular in the upwind side of the disc where the external pressure is significantly increased.

The lag gradient of NGC 4330 can be estimated by measuring the change in velocity along the z axis using the PVDs parallel to the major axis. Once corrected from the morphological inclination of the galaxy (89.5°), the envelopes of the different PVDs shown in Fig. 8 become the RCs given in Fig. 10. The RC with the highest velocity amplitude is obtained in the plane of the disc (z = 0 arcsec). To the south-west, at the effective radius, the velocity decreases by ⟨ΔV⟩≃9 km s⁻¹ each Δz = 1.8 arcsec (∼0.14 kpc) from the kinematic major axis up to z = −5.5 arcsec (∼0.4 kpc), where, because the ram pressure stripping compresses the gas, there is no change in velocity along the z axis up to -11 arcsec (∼0.9 kpc), where the most diffuse gas is. To derive the velocity gradient along the z axis of NGC 4330 we used the fact that the asymptotic velocity V₀ of the relation (B.4) is one of the free parameters of the fitting to the terminal velocity of each PVD, so that, from z = 0 to −5.5 arcsec (∼0.4 kpc) there is a lag gradient of ΔV/Δz ∼ 64 ± 9 km s⁻¹ kpc⁻¹. The rotation in the north along the z axis has been modified by the wind, it is slower by ΔV = 49 km s⁻¹ with respect to the one measured at z = 0 arcsec. The difference in velocity between the northern RCs is ΔV = 16 km s⁻¹. Although the rotation velocity decreases with the distance from the stellar disc, as indeed expected, it is not possible to determine a constant lag gradient in this direction.

Fig. 10. Envelopes corrected from the morphological inclination (89.5°) of the different PVDs shown in Fig. 8. The vertical blue dotted line indicates the effective radius Re(i).

6.3. PVDs along the minor axis

The ram pressure stripping event, which is acting almost face-on in NGC 4330, has reduced the vertical height of the ionised gas disc in the northern direction (upwind). We thus study the PVDs along the kinematic minor axis to determine whether galactic fountains or inflows have been suppressed by the external pressure (e.g., Heald et al. 2006, 2007; Rosado et al. 2013). Figure 9 shows the PVDs perpendicular to the kinematic major axis obtained from the PUMA dataset. The pseudo-slits are located in specific places to highlight the main features of the disc and its asymmetries, and they cross the regions labelled from A to M identified in Sect. 2. We also compare the PVDs from the observational data with the PVDs at the same position obtained from that data cube model computed with the ITM (see Sect. 4).

In general, these PVDs show that the brightest Hα emitting gas is compressed to the galactic mid plane and that it follows mostly the rotation predicted by the cylindrical model. Upwind (north-east) there is no Hα emitting gas beyond z ≃ 5 arcsec (0.4 kpc). Only near the midplane of the disc, downwind (south-west) the gas decreases its brightness along the z−axis following the exponential structure used to construct the model. Here the velocity of the most diffuse gas generally increases with z all over the galactic disc. This increase in velocity with z is not reproduced by the ITM model.

In region A, at 65.5 arcsec (5.2 kpc) from the kinematic centre, the ionised gas has exactly the same velocity of the cylindrical model. At 31.7 arcsec (2.5 kpc), the low radial velocity observed for the gas in region D in the velocity field is now evident as long as its corresponding position on the disc follows the model. We only detected some diffuse gas to the north-east along the z−axis inside region F at a radius of 4.1 arcsec (0.3 kpc) from the kinematic centre, where the diffuse gas rotates a little faster than the model. The pseudo-slit at the kinematic centre shows that the gas in region E is completely compressed to the galactic mid plane. To the south-western edge of the disc, the most diffuse gas inside regions G, H, J, and K rotates slower than what expected from the model, in agreement with the residual velocity map. Finally, we have observed that the gas inside region L (-121.5 and -130.9 arcsec, 9.7 and 10.5 kpc, respectively) has a high velocity dispersion, radial and residual velocity likely because the gas is moving with non-circular motions. Figure 9 also shows that the diffuse gas in region L has a velocity higher by ∼200 km s⁻¹ than the one expected from the model.

Instead of a heart-shaped profile, which would indicate the presence of an inflow in the galaxy (e.g., Fraternali & Binney 2006; Rosado et al. 2013; Zschaechner et al. 2015), the PVDs suggest features similar to the intermediate- and high-velocity clouds of the Milky Way (e.g., Wakker & van Woerden 1997; Fraternali & Binney 2006) whose origin is attributed to the interaction with the Magellanic Clouds. This suggests that the kinematics of the extraplanar diffuse ionised gas of NGC 4330 is likely dominated by the ram-pressure stripping process.

7. The dynamical model

We present here the kinematical results of the simulation of NGC 4330 presented in Vollmer et al. (2021). They have been obtained using a N-body code, described in Vollmer et al. (2001), which consists of a non-collisional component plus a collisional one. The non-collisional component of 81 920 particles simulates the stellar disc, the stellar bulge and the dark halo. The collisional component simulates the ISM as an ensemble of 20 000 gas clouds that have inelastic collisions (sticky particles) that might coalesce or fragment, which evolve in the gravitational potential of the galaxy. A scheme for star formation was implemented, in which stars are formed during cloud collisions to further evolve as non-collisional particles. The star formation rate is proportional to the cloud collision rate. Ram pressure has been included as an additional acceleration on the sticky particles that are not protected by other particles. Once stripped out from the galactic plane, the warm gas clouds become diffuse if their density falls below a critical density (Vollmer et al. 2021). The model consists of two components: the HII regions ionised by young massive stars with ages less than 20 Myr and the diffuse gas ionised by the stellar UV radiation or by strong shocks induced by the ram pressure stripping process. As discussed in Vollmer et al. (2021), collisional ionisation through the thermal electrons of the ambient ICM that confines the filaments is the most probable ionisation mechanism in the extraplanar ionised gas filaments detected in the VESTIGE narrow-band image.

In Figs. C.8 and 11, we compare the observed Hα data to the simulated ones. For each figure we provide three panels, the top one corresponds to the data whereas the middle and bottom ones to the models. Figures 11 (panels b and e) and C.8 (panel b) give the lowest surface brightness achieved by the model while Figs. 11 (panels c and f) and C.8 (panel c) show the model at a depth comparable to the sensitivity of our observations. The plotted Hα emission map is the one extracted from the model for the snapshots of highest goodness for the diffuse gas. This corresponds to a peak of ram pressure occurred 140 Myr ago. Similarly to the deep Hα observations presented in Fossati et al. (2018), the model shows extra-planar linear filaments in the downwind region of the galactic disc, which correspond to the low-surface-brightness filaments located in the south-west tail and close to the downturn in the north-east (Vollmer et al. 2021). These filaments indicate the direction of the motion on the plane of the sky of the galaxy within the cluster.

Fig. 11. GHASP (panels a and d) and simulated (panels b, c, e and f) Hα velocity (left column) and velocity dispersion (right column) fields of NGC 4330. The simulated Hα velocity and velocity dispersion fields are given down to the faintest levels used in the middle panel of Fig. C.8 to show the kinematical properties of the low column density ionised gas reached by the simulations (panels b and e) and at a higher surface brightness level corresponding to the sensitivity of the FP observations (Σ(Hα)≃ 2 × 10−17 erg s−1 cm−2 arcsec−2, panel c). The contours on the velocity fields indicate the Hα surface brightness Σ(Hα)∼ 2×10−18 erg s−1 cm−2 arcsec−2 taken from the VESTIGE data (Fossati et al. 2018).

The agreement between the gas distribution in the observations and in the simulations is overall very satisfactory. The brightest Hα knots detected in the FP observation (panel a) correspond to the brightest Hα regions in the simulation (panels b and c). The diffuse Hα emission in the model, given (in log scale) in panel b, is far too faint to be detected in the FP observations, but it matches fairly well the diffuse emission observed using the Hα narrow-band filter by Fossati et al. (2018). While the tuning is very good at the north-east edge of the disc, the simulated south-west tail has a larger curvature than that observed in the diffuse ionised gas detected by VESTIGE. This difference might be due to the fact that, given the orientation of NGC 4330 and its deduced orbit within the cluster, the north-east edge of the disc is the region of the galaxy that first suffered the stripping process.

As suggested by the simulations, all the diffuse gas component (mainly HI gas) has been already stripped, while that in the south-west direction, which was perturbed later, is still partly present, as suggested by the diffuse HI tail observed by Chung et al. (2007). This effect is probably combined with an asymmetric distribution of the star-forming regions over the disc of the galaxy, possibly formed after the external perturbation, with more compact HII regions along the south-west edge of the disc. Given their compactness and high density, these regions are less subject to the external perturbation than the diffuse gas. The difference between observations and simulations might thus result from a combined effect of sensitivity, gas density, and time elapsed since the beginning of the perturbation. Indeed, as shown in Fig. 19 of Abramson et al. (2011), the south-west curvature of the gas is more pronounced in the cold and diffuse atomic phase (the first one stripped) than in the compact Hα, more resistant to the external pressure.

The Hα LoS velocity fields are shown in panel a of Fig. 11. In order to compare the model to the data, a systemic velocity of 1551 km s⁻¹ has been added to the model. The velocity amplitudes are the same for the three panels. Hα isocontours from Fossati et al. (2018), corresponding to a surface brightness of Σ(Hα)∼ 2 × 10⁻¹⁸ erg s⁻¹ cm⁻² arcsec⁻², have been displayed in order to facilitate the comparison between the panels. We observe a very good agreement between the simulation and the observations. Where the simulations can be compared to the data, the amplitudes of the velocity fields are consistent. We note that if we include the low-surface-brightness tails, unfortunately not reached by the shallow FP observations, the velocity range in the simulations is even larger.

On the other hand, the HI velocity field (Fig. C.5, bottom panel) seems to show that the LoS velocities increase at the edge of the tail, although this rising is not very sharp, probably because of the lack of resolution in the HI data. The same excellent agreement is present also between the observed and simulated unfolded RCs shown in Fig. C.7 on which the simulated RC measured along the major axis of the galaxy has been compared to the GHASP observed and modelled RCs. The RC extracted from the simulations well predict the observed asymmetries between the receding (north-east) and approaching (south-west) sides of the disc up to ∼70 arcsec (∼5.6 kpc). Since the RC extracted from the simulations is taken on the major axis, we do not have any data outside this radius because the ionised gas emission is bent to the south-west.

The Hα LoS velocity dispersion fields are compared in the right panels of Fig. 11 using the same velocity range in the three panels. The mean velocity dispersion in the data is ∼20.5 km s⁻¹, with a standard deviation of ∼6.7 km s⁻¹, while in the model those quantities are respectively equal to ∼10.6 and ∼6.7 km s⁻¹. The smaller mean velocity dispersion in the model with respect to the observations (which have an instrumental limiting resolution of ∼13 km s⁻¹), is due to the fact that the model estimate is strongly weighted by the diffuse gas component. This component is not reached by the observations. If we consider only the area of the model corresponding to the one detected in the data, the mean model velocity dispersion becomes ∼17.5 km s⁻¹, with a standard deviation of ∼6.4 km s⁻¹, and matches fairly well the observations. Striking is the similarity between models and observations around RA = 12^h23^m15.5^s, Dec = 11° 21′42″, where the velocity dispersion increases by one sigma with respect to other regions along the galaxy disc up to 28–45 km s⁻¹ in the data and 30–40 km s⁻¹ in the model (in red in both cases).

8. Discussion and conclusion

Located in the nearby Virgo cluster, where 1 arcsec corresponds to ∼80 pc, the edge-on late-type NGC 4330 is the perfect candidate for studying the effects of ram pressure stripping in the direction perpendicular to the disc plane of a perturbed system with an angular and spectral resolution unreachable elsewhere. The presence of extended tails in the atomic (Chung et al. 2007) and ionised (Fossati et al. 2018) gas and in the dust component (Longobardi et al. 2020) observed in this object testify to an ongoing stripping process that is able to remove the different phases of the ISM in the outer regions and produce a truncated disc.

Using two different and independent sets of FP observations with a spectral resolution of up to R ∼ 10 000, combined with our own tuned hydrodynamic simulations (Vollmer et al. 2021), we studied the kinematics of the ionised gas. Although the diffuse emission of the ionised gas in the stripped tails, which has a surface brightness of Σ(Hα) ≃ 3 × 10⁻¹⁸ erg s⁻¹ cm⁻² arcsec⁻², has not been reached by the present observations, limited to Σ(Hα) ≃ 10⁻¹⁷ erg s⁻¹ cm⁻² arcsec⁻², our data allowed us to study the kinematics of the gas along the disc and in a few HII regions located outside the disc plane and probably formed within the stripped gas after the interaction.

The Fourier transform analysis of the deep NGVS i-band image reveals that the B₄ parameter is negative, indicating that the galaxy disc has a boxy shape outside a radius of ≃60 arcsec (≃4.8 kpc). This boxy structure is present mostly in the south-west direction, where the tail of stripped material is more pronounced. Simulations suggests that a thicker disc could result from the perturbation of the gravitational potential well due to the displacement of the gaseous component from the disc plane after a ram pressure stripping event (Farouki & Shapiro 1980; Clarke et al. 2017; Safarzadeh & Scannapieco 2017; Steyrleithner et al. 2020). The boxy shape would thus result from the induced perturbation on the stellar orbits. The datasets used in this work, which are sensitive only to the kinematics of the ionised gas, do not allow us to test this hypothesis. However, in the central regions, the kinematic major axis matches with the dust lane position angles, indicating that the perturbation probably had only a moderate effect on the dynamics of the inner gas.

It is thus conceivable that the dynamics of the stars are also moderately perturbed. Indeed, within the inner stellar disc (r ≲ 2 kpc), the velocity field of the galaxy dominated by the emission of HII regions is characterised by isovelocities parallel to the minor axis. The analysis of the RCs and of the PVDs consistently indicates a solid-body rotation. Outside this radius the isovelocities are bent, becoming parallel to the major axis, and show a low velocity gradient in the south-western edges of the disc, indicating non-circular motions.

The FP data allowed us to detect several ionised gas features formed after the dynamical interaction of the galaxy with the surrounding ICM (ram pressure) and located at the edges or outside the stellar disc. These are i) a hook-like structure at the north-east of the galaxy in the region that first came into contact with the ICM, and ii) several bright and extended regions formed along the downstream direction of the wind, principally located in the south-west of the galaxy, including a few HII regions well outside the disc and located along the filaments of stripped gas.

The analysis of the velocity field, of the velocity dispersion, of the RCs (derived using the ETM necessary to correct for beam smearing effects due to the velocity superposition along the LoS), and of the PVDs consistently indicates that these regions do not follow the solid-body rotation but have peculiar velocities indicative of streaming motions out of the plane of the galaxy. We also observe a small decrease in the rotation of the stripped gas with increasing distance from the galaxy disc.

This observational evidence can be explained considering that the stripping process is mainly perpendicular to the disc plane, as indeed suggested by our hydrodynamic simulations. Under this geometrical configuration, the gas is expected to keep its rotation with only a partial and gradual loss of angular momentum. Furthermore, the similarity of the Hα, CO, and HI velocity fields and RCs, particularly in the inner disc, indicates that here ram pressure stripping significantly affects the gas kinematics. Since in this scenario the acceleration caused by ram pressure is inversely proportional to the stripped gas surface density, the observed similarity in the kinematical perturbations of the different gas components (ionised, molecular, atomic) is expected given their comparable surface density.

Overall, the kinematical properties of the galaxy and of all these extraplanar features formed after the ram pressure stripping event are well reproduced by our simulations, at least up to the limited sensitivity of our FP observations. Nevertheless, it is hard to drive general conclusions regarding the effects of ram pressure stripping on the kinematical properties of galaxies and of the stripped material from the analysis of a single object.

The present results, however, combined with those obtained from the analysis of the kinematical properties of a few other well-resolved galaxies with available IFU data (e.g., Chemin et al. 2006; Merluzzi et al. 2013; Fumagalli et al. 2014; Consolandi et al. 2017; Boselli et al. 2022), consistently suggest that the gas removed during a ram pressure stripping event tends to keep its rotation, albeit with a possible loss of angular momentum. The perturbation can affect the kinematics of the gas even inside the disc of the perturbed system (Boselli et al. 2021). Although we still lack direct evidence, indirect observations and simulations consistently indicate that a ram pressure stripping event can also perturb, although to a lower extent, the kinematics of the stellar component. Deep IFU spectroscopic observations with a sufficient spectral resolution of perturbed systems will be required to probe this scenario.

Acknowledgments

Based on observations collected at the Observatorio Astronómico Nacional at San Pedro Mártir, Baja California, México (OAN -SPM). We thank the daytime and night support staff at the OAN-SPM for facilitating and helping obtain our observations. Also, based on observations taken with the GHASP and MISTRAL spectrographs at the Observatoire de Haute Provence (OHP, France), operated by the French CNRS. The authors warmly thank Olivier Boissin from LAM and the OHP team for its technical assistance before and during the observations, namely the night team: Jean Balcaen, Stéphane Favard, Jean-Pierre Troncin, Didier Gravallon and the day team led by François Moreau as well as Dr. Auguste Le Van Suu, the Head of Observatoire de Haute Provence-Institut Pythéas. We are grateful to the whole CFHT team who assisted us in the preparation and in the execution of the observations and in the calibration and data reduction: Todd Burdullis, Daniel Devost, Bill Mahoney, Nadine Manset, Andreea Petric, Simon Prunet, Kanoa Withington. We acknowledge financial support from “Programme National de Cosmologie and Galaxies” (PNCG) funded by CNRS/INSU-IN2P3-INP, CEA and CNES, France, and from “Projet International de Coopération Scientifique” (PICS) with Canada funded by the CNRS, France. This research has made use of the NASA/IPAC Extragalactic Database (NED), which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration and of the GOLDMine database (http://goldmine.mib.infn.it/) (Gavazzi et al. 2003). M. M. S. warmly thanks the Mexican National Council on Science and Technology (CONACyT) who found her through the program “Becas CONACyT al Extranjero 2017”, CVU 666085, and to the Secretariat of Public Education (SEP) of the Mexican Government through the scholarship “Becas Complemento de Apoyo al Posgrado ciclo 2018-2019”. M. R. thanks also the grants IN109919 of DGAPA-UNAM and CY-253085 and CF-86367 of CONACyT. M.B. acknowledges FONDECYT regular grant 1211000. A. L. is supported by Fondazione Cariplo, grant No. 2018-2329. M.F. acknowledges funding from the European Research Council (ERC) (grant agreement No. 757535).

References

Appendix A: Fabry-Perot observations

A.1. PUMA observations

Three-dimensional FP spectroscopic observations of NGC 4330 were gathered using PUMA at the 2.1 m telescope at the Observatorio Astronómico Nacional in San Pedro Mártir, Baja California, Mexico (OAN-SPM; Rosado et al. 1995). The PUMA focal reducer hosts an FP interferometer with a field of view of ∼10 × 10 arcmin². The camera is a 512 × 512 CCD detector³ with a pixel scale of ∼1.27 × 1.27 arcsec².

NGC 4330 was observed during dark time, in February 2017⁴ (see Table 1), as part of the FP survey of the Herschel Reference Survey (Boselli et al. 2010; Gómez-López et al. 2019). To increase the S/N, the pixels have been electronically re-binned by 4 × 4, leading to a final spatial sampling of ∼ 1.3 arcsec. The galaxy was observed during poor seeing conditions (FWHM ∼2.9±0.3 arcsec) but excellent transparency with a total exposure time of 96 min, with 120 seconds per channel. The spectral domain has been selected using a FWHM = 90 Å wide filter centred at 6607 Å. The spectroscopic calibration of the data was made using the narrow Ne line at 6598.95 Å in the same interference filter. The choice of this line is optimal to minimise the phase shift effects due to the coating of the interferometer since it is very close to the redshifted Hα emission line of the galaxy (∼6597 Å). The theoretical PUMA effective finesse is F∼ 24 but the one measured using the narrow Ne spectral line was F_e ∼ 19.9. The FP interference order at the mean wavelength of the observation is p ∼ 328, which gives a resolution R = pF_e ∼ 6533. The free spectral range (FSR) of the instrument at the redshifted wavelength of the galaxy, ∼20.04 Å (911 km s⁻¹), was scanned through 48 channels, corresponding to a spectral sampling of ∼ 0.42 Å (19 km s⁻¹). The instrumental and observational parameters are listed in Table A.1.

Standard corrections were applied to the CCD images. The astrometric calibration of the data has been done using the KOORDS task of the KARMA⁵ package (Gooch 1996) and WCS (Word Coordinate System) textscIRAF⁶ task. The FP data were reduced and analysed using the ADOCHW⁷ software, the IDL-based COMPUTEEVERYTHING⁸, REDUCWIZARD interface ⁹, and our own PYTHON scripts.

The filter includes the Hα line and the [NII] lines at 6548.03 and 6583.48 Å (see Fig. C.1). Since the FSR of the interferometer is ∼ 4.5 times smaller than the filter width, the [NII] lines are present in the data cube but are respectively observed at the interference order p + 1 and p − 1 with respect to Hα, which is observed at order p. Furthermore, the apparent separation between Hα and [NII]6548 lines is ∼+5.0 Å and the one with the [NII]6583 is ∼+0.9 Å, as illustrated in the bottom panel of Fig. C.1.

Since the spectral resolution at Hα redshifted is ∼1.01 Å, the [NII]6548 line is resolved from the Hα line while the [NII]6583 line partially overlaps with Hα, affecting the barycentre position and the line width measurements of the Hα line. In order to correct for this effect, we make a model that takes into account (i) a constant line ratio of [NII]6583/Hα = 0.35, as measured using the VLT/FORS long slit spectrum published in Fossati et al. (2018), (ii) the known wavelength separation of ∼20.73 Å (∼ 942 km s⁻¹) between the redshifted Hα and [NII]6583 lines, and (iii) the Hα₀ FSR at rest ∼19.8Å (∼ 905 km s⁻¹). Considering that the dispersion of the instrumental line spread function (LSF) is σ_LSF∼ 0.43 Å (∼ 19.5 km s⁻¹), the thermal broadening σ_T∼ 0.18 Å (∼9.1 km s⁻¹), and the intrinsic mean velocity dispersion σ∼ 0.44 Å (∼ 20.0 km s⁻¹), the model indicates that the observed Hα barycentre and the velocity dispersion should be corrected by -0.20 Å (-8.9 km s⁻¹) and -0.12 Å (-5.3 km s⁻¹), respectively.

We computed a parabolic phase map from the calibration cube in order to obtain the reference wavelength for the line profile observed inside each pixel. This phase map provides the shift that has to be applied in the spectral dimension to every pixel of the interferogram cube to bring all channels to the same wavelength. The wavelength-sorted data cube is then created by applying the phase map correction to the interferogram data cube. A detailed explanation about the data reduction process for data obtained at OAN-SPM with PUMA can be found in Fuentes-Carrera et al. (2004).

In order to work with profiles and maps with a similar S/N and at the same time optimise the angular resolution, we applied a Voronoi tessellation with a constant S/N = 9 to the data cube Cappellari & Copin (2003), Daigle et al. (2006a,b). The brightest HII regions are not affected by this smoothing procedure because their S/N per pixel is larger than 9 and are thus displayed at full resolution.

The resulting data cubes are used to extract the monochromatic, continuum, LoS radial velocity and velocity dispersion maps of the galaxy. The continuum map is computed considering the average of the 3 lowest intensities of the 48 channels of the cube. For the monochromatic image, the intensity of the Hα line is obtained by integrating the flux of the line profile for each pixel. Following Daigle et al. (2006a), we calculate the radial velocity in each pixel by measuring the barycentre of the profile of the Hα line. This barycentre is computed as the intensity weighted centroid of the spectral bins falling into the emission line boundaries, with the continuum subtracted. Radial velocities and velocity dispersions are extracted using a single emission-line detection algorithm. In the case where more than one velocity component is present in the spectrum, only the strongest emission line is taken into account. When two emission lines are spectrally close and have comparable amplitudes, they might be taken as a single one with a larger velocity dispersion.

Assuming that all the Hα emission-profiles are described by Gaussian functions, from the emission profiles width map the velocity dispersion σ was estimated as

$$\begin{array}{c}\hfill \sigma ={({\sigma }_{\mathrm{LoS}}^2-{\sigma }_{\mathrm{LSF}}^2-{\sigma }_{\mathrm{T}}^2)}^{1/2},\end{array}$$(A.1)

where σ_LoS is the observed LoS velocity dispersion directly measured from the data cube corrected (i) from the instrumental LSF (e.g. Valdez-Gutiérrez et al. 2001; Rosado et al. 2013; Cárdenas-Martínez & Fuentes-Carrera 2018), which is modelled using the instrumental finesse F and the FSR by σ_LSF = FSR/ $(\mathrm{F}\times 2\sqrt{2\ln 2})\sim$19.5  km s⁻¹ (see Table A.1) and (ii) from the thermal broadening σ_T = (kT_e/m_H)^1/2 = 9.1   km s⁻¹, assuming an electronic temperature of HII regions of T_e = 10⁴   K (Osterbrock & Ferland 2006).

A.2. GHASP observations

An independent set of FP spectroscopic observations was taken using the GHASP instrument on the 1.93 m telescope at the OHP in spring 2021 (northern emisphere). These observations were gathered using the same configuration extensively described in Gómez-López et al. (2019). The GHASP focal reducer is equipped with an FP with a field of view of 5.8 × 5.8 arcmin² and a 512 × 512 Imaging Photon Counting System¹⁰ with a pixel scale of 0.68 × 0.68 arcsec² (Gach et al. 2002). The observations were taken during poor seeing conditions (2.5 < FWHM < 4.0 arcsec), with medium transparent sky, and with a total integration time of 645 minutes. The FSR of the FP at the redshifted wavelength of the galaxy (379 km s⁻¹) was scanned through 32 channels, with a typical spectral resolution of R∼ 9704 at Hα redshifted. NGC 4330 was observed using two different narrow-band interference filters. Indeed, the mean redshifted Hα wavelength of the galaxy is ∼ 6597 Å, which unfortunately falls just in the middle of two available interference filters centred on 6591 Å and 6601 Å. Both filters have the same FWHM ∼15 Å useful to overcome [NII] line contamination. To cover the whole velocity field of the galaxy, the first filter was used for the blueshifted part of the galaxy, the second one for the redshifted side. This does not exactly double the observing time because approximatively two-thirds of the velocity amplitude was covered by the two filters that overlap by ∼5 Å. Similarly to the PUMA observation, the wavelength calibrations were done using the NeI emission line at λ = 6598.95 Å. The data reduction procedure adopted to reduce the GHASP data have been extensively described in Epinat et al. (2008b) and Gómez-López et al. (2019). The LoS velocity dispersion has been corrected from the LSF and thermal broadening using the method and the relations described in the previous section (Sect. A.1). The main observational parameters of the FP observations are summarised in Table A.1.

We decided not to combine the FP observations obtained with GHASP and PUMA because of their different spatial and spectral resolutions. On the contrary, we decided to keep them separated to i) check the consistency of the results in particular in the low-surface-brightness regions where the signal-to-noise is limited, ii) take benefit of the slightly better angular resolution of the PUMA data for the study of the extraplanar gas, and iii) of the higher spectral resolution of the GHASP data for the study of the velocity dispersion of the gas, and finally iv) use the deeper GHASP data for the comparison with the output of the simulations. We notice that for these purposes the observations were analysed differently: using a Voronoi adaptive spatial smoothing for the PUMA data in order not to thicken the almost edge-on disc of the galaxy, while using a Gaussian smoothing for the GHASP observations to detect the weakest structures.

The different data cubes and maps presented in the next sections were thresholded using the VESTIGE image; all pixels below a surface brightness limit of 1 and 3.5 × 10⁻¹⁷ erg s⁻¹ cm⁻² arcsec⁻², for GHASP and PUMA, respectively, have been masked. Nevertheless, because of the night sky lines contamination, for the GHASP data, even above this threshold limit, it was not possible to calculate a consistent line flux nor velocity, and velocity dispersion for each pixel in the regions labelled B, C, D, H, J (western tail only) and M. For each region, the profiles were summed to give average values, which in turn lead to average flux, LoS, and velocity dispersions.

Appendix B: Derivation of the rotation curve

B.1. Position velocity diagrams

In order to study the velocity distribution as a function of the Hα surface brightness, we built PVDs from the FP data cube after having subtracted the stellar continuum (e.g. Epinat et al. 2008b; Rosado et al. 2013). Using a pseudo-slit of 2 pixels in width (i.e. ∼2.6/1.4 arcsec for PUMA/GHASP, close to or below the seeing) centred at the photometric centre and aligned to the major axis defined in Sect. 3, we extracted the emission intensity projected along the major axis of the galaxy after rotating the data cube in order to align the position angle of the galaxy with the north-south direction (see Figure 8). We modified the PA and the centre of the pseudo-slit in order to obtain the most symmetric PVD in intensity, the final values correspond to the kinematic major axis (PA_kin) and kinematic centre (see Table B.1).

B.2. The intensity-peak method

The IPM is based either on a barycentric measurement of the line profiles for each pixel or by fitting it with a theoretical function such as a Gaussian (e.g. Sofue & Rubin 2001; Rosado et al. 2013). With this method, the tangential velocity in the galactic plane, V_rot, is computed under the assumption that the galaxy has an infinitesimally flat disc, inclined by an angle i with respect to the sky plane and dominated by rotational motions around an axis perpendicular to the galactic plane. Assuming an axial symmetry about the galactic centre and considering that the velocity field is dominated by circular motions, we can measure the velocity of the pixels that are at the same distance from the kinematic centre along the major axis by dividing the galaxy on the sky-plane into rings. The rotation velocity V_rot(r) for a ring of radius r is given as the azimuthal average over the ring of the deprojected velocities:

$$\begin{array}{c}\hfill V_{\mathit{rot}}(r,\theta )=\frac{V_{\mathrm{LoS}}(r,\theta )-V_{\mathit{sys}}}{cos\theta sini},\end{array}$$(B.1)

where r and θ are the radial and angular coordinates in the plane of the galaxy (Mihalas & Binney 1981) and V_LoS(r,  θ) is the LoS velocity extracted from each pixel of the FP velocity field (e.g. Fuentes-Carrera et al. 2004; Epinat et al. 2008a; Cárdenas-Martínez & Fuentes-Carrera 2018; Sardaneta et al. 2020).

This method has been used for the PUMA and GHASP datasets; the results are shown by the black filled and open squares in Fig. 3. In both cases, the RCs were computed considering all the pixels located within an angular sector of ±5° from the kinematic major axis, V_rot have been computed within circular rings of 6.5 arcsecond a galaxy kinematical inclination i = 88°. Both sides of each RC show a slowly rising solid-body shape despite the fact that both sides do not match outside the inner radius r₀∼30 arcsec (∼2.4 kpc), but the increasing trend of the curve continue up to ∼95 arcsec (∼8 kpc). We have indicated r₀ by dashed lines in all the maps and PVDs all over the paper. Beyond this radius, differences are observed between the receding and approaching sides. Discrepant points are associated with several regions previously identified in Fig. 2: region A on the hook-like structure and region J located at about the same galactocentric distance show a velocity discrepancy of ∼110 km s⁻¹ in the GHASP RC and ∼50 km s⁻¹ in the PUMA RC; regions K and L in the south-west tail do not have a counterpart on the other side of the galaxy. The deepest GHASP data show velocities fairly constant (V_rot ∼ 122 ± 10 km s⁻¹) in the south-west tail, consistent with a streaming motion.

B.3. The envelope-tracing method

The ETM uses the PVD traced along the major axis of the galaxy to simulate a pseudo-slit that has a width equal to two pixels (see Sect. B.1) and adopts all the kinematic parameters (centre, position angle, and systemic velocity; see Table B.1) derived from it, with the exception of the inclination, which cannot be determined using this method. Following Sofue et al. (1997, 1999a), García-Ruiz et al. (2002), the rotation velocity is defined as

$$\begin{array}{c}\hfill V_{\mathrm{rot}}sini=(V_{\mathrm{env}}-V_{\mathrm{sys}})-{({\sigma }_{\mathrm{T}}^2+{\sigma }_{\mathrm{LSF}}^2)}^{1/2},\end{array}$$(B.2)

where V_sys is the systemic velocity, i is the galaxy inclination with respect to the sky plane, σ_LSF is the instrumental velocity dispersion, σ_T is the thermal broadening of the gas, and V_env is the terminal velocity characterised by the intensity of the envelope of an observed PVD along the line profile¹¹. The intensity of the envelope is defined as

$$\begin{array}{c}\hfill I_{\mathrm{env}}={[{(\eta I_{\max })}^2+{(I_{\min })}^2]}^{1/2},\end{array}$$(B.3)

where I_max is the maximum intensity in the line profile, I_min is a minimum value of the contour, typically taken to be 3 × rms noise in the PVD, and η is a constant that determines the fraction of maximum intensity to set the envelope intensity, which is usually taken in the range 0.2 − 0.5 (see references above). An arc-tangent model reproduces fairly well the shape of the envelope of the PVD with the smallest number of arguments and provides an adequate match to most RCs (Courteau 1997). This function is given by

$$\begin{array}{c}\hfill V_{\mathrm{env}}(r)=V_{\mathrm{sys}}+V_0{tan}^{-1}\left(\frac{r-r_0}{r_t}\right)+C,\end{array}$$(B.4)

where V_sys, V₀, r₀, r_t and C are five free parameters: V_sys is the systemic velocity, V₀ is the velocity amplitude with respect to the velocity C at r = r₀, r₀ and r_t are two characteristic radii that shape the curve marking the radius between the rising and flat part of the RC and the sharpness of this transition, respectively. r₀ is indeed the inflection point of the function visible on the purple curve on Fig. 3, and represented by two dashed lines on the Hα velocity fields (Fig. 5). This arc-tangent model has also been used by other authors Weiner et al. (2006), Drew et al. (2018), Zhao et al. (2021). Because the PVD is traced along the major axis, relation (B.4) can be transformed into a RC using relation (B.2). Therefore, to apply the ETM, we determined the terminal velocity V_env by fitting the arc-tangent model to both the maximum I_max and minimum I_min intensity line profiles. We corrected the velocities for σ_LSF for both instruments (see Table B.1) and for σ_T = 9.1 km s⁻¹ to account for the dispersion due to the gas temperature along the LoS (see Sect. 2). The results of the ETM are shown in Fig. 3. The rotation velocities and its best fit are given on the upper panels a and c while the PVDs from which the rotation velocities have been computed are shown on the lowest panel b and d. Using Eq. (B.3), we computed the intensity of the envelope (I_env) for η=[0.2-0.5] and then we applied relation (B.2). The shaded area corresponds to η=[0.2-0.5] and the purple curve to η = 0.3. For η = 0.3, the maximum coefficient of determination is R² = 0.8 and was obtained with the parameters given in Table B.1.

B.4. The iteration method

The IPM and the ETM are limited by the observational parameters such as the seeing, the spectral resolution, the slit width (or pseudo-slit width in case of two-dimensional velocity fields). They also depend on the physical properties of the galaxy such as its inclination, the velocity dispersion of the gas and its distribution within the disc. To overcome these sources of uncertainties, Takamiya & Sofue (2000) proposed a different model, the ITM, to generate a model data cube using specific radial and vertical density profiles of the ionised gas (e.g. exponentially declining), RCs, velocity dispersions, and viewing angles. First, to provide the initial conditions from which a first model data cube is computed, we used a RC derived from a PVD. Second, the ETM is used to derive a new RC from the first model data cube, which is compared to the data. Third, the differences between these two RCs are added to the input RC to generate a new model. This procedure is repeated until the difference between the RCs derived from the data and the final model satisfies a convergence criterion or until a maximum number of iterations has been performed (e.g. Takamiya & Sofue 2000; Sofue & Rubin 2001; Heald et al. 2006, 2007). As suggested by the intensity-weighted and the envelope-tracing RCs of NGC 4330, the inner regions of the disc within R₂₃ show a solid-body shape (see Fig. 3 and 5).

Therefore, following Swaters et al. (1997), Heald et al. (2006), Rosado et al. (2013) and Moiseev (2014), we applied to NGC 4330 the classically called cylindrical model. This model assumes, as a first grade approximation, that the galaxy has the shape of a cylinder that rotates as a rigid body around the z axis perpendicular to the x, y plane, which defines the disc. In order to compute the model data cube, we proposed an initial RC that describes the rotation of the cylindrical shaped rigid body based on the arc-tangent model of a RC (Eq. B.4). Using the non-linear least-squares minimisation of the library LMFIT from PYTHON, we computed a model data cube for each FP dataset whose emission profiles were represented by Gaussian functions. The velocity dispersion of the emission peaks was fixed for all profiles to the medium value of the velocity dispersion of the galaxy ⟨σ_LoS⟩. For the amplitude of the Gaussian functions, we assumed that the gas has an exponentially declining distribution (Mo et al. 2010; Takamiya & Sofue 2002) within the disc of the galaxy along the z direction, B(z) = B₀exp(−|z|/z_d), and we measured the vertical disc length z_d = 6 arcsec using the Hα surface brightness distribution. We used the value of the kinematic centre of each dataset as the systemic velocity of the galaxy (V_sys) and the kinematic PA obtained from the PVDs. This implies that the free parameters of the model are: the transition radius between the rising and flat part of the RC (r_t), the asymptotic velocity (V₀), and the inclination (i) of the galaxy.

Appendix C: Additional figures

Fig. C.1. Model of the PUMA spectrum. The Hα plus the two [NII]6548, 6583 lines, transmitted by the filter, are represented by Gaussian functions of width σLSF (see Table 1). Normalised units to the Hα line are used for the y axis. The relative intensity between the Hα and the [NII] lines are given in Sect. A.1. The distance between the two green dashed vertical lines represents the FSR ( ∼ 20 Å). Top panel: Redshifted Hα and [NII]6548, 6583 lines are located at their right position. The distance between the two grey dashed vertical lines represents the filter width (FWHM ∼ 90 Å). Bottom panel: Redshifted lines located at their apparent position due to the FSR recovering. The [NII]6548 line is shifted by plus one FSR, and the [NII]6583 is shifted by minus one FSR. The green curve represents the sum of the Hα plus [NII]6583 lines, whereas the red one represents the Gaussian that best fits the sum of the two lines (best fit method).

Fig. C.2. Long-slit spectra obtained with MISTRAL at the 1.93 m OHP telescope. The Hα lines of region D (red spectrum) and the position of the orthogonal projection of region D on the major axis of the galaxy (blue spectrum) are fitted by Gaussian functions, with intensities in normalised units. The two vertical dashed lines on Hα indicate the measured velocity of region D (1458±46 km s−1) and on the disc (1629±46 km s−1).

Fig. C.3. GHASP velocity dispersion properties of NGC 4330. Panel (a): Intensity histogram. Panel (b): Loci, in the galaxy, of the different areas shown in panel (e) with different colours (red, green, blue, and light blue), on top of which Hα isocontours have been superimposed. Panel (c): Velocity dispersion versus emission intensity diagnostic diagram. Each black dot represents a pixel, and the blue dots represent the mean of the intensity in each velocity dispersion bin. Panel (d): Velocity dispersion histogram. Panel (e): Pixel density map of the velocity dispersion versus intensity diagram. The colour bar represents pixel density. The dashed line in plots (c) and (d) represents the LSF σLSF.

Fig. C.4. HI emission data of NGC 4330 from Chung et al. (2009). Panel (a): Cold gas velocity field on which the positions of pseudo-slits are superimposed (∼28 arcsec wide) parallel to the kinematic major axis and the ellipse fitted to the isophote of surface brightness μ(i) = 23.5  mag arcsec−2 of the NGVS i-band image, which traces the stellar disc. Panels (b) and (d): Outermost contour of the PVDs. Panel (b) corresponds to the three pseudo-slit positions drawn in panel (a), and panel (d) corresponds to the central pseudo-slit position. The red and blue circles indicate the position of the minimum intensity of the receding and approaching side, respectively, while the dashed line indicates the maximum intensity. Both entities are used to estimate the terminal velocity and to compute the RC using the ETM method. Panel (c): RC computed using the ETM, where the terminal velocity is represented by empty blue circles for the receding side and the filled red circles for the approaching side. The black solid line indicates the best fit to these points using η = 0.3, the shaded area shows the rotation velocity amplitude for η = [0.2, 0.5], the purple dashed line indicates the RC computed from the PVD maximum intensity, and the red and blue vertical lines indicate the photometric and effective radius, respectively.

Fig. C.5. Comparison between the GHASP Hα (upper panel) and the HI (lower panel, Chung et al. 2009) velocity fields of NGC 4330. The inner and outer black contours are at Σ(Hα) = 3 × 10−18 and 1.25 × 10−16 erg s−1 cm−2 arcsec−2, respectively, and the red contour Σ(Hα) = 10−17 erg s−1 cm−2 arcsec−2. The black diamond, the cross, and the star indicate the position of the CO, Hα, and HI emission kinematic centres, respectively. To align the maps horizontally, both are rotated i = 59°, the photometric PA of the galaxy. The grey ellipses encapsulated in a box located in the upper-right corner of each panel represent the sizes of the Hα mean seeing disc (upper panel) and of the HI beam (upper panel).

Fig. C.6. Comparison between the 12CO (upper panel; Lee et al. 2017) and the PUMA Hα (lower panel) velocity fields of NGC 4330. The inner black contours and the red outer contours are at Σ(Hα) = 10−17 and 1.25×10−16 erg s−1 cm−2 arcsec−2, respectively. The black diamond and the cross indicate the position of the CO and Hα emission kinematic centres, respectively. To align the maps horizontally, both are rotated i = 59°, the photometric PA of the galaxy. The grey ellipses located in the upper-right corner of each panel represent the sizes of the CO beam (upper panel) and of the Hα mean seeing disc (lower panel).

Fig. C.7. Comparison between observed and simulated unfolded RCs. The simulated velocity field shown in panel (c) of Fig. 11 has been used to measure the RC along the major axis of the galaxy; the receding and approaching sides are respectively represented by filled orange and opened cyan stars (⋆) and superimposed onto the GHASP observed and modelled RCs already displayed in panel (c) of Fig. 3, using the same symbols, colours, and lines. The blue vertical lines indicate the position of the effective radius.

Fig. C.8. GHASP (panel (a)) and simulated Hα surface brightness (panels (b) and (c)) of NGC 4330, with scales identical to those used in Fig. 2 (d), i.e. with a threshold limit Σ(Hα) ≃ 10−17 erg s−1 cm−2 arcsec−2. The low level cuts used to display the simulated Hα images have been chosen to show the low column density ionised gas reached by the simulations (panel (b)) and at higher surface brightness levels (Σ(Hα) ≃ 2×10−17 erg s−1 cm−2 arcsec−2; panel (c)), matching the typical surface brightness limit of VESTIGE and FP observations, respectively.

All Tables

All Figures

Fig. 1. Photometric properties of NGC 4330. Left panels: NGVS i-band image at low contrast (upper panel) to see the distribution of the absorbing dust within the disc, and at high contrast in logarithmic scale (lower panel) to see the boxy shape of the disc at its edges. The red and blue ellipses indicate the i-band mean surface brightness μ(i) = 23.5 and the effective surface brightness μe = 21.5 mag arcsec−2, respectively. The green bullet and the red cross indicate the position of the photometric and kinematical centres, respectively. Right panels, from top to bottom: radial variation in the i-band surface brightness (μ), position angle (PA), ellipticity (ε), morphological inclination (iph), and the fourth degree cosine coefficient of the Fourier series (B4) parameters. The red and blue vertical lines correspond to the ellipses drawn in the upper left panel and indicate the radius R23.5(i)∼145.0 arcsec (∼11.6 kpc) and the effective radius Re(i)∼66.5 arcsec (∼5.3 kpc), respectively.

In the text

Fig. 2. Images of NGC 4330. Panel a: full resolution VESTIGE continuum-subtracted Hα narrow-band image of NGC 4330 in linear scale. Panel b: NGVS i-band image. Panels c,d: PUMA and GHASP FP Hα monochromatic images, respectively. Panel e: NGVS i-band image. The dashed black ellipse and the inner ellipse indicate the surface brightness, μ(i) = 23.5, and the effective surface brightness, μe(i) = 21.5 mag arcsec−2, respectively, and the green bullet indicates the photometric centre (see Fig. 1). Panel f: NGVS i-band image with marks on the zones of interest labelled by letters (A to M) and two parallel dashed lines at r0 ∼ 30 arcsec (∼2.4 kpc), which divide the galactic disc into three main regions: north-eastern, central, and south-western (see Sect. 4). The cyan, red, and blue contours in panels a to d show the VESTIGE surface brightness Σ(Hα) = 0.3, 1.0, and 12.5 × 10−17 erg s−1 cm−2 arcsec−2, respectively. The red contours from panels e and f show the VESTIGE surface brightness Σ(Hα) = 10−17 erg s−1 cm−2 arcsec−2. The blue contours in panel e indicate the surface brightness limit of PUMA, Σ(Hα)∼3.5 × 10−17 erg s−1 cm−2 arcsec−2, while those in panel f show the surface brightness limit of GHASP, Σ(Hα)∼10−17 erg s−1 cm−2 arcsec−2.

In the text

Fig. 3. PUMA (panels a and b) and GHASP (panels c and d) RCs (panels a and c) and PVDs (panels b and d) of NGC 4330. Panels a,c: RCs derived using three different methods: the IPM, the ETM, and the ITM (see Appendix B.3 for a further explanation on the methodology). In the case of the IPM (black squares), an azimuthal sector of ±5° in the sky plane has been taken around the major axis to compute the rotation velocities; filled and empty black squares indicate the receding and the approaching sides, respectively. The upper-left insert shows the azimuthal sector and the location of the regions pointed on the RC. For the ETM, the rotation velocities have been derived from the PVD along the major axis (z = 0), shown in panels b and d. The terminal velocities are represented by empty blue circles for the receding side and filled red circles for the approaching side; the pink line corresponds to the best fit to these points using η = 0.3, and the shaded area shows the rotation velocity amplitude for η = [0.2, 0.5] (see relation (B.3)). The ITM based on a cylindrical model (described in Appendix B.4) is represented by the black line. Panels b,d: PVDs showing the LoS velocities computed with the ETM from which the rotation velocities are computed in the upper panels. On top of the PVD, the filled red and open blue circles indicate the extremal LoS velocities and the pink curve the arctan model fitting those velocities. The red and blue vertical dotted lines in panels a and c show the photometric and effective radius, and the grey dotted lines in the lowest panels show the position of r0 (see Sect. 4).

In the text

Fig. 4. RCs of NGC 4330 derived using the 12CO(2 − 1) at λ = 1.3 mm data from Lee et al. (2017) (green), the HI (λ = 21 cm) data from Chung et al. (2009) (purple), and the FP Hα data from PUMA (black) and GHASP (blue). The RCs from CO and HI data were computed using the ETM. The solid lines correspond to the best fit to the points of the terminal velocity using η = 0.3, and the shaded areas show the rotation velocity amplitude for η = [0.2, 0.5] (see relation (B.3)); the dashed lines show the RC computed from the PVD maximum intensity. For comparison, the Hα RCs shown in Fig. 3 are also reported on this plot. The green, purple, and red vertical lines indicate the CO extension (∼47 arcsec, ∼3.8 kpc), the HI radius (∼141 arcsec, ∼11.3 kpc), and the outermost photometric radius (R23.5 ∼ 145 arcsec, ∼11.6 kpc), respectively. The top-right insert indicates the angular resolution of each dataset.

In the text

Fig. 5. Velocity fields of the ionised gas of NGC 4330 (left PUMA, right GHASP). The general features of the main disc regions (north-east, central, and south-west), identified by two parallel dashed lines, and the regions with particular features, identified with letters A to H, are discussed in the text. The black bullet and the cross indicate the position of the photometric and kinematical centre, respectively. The black arrow shows the wind direction as derived from the models of Vollmer et al. (2012, 2021).

In the text

Fig. 6. Velocity dispersion fields of NGC 4330 (left PUMA, right GHASP). The general features of the main disc regions (north-east, central, and south-west), identified by two parallel dashed lines, and the regions with particular features, identified with letters A to H, are discussed in the text. The black bullet indicates the position of the photometric centre. The black arrow shows the wind direction as derived from the models of Vollmer et al. (2012, 2021).

In the text

Fig. 7. Residual velocity fields obtained after subtracting the velocity field computed from the ITM model data cube (left PUMA, right GHASP). The model, based on a solid-body rotation, is shown in the inset. The purple parallel lines show the isovelocity contours of the edge-on rotating disc model. The bullet and the cross indicate the position of the photometric and kinematical centre, respectively. The arrow indicates the wind direction as derived from the models of Vollmer et al. (2012, 2021).

In the text

Fig. 8. PVDs parallel to the kinematic major axis extracted from the PUMA dataset. The distance of the pseudo-slit from the disc plane along the z axis (positive to the north) is given in the upper-left corner of each panel. The positions of all the pseudo-slits are shown on the Hα monochromatic map in the right panel. From the top at +3.7 to −5.5 arcsec (from ∼0.3 to ∼0.4 kpc), the intensity of the contour level is at 1.5 σ, and from −7.4 to −11 arcsec (from ∼0.6 to ∼0.9 kpc) the intensity level is at 1.2σ. In the left panels, coloured lines overlaid on the PVD represent the best fit of the terminal velocity (Relation (B.3) with η = 0.3) according to the arc-tangent model (Eq. (B.4)). There is a correspondence in colour between the terminal velocities and the pseudo-slits plotted on the monochromatic map in the right panel. In the middle panels, the PVDs derived from the data cube model, which was built in the context of the ITM to compute the RC (see Appendix B.4), are overlaid on the observed PVD in grey. All the PVDs extracted from the data cube model are at 8% of the maximum intensity of the model. The two parallel dashed lines are at r0 = ±30 arcsec (±2.4 kpc) and indicate the three main disc regions: north-east, central, and south-west. The different features of interest are indicated with letters A to H.

In the text

Fig. 9. PVDs parallel to the kinematic minor axis obtained from the PUMA dataset. The locations of the pseudo-slits are in the upper-left corner of each panel, and they are traced on the Hα monochromatic map in the left panel with colour correspondence. The distance increases from the south to the north of the disc, with zero at the kinematic centre. The black dotted ellipse on the Hα monochromatic map indicates the i-band μe(i) = 23 mag arcsec2 isophote. The intensity of the contour levels is 1.3, 1.4, and 1.6 σ from −45.5 to +63 arcsec (∼−3.6 to 5 kpc) and 1, 1.1, and 1.25σ from −130.9 to −52 arcsec (∼−10.5 to −4.2 kpc). The PVDs extracted from the data cube model (black contours) are at 1, 10, and 50% of the maximum intensity of the model.

In the text

	Fig. 10. Envelopes corrected from the morphological inclination (89.5°) of the different PVDs shown in Fig. 8. The vertical blue dotted line indicates the effective radius Re(i).
In the text

Fig. 11. GHASP (panels a and d) and simulated (panels b, c, e and f) Hα velocity (left column) and velocity dispersion (right column) fields of NGC 4330. The simulated Hα velocity and velocity dispersion fields are given down to the faintest levels used in the middle panel of Fig. C.8 to show the kinematical properties of the low column density ionised gas reached by the simulations (panels b and e) and at a higher surface brightness level corresponding to the sensitivity of the FP observations (Σ(Hα)≃ 2 × 10−17 erg s−1 cm−2 arcsec−2, panel c). The contours on the velocity fields indicate the Hα surface brightness Σ(Hα)∼ 2×10−18 erg s−1 cm−2 arcsec−2 taken from the VESTIGE data (Fossati et al. 2018).

In the text

Fig. C.1. Model of the PUMA spectrum. The Hα plus the two [NII]6548, 6583 lines, transmitted by the filter, are represented by Gaussian functions of width σLSF (see Table 1). Normalised units to the Hα line are used for the y axis. The relative intensity between the Hα and the [NII] lines are given in Sect. A.1. The distance between the two green dashed vertical lines represents the FSR ( ∼ 20 Å). Top panel: Redshifted Hα and [NII]6548, 6583 lines are located at their right position. The distance between the two grey dashed vertical lines represents the filter width (FWHM ∼ 90 Å). Bottom panel: Redshifted lines located at their apparent position due to the FSR recovering. The [NII]6548 line is shifted by plus one FSR, and the [NII]6583 is shifted by minus one FSR. The green curve represents the sum of the Hα plus [NII]6583 lines, whereas the red one represents the Gaussian that best fits the sum of the two lines (best fit method).

In the text

Fig. C.2. Long-slit spectra obtained with MISTRAL at the 1.93 m OHP telescope. The Hα lines of region D (red spectrum) and the position of the orthogonal projection of region D on the major axis of the galaxy (blue spectrum) are fitted by Gaussian functions, with intensities in normalised units. The two vertical dashed lines on Hα indicate the measured velocity of region D (1458±46 km s−1) and on the disc (1629±46 km s−1).

In the text

Fig. C.3. GHASP velocity dispersion properties of NGC 4330. Panel (a): Intensity histogram. Panel (b): Loci, in the galaxy, of the different areas shown in panel (e) with different colours (red, green, blue, and light blue), on top of which Hα isocontours have been superimposed. Panel (c): Velocity dispersion versus emission intensity diagnostic diagram. Each black dot represents a pixel, and the blue dots represent the mean of the intensity in each velocity dispersion bin. Panel (d): Velocity dispersion histogram. Panel (e): Pixel density map of the velocity dispersion versus intensity diagram. The colour bar represents pixel density. The dashed line in plots (c) and (d) represents the LSF σLSF.

In the text

Fig. C.4. HI emission data of NGC 4330 from Chung et al. (2009). Panel (a): Cold gas velocity field on which the positions of pseudo-slits are superimposed (∼28 arcsec wide) parallel to the kinematic major axis and the ellipse fitted to the isophote of surface brightness μ(i) = 23.5  mag arcsec−2 of the NGVS i-band image, which traces the stellar disc. Panels (b) and (d): Outermost contour of the PVDs. Panel (b) corresponds to the three pseudo-slit positions drawn in panel (a), and panel (d) corresponds to the central pseudo-slit position. The red and blue circles indicate the position of the minimum intensity of the receding and approaching side, respectively, while the dashed line indicates the maximum intensity. Both entities are used to estimate the terminal velocity and to compute the RC using the ETM method. Panel (c): RC computed using the ETM, where the terminal velocity is represented by empty blue circles for the receding side and the filled red circles for the approaching side. The black solid line indicates the best fit to these points using η = 0.3, the shaded area shows the rotation velocity amplitude for η = [0.2, 0.5], the purple dashed line indicates the RC computed from the PVD maximum intensity, and the red and blue vertical lines indicate the photometric and effective radius, respectively.

In the text

Fig. C.5. Comparison between the GHASP Hα (upper panel) and the HI (lower panel, Chung et al. 2009) velocity fields of NGC 4330. The inner and outer black contours are at Σ(Hα) = 3 × 10−18 and 1.25 × 10−16 erg s−1 cm−2 arcsec−2, respectively, and the red contour Σ(Hα) = 10−17 erg s−1 cm−2 arcsec−2. The black diamond, the cross, and the star indicate the position of the CO, Hα, and HI emission kinematic centres, respectively. To align the maps horizontally, both are rotated i = 59°, the photometric PA of the galaxy. The grey ellipses encapsulated in a box located in the upper-right corner of each panel represent the sizes of the Hα mean seeing disc (upper panel) and of the HI beam (upper panel).

In the text

Fig. C.6. Comparison between the 12CO (upper panel; Lee et al. 2017) and the PUMA Hα (lower panel) velocity fields of NGC 4330. The inner black contours and the red outer contours are at Σ(Hα) = 10−17 and 1.25×10−16 erg s−1 cm−2 arcsec−2, respectively. The black diamond and the cross indicate the position of the CO and Hα emission kinematic centres, respectively. To align the maps horizontally, both are rotated i = 59°, the photometric PA of the galaxy. The grey ellipses located in the upper-right corner of each panel represent the sizes of the CO beam (upper panel) and of the Hα mean seeing disc (lower panel).

In the text

Fig. C.7. Comparison between observed and simulated unfolded RCs. The simulated velocity field shown in panel (c) of Fig. 11 has been used to measure the RC along the major axis of the galaxy; the receding and approaching sides are respectively represented by filled orange and opened cyan stars (⋆) and superimposed onto the GHASP observed and modelled RCs already displayed in panel (c) of Fig. 3, using the same symbols, colours, and lines. The blue vertical lines indicate the position of the effective radius.

In the text

Fig. C.8. GHASP (panel (a)) and simulated Hα surface brightness (panels (b) and (c)) of NGC 4330, with scales identical to those used in Fig. 2 (d), i.e. with a threshold limit Σ(Hα) ≃ 10−17 erg s−1 cm−2 arcsec−2. The low level cuts used to display the simulated Hα images have been chosen to show the low column density ionised gas reached by the simulations (panel (b)) and at higher surface brightness levels (Σ(Hα) ≃ 2×10−17 erg s−1 cm−2 arcsec−2; panel (c)), matching the typical surface brightness limit of VESTIGE and FP observations, respectively.

In the text