# A Virgo Environmental Survey Tracing Ionised Gas Emission (VESTIGE)

M. Sardaneta, P. Amram, A. Boselli, B. Vollmer, M. Rosado, M. Sánchez-Cruces, A. Longobardi, C. Adami, M. Fossati, B. Epinat, et al.

► To cite this version:

M. Sardaneta, P. Amram, A. Boselli, B. Vollmer, M. Rosado, et al., A Virgo Environmental Survey Tracing Ionised Gas Emission (VESTIGE). Astronomy and Astrophysics - A&A, 2022, 659, pp.A45. 10.1051/0004-6361/202142270. [hal-03607810](https://hal.science/hal-03607810)

HAL Id: hal-03607810

https://hal.science/hal-03607810

Submitted on 17 Mar 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# A Virgo Environmental Survey Tracing Ionised Gas Emission
(VESTIGE)

## XI. Two-dimensional Hα kinematics of the edge-on ram pressure stripped galaxy NGC 4330*,**

M. M. Sardaneta<fnref n="1" />, P. Amram<fnref n="1" />, A. Boselli<fnref n="1" />, B. Vollmer<fnref n="2" />, M. Rosado<fnref n="3" />, M. Sánchez-Cruces<fnref n="1" />,<fnref n="3" />, A. Longobardi<fnref n="4" />,
C. Adami<fnref n="1" />, M. Fossati<fnref n="4" />, B. Epinat<fnref n="1" />, M. Boquien<fnref n="5" />, P. Côté<fnref n="6" />, G. Hensler<fnref n="7" />, Junais<fnref n="1" />, H. Plana<fnref n="8" />, J. C. Cuillandre<fnref n="9" />,
L. Ferrarese<fnref n="6" />, J. L. Gach<fnref n="1" />, J. A. Gomez-Lopez<fnref n="1" />, S. Gwyn<fnref n="6" />, and G. Trinchieri<fnref n="10" />

Received 21 September 2021 / Accepted 17 December 2021

## ABSTRACT

Using the VESTIGE survey, a deep narrow-band Hα imaging survey of the Virgo cluster carried out at the CFHT with MegaCam, we discovered a long and diffuse tail of ionised gas in the edge-on late-type galaxy NGC 4330. This peculiar feature indicates an ongoing ram pressure stripping event able to remove the gas in the outer region of the disc. Tuned hydrodynamic simulations suggest that the ram pressure stripping event is occurring almost face-on, making NGC 4330 the ideal candidate for studying the effects of the perturbation in the direction perpendicular to the disc plane. We present here two new independent sets of Fabry-Perot observations ($R \simeq 10\,000$) with the purpose of understanding the effects of the ram pressure stripping process on the kinematics of the ionised gas. Despite their limited sensitivity to the diffuse gas emission, the data allowed us to measure the velocity and the velocity dispersion fields over the galaxy disc and in several features at the edges or outside the stellar disc formed after the ram pressure stripping event. We constructed the position-velocity diagrams and the rotation curves of the galaxy using three different techniques. The data show, consistent with the hydrodynamic simulations, that the galaxy has an inner solid-body rotation up to ~2.4 kpc, with non-circular streaming motions outside this radius and in the several external features formed during the interaction of the galaxy with the surrounding intracluster medium. The data also indicate a decrease in the rotational velocity of the gas with increasing distance from the galaxy disc along the tails, suggesting a gradual but not linear loss of angular momentum in the stripped gas. Consistent with a ram pressure stripping scenario, the i-band image shows a boxy shape at the south-west edge of the disc, where the stellar orbits might have been perturbed by the modification of the gravitational potential well of the galaxy due to the displacement of the gas in the z direction.

**Key words.** galaxies: clusters: general – galaxies: clusters: individual: Virgo – galaxies: evolution – galaxies: interactions – galaxies: kinematics and dynamics

<footnote n="1">¹ Aix-Marseille Univ., CNRS, CNES, LAM, Marseille, France
e-mail: minerva.munoz@lam.fr, philippe.amram@lam.fr, alessandro.boselli@lam.fr</footnote>

<footnote n="2">² Université de Strasbourg, CNRS, Observatoire astronomique de Strasbourg, UMR 7550, 67000 Strasbourg, France</footnote>

<footnote n="3">³ Instituto de Astronomía, Universidad Nacional Autónoma de México, Apartado Postal 70-264, Coyoacán, Ciudad de México 04510, Mexico</footnote>

<footnote n="4">⁴ Università di Milano-Bicocca, piazza della scienza 3, 20100 Milano, Italy</footnote>

<footnote n="5">⁵ Centro de Astronomía (CITEVA), Universidad de Antofagasta, Avenida Angamos 601, Antofagasta, Chile</footnote>

<footnote n="6">⁶ National Research Council of Canada, Herzberg Astronomy and Astrophysics, 5071 West Saanich Road, Victoria, BC V9E 2E7, Canada</footnote>

<footnote n="7">⁷ Department of Astrophysics, University of Vienna, Türkenschanzstrasse 17, 1180 Vienna, Austria</footnote>

<footnote n="8">⁸ Laboratório de Astrofísica Teórica e Observacional, Universidade Estadual de Santa Cruz, 45650-000 Ilhéus-BA, Brazil</footnote>

<footnote n="9">⁹ AIM, CEA, CNRS, Université Paris-Saclay, Université Paris Diderot, Sorbonne Paris Cité, Observatoire de Paris, PSL University, 91191 Gif-sur-Yvette Cedex, France</footnote>

<footnote n="10">¹⁰ INAF – Osservatorio Astronomico di Brera, via Brera 28, 20159 Milano, Italy</footnote>

<footnote>* Reduced datacubes are also available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsarc.u-strasbg.fr/viz-bin/cat/J/A+A/659/A45</footnote>

<footnote>** Based on observations obtained: (1) at the Observatoire de Haute Provence (OHP; France), operated by the French Centre National de la Recherche Scientifique (CNRS); (2) at the Observatorio Astronómico Nacional in San Pedro Mártir, Baja California, Mexico (OAN-SPM); and with (3) MegaPrime/MegaCam, a joint project of CEA/DAPNIA and the Canada-French-Hawaii Telescope (CFHT), which is operated by the National Research Council (NRC) of Canada, the Institut National des Sciences de l'Univers of the CNRS of France, and the University of Hawaii.</footnote>
# 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 \sim 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 \approx 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 $\sigma \sim 13 \text{ km s}^{-1}$. 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 ($\lambda_c = 6591$ Å; $\Delta\lambda = 106$ Å), which includes the Hα line and the two [NII] lines at $\lambda 6548$ Å and $\lambda 6584$ Å<fnref n="1" />. 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<fnref n="2" />. These deep observations have a typical sensitivity of $f(H\alpha) \approx 2 \times 10^{-17} \text{ erg s}^{-1} \text{ cm}^{-2} (5\sigma)$ for point sources and $\Sigma(H\alpha) \approx 10^{-18} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$ ($1\sigma$ 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<fnref n="2" /> coupled with a 512×512 charge-coupled device (CCD) camera, with pixel size of 1.27×1.27 arcsec<fnref n="2" /> and has a spectral resolution R ~ 6500. GHASP has a field of view of 5.8×5.8 arcmin<fnref n="2" />, is coupled with a 512×512 Image Photon Counting System (IPCS) with pixel size of 0.68×0.68 arcsec<fnref n="2" />, 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⁻<fnref n="1" /> 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⁻<fnref n="1" /> and an angular resolution of full width at half maximum (FWHM) of 26.4×24.0 arcsec<fnref n="2" />. ¹²CO(2–1) observations obtained with the Submillimeter Array (SMA), with a synthesised beam of 6.35×4.47 arcsec<fnref n="2" /> and a spectral resolution of 5 km s⁻<fnref n="1" />, 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 ($\mu_i = 27.4$ mag arcsec⁻<fnref n="2" />), 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).

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:

$$ \mu(\theta) = \mu_0 + \sum_{n=1,\infty} (A_n \sin n\theta + B_n \cos n\theta), \quad (1) $$

where $\mu_0$ 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_1, B_1, A_2$, and $B_2$ indicate errors in the fitting

<footnote n="1">¹ Hereafter we refer to the Hα+[NII] band simply as Hα, unless otherwise stated.</footnote>

<footnote n="2">² The total integration time in the r-band filter reported in Fossati et al. (2018) includes archival data not considered here, for the sake of consistency with the rest of the VESTIGE survey.</footnote>
<figure><img src="image_2.png" /><figcaption>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 $\mu(i) = 23.5$ and the effective surface brightness $\mu_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 ($\mu$), position angle (PA), ellipticity ($\epsilon$), morphological inclination ($i_{ph}$), and the fourth degree cosine coefficient of the Fourier series ($B_4$) parameters. The red and blue vertical lines correspond to the ellipses drawn in the upper left panel and indicate the radius $R_{23.5}(i) \sim 145.0$ arcsec ($\sim 11.6$ kpc) and the effective radius $R_e(i) \sim 66.5$ arcsec ($\sim 5.3$ kpc), respectively.</figcaption></figure>

procedure. The terms $A_3$ and $B_3$ give the ‘egg-shaped’ isophote but the most interesting is $B_4$. If $B_4 > 0$, the galaxy has discy-shape isophotes, while if $B_4 < 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 ($\epsilon$), 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) \sim 145.0$ arcsec ($\sim 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) \sim 140.3$ arcsec, $\sim 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) \sim 66.5$ arcsec ($\sim 5$ kpc) and the effective surface brightness $\mu_e(i) \sim 21.5$ mag arcsec$^{-2}$, consistently with those given in Cortese et al. (2012a). The photometric centre is located at RA(J2000) = 12h23m17.26s, Dec = +11°22′07″.9, with an

uncertainty of $\sim 0.9$ arcsec ($\sim 76$ pc). The sudden jump of the parameters observed in Fig. 1 at a radius of $\sim 5$ arcsec ($\sim 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

$$ \cos i_{ph} = \frac{(1 - \epsilon)^2 - q_0^2}{1 - q_0^2}, \qquad (2) $$

where $q_0 = 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 \pm 1°$ for $R(i) \ge 75$ arcsec ($\sim 6$ kpc) ($\mu_i \ge 22.75$ mag arcsec$^{-2}$). The $B_4$ 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 stel-
lar 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 particu-
larly 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 trun-
cated with respect to the stellar disc, indicating that the activ-
ity 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 north-
wards 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 monochro-
matic 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).

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 gen-
erally 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 method-
ology 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 tran-
sition, the RCs computed with the ETM method are flat up to
the photometric radius R<sub>23.5</sub> ~ 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 veloc-
ity obtained for each dataset with the ETM differs by 8 km s<sup>-1</sup>.
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.

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 diame-
ter 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_0$ ~ 30 arcsec (~2.4 kpc) and slightly before the photometric
radius.

The IPM allows us to analyse the rotation on both the reced-
ing 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 with-
draw 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_{\text{rot}}^{\text{max}} = 158 \pm 12$ km s$^{-1}$, the uncertainty being estimated from
their standard deviation.
<figure><img src="image_3.png" /></figure>

**Fig. 2.** Images of NGC 4330. **Panel a:** full resolution VESTIGE continuum-subtracted $H\alpha$ narrow-band image of NGC 4330 in linear scale. **Panel b:** NGVS i-band image. **Panels c,d:** PUMA and GHASP FP $H\alpha$ monochromatic images, respectively. **Panel e:** NGVS i-band image. The dashed black ellipse and the inner ellipse indicate the surface brightness, $\mu(i) = 23.5$, and the effective surface brightness, $\mu_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 $r_0 \sim 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 $\Sigma(H\alpha) = 0.3, 1.0$, and $12.5 \times 10^{-17}$ erg s$^{-1}$ cm$^{-2}$ arcsec$^{-2}$, respectively. The red contours from panels e and f show the VESTIGE surface brightness $\Sigma(H\alpha) = 10^{-17}$ erg s$^{-1}$ cm$^{-2}$ arcsec$^{-2}$. The blue contours in panel e indicate the surface brightness limit of PUMA, $\Sigma(H\alpha) \sim 3.5 \times 10^{-17}$ erg s$^{-1}$ cm$^{-2}$ arcsec$^{-2}$, while those in panel f show the surface brightness limit of GHASP, $\Sigma(H\alpha) \sim 10^{-17}$ erg s$^{-1}$ cm$^{-2}$ arcsec$^{-2}$.
<figure><img src="image_4.png" /></figure>

**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 r₀ (see Sect. 4).

## 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_{\text{rot}}^2(R)/G$, where G is the universal gravity constant, $V_{\text{rot}}$ the flat rotation velocity, and $\alpha$ a parameter spanning between the lowest value of $\alpha = 0.6$ for a galaxy dominated by a flat disc component and the higher value of $\alpha = 1.0$ if the spherical halo components dominate. Inside the radius $R_{23.5}(i$-band) ~ 11.6 kpc, using

$V_{\text{rot}}^{\text{max}} = 158 \pm 12 \text{ km s}^{-1}$ and using $\alpha = 0.8 \pm 0.2$ to be conservative, the dynamical mass of NGC 4330 within $R_{23.5}(i$-band) is $M_{\text{dyn}}(R_{23.5,i}$-band) = 4.4 \pm 0.6 \times 10^{10} M_{\odot}$. For the same maximum rotation velocity, the mass within the HI radius $R_{\text{HI}} = 11.3$ kpc is $M_{\text{dyn}}(R_{\text{HI}}) = 4.9 \pm 1.2 \times 10^{10} M_{\odot}$. These estimates can be compared to the stellar mass of the galaxy, $M_{\text{star}} = 5 \times 10^9 M_{\odot}$ (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$\alpha$ RC, except maybe within the first 2 kpc where the HI velocities have probably been overestimated (see Fig. 4). The comparison between the H$\alpha$, 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).
<table><caption>**Table 1.** General properties of NGC 4330.</caption><thead><tr><th>Variable</th><th>Value</th><th>Ref.</th></tr></thead><tbody><tr><td>Alternative names <sup>(1)</sup></td><td>HRS 124</td><td>a</td></tr><tr><td>α (J2000) <sup>(2)</sup></td><td>12<sup>h</sup>23<sup>m</sup>17.10<sup>s</sup></td><td>b</td></tr><tr><td>δ (J2000) <sup>(3)</sup></td><td>+11°22′5″.7</td><td>b</td></tr><tr><td>Type</td><td>Scd?</td><td>a</td></tr><tr><td>Distance</td><td>16.5 Mpc</td><td>c</td></tr><tr><td>Dist. from M87 <sup>(4)</sup></td><td>2.1°/604.8 kpc</td><td>TW</td></tr><tr><td>D<sub>25</sub></td><td>5.86 arcmin/28.1 kpc</td><td>d</td></tr><tr><td>i<sub>D25</sub></td><td>11.925 mag</td><td>d</td></tr><tr><td>V<sub>heliocentric</sub></td><td>1551 km s<sup>-1</sup></td><td>TW</td></tr><tr><td>V<sub>max</sub><br>rot</td><td>140 km s<sup>-1</sup></td><td>e</td></tr><tr><td>PA</td><td>59 ± 1°</td><td>TW</td></tr><tr><td>Morphological Inclination</td><td>89.5°</td><td>TW</td></tr><tr><td>M<sub>star</sub></td><td>5 × 10<sup>9</sup> M<sub>☉</sub></td><td>f</td></tr><tr><td>M(HI)</td><td>4.73 ± 1.04 × 10<sup>8</sup> M<sub>☉</sub></td><td>e</td></tr><tr><td>M(H<sub>2</sub>)</td><td>1.27 ± 0.05 × 10<sup>8</sup> M<sub>☉</sub></td><td>g</td></tr></tbody></table>

**Notes.** <sup>(1)</sup>UGC 7456, VCC 630, IRAS 12207+1138. <sup>(2)</sup>Right ascension and <sup>(3)</sup>declination of the galaxy. <sup>(4)</sup>Projected distance for M 87.

**References.** TW: This work; a: NED; b: NGVS, Ferrarese et al. (2012); c: Gavazzi et al. (1999), Mei et al. (2007); d: Cortese et al. (2012a); e: Chung et al. (2009); f: Boselli et al. (2015); g: Lee et al. (2017).

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 <sup>12</sup>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.

For the HI RC, the ETM fit gave a coefficient of determination R² = 0.97, with r₀ = 18 ± 2 arcsecond rᵢ = 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ₘₐₓ (in setting η = Iₘₐₓ = 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

<figure><img src="image_5.png" /></figure>

**Fig. 4.** RCs of NGC 4330 derived using the <sup>12</sup>CO(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 (R<sub>23.5</sub> ~ 145 arcsec, ~11.6 kpc), respectively. The top-right insert indicates the angular resolution of each dataset.

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ᵢ = 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
<figure><img src="image_6.png" /><figcaption>**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).</figcaption></figure>

<table><caption>**Table 2.** Velocities of peculiar regions derived from the PUMA and GHASP datasets.</caption><thead><tr><th rowspan="2">Region</th><th colspan="3">PUMA</th><th colspan="3">GHASP</th></tr><tr><th>V<sub>rad</sub><sup>(1)</sup></th><th>⟨σ⟩<sup>(2)</sup></th><th>V<sub>res</sub><sup>(3)</sup></th><th>V<sub>rad</sub><sup>(1)</sup></th><th>⟨σ⟩<sup>(2)</sup></th><th>V<sub>res</sub><sup>(3)</sup></th></tr></thead><tbody><tr><td>A</td><td>1659</td><td>35</td><td>-8</td><td>1673</td><td>19</td><td>8</td></tr><tr><td>B</td><td>...</td><td>...</td><td>...</td><td>1564</td><td>9</td><td>-41</td></tr><tr><td>C</td><td>...</td><td>...</td><td>...</td><td>1601</td><td>9</td><td>-21</td></tr><tr><td>D</td><td>1580</td><td>37</td><td>-54</td><td>1585</td><td>11</td><td>34</td></tr><tr><td>E</td><td>1562</td><td>38</td><td>20</td><td>1579</td><td>19</td><td>2</td></tr><tr><td>F</td><td>1543</td><td>38</td><td>21</td><td>1533</td><td>20</td><td>-2</td></tr><tr><td>G</td><td>1477</td><td>53</td><td>48</td><td>1493</td><td>19</td><td>22</td></tr><tr><td>H</td><td>1548</td><td>37</td><td>109</td><td>1496</td><td>21</td><td>51</td></tr><tr><td>I</td><td>1460</td><td>37</td><td>13</td><td>1482</td><td>25</td><td>3</td></tr><tr><td>J</td><td>1435</td><td>33</td><td>19</td><td>1542</td><td>14</td><td>116</td></tr><tr><td>K</td><td>1412</td><td>45</td><td>9</td><td>1480</td><td>15</td><td>72</td></tr><tr><td>L</td><td>1450</td><td>38</td><td>49</td><td>1501</td><td>21</td><td>99</td></tr><tr><td>M</td><td>...</td><td>...</td><td>...</td><td>1541</td><td>18</td><td>139</td></tr></tbody></table>

**Notes.** <sup>(1)</sup>Mean LoS velocity in km s<sup>-1</sup>. <sup>(2)</sup>Mean velocity dispersion in km s<sup>-1</sup>. <sup>(3)</sup>Mean residual velocity in km s<sup>-1</sup>.

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 south-ern 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.

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<sup>-1</sup>, which is ~25 km s<sup>-1</sup> 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<sup>-1</sup>) 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<sup>-1</sup> 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<sup>-1</sup>), 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 isov velocities 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\alpha$ line is marginally detected with a signal-to-noise ratio (S/N) ~2 and a mean LoS velocity of $1458 \pm 92$ km s$^{-1}$, only marginally consistent with the FP measurements of $1580 \pm 30$ km s$^{-1}$. Figure C.2 also shows the corresponding MISTRAL spectrum obtained on the galaxy disc north of region D, where the average velocity is $1629 \pm 46$ km s$^{-1}$, in very good agreement with the FP data. As for the other extrapolar 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).

At first order, the velocity dispersion of the gas is fairly constant over the disc of the galaxy, with $\sigma \sim 35 \pm 10$ km s$^{-1}$ and $\sigma \sim 21 \pm 7$ km s$^{-1}$ 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 ($\sigma_{LSF} = 19.5$ and $13.1$ km s$^{-1}$, 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\alpha$ 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\alpha$ 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$^{-1}$ and $26$ km s$^{-1}$.

## 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.

The residual velocity maps given in Fig. 7 do not show any systematic pattern in the galactic disc and have a velocity dispersion of $\sigma_{res} \sim 19$ km s$^{-1}$ for the GHASP dataset and $\sigma_{res} \sim 24$ km s$^{-1}$ for the PUMA dataset. The residual velocity is $\sim 0$ km s$^{-1}$, suggesting that the adopted cylindrical model (described in Appendix B.4) reproduces the axisymmetric
<figure><img src="image_7.png" /><figcaption>**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).</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**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).</figcaption></figure>

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\alpha$ 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).

### 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 kinematic behaviour of the major axis PVD.

Some low-surface-brightness plumes of $H\alpha$ 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\alpha$ 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.
<figure><img src="image_9.png" /><figcaption>**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 r₀ = ±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.</figcaption></figure>
<figure><img src="image_10.png" /><figcaption>**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 *μ*<sub>e</sub>(i) = 23 mag arcsec<sup>2</sup> 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.</figcaption></figure>

The PVD along the major axis presents a kinematic symme-
try 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<sup>-1</sup> to the south-west, while to the
north-east the PVDs show velocities that change by ~±25 km s<sup>-1</sup>
perpendicular to the disc. To the north-east, the gas in region A
has a radial velocity of ~200 km s<sup>-1</sup>. 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 pres-
ence 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;
<figure><img src="image_11.png" /><figcaption>**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 $R_e(i)$.</figcaption></figure>

Levy et al. 2019) and radio (e.g., Miller & Veilleux 2003; Fraternali & Binney 2006; Zschaechner & Rand 2015) wave-lengths, 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 ($\Delta V$) ≈ 9 km s⁻¹ each $\Delta 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_0$ 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 $\Delta V/\Delta 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 $\Delta V$ = 49 km s⁻¹ with respect to the one measured at $z=0$ arcsec. The difference in velocity between the northern RCs is $\Delta 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.

## 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\alpha$ 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\alpha$ emitting gas beyond $z \approx 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
<figure><img src="image_12.png" /><figcaption>**Fig. 11.** GHASP (panels a and d) and simulated (panels b, c, e and f) $H\alpha$ velocity (left column) and velocity dispersion (right column) fields of NGC 4330. The simulated $H\alpha$ 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 ($\Sigma(H\alpha) \approx 2 \times 10^{-17} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$, panel c). The contours on the velocity fields indicate the $H\alpha$ surface brightness $\Sigma(H\alpha) \sim 2 \times 10^{-18} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$ taken from the VESTIGE data (Fossati et al. 2018).</figcaption></figure>

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\alpha$ 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\alpha$ 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\alpha$ 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.

The agreement between the gas distribution in the observations and in the simulations is overall very satisfactory. The brightest $H\alpha$ knots detected in the FP observation (panel a) correspond to the brightest $H\alpha$ regions in the simulation (panels b and c). The diffuse $H\alpha$ 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\alpha$ 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\alpha$, more resistant to the external pressure.

The $H\alpha$ 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\alpha$ isocontours from Fossati et al. (2018), corresponding to a surface brightness of $\Sigma(H\alpha) \sim 2 \times 10^{-18} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$, 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\alpha$ 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^\circ21'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 \sim 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 $\Sigma(H\alpha) \approx 3 \times 10^{-18} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$, has not been reached by the present observations, limited to $\Sigma(H\alpha) \approx 10^{-17} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$, 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_4$ 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 suggest 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 \le 2$ kpc), the velocity field of the galaxy dominated by the emis-
sion of HII regions is characterised by isovelocities parallel to
the minor axis. The analysis of the RCs and of the PVDs con-
sistently 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 sur-
rounding 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 out-
side 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 angu-
lar momentum. Furthermore, the similarity of the Hα, CO, and
HI velocity fields and RCs, particularly in the inner disc, indi-
cates 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 strip-
ping 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; Con-
solandi 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 momen-
tum. 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.

**Acknowledgements.** 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

Abramson, A., Kenney, J. D. P., Crowl, H. H., et al. 2011, AJ, 141, 164

Adami, C., Basa, S., Brunel, J. C., et al. 2018, in SF2A-2018: Proceedings of the Annual meeting of the French Society of Astronomy and Astrophysics, eds. P. Di Matteo, F. Billebaud, F. Herpin, et al., 357

Amram, P., Le Coarer, E., Marcelin, M., et al. 1992, A&AS, 94, 175

Barbosa, C. E., Mendes de Oliveira, C., Amram, P., et al. 2015, MNRAS, 453, 2965

Bellhouse, C., Jaffé, Y. L., Hau, G. K. T., et al. 2017, ApJ, 844, 49

Bellhouse, C., Jaffé, Y. L., McGee, S. L., et al. 2019, MNRAS, 485, 1157

Biziyaev, D., Walterbos, R. A. M., Yoachim, P., et al. 2017, ApJ, 839, 87

Boissier, S., Boselli, A., Duc, P. A., et al. 2012, A&A, 545, A142

Boselli, A., & Gavazzi, G. 2006, PASP, 118, 517

Boselli, A., & Gavazzi, G. 2014, A&ARv, 22, 74

Boselli, A., Boissier, S., Cortese, L., et al. 2006, ApJ, 651, 811

Boselli, A., Boissier, S., Cortese, L., & Gavazzi, G. 2008a, ApJ, 674, 742

Boselli, A., Boissier, S., Cortese, L., & Gavazzi, G. 2008b, A&A, 489, 1015

Boselli, A., Eales, S., Cortese, L., et al. 2010, PASP, 122, 261

Boselli, A., Cortese, L., & Boquien, M. 2014a, A&A, 564, A65

Boselli, A., Cortese, L., Boquien, M., et al. 2014b, A&A, 564, A66

Boselli, A., Fossati, M., Gavazzi, G., et al. 2015, A&A, 579, A102

Boselli, A., Cuillandre, J. C., Fossati, M., et al. 2016a, A&A, 587, A68

Boselli, A., Roehlly, Y., Fossati, M., et al. 2016b, A&A, 596, A11

Boselli, A., Fossati, M., Cuillandre, J. C., et al. 2018a, A&A, 615, A114

Boselli, A., Fossati, M., Ferrarese, L., & Boissier, S. 2018b, A&A, 614, A56

Boselli, A., Fossati, M., Consolandi, G., et al. 2018c, A&A, 620, A164

Boselli, A., Fossati, M., Longobardi, A., et al. 2019, A&A, 623, A52

Boselli, A., Lupi, A., Epinat, B., et al. 2021, A&A, 646, A139

Boselli, A., Fossati, M., Longobardi, A., et al. 2022, A&A, 659, A46

Byrd, G., & Valtonen, M. 1990, ApJ, 350, 89

Cappellari, M., & Copin, Y. 2003, MNRAS, 342, 345

Cárdenas-Martínez, N., & Fuentes-Carrera, I. 2018, ApJ, 868, 141

Carvalho, M. S., & Plana, H. 2018, MNRAS, 481, 122

Cayatte, V., van Gorkom, J. H., Balkowski, C., & Kotanyi, C. 1990, AJ, 100, 604

Chemin, L., Balkowski, C., Cayatte, V., et al. 2006, MNRAS, 366, 812

Chung, A., van Gorkom, J. H., Kenney, J. D. P., & Vollmer, B. 2007, ApJ, 659, L115

Chung, A., van Gorkom, J. H., Kenney, J. D. P., Crowl, H., & Vollmer, B. 2009, AJ, 138, 1741

Clarke, A. J., Debattista, V. P., Roškar, R., & Quinn, T. 2017, MNRAS, 465, L79

Consolandi, G., Gavazzi, G., Fossati, M., et al. 2017, A&A, 606, A83

Cortese, L., Davies, J. I., Pohlen, M., et al. 2010, A&A, 518, L49

Cortese, L., Ciesla, L., Boselli, A., et al. 2012a, A&A, 540, A52
Cortese, L., Boissier, S., Boselli, A., et al. 2012b, A&A, 544, A101

Cortese, L., Catinella, B., & Smith, R. 2021, PASA, 38

Courteau, S. 1997, AJ, 114, 2402

Cowie, L. L., & Songaila, A. 1977, Nature, 266, 501

Daigle, O., Carignan, C., Amram, P., et al. 2006a, MNRAS, 367, 469

Daigle, O., Carignan, C., Hernandez, O., Chemin, L., & Amram, P. 2006b, MNRAS, 368, 1016

Dressler, A. 1980, ApJ, 236, 351

Drew, P. M., Casey, C. M., Burnham, A. D., et al. 2018, ApJ, 869, 58

Epinat, B., Amram, P., & Marcelin, M. 2008a, MNRAS, 390, 466

Epinat, B., Amram, P., Marcelin, M., et al. 2008b, MNRAS, 388, 500

Erroz-Ferrer, S., Knapen, J. H., Font, J., et al. 2012, MNRAS, 427, 2938

Farouki, R., & Shapiro, S. L. 1980, ApJ, 241, 928

Ferrarese, L., Côté, P., Cuillandre, J.-C., et al. 2012, ApJS, 200, 4

Fossati, M., Fumagalli, M., Boselli, A., et al. 2016, MNRAS, 455, 2028

Fossati, M., Mendel, J. T., Boselli, A., et al. 2018, A&A, 614, A57

Fraternali, F., & Binney, J. J. 2006, MNRAS, 366, 449

Fuentes-Carrera, I., Rosado, M., Amram, P., et al. 2004, A&A, 415, 451

Fumagalli, M., Krumholz, M. R., Prochaska, J. X., Gavazzi, G., & Boselli, A. 2009, ApJ, 697, 1811

Fumagalli, M., Fossati, M., Hau, G. K. T., et al. 2014, MNRAS, 445, 4335

Gach, J.-L., Hernandez, O., Boulesteix, J., et al. 2002, PASP, 114, 1043

García-Ruiz, I., Sancisi, R., & Kuijken, K. 2002, A&A, 394, 769

Gavazzi, G., Catinella, B., Carrasco, L., Boselli, A., & Contursi, A. 1998, AJ, 115, 1745

Gavazzi, G., Boselli, A., Scodeggio, M., Pierini, D., & Belsole, E. 1999, MNRAS, 304, 595

Gavazzi, G., Boselli, A., Mayer, L., et al. 2001, ApJ, 563, L23

Gavazzi, G., Bonfanti, C., Sanvito, G., Boselli, A., & Scodeggio, M. 2002a, ApJ, 576, 135

Gavazzi, G., Boselli, A., Pedotti, P., Gallazzi, A., & Carrasco, L. 2002b, A&A, 396, 449

Gavazzi, G., Boselli, A., Donati, A., Franzetti, P., & Scodeggio, M. 2003, A&A, 400, 451

Gavazzi, G., Boselli, A., van Driel, W., & O'Neil, K. 2005, A&A, 429, 439

Gavazzi, G., Boselli, A., Cortese, L., et al. 2006, A&A, 446, 839

Gavazzi, G., Fumagalli, M., Cucciati, O., & Boselli, A. 2010, A&A, 517, A73

Gavazzi, G., Fumagalli, M., Fossati, M., et al. 2013, A&A, 553, A89

Gómez, P. L., Nichol, R. C., Miller, C. J., et al. 2003, ApJ, 584, 210

Gómez-López, J. A., Amram, P., Epinat, B., et al. 2019, A&A, 631, A71

Gooch, R. 1996, in Astronomical Data Analysis Software and Systems V, eds. G. H. Jacoby, & J. Barnes, ASP Conf. Ser., 101, 80

Gullieuszik, M., Poggianti, B. M., Moretti, A., et al. 2017, ApJ, 846, 27

Gunn, J. E., Gott, J., & Richard, I. 1972, ApJ, 176, 1

Gwyn, S. D. J. 2008, PASP, 120, 212

Haynes, M. P., & Giovanelli, R. 1984, AJ, 89, 758

Heald, G. H., Rand, R. J., Benjamin, R. A., Collins, J. A., & Bland-Hawthorn, J. 2006, ApJ, 636, 181

Heald, G. H., Rand, R. J., Benjamin, R. A., & Bershady, M. A. 2007, ApJ, 663, 933

Jáchym, P., Combes, F., Cortese, L., Sun, M., & Kenney, J. D. P. 2014, ApJ, 792, 11

Jáchym, P., Sun, M., Kenney, J. D. P., et al. 2017, ApJ, 839, 114

Jáchym, P., Kenney, J. D. P., Sun, M., et al. 2019, ApJ, 883, 145

Jedrzejewski, R. I. 1987, MNRAS, 226, 747

Junais, S., Boissier, A., Boselli, M., et al. 2021, A&A, 650, A99

Kalberla, P. M. W., Kerp, J., Dedes, L., & Haud, U. 2014, ApJ, 794, 90

Korsaga, M., Epinat, B., Amram, P., et al. 2019, MNRAS, 490, 2977

Kronberger, T., Kapferer, W., Ferrari, C., Unterguggenberger, S., & Schindler, S. 2008, A&A, 481, 337

Larson, R. B., Tinsley, B. M., & Caldwell, C. N. 1980, ApJ, 237, 692

Lee, B., Chung, A., Tonnesen, S., et al. 2017, MNRAS, 466, 1382

Lequeux, J. 1983, A&A, 125, 394

Levy, R. C., Bolatto, A. D., Sánchez, S. F., et al. 2019, ApJ, 882, 84

Lewis, I., Balogh, M., De Propris, R., et al. 2002, MNRAS, 334, 673

Longobardi, A., Boselli, A., Fossati, M., et al. 2020, A&A, 644, A161

Mei, S., Blakeslee, J. P., Côté, P., et al. 2007, ApJ, 655, 144

Merluzzi, P., Busarello, G., Dopita, M. A., et al. 2013, MNRAS, 429, 1747

Merritt, D. 1983, ApJ, 264, 24

Mihalas, D., & Binney, J. 1981, Galactic astronomy: Structure and Kinematics, 2nd edn. (W H Freeman & Co.)

Miller, S. T., & Veilleux, S. 2003, ApJ, 592, 79

Mo, H., van den Bosch, F. C., & White, S. 2010, Galaxy Formation and Evolution (Cambridge University Press)

Moiseev, A. V. 2014, Astrophys. Bull., 69, 1

Moiseev, A. V., & Lozinskaya, T. A. 2012, MNRAS, 423, 1831

Moore, B., Lake, G., & Katz, N. 1998, ApJ, 495, 139

Moretti, A., Paladino, R., Poggianti, B. M., et al. 2020, ApJ, 897, L30

Munoz-Tunon, C., Tenorio-Tagle, G., Castaneda, H. O., & Terlevich, R. 1996, AJ, 112, 1636

Osterbrock, D. E., Fulbright, J. P., Martel, A. R., et al. 1996, PASP, 108, 277

Osterbrock, D. E., & Ferland, G. J. 2006, Astrophysics of Gaseous Nebulae and Active Galactic Nuclei (University Science Books)

Peng, Y.-J., Lilly, S. J., Kovač, K., et al. 2010, ApJ, 721, 193

Poggianti, B. M., Gullieuszik, M., Tonnesen, S., et al. 2019, MNRAS, 482, 4466

Quilis, V., Moore, B., & Bower, R. 2000, Science, 288, 1617

Roediger, E., & Hensler, G. 2005, A&A, 433, 875

Rosado, M., Langarica, R., Bernal, A., et al. 1995, Rev. Mex. Astron. Astrofis. Conf. Ser., 3, 263

Rosado, M., Gabbasov, R. F., Repetto, P., et al. 2013, AJ, 145, 135

Safarzadeh, M., & Scannapieco, E. 2017, ApJ, 850, 99

Sancisi, R., & Allen, R. J. 1979, A&A, 74, 73

Sardaneta, M. M., Rosado, M., & Sánchez-Cruces, M. 2020, Rev. Mex. Astron. Astrofis., 56, 71

Sheen, Y.-K., Smith, R., Jaffé, Y., et al. 2017, ApJ, 840, L7

Sofue, Y., & Rubin, V. 2001, ARA&A, 39, 137

Sofue, Y., Tutui, Y., Honma, M., & Tomita, A. 1997, AJ, 114, 2428

Sofue, Y., Tomita, A., Honma, M., & Tutui, Y. 1999a, PASJ, 51, 737

Sofue, Y., Tutui, Y., Honma, M., et al. 1999b, ApJ, 523, 136

Solanes, J. M., Manrique, A., García-Gómez, C., et al. 2001, ApJ, 548, 97

Spector, O., Finkelman, I., & Brosch, N. 2012, MNRAS, 419, 2156

Steyrleithner, P., Hensler, G., & Boselli, A. 2020, MNRAS, 494, 1114

Sun, M., Donahue, M., & Voit, G. M. 2007, ApJ, 671, 190

Swaters, R. A., Sancisi, R., & van der Hulst, J. M. 1997, ApJ, 491, 140

Takamiya, T., & Sofue, Y. 2000, ApJ, 534, 670

Takamiya, T., & Sofue, Y. 2002, ApJ, 576, L15

Tonnesen, S., & Bryan, G. L. 2009, ApJ, 694, 789

Valdez-Gutiérrez, M., Rosado, M., Georgiev, L., Borissova, J., & Kurtev, R. 2001, A&A, 366, 35

Vollmer, B., Cayatte, V., Balkowski, C., & Duschl, W. J. 2001, ApJ, 561, 708

Vollmer, B., Beck, R., Kenney, J. D. P., & van Gorkom, J. H. 2004, AJ, 127, 3375

Vollmer, B., Soida, M., Otmianowska-Mazur, K., et al. 2006, A&A, 453, 883

Vollmer, B., Braine, J., Pappalardo, C., & Hily-Blant, P. 2008, A&A, 491, 455

Vollmer, B., Soida, M., Braine, J., et al. 2012, A&A, 537, A143

Vollmer, B., Fossati, M., Boselli, A., et al. 2021, A&A, 645, A121

Wakker, B. P., & van Woerden, H. 1997, ARA&A, 35, 217

Warner, P. J., Wright, M. C. H., & Baldwin, J. E. 1973, MNRAS, 163, 163

Weiner, B. J., Willmer, C. N. A., Faber, S. M., et al. 2006, ApJ, 653, 1027

Yagi, M., Yoshida, M., Komiyama, Y., et al. 2010, AJ, 140, 1814

Zhao, Q., Sun, L., Shen, L., et al. 2021, ApJ, 913, 111

Zschaechner, L. K., & Rand, R. J. 2015, ApJ, 808, 153

Zschaechner, L. K., Rand, R. J., & Walterbos, R. 2015, ApJ, 799, 61
# 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 Observatory 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<fnref n="3" /> 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<sub>e</sub> ~ 19.9. The FP interference order at the mean wavelength of the observation is p ~ 328, which gives a resolution R = p F<sub>e</sub> ~ 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<fnref n="5" /> package (Gooch 1996) and WCS (Word Coordinate System) textscIRAF<fnref n="6" /> task. The FP data were reduced and analysed using the ADOCHW<fnref n="7" /> software, the IDL-based COMPUTEVERYTHING<fnref n="8" />, REDUCWIZARD interface<fnref n="9" />, 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 σ<sub>LSF</sub> ~ 0.43 Å (~19.5 km s⁻¹), the thermal broadening σ<sub>T</sub> ~ 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

$$ \sigma = (\sigma_{\text{LoS}}^2 - \sigma_{\text{LSF}}^2 - \sigma_{\text{T}}^2)^{1/2}, \qquad (A.1) $$

where σ<sub>LoS</sub> 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 σ<sub>LSF</sub> = FSR/

<footnote n="3">³ The PUMA CCD has in fact 2048 × 2048 px², which are electronically binned 4×4 to increase the S/N per pixel by a factor of 16 and to take into account the mean seeing of the site.</footnote>

<footnote n="4">⁴ Additional deeper PUMA observations of NGC 4330 could not be made due to the lockdown in spring 2020 (northern emisphere).</footnote>

<footnote n="5">⁵ https://www.atnf.csiro.au/computing/software/karma/</footnote>

<footnote n="6">⁶ ‘Image Reduction and Analysis Facility’, http://iraf.noao.edu/; IRAF is distributed by the National Optical Astronomy Observatories, which are operated by the Association of Universities for Research in Astronomy, Inc., under cooperative agreement with the National Science Foundation.</footnote>

<footnote n="7">⁷ ADHOCw, http://cesam.lam.fr/fabryperot/index/softwares, developed by J. Boulesteix.</footnote>

<footnote n="8">⁸ https://projets.lam.fr/projects/computeeverything</footnote>

<footnote n="9">⁹ https://projets.lam.fr/projects/fpreducwizard</footnote>
<table><caption>**Table A.1.** Set-up of the FP observations.</caption><thead><tr><th>Parameter</th><th>OAN-SPM</th><th>OHP</th></tr></thead><tbody><tr><td>Telescope</td><td>2.1 m</td><td>1.93 m</td></tr><tr><td>Instrument</td><td>PUMA<sup>(1)</sup></td><td>GHASP</td></tr><tr><td>Seeing</td><td>2.9 arcsec</td><td>2.5 - 4.5 arcsec</td></tr><tr><td>Galaxy Hα wavelength range</td><td>6593 - 6600 Å</td><td></td></tr><tr><td>velocity range</td><td>1380 - 1672 km s<sup>-1</sup></td><td></td></tr><tr><td>Observation date</td><td>2017-02-23</td><td>2021-01-14 to 18<br>2021-02-11</td></tr><tr><td>Number of cycles</td><td>1</td><td>80</td></tr><tr><td>Total exposure time</td><td>96 min</td><td>427 min</td></tr><tr><td>Scanning FP interferometer</td><td>ICOS - ET-50<sup>(2)</sup></td><td></td></tr><tr><td>Calibration (Neon line)</td><td>6598.95 Å</td><td></td></tr><tr><td>Finesse (F)</td><td>19.9</td><td>12.3±0.3<sup>(3)</sup></td></tr><tr><td>p<sup>(4)</sup> at Hα<sup>(5)</sup>=6562.8 (Hα<sup>(5)</sup> = 6596.8Å)</td><td>330 (328)</td><td>793 (789)</td></tr><tr><td>Resolution power (pF) at Hα<sup>(5)</sup></td><td>6533</td><td>9704±25</td></tr><tr><td>Free spectral range at Hα<sup>(5)</sup></td><td>20.05 Å / 911 km s<sup>-1</sup></td><td>8.35 Å / 379 km s<sup>-1</sup></td></tr><tr><td>Spectral resolution at Hα<sup>(5)</sup></td><td>1.01 Å / 45.8 km s<sup>-1</sup></td><td>0.68 Å / 30.9 km s<sup>-1</sup></td></tr><tr><td>Spectral sampling at Hα<sup>(5)</sup></td><td>0.42 Å / 19.0 km s<sup>-1</sup></td><td>0.26 Å / 11.9 km s<sup>-1</sup></td></tr><tr><td>Line spread function (σ<sub>LSF</sub>)</td><td>0.43 Å / 19.5 km s<sup>-1</sup></td><td>0.29 Å / 13.1 km s<sup>-1</sup></td></tr><tr><td>Scanning steps</td><td>48</td><td>32</td></tr><tr><td>Exposure time per channel</td><td>120 s</td><td>10 s</td></tr><tr><td>Detector Name</td><td>Spectral 2<sup>(6)</sup></td><td>IPCS (GaAs)</td></tr><tr><td>Detector size</td><td>512 × 512 pixels<sup>2</sup></td><td>512 × 512 pixels<sup>2</sup></td></tr><tr><td>Field-of-view</td><td>10 × 10 arcmin<sup>2</sup></td><td>5.8 × 5.8 arcmin<sup>2</sup></td></tr><tr><td>Image scale</td><td>1.27<sup>(7)</sup> arcsec pix<sup>-1</sup></td><td>0.68<sup>(8)</sup> arcsec pix<sup>-1</sup></td></tr></tbody></table>

Notes. (1) https://www.astrossp.unam.mx/instrumentos/interferometria/puma/docs/manualespuma/puma.html. (2) *Queensgate Instruments*, which is now called *IC Optical Systems*. (3) Mean and standard deviation of the finesse (F) measured from eight different calibration cubes. (4) p is the interference order. (5) $Hα_0$ is the wavelength of the Hα-line at rest, and $Hα$ the wavelength at the galaxy heliocentric LoS velocity (see Table 1). (6) https://www.astrossp.unam.mx/en/users/ccd-s. (7) After 4 × 4 electronic binning on the CCD. (8) After 2 × 2 electronic binning on the IPCS.

($F \times 2\sqrt{2\ln2}$) ~19.5 km s⁻¹ (see Table A.1) and (ii) from the thermal broadening $\sigma_T = (kT_e/m_H)^{1/2} = 9.1$ km s⁻¹, assuming an electronic temperature of HII regions of $T_e = 10^4$ 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 hemisphere). 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<fnref n="10" /> 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 \sim 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 approximately 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 $\lambda = 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

<footnote n="10">¹⁰ The GHASP image photon counting system (IPCS) has in fact 1024 × 1024 px², which are electronically binned 2×2 to increase by a factor of four the frame readout frequency and to take into account the mean seeing of the site. We remember that for an IPCS the S/N does not increase with the pixel binning.</footnote>
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 \times 10^{-17}$ erg s$^{-1}$ cm$^{-2}$ arcsec$^{-2}$, 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\alpha$ 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:

$$V_{rot}(r, \theta) = \frac{V_{LoS}(r, \theta) - V_{sys}}{\cos \theta \sin i}, \qquad (B.1)$$

where $r$ and $\theta$ are the radial and angular coordinates in the plane of the galaxy (Mihalas & Binney 1981) and $V_{LoS}(r, \theta)$ 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 $\pm 5^\circ$ from the kinematic major axis, $V_{rot}$ have been computed within circular rings of 6.5 arcsecond a galaxy kinematical inclination

$i = 88^\circ$. 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_0 \sim 30$ arcsec (~2.4 kpc), but the increasing trend of the curve continue up to ~95 arcsec (~8 kpc). We have indicated $r_0$ 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$^{-1}$ in the GHASP RC and ~50 km s$^{-1}$ 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} \sim 122 \pm 10$ km s$^{-1}$) 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

$$V_{rot} \sin i = (V_{env} - V_{sys}) - (\sigma_{T}^{2} + \sigma_{LSF}^{2})^{1/2}, \qquad (B.2)$$

where $V_{sys}$ is the systemic velocity, $i$ is the galaxy inclination with respect to the sky plane, $\sigma_{LSF}$ is the instrumental velocity dispersion, $\sigma_{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$^{11}$. The intensity of the envelope is defined as

$$I_{env} = [(\eta I_{max})^2 + (I_{min})^2]^{1/2}, \qquad (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 \times rms$ noise in the PVD, and $\eta$ 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

$$V_{env}(r) = V_{sys} + V_0 \tan^{-1} \left( \frac{r - r_0}{r_t} \right) + C, \qquad (B.4)$$

where $V_{sys}$, $V_0$, $r_0$, $r_t$ and $C$ are five free parameters: $V_{sys}$ is the systemic velocity, $V_0$ is the velocity amplitude with respect to the velocity $C$ at $r = r_0$, $r_0$ 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_0$ 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\alpha$ velocity fields (Fig. 5). This arc-tangent model has also been used by other authors Weiner et al. (2006), Drew et al. (2018),

<footnote n="11">11 Sofue & Rubin (2001) provide a slightly different expression that corrects for instrumental and thermal broadening after deprojecting from inclination. Despite their work often being cited as reference since then, we believe this was a typo since the term that needs to be corrected for broadening is the observable, $V_{env} - V_{sys}$, as done in the previous work of Sofue et al. (1997).</footnote>
<table><caption>**Table B.1.** Photometric and kinematic parameters.</caption><thead><tr><th rowspan="2">Parameter</th><th rowspan="2">Photometric</th><th colspan="2">Kinematic</th></tr><tr><th>PUMA</th><th>GHASP</th></tr></thead><tbody><tr><td>C<sub>α</sub> (J2000)<sup>1</sup></td><td>12<sup>h</sup>23<sup>m</sup>17.1<sup>s</sup></td><td>12<sup>h</sup>23<sup>m</sup>17.0<sup>s</sup></td><td>12<sup>h</sup>23<sup>m</sup>17.0<sup>s</sup></td></tr><tr><td>C<sub>δ</sub> (J2000)<sup>1</sup></td><td>+11° 22′ 5.7″</td><td>+11° 22′ 06.7″</td><td>+11° 22′ 06.7″</td></tr><tr><td>PA (deg)<sup>2</sup></td><td>59.0 ± 0.1</td><td>56 ± 1</td><td>59 ± 1</td></tr><tr><td>i (deg)<sup>3</sup></td><td>89.5 ± 0.1</td><td>88 ± 2</td><td>88 ± 2</td></tr><tr><td>V<sub>sys</sub> (km s<sup>-1</sup>)<sup>4</sup></td><td>...</td><td>1569 ± 2</td><td>1563 ± 2</td></tr><tr><td>V<sup>ETM</sup><sub>rot</sub> (km s<sup>-1</sup>)<sup>5</sup></td><td>...</td><td>164 ± 11</td><td>158 ± 8</td></tr><tr><td>V<sup>ITM</sup><sub>rot</sub> (km s<sup>-1</sup>)<sup>6</sup></td><td>...</td><td>162 ± 8</td><td>150 ± 21</td></tr><tr><td>⟨σ<sub>obs</sub>⟩ (km s<sup>-1</sup>)<sup>7</sup></td><td>...</td><td>35 ± 10</td><td>21 ± 7</td></tr><tr><td>σ<sub>LSF</sub> (km s<sup>-1</sup>)<sup>8</sup></td><td>...</td><td>19.5</td><td>13.1</td></tr><tr><td>V<sup>ETM</sup><sub>0</sub> (km s<sup>-1</sup>)<sup>9</sup></td><td>...</td><td>46 ± 8</td><td>34 ± 5</td></tr><tr><td>r<sup>ETM</sup><sub>0</sub> (arcsec)<sup>10</sup></td><td>...</td><td>27 ± 2</td><td>28 ± 2</td></tr><tr><td>r<sup>ETM</sup><sub>t</sub> (arcsec)<sup>11</sup></td><td>...</td><td>17 ± 4</td><td>15 ± 4</td></tr><tr><td>C<sup>ETM</sup><sub>0</sub> (km s<sup>-1</sup>)<sup>12</sup></td><td>...</td><td>98 ± 3</td><td>100 ± 2</td></tr><tr><td>V<sup>ITM</sup><sub>0</sub> (km s<sup>-1</sup>)<sup>13</sup></td><td>...</td><td>107 ± 5</td><td>110 ± 15</td></tr><tr><td>r<sup>ITM</sup><sub>t</sub> (arcsec)<sup>14</sup></td><td>...</td><td>8 ± 3</td><td>20 ± 1</td></tr></tbody></table>

**Notes.** (1): Right ascension and declination of the galaxy centre. The photometric and the kinematic centre differ by less than half the seeing disc. (2): Position angle of the major axis. (3): Column 2: photometric inclination computed with the ELLIPSE task from IRAF. Columns 3 and 4: kinematic inclination computed with the ITM. (4): Systemic velocity. (5): Rotation velocity measured with the RC computed with the ETM at R<sub>23.5</sub>. (6): Rotation velocity measured with the RC computed with the ITM at R<sub>23.5</sub>. (7): Mean velocity dispersion. (8): Instrumental velocity dispersion (also called the LSF). (9) to (14): Best fit parameters of equation B.4 to compute the RC with the ETM (see text).

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_{\text{env}}$ by fitting the arc-tangent model to both the maximum $I_{\text{max}}$ and minimum $I_{\text{min}}$ intensity line profiles. We corrected the velocities for $\sigma_{\text{LSF}}$ for both instruments (see Table B.1) and for $\sigma_{\text{T}} = 9.1 \text{ km s}^{-1}$ 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_{\text{env}}$) for $\eta=[0.2-0.5]$ and then we applied relation (B.2). The shaded area corresponds to $\eta=[0.2-0.5]$ and the purple curve to $\eta=0.3$. For $\eta=0.3$, the maximum coefficient of determination is $R^2 = 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<sub>23</sub> 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 $\langle\sigma_{\text{LoS}}\rangle$. 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_0 \exp(-|z|/z_d)$, and we measured the vertical disc length $z_d = 6 \text{ arcsec}$ using the H$\alpha$ surface brightness distribution. We used the value of the kinematic centre of each dataset as the systemic velocity of the galaxy ($V_{\text{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_0$), and the inclination ($i$) of the galaxy.
Appendix C: Additional figures

<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /></figure>

**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 σ<sub>LSF</sub> (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).

<figure><img src="image_15.png" /></figure>

**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<sup>-1</sup>) and on the disc (1629±46 km s<sup>-1</sup>).
<figure><img src="image_16.png" /></figure>

**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\alpha$ 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 $\sigma_{LSF}$.

<figure><img src="image_17.png" /></figure>

**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 $\mu(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 $\eta = 0.3$, the shaded area shows the rotation velocity amplitude for $\eta = [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.
<figure><img src="image_18.png" /></figure>

**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 $\Sigma(H\alpha) = 3 \times 10^{-18}$ and $1.25 \times 10^{-16}$ erg s⁻¹ cm⁻² arcsec⁻², respectively, and the red contour $\Sigma(H\alpha) = 10^{-17}$ erg s⁻¹ cm⁻² arcsec⁻². 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^\circ$, 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).

<figure><img src="image_19.png" /></figure>

**Fig. C.6.** Comparison between the $^{12}$CO (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 $\Sigma(H\alpha) = 10^{-17}$ and $1.25 \times 10^{-16}$ erg s⁻¹ cm⁻² arcsec⁻², 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^\circ$, 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).

<figure><img src="image_20.png" /></figure>

**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.
<figure><img src="image_21.png" /></figure>

**Fig. C.8.** GHASP (panel (a)) and simulated $H\alpha$ 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 $\Sigma(H\alpha) \approx 10^{-17} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$. The low level cuts used to display the simulated $H\alpha$ images have been chosen to show the low column density ionised gas reached by the simulations (panel (b)) and at higher surface brightness levels ($\Sigma(H\alpha) \approx 2 \times 10^{-17} \text{ erg s}^{-1} \text{ cm}^{-2} \text{ arcsec}^{-2}$; panel (c)), matching the typical surface brightness limit of VESTIGE and FP observations, respectively.