NGC628 with SITELLE : I. Imaging Spectroscopy of 4285 HII region candidates

TL;DR

This paper introduces a new imaging spectroscopy technique applied to NGC628, identifying and analyzing over 4000 HII regions to study their physical properties, luminosity functions, and ionization conditions.

Contribution

The paper presents a novel method for extracting detailed parameters from HII regions using SITELLE imaging spectroscopy, creating a comprehensive catalog for galaxy NGC628.

Findings

HII region luminosity function slope of 1.12±0.03

Most detections are confirmed as HII regions via BPT diagrams

Complex ionization variations observed within HII regions

Abstract

This is the first paper of a series dedicated to nebular physics and the chemical evolution of nearby galaxies by investigating large samples of HII regions with the CFHT imaging spectrograph SITELLE. We present a technique adapted to imaging spectroscopy to identify and extract parameters from 4285 HII region candidates found in the disc of NGC628. Using both the spatial and spectral capabilities of SITELLE, our technique enables the extraction of the position, dust extinction, velocity, Halpha profile, diffuse ionised gas (DIG) background, luminosity, size, morphological type, and the emission line fluxes for individual spaxels and the integrated spectrum for each region. We have produced a well-sampled HII region luminosity function and studied its variation with galactocentric radius and level of the DIG background. We found a slope of 1.12+/-0.03 with no evidence of a break at high luminosity. Based on the width of the region profile, bright regions are rather compact, while faint regions are seen over a wide range of sizes. The radius function reveals a slope of 1.81+/-0.02. BPT diagrams of the individual spaxels and integrated line ratios confirm that most detections are HII regions. Also, maps of the line ratios show complex variations of the ionisation conditions within HII regions. All this information is compiled in a new catalogue for HII regions. The objective of this database is to provide a complete sample which will be used to study the whole parameter space covered by the physical conditions in active star-forming regions.

Tables

Table 1. Table 1: Observing Parameters

Filter	SN1	SN2	SN3
Spectral range [Å]	3640-3850	4840-5120	6480-6860
Observing date	2016/01/13	2016/01/11	2015/08/09
Exposure time/step [s]	73.0	103.3	25.0
Number of steps	105	138	323
Total exposure time [h]	2.13	3.96	2.24
Mean resolution R	$\sim$ 600	$\sim$ 600	$\sim$ 1800
Image Quality [arcsec]	1.28	0.96	0.90

Table 2. Table 2: Dispersion thresholds for the morphological classification. The dispersion thresholds are given in % of the region peak value.

Category	Centrea	Sigmaa	Outer Limita	Regions
Symmetrical	$\le$ 6 %	$\le$ 6 %	none	448
Asymmetrical	$\le$ 6 %	$>$ 6 %	$\le$ 3 %	772
Transient	$\le$ 6 %	$>$ 6 %	$>$ 3 %	1785
Diffuse	$\ge$ 6 %	$>$ 6 %	none	1280

Table 3. Table 3: Catalogue Parameters

ID	RA	DEC	$R_G$	Log LHα	log DIGSB	Cat	$I_{0profile}$	Aprofile	$\sigma$profile	$\alpha$profile	R2profile	Size	E(B$-$V)	E(B$-$V)err
#	–	–	kpc	–	erg s-1 cm-2	#	erg s-1	erg s-1	pc	%	–	pc	mag	$\pm$
					spaxel-1		cm-2	cm-2
2430	24.1818	15.7876	1.39	37.20	-17.26	1	-15.89	-14.22	39.34	0.72	0.96	46.0	0.408	0.172
2431	24.1375	15.7867	5.48	37.21	-16.85	3	-15.87	-13.92	64.66	0.50	0.95	64.3	0.227	0.048
2432	24.1368	15.7868	4.60	36.49	-17.65	2	-16.92	-14.84	84.41	0.44	0.51	84.2	1.265	0.090
2433	24.2040	15.7882	3.89	37.10	-16.91	3	-16.33	-14.34	77.23	0.60	0.62	97.0	0.347	0.147
2434	24.1891	15.7879	0.83	37.87	-17.06	3	-15.38	-13.61	43.87	0.50	0.96	43.7	0.356	0.053
2435	24.1773	15.7876	5.64	37.70	-16.89	2	-15.49	-13.54	65.68	0.50	0.95	83.0	0.286	0.029
2436	24.2416	15.7892	10.28	37.22	-17.36	1	-15.95	-14.21	40.61	0.50	0.97	48.0	0.308	0.057
2437	24.2121	15.7885	5.82	38.63	-16.72	1	-14.46	-12.79	34.29	0.50	0.99	35.0	0.349	0.026
2438	24.2021	15.7883	4.31	37.35	-17.81	3	-16.02	-14.25	59.17	1.00	0.87	75.0	0.479	0.275
2439	24.1997	15.7881	4.06	37.06	-17.34	3	-16.40	-14.21	114.01	0.48	0.30	140.0	0.000	0.091
2440	24.1990	15.7881	2.52	37.04	-17.45	3	-16.09	-14.39	44.32	0.77	0.84	54.00	0.261	0.188
2441	24.1959	15.7881	1.68	37.28	-17.28	3	-16.23	-14.08	102.78	0.50	0.44	127.00	0.067	0.168
2442	24.1840	15.7878	0.67	–	–	4	-16.61	-14.70	60.00	0.50	0.00	76.00	0.607	0.994
2443	24.1736	15.7876	0.79	37.45	-16.89	3	-15.64	-13.82	48.25	0.50	0.92	60.00	0.448	0.080
2444	24.1702	15.7876	0.91	37.06	-17.26	4	-16.51	-14.48	79.00	0.50	0.48	99.00	0.000	0.212
2445	24.1630	15.7874	1.68	36.67	-17.51	4	-16.47	-14.48	72.77	0.50	0.86	72.45	1.154	0.708
…	…	…	…	…	…	…	…	…	…	…	…	…	…	…

log ([NII]/H$\alpha$)	log ([NII]/H$\alpha$)err	log ([NII]/H$\alpha$) SNRcross	log ([SII]/H$\alpha$)	log ([SII]/H$\alpha$)err	log ([SII]/H$\alpha$) SNRcross	…∗
–	$\pm$	–	–	$\pm$	–	…
-0.450	0.069	20.919	-0.657	0.114	20.919	…
-0.629	0.023	56.179	-0.615	0.042	56.179	…
-0.898	0.162	8.021	-0.413	0.329	8.021	…
-0.487	0.048	29.180	-0.452	0.086	29.180	…
-0.597	0.024	54.618	-0.768	0.042	54.618	…
-0.618	0.014	88.427	-0.627	0.026	88.427	…
-0.769	0.023	52.319	-0.504	0.047	52.319	…
-0.586	0.017	77.199	-0.696	0.030	77.199	…
-0.592	0.069	19.538	-0.443	0.130	19.538	…
-0.493	0.056	25.023	-0.396	0.104	25.023	…
-0.698	0.067	19.171	-0.444	0.133	19.171	…
-0.610	0.061	21.629	-0.505	0.117	21.629	…
-0.330	0.302	5.915	0.034	0.631	6.344	…
-0.461	0.035	40.596	-0.801	0.056	40.596	…
-0.109	0.182	10.803	-0.213	0.274	10.803	…
-0.243	0.178	9.849	-0.201	0.400	9.849	…
…	…	…	…	…	…
*All line ratios listed in Section  6.2 are included in the catalogue.
Note: More parameters can be provided upon request.

Equations

$$E(B-V)=\frac{2.5}{1.07}\log \left(\frac{F_{H\alpha }/F_{H\beta obs}}{2.87}\right),$$

$$A_V=-0.083R_G+\mathrm{1.16.}$$

$$f(x,A,x_0,\sigma ,\alpha )=\frac{1-\alpha }{\sqrt{2\pi }{\sigma }_g}A{e}^{\frac{-(x-x_0{)}^2}{2{\sigma }_g^2}}+\frac{\alpha A}{\pi }\frac{\sigma }{(x-x_0{)}^2+{\sigma }^2},$$

$$R_S^3=\frac{3Q_H}{4\pi ϵ{n_e}^2{\alpha }_B(H,T_e)},$$

$$\begin{array}{rl}\log (Q_H)=& 5.0754+42.034x\\ & -15.0797x^2+2.4439x^3-0.1474x^4,\end{array}$$

$$\begin{array}{rl}\log ([\mathrm{NII}]/H\alpha )=& -0.6841-0.651x-0.9026x^2\\ & -0.4249x^3-0.06401x^4,\end{array}$$

Full text

galaxies: HII regions, galaxies: diffused ionized gas, galaxies: nebular abundances, catalogues: HII regions, instrumentation: IFTS, IFTS : SITELLE

††pagerange: NGC628 with SITELLE : I. Imaging Spectroscopy of 4285 HII region candidates.–51††pubyear: 2018

1 INTRODUCTION

The study of the ionized gas components in galaxies, for instance the population of HII regions as well as the diffuse ionized gas (DIG, also referred as the diffuse warm ionized medium WIM), helps gather important clues on the mechanisms that drives star formation and galaxy evolution. The gas emission lines are sensitive to the physical conditions in the interstellar medium (ISM) and their relative intensities probe underlying mechanisms involved in ionisation processes. One predominant example related to the ISM ionisation is the recent star formation sites where massive stars act as the main source of ionizing photons. Gas ionisation can also be triggered by other sources such as the energetic emission released by an active galactic nucleus (e.g. Ho 2008; Davies et al. 2014), the global evolved stars contribution to the heating of the ISM (e.g. Flores-Fajardo et al. 2011), different sources of small- and large-scale shocks in the ISM (e.g. Allen et al. 2008), etc.

The gas strong emission lines can be observed at cosmological distances, allowing the study of various star-forming environments through time. Establishing reliable HII region and DIG physical characteristics are therefore key elements in probing galaxy’s chemo-dynamical evolution and improving star formation models. The extraction of these physical characteristics often involves photoionisation codes, like MAPPINGS (Sutherland & Dopita 1993; Nicholls et al. 2013) and CLOUDY (Ferland et al. 1998). These codes are used to model the complex interaction between the ionizing spectrum of different sources and the ISM. To constrain photoionisation models, one needs spatially resolved star-forming region spectra including a set of different emission lines.

From a historical perspective, imagers and Fabry-Perot instruments were often used to observe simultaneously a large number of HII regions in nearby galaxies. Imagers combined with different narrow-band filters allow us to target a few important emission lines and continuum windows (used to isolate the emission line flux from the presence of the underlying stellar and nebular continuum; Mollá et al. 2009, Sánchez et al. 2011). Nevertheless, imaging techniques are limited by the number of emission lines accessible and are also affected by problems related to the spectral shift of the lines due to nebular and galactic dynamics, which prevent the accurate measurement of the lines (as mentioned byYoung et al. 1996, Hoopes et al. 2001, and Kennicutt et al. 2008). High-resolution spectra produced from Fabry-Perot interferometers can provide additional information on stellar and gas dynamics. However, as their spectral range is limited, Fabry-Perot interferometers were mostly used for the observation of one emission or absorption line (e.g. Fathi et al. 2007 and Epinat et al. 2008). For this reason, most of the complex spectral analysis of extragalactic HII regions have been performed using slit spectrograph data, as they enable accurate measurements of multiple emission lines over a wide spectral window (e.g. Berg et al. 2015). More recently, integral field spectrographs (IFS) using dispersive techniques offered a breakthrough by giving simultaneously spectral and spatial information on extragalactic HII regions (e.g. Kreckel et al. 2016). Still, the incomplete spatial coverage of these instruments is not convenient for a detailed study of extragalactic HII regions and to extract some of their intrinsic properties (e.g. total luminosity and morphology).

Imaging Fourier Transform Spectrograph (IFTS) have advantages over these other intruments. The multiplex advantage has to be mentioned first as multiple scans, or spectra, can be performed at a faster rate than most dispersive spectrographs. Secondly, the throughput of an IFTS is significantly higher than slit spectrograph while nothing limits the amount of light getting in the cameras and fewer mirrors are required. Finally, the sampling of the scans can be adjusted to reach a desired spectral resolution. Nevertheless, for an IFST working in the visible, constraining the bandpass using narrow- to wide-band filters is necessary to reach an appropriate spectral resolution in a reasonable amount of time (a few hours) and limit the photon noise. All things considered, the large field-of-view resulting from the multiplex advantage, the high throughput, and the complete spatial coverage make the IFST an ideal instrument to study extended object, and more particularly here, large nearby galaxies.

Some fundamental problems have been outlined from the analysis of extragalactic HII and DIG regions in nearby galaxies with integral field spectroscopy using disperse units:

The lack of spatial resolution is preventing us from resolving individual HII regions and is imposing a great limitation on the scale of the physical processes studied. Even when analyzing the large scale properties of galaxies, it has been shown that the gas abundance [O/H] gradient (derived from emission line diagnostics) is affected by a low spatial resolution, resulting in an artificial flattening of the gradients with respect to the redshift (Yuan et al. 2013; Mast et al. 2014). In addition, the lack of spatial resolution induces an attenuation on the extinction correction factor obtained from the Balmer decrement in integrated spectra as the extinction is biased toward the brightest and less dusty regions (Vale Asari et al. 2014). As a result, the derived star formation rate (SFR) and stellar mass content of HII regions are only upper limits.

• 2)

A variable contribution from the superimposed DIG emission is often pointed out as a source of uncertainties when trying to extract the HII region properties (Blanc et al. 2009; Kreckel et al. 2016; Sanders et al. 2017).

• 3)

A degeneracy effect between the age of the HII regions and the under-sampling of the massive end of the initial mass function (IMF) precludes proper modeling of the ionized gas (Villaverde et al. 2010). Without resolving the stellar content, additional information on the ionisation structure surrounding the stars is required using photoionisation models to determine the precise stellar content of the regions.

• 4)

Most of the diagnostic tools proposed to gather the abundances and other ionized gas properties are designed for integrated spectra (i.e. at low spatial resolution) and are calibrated using a small, sometimes biased, sample of HII regions with well-defined properties. They are therefore limited tools, not representative of a large variety of environments within galaxies, and are not suited for a desired detailed analysis of resolved HII regions.

To overcome these numerous issues, we have conducted new observations of the nearby spiral galaxy NGC 628 with SITELLE, the Canada-France-Hawaii IFTS. SITELLE allows a high-spatial resolution (35 pc) over the whole galaxy disc (FOV: 11′$\times$11′) with full coverage (100% filling factor). NGC 628 is a well-known spiral galaxy seen almost face-on. It has been widely observed by others with various instruments. For instance, Rosales-Ortega et al. (2011) and Sánchez et al. (2011) performed a spectral analysis and studied the global properties (abundance gradients) of NGC 628 using a mosaic of the PPAK IFS ($\sim$6′ circular aperture) with individual fiber spectra and 94 selected apertures. Berg et al. (2015) used the Multi Object Double Spectrograph (MODS), as part of the CHAOS project, to target 62 bright HII regions over the galaxy disc and investigated the abundance dispersion using direct temperature measurements from the auroral line. Grasha et al. (2015) used multiple Hubble Space Telescope broad band UV and optical images over two 2.7′$\times$2.7′ fields in the disc of NGC 628 to identify 1392 young stellar cluster candidates ($<$ 100 Myr) and studied their size, mass, and age, and the concentration of the stars embedded in the clusters. More recently, with the Multi-Unit Spectroscopic Explorer (MUSE), Kreckel et al. (2016) sampled 391 HII regions at the spatial resolution of 35 pc, analyzed the DIG contribution, and studied the differences between the HII region properties in the arms and inter-arm regions. These previous studies have targeted different science goals and their contribution to the understanding of star formation in galaxies is notable. Complementary to these studies, we believe that SITELLE has the potential to overcome several of the problems related to sample selection biases, the lack of spatial resolution, the contribution of the DIG, as well as to explore line ratio variations in the ISM with greater details.

The spectral datacubes obtained for this study with SITELLE triggered the development of new analysis tools tailored to its unique capabilities. In this paper, we present our approach for the systematic analysis of the emission line spectra of a large sample of regions observed with a high-spatial resolution. Our results are meticulously compiled in a general catalogue for the HII regions which is made available to the entire astronomical community. Our immediate goal is to demonstrate the impact of the spatial resolution, mainly through the evaluation of the DIG background and the effect of sampling on characterizing the HII regions.

This paper is the first of a series where the systematic study of many nearby galaxies will be done with SITELLE to enlarge the domain of HII region physical parameters analysed. All the data gathered for the HII regions will be added to the catalogue. It aims to be the largest and most detailed spectroscopic study of star-forming region properties in different galactic environments, and represents an important statistical tool to study: small-scale processes in chemical enrichment and mixing; the balance between ionizing photons absorbed by the gas, the UV flux observed from a young population, and the overall extinction; the DIG in relation with HII regions; the star formation efficiency and rate directly from the gas and stellar population characteristics; and the massive stars IMF in different environments. This catalogue is tailored for the development of new diagnostics with photoionisation codes and for the study of HII regions and their content at different spatial resolution.

Sections 2 and 3 describe the observations and the technique used to properly measure the gas emission lines. Section 4 focuses on the emission region definition, while Section 5 presents the luminosity function, morphology, and size of the regions. Section 6 presents the emission line ratios measured in individual spaxel within the emission regions, as well as global ratios for each region for a comparison with the literature. The results, along with the database content, are summarized in Section 7.

2 OBSERVATIONS AND DATA REDUCTION

NGC 628 was observed with SITELLE (Spectro-Imageur à Transformée de Fourier pour l’Étude en Long et en Large des raies d’Émission) at the Canada-France-Hawaii Telescope during the instrument commissioning and science verification runs in August 2015 and January 2016. Table 1 gives the main observing parameters. SITELLE is a Michelson interferometer inserted, with an inclination, into the collimated beam of an astronomical camera system (Drissen et al. 2014). Its has a large field-of-view (FOV: 11′$\times$11′) with complete spatial coverage, a high-resolution power up to R $\simeq$10 000, and a broad wavelength range from 3500 to 9000 Å  with an excellent efficiency all the way to the blue part of the spectrum. SITELLE is equipped with two E2V detectors (2048$\times$2064 pixels, resulting in a mean pixel scale of 0.321*′′$\times$0.321′′*). Filters are necessary with SITELLE in order to reduce the noise in a selected bandpass. A raw datacube from SITELLE is composed of a sequence of images (one for each position of the movable mirror) that, when taken together, results in more than 4 million interferograms (one interferogram for each pixel of the detectors). After the data reduction, SITELLE produces 4 million spectra in a single datacube for a selected filter.

For this project, three filters have been used : SN1, SN2, and SN3; Table 1 gives the spectral range and the resolution adopted for each filter. These configurations enabled the measurement of multiple strong emission lines: [OII]$\lambda$3727, H$\beta$ (4861 Å), [OIII]$\lambda$4959, [OIII]$\lambda$5007, [NII]$\lambda$6548, H$\alpha$ (6563 Å), [NII]$\lambda$6583, HeI$\lambda$6678, [SII]$\lambda$6716, and [SII]$\lambda$6731. Figure  1 shows the integrated spectrum in all three filters for one HII region. Figure  2 shows the deep image produced by adding, for each pixel, the whole signal from the three filters together with an enhanced contribution of the H$\alpha$ intensity map extracted from the line fitting procedure (as described in $\S$  3). SITELLE’s FOV was slightly offset from the galactic centre of NGC 628 in order to cover a vast portion of the disc as well as the extended northern spiral arm; the FOV was centered at RA 01h36m46.2s and DEC +15$\degree$48′42.8*′′*.

The data reduction was performed using the fully-parallelized reduction software ORBS (data-release 1) specifically developed for SITELLE (Martin 2015). The data reduction process requires additional calibrations from the instrument. Typical flat (using a white light source) were acquired during each night and a pre-scan section of the CCD was preserved to measure the bias contribution in each exposure. The flux calibration was performed using a datacube of the standard star GD 71 secured with each filter to define the global transmission function (taking into account the minimum airmass value for the standard star and the galaxy, which is 1.002 for both objects). The flux calibrated spectrum of NGC 628 is in good agreement with those from other studies (as shown in Martin & Drissen 2017). The spectral calibration was done using a Helium-Neon laser datacube acquired during daytime. For the SN3 filter (at higher R), greater precision on the spectral calibration was obtained using the centroid position of sky lines ($\S$  3.1). Additionally, three flat datacubes were acquired using a white light source (one for each filter) in order to apply high-order phase correction in the data reduction process (Martin et al. 2018). To match the image quality obtained for the different datacubes (see Tab. 1) and attenuate the impact of a small misalignment of the datacubes, we applied a Gaussian convolution (with a FWHM $=$ 3 pixels) on all the images within the datacubes. We did not correct for the Milky Way extinction knowing that it is very small (E(B$-$V) = 0.062 according to NED111The NASA/IPAC Extragalactic Database (NED) is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration.).

3 MEASUREMENTS

A few manipulations are required before proceeding with the accurate measurements of the emission lines. These steps include the refinement of the spectral calibration for the SN3 datacube, sky subtraction, corrections for the astrometry and image distortion, and subtraction of the galaxy stellar contributions along the line of sight of each spaxel. Fitting of the emission lines is then performed. These operations are explained in the following subsections.

3.1 Spectral Calibration Refinement

The spectral resolution of the SN3 filter was pushed further on, compared to the other filters (see Tab. 1), to allow us to fully separate the [SII] doublet, but also to enable the use of the strong H$\alpha$ emission line to extract the ionized gas kinematics. At this resolution, higher precision for the spectral calibration of the SN3 datacube can be achieved by fitting, over the FOV, an isolated and strong sky line available in the filter (a similar technique than the one used by Martin et al. 2016 and Shara et al. 2016; this technique is now a standard procedure and is well documented in Martin et al. 2018). Since the Helium-Neon laser datacube used for the spectral calibration was obtained with the telescope at zenith, some telescope and instrument flexures affecting the calibration may occur while pointing at the target. These flexures are then corrected by using sky-line centroid positions observed simultaneously. For this task, we selected the OH sky line at 6498.729 Å  because of its proximity to the H$\alpha$ line while it remains uncontaminated by emission or absorption lines from other sources. The sky line signal is enhanced by summing its flux over boxes of 10$\times$10 spaxels. To get the spectral calibration correction, the sky line was fitted (see the description of line measurement technique in $\S$  3.5) in order to obtain its centroid shift with respect to the expected rest frame wavelength. A constant shift of 83.2 km s*-1* was determined over the FOV and considered for the correction of the velocity map (as expected according to Martin et al. 2016, 2018).

3.2 Sky Subtraction

For each datacube, sky subtraction was performed using the median sky spectrum extracted from a region located as far as possible from the galaxy disc. A total of 40 000 spaxels (200$\times$200 spaxels) were used centred at the position RA 01h37m04.2s and DEC +15°49′23.3*′′*. The sky spectrum in the SN3 datacube is shown in Figure  3.

3.3 Astrometry and Distortion Corrections

The astrometric solution for each datacube was calculated using bright stars in the FOV of the deep image and applied using ccmap and ccsetwcs (STSDAS/IRAF222IRAF 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.). We selected a TNX astrometric solution with a polynomial degree of 4. The accuracy of the solutions is below 0.324*′′* ($\sim$ 1 pixel) for most regions (all except a small group of $\sim$ 20 regions located in the SE corner of the FOV) detected in Section  4 for all filters. Ultimately distortions were corrected using remap (WCStools) and all line flux maps (after the line fitting procedure described in $\S$  3.5) were adjusted to correspond to the astrometric solution of the SN3 deep image. This step is crucial to measure reliable line ratios extracted from lines lying in different filters.

3.4 Subtraction of the Stellar Population Contribution

Stellar populations in galaxy discs are responsible for continuum emission and absorption features that add up to the spectrum of the ionized gas regions. In the visible domain, emission lines (and most importantly Balmer lines) are then a combination of a stellar absorption component (from the older galaxy stellar populations and, if the case, from the young stellar population embedded within an HII region) and the emission line from the ionized gas itself. Depending on the galaxy mass and type, and also on where an HII region is located (e.g. far out in the disc or near the galaxy centre), the measurement of its emission lines can be significantly affected by the disc and bulge older stellar populations (Rousseau-Nepton 2017). Although Gonzalez Delgado et al. (2005) indicate, for example, that the equivalent width of the H$\beta$ absorption line can reach 5 Å for a young population below 10 Myr, the young cluster absorption features may also be important after its superposition with the old population, depending on the spatial resolution and morphology of the young cluster relative to the emission region.

We first created a reference spectrum for the old disc/bulge stellar population for each datacube considering a small circular aperture of 5.5*′′* centered on the galaxy continuum emission peak (emission line regions are not detected in this aperture). The choice of this central aperture (compared to specific regions in the galaxy disc) allows us to minimise the noise in the reference spectrum. This reference spectrum also offers the advantage of containing the exact information about the instrumental lines profile (i.e. a sinc profile) as compared to a stellar population model. A comparison of spectra created at different galactocentric radii does not show a significant variation in the absorption lines strength and width when the noise is taken into account. This is consistent with the study of the stellar population characteristics and the velocity dispersion done by Rousseau-Nepton (2017) using long-slit spectroscopy data (between 3500 and 7500 Å) across the disc and nucleus of NGC 628. It supports the idea that our reference spectrum is a good representation of the old galaxy stellar contribution all over the disc, considering an adjustment of its velocity, using a preliminary interpolated velocity map, and an adjustment of its continuum level. For this last step we decided to simply match the continuum level of the old population with the continuum level of each spaxel studied, assuming that it will be a sufficient correction for the young population contribution, if present. For NGC 628, SITELLE gives us a spatial resolution of 35 pc, which allows us to see, with limited details, the spatial distribution of the ionized gas relative to the young stellar clusters: a young cluster (considering variation in the continuum level) may sometimes be seen sitting in the centre of an emission region, or it is not seen at all, or it is sitting at the edge of it. Not to mention as well that a young cluster and its related emission region may not share the same extinction; this is a complex problem that we will address in another paper. Nevertheless, here we ran some tests to see what would be the error on the emission line flux if we had specifically considered young population spectra. To do so, we combined the average old population spectrum with young (6 and 10 Myr) population models from Starburst99 (Leitherer et al. 1999) in various proportions. For example, for HII regions at a large galactocentric radius (we considered regions at 6 and 14 kpc), if we suppose that the young population contributes to 100% of the continuum, relative to the case where 100% of the continuum is associated with the old population, we see a maximum variation of 5% and 1% in the H$\beta$ and H$\alpha$ flux, respectively, if the region has weak emission lines and presents a high continuum level. We find a similar variation of the line flux for HII regions at a small galactic radius if we consider that the young population accounts for 50% of the continuum, relative to the case where 100% of the continuum is associated with the old population. These are extreme examples for the integrated flux from HII regions, considering that while the old population is present in all the spaxels where the emission lines are measured, it is probably not the case for the young population clusters. But it is more realistic when studying individual spaxels. Nevertheless, the flux variations measured in our tests are small and we are therefore confident that our flux measurements are good when considering the old galaxy stellar population with its continuum level scaled with the one of the individual spaxels or the integrated regions.

Figures 4 and 5 show the effect of the subtraction of the stellar population reference spectrum for the H$\beta$ (SN2) and H$\alpha$ (SN3) emission line of one spaxel selected in the centre of a faint HII region, and for the H$\alpha$ (SN3) emission line in a DIG region. The stellar population reference spectrum is superimposed to the spectrum of these regions in the figures. The final spectrum obtained after the subtraction is also shown. It is important to realize that the oscillations around the emission lines are mostly the signature of the sinc instrumental profile. These oscillations naturally display large negative side-lobes on each side of the line centre, which must not be confused with the shape of an absorption component. These figures clearly show that it is important to subtract the stellar population prior to measure the emission lines. For weak emission regions, an emission component in the lines, especially for H$\beta$, may suddenly show up after the subtraction of the stellar population spectrum.

3.5 Line Fitting Procedure and the Final Velocity Map

For each spaxel in all three datacubes, the emission lines are fitted simultaneously using the extraction software ORCS (Martin et al. 2015333https://sourceforge.net/projects/orb-orcs/). This software was specifically created for SITELLE and therefore uses sinc-shaped line profiles. ORCS output map for each line includes the maximum intensity $I_{0}$, FWHM, velocity, continuum level, all the corresponding uncertainties, the SNR, as well as the standard deviation $\chi^{2}$ of the fit. These parameters have been used to get the flux and its uncertainty map for all lines (the latest version of ORCS automatically gives all these quantities). Line flux uncertainties were here calculated using the relation for a non-resolved line (i.e. below the Nyquist criteria) in the case of an FTS instrument according to Lenz & Ayres (1992); they consist of the product of the line flux with its ratio $\chi^{2}$/$I_{0}$. The uncertainty for line flux ratios were then put equal to half of the difference between the maximum and minimum ratio values obtained considering the line flux uncertainties.

Examples of the line fits are shown in Figure  6 for all three filters. The fluctuation seen in the residual underneath H$\alpha$ is caused by the subtraction of the galaxy stellar contribution which, as we discussed in the previous section, may not be perfectly well adapted (in content and velocity dispersion). But, as seen here, the impact on the line fitting is negligible ($<$ 1$\,\%$). Appendix A contains the resulting flux maps for all the lines (Fig. 27 and 28).

Appendix A also presents the final velocity map (Fig. 29) obtained by fitting simultaneously all the SN3 emission lines with ORCS. As explained in the appendix, the velocity precision achieved here is often better than 30 km s*-1*, i.e. the limit when SNR${}_{{\rm H}\alpha}=3$.

4 REGION DEFINITION

Although one can study individual spaxel within an ionized gas region in a nearby galaxy like NGC 628 with SITELLE, some properties must be investigated using their integrated flux (e.g. the region luminosity function) and a description of their morphology (e.g. their size). It is therefore necessary to distinguish which pixels belong to a given region. A new code is presented here for identifying the position of the ionizing source associated to a region and its boundaries. Compared to other existing codes, like those based on the H$\alpha$-continuum-subtracted images and a standard-fixed-threshold-photometry method, such as the percent-of-peak-photometry method (McCall et al. 1996), or HIIphot (Thilker et al. 2000, 2002), this new code gives us flexibility in defining regions and their DIG background.

Our method proceeds with multiple steps : 1)  the identification of the emission peaks, 2)  the determination of the zone of influence around each emission peak, and 3)  the definition of the outer limit of a region and its DIG background. These operations are done prior to the extinction correction, to avoid introducing unnecessary noise. The different steps of our method are described in the following subsections.

4.1 Definition of the Peaks

The first step involves adopting a definition for the emission peaks and identifying their positions. We use a combination of the H$\alpha$, H$\beta$, and [OIII]$\lambda$4959,5007 continuum subtracted images. This choice is motivated by the increased SNR of the combined image. Moreover, as found from our high-spatial resolution image, the maximal intensity of these three lines coincides with the centre of the emission regions. An emission peak is identified for a pixel when: 1)  the pixel intensity is greater than the intensity of at least 5 immediate surrounding pixels (i.e. all these pixels must be in contact with the maximum intensity pixel), and 2)  the total intensity in a 3$\times$3 pixels box centred on the emission peak is above the adopted detection threshold. This threshold is fixed by the 3$\sigma$ noise level of the H$\alpha$+H$\beta$+[OIII] flux map (a slightly different value is adopted for each of the CCD quadrant according to their individual noise level). Furthermore, if two emission peaks are separated by a distance smaller than the image quality of the observation, only the brightest peak is preserved. For each so-defined emission peak, two values of its mean flux are evaluated considering a circular aperture with a radius of 1 and 2 pixels centered on the emission peak. The peak is rejected when the mean flux does not significantly decrease between the small and large aperture (i.e when a drop less than $\sim$ 2$\%$ between adjacent pixels is seen for the bright regions, with a flux FHα $\simeq$1$\times$10*-15* erg s*-1* cm*-2*, and when a drop less than $\sim$ 5$\%$ is seen for the faint regions). This criterion allows us to remove false identifications inside a bright ionized gas region, that are likely to be associated with small surface brightness variations, while eliminating false detections due to random noise fluctuations. To find the accurate position of the centroid, a 2D Gaussian model is fitted on the remaining emission peaks using the IDL procedure MPFIT2DPEAK. A total of 4285 emission peaks is finally detected for NGC 628 with our extraction technique. After a visual verification of a representative group of emission peaks, we estimate that only a small fraction, $\sim$ 0.5$\%$, of these detections are false detections (not an ionized gas region). It is also important to note that most of the remaining regions are HII regions while a few percentage may be pure DIG regions, supernova remnants, and planetary nebulae, as it will be demonstrated from their location within the BPT diagrams in Section 7.

4.2 Zone of Influence

Figure  7 shows the emission peaks and their zone of influence for a section of the northern arm and for four reference regions selected in the galaxy disc. We name the area that is influenced by the ionizing photons from an emission peak, the zone of influence. It is defined by studying the distance between each pixel and all its surrounding emission peaks. In order to consider all the possible sizes for the zone of influence, we looked for emission peaks within a maximum distance of 425 pc around each pixels (corresponding to 30 pixels). This distance is greater than the maximal outer limit distance found for the regions of NGC 628, as seen in the next paragraph. Most of the time a pixel is then associated with the closest emission peak. In the case of multiple emission peaks at the same distance (considering a distance uncertainty given by the image quality), the pixel is then associated with the brightest emission peak. The most distant pixels associated with an emission peak define the boundary of the zone of influence for this ionizing source.

Figure  8 shows an example, for a typical region, of the H$\alpha$ flux of all the individual pixels included inside the zone of influence as a function of their distance from the emission peak centroid. Although the average H$\alpha$ flux profile for most regions presents a smooth drop, the dispersion around the mean profile can, in some cases, be large (see § 5.2). This dispersion is explained by the asymmetry of the emission region, inhomogeneous DIG emission, and also crowding and density effects (in these last cases, pixels are not associated with the right emission peak).

4.3 Outer Limit of a Region and the DIG Background

Within the zone of influence of an ionizing source, the intensity profile changes rapidly (as shown in Fig.  8). Figure  9 shows the total flux profile for the region presented in Figure 8. The total flux profile is calculated by summing the flux of the pixels within a circular annulus 25 pc thick centered on the emission peak. This fixed 25 pc annulus, corresponding to $\sim$2 pixels, allows us to properly sample the profile (e.g. its slope) with a sufficient amount of pixels. The outer limit of the region is then defined at a distance where the slope of the total flux profile decreases by less than 2 % within the zone of influence. The DIG background surface brightness (DIGSB) is given by the median intensity value for all the pixels included in an annulus 50 pc thick centered on the position of the outer limit. The zone of influence and the DIG background are indicated for the example in Figure 8. The 2 % threshold was selected after a few iterations: a larger value overestimates the DIG background for numerous HII regions, whereas a smaller value is not sensitive to variations of the slope. With these constrains, an outer limits is not found for only three regions; in these cases, the outer limit is simply set at the radius of the zone of influence. For our sample of regions, the average outer limit is 71 pc (with a standard deviation of 33 pc) and the smallest value is 21 pc. In Figure  9, the effect of the DIG background subtraction on the total intensity profile is also shown.

5 MORPHOLOGY, SIZE, AND LUMINOSITY FUNCTION

5.1 Extinction Correction

To evaluate the colour excess ${E(B-V)}$, we used the following equation:

[TABLE]

with ${A_{\lambda}=2.5log{I_{\lambda}}/{I_{\lambda 0}}}$, and $E(\beta-\alpha)$/$E(B-V)$ = 1.07 (extracted from the extinction curve of Cardelli et al. 1989), and where $F_{H_{\alpha}}/F_{H_{\beta}\,obs}$ is the observed flux ratio of the two main Balmer lines. We adopted a theoretical $F_{H_{\alpha}}/F_{H_{\beta}}$ ratio of 2.87 (as expected for the Case B of Osterbrock 1989, at a temperature of 10 000 K). As discussed at the end of this subsection, the value selected here are for an average HII region.

Using the previous equation and the relation ${A_{V}=3.1\,E(B-V)}$, we obtained a pixel-by-pixel SITELLE extinction map for NGC 628. This map is shown in Figure  10, where only the spaxels with a 3$\sigma$ detection for both lines have been considered. Figure  11 presents the extinction as a function of the galactocentric radius ($R_{G}$) for all the spaxels (the observed radii were deprojected using a galaxy inclination of 21$\degree$ from NED). The extinction is clearly decreasing as we look further out in the disc. Fitting a linear function through the data between the center and a galactocentric radius of 15 kpc, we obtained the following relation:

[TABLE]

The highest density points in Figure  11 (in green) show an extinction value in good agreement with Berg et al. (2015), where $A_{V}\simeq 0$ to 1.28 mag for a sample of 45 HII regions. Our mean value, considering only the inner part of the galaxy ($R_{G}\leq 9$ kpc), is 0.86$\pm$0.81 mag. This value is very close to the mean value of 0.79 obtained by Cedrés et al. (2012) derived from a sample of 209 HII regions observed with narrow-band imaging. Within the uncertainties, it is also in agreement with the result of Sánchez et al. (2011), 1.24$\pm$ 0.76 mag, gathered using an IFS over the same area of the galaxy disc.

A close examination of the extinction map as displayed in Figure  10 suggests significant fluctuations in the extinction within individual regions. Furthermore, $\sim$ 10$\%$ of the spaxels have an observed ratio $F_{H\alpha}$/$F_{H\beta}$ below the theoretical value (this is seen in Fig.  11 for AV $<$ 0). Some tests were done to investigate the cause of these variations and the low flux ratios. First, the subtraction of the stellar population reference spectrum was considered. From our tests, using a combination of young burst (6 to 10 Myr) with the old stellar population in various proportions of the continuum level (§ 3.4), we get a maximum effect of $-$0.2 mag for AV, which is of interest to explain some of the variations seen in the extinction map of individual regions (Fig. 10), but below the dispersion seen in Figure 11. Next, by combining many spaxels within a region, we used various aperture sizes to see the impact of atmospheric refraction on the measured flux ratio. The accuracy of the flux calibration and the alignment of the datacubes were also tested. Although these effects can increase the flux uncertainty, we find that they cannot be held responsible for all the spaxels with a Balmer line flux ratio below the theoretical value. A simple explanation often proposed for low Balmer line flux ratios invokes the low SNR of the H$\beta$ line. But this is not sufficient here as we estimate that only $\sim$ 4 $\%$ of the spaxels with a ratio $F_{H\alpha}$/$F_{H\beta}$ below the theoretical value have SNRHβ $<$ 10 (including the effect associated to atmospheric refraction, flux calibration, and alignment of the datacubes). Therefore, some of the low Balmer line flux ratios, along with the extinction variation seen inside a region itself, could be real and possibly related to the change of the physical parameters (e.g. the electron temperature and density) within that region. For an HII region, the theoretical Balmer line flux ratio can reach 2.75 at 20 000 K according to Osterbrock (1989). Depending on the morphology of the HII regions, the dust backscattering (O’Dell et al. 2017) can also contribute to the observed variation of the flux ratio. In regions highly dominated by the DIG, the electron temperature and density would be very different than in an HII region (Hafner et al. 1999; Elwert & Dettmar 2005). SITELLE’s data is therefore revealing the detailed physical conditions within emission regions and their close surroundings.

Using the extinction law of Cardelli et al. (1989) and considering that our emission regions are mainly HII regions, we corrected all our flux maps with the extinction calculated fror each spaxel or with the relation given by Equation 2 for spaxels where no precise extinction was measured (for AV $<$ 0 and with SNRHα or SNRHβ $<$ 3). As the extinction correction is subject to uncertainties, line ratios (as presented in $\S$ 6) built using lines that are far away from each other (i.e. lines from different filters) also have greater uncertainties.

5.2 Morphology of the Emission Regions

Constraining the morphology and the size of the emission regions requires modeling their intensity profiles in more details. The morphological classification of HII regions is often performed by eye, but with a large sample of candidates as obtained here, an automated approach was considered. At the spatial resolution of our survey ($\sim$ 0.8*′′* or 35 pc), we were mostly interested in obtaining the general shape of the regions and in separating compact from extended regions.

To establish the general shape of the emission regions, we divided them into four categories: symmetrical, asymmetrical, transient, and diffuse. Figure  12 presents an example for each category. Before explaining our criteria for the classification, we can see that the symmetrical region looks almost spherical on the H$\alpha$ image and presents very low dispersion in its luminosity profile, as shown in Figure  13. The ionizing source (i.e. possibly a young stellar population here) of this region is likely to be concentrated in the centre with the surrounding gas equally distributed. The asymmetrical region shows an important emission peak but also a large variety of shapes diverging from the spherical case. In this category, the ionizing source is expected to be either more dispersed or to encompass multiple components, and/or the surrounding gas is not equally distributed in all directions, and/or emitting pixels have been associated with the wrong peak due to a crowding effect. The diffuse region is less peaked and its intensity profile presents even more dispersion; a fraction of the regions in this category may be more easily related to DIG regions. The dispersion seen for transient regions impedes our ability to distinguish between an asymmetrical region and a region with diffuse structures such as large filaments of ionized gas.

The HII region morphologies are fully quantified using the intensity profile. We selected a pseudo-Voigt function to fit the profiles (a linear combination of a Gaussian and a Lorentzian profile) using the Nelder-Mead method (Nelder & Mead 1965). The pseudo-Voigt profile is defined by the equation:

[TABLE]

where the constant $A$ is the amplitude, $\sigma_{g}=\sigma/\sqrt{2\,\mathrm{ln}2}$, and 2$\sigma$ = FWHM. The constant $\alpha$ is a fraction revealing the relative importance of the Gaussian and the Lorentzian terms to the profile. The centroid position ($x_{0}$) was fixed to 0 as it is already constrained by our peak identification procedure. The use of a pseudo-Voigt function enables us to reproduce the observed H$\alpha$ flux profile very well in the case of a peaked emission with a spherical morphology (the pseudo-Voigt profile is shown in Fig.  13 for the example selected). The fraction $\alpha$ varies significantly from one region to another, and allows us to recover the outer part of the profile more accurately than by using a simple Gaussian profile. A correlation coefficient (R2) for the fit was also evaluated for each region by considering the whole profile from the centre to the outer boundary. Only 13$\%$ of the regions have a correlation coefficient R2 $<$ 0.5 and 7$\%$ have a correlation coefficient R2 $<$ 0.3.

For each emission region, we used the pseudo-Voigt H$\alpha$ intensity profile to evaluate the average pixel dispersion as a function of the distance over three 15 pc annuli: in the very centre, at the position of 1$\sigma$, and at the outer limit. We also evaluated the average pixel dispersion within these same annuli by using the mean value instead of the fitted profile. This second method was required for specific regions where the pseudo-Voigt profile could not properly fit the peak or the background intensity. The smallest dispersion values, evaluated using either the mean or the pseudo-Voigt fitted profile, were considered to assign a morphological category to an emission region. The dispersion thresholds for all four categories were first selected arbitrarily (by eye) to finally follow the rules listed in Table 2.

Figure  14 shows, for each morphological category individually, the integrated H$\alpha$ luminosity of the emission regions as a function of the $\sigma$ value of their pseudo-Voigt fitted profile. The symmetrical morphology category contains on average regions that are bright and more compact (with an outer limit between 21 and 173 pc), while the diffuse category contains regions that are among the faintest ones and among the most extended ones (covering the whole domain of radius up to the limit of 300 pc). Asymmetric regions, which may contain multiple sources, are among the brightest ones.

Figure  15 shows the integrated H$\alpha$ luminosity of the emission regions as a function of the $\sigma$ value of their pseudo-Voigt fitted profile and their corresponding correlation coefficient. As suspected for the diffuse category, low-luminosity and extended regions have a weaker correlation coefficient. Furthermore, Figures 14 and 15 clearly show a lack of extended-high luminosity regions and reveal a possible maximum envelope for a correlation between the compactness and integrated H$\alpha$ luminosity of the regions. Multiple mechanisms could explain this behavior. Firstly, considering that our sample of emitting regions is dominated by HII regions (this is well supported by the BPT diagrams presented in Section 6), the most massive interstellar clouds are the progenitors of the most luminous regions and these clouds are known to be denser than the smaller ones (almost a linear function between the mass of the clouds and the ISM density; Morrison MacLachlan 2015). Due to the global equilibrium between the ISM density and the ionizing flux of the source, the high density of these clouds could explain why the massive regions are smaller in size than the low mass, less dense, regions. Secondly, the number of ionizing photons (QH) emitted by massive stars in the clusters is well constraint by the IMF while more dispersion is observed in low mass clusters where the ionizing stars content reflects a more stochastic sampling of the IMF (Morrison MacLachlan 2015). Additionally, the very bright HII regions observed here (log LHα $>$ 38.5) could be density bounded rather than ionisation bounded. This would also explain that bright HII regions have on average smaller $\sigma$ value compared to fainter regions. In the same way, the fainter regions observed (between log LHα = 37 and 38.5) can clearly be more extended; ionisation bounded regions in a lower density medium could produce this result.

5.3 Size of the Emission Regions

The size of HII regions can be compared to spherical models of uniform density (i.e. Strömgren spheres; Osterbrock & Ferland 2006). This spherical model was chosen here for its simplicity but is certainly not representative of all the observed regions (Wood et al. 2005). The radius of a Strömgren sphere is given by :

[TABLE]

where $Q_{H}$ is the number of ionizing photons (with E $>$ 13.6 eV), $\epsilon$ is the volume filling factor of the ionized gas, $n_{e}$ is the electron density, and $\alpha_{B}$($H$, $T_{e}$) is the recombination coefficient of the hydrogen atom for a given electron temperature ($T_{e}$). Typical values of $\epsilon$ for large HII regions are between 10*-3* and 10*-4* (according to Giammanco et al. 2004 and Cedrés et al. 2013) and $\alpha_{B}$($H$, 10 000 K)\,\,=\,\,2.59$$\times10*-13cm3* s*-1* (Osterbrock 1989). The number of ionizing photons was estimated using the relation given by Morrison MacLachlan (2015):

[TABLE]

where $x$ is a function of the cluster stellar mass : $x\,=\,$log$(M_{cl}/M_{\odot})$. The luminosity of a region was calculated using the equation of Osterbrock & Ferland (2006), $L_{{\rm H}\alpha}=1.36\times 10^{-12}Q_{H}$. It is important to realize that with the Strömgren sphere model selected here, the regions are ionisation bounded and their luminosities are therefore upper limits.

Leaving aside emission regions with a poor quality H$\alpha$ flux profile (i.e. with a correlation coefficient R${}^{2}<0.6$), we estimate the radius $R_{HII}$ of the 3466 remaining emission regions, considering that they are mainly HII regions. To convert the pseudo-Voigt profile half-width $\sigma$ into a radius, the uniform spherical models intensity profiles were convoluted with a Gaussian to account for differences in image quality and manipulations during data reduction. Then, the H$\alpha$ intensity profile half-width $\sigma_{g}$ of the models were obtained by fitting a Gaussian. A relation between $\sigma_{g}$ and $R_{S}$ was derived and transposed to our measurements (see the insurt in Fig. 16) to get the corresponding radius $R_{HII}$.

Figure  16 shows the HII region integrated H$\alpha$ luminosity as a function of their radius $R_{HII}$. Theoretical curves for Strömgren spheres with various stellar cluster masses, gas electron temperatures and densities, and volume filling factors are also shown on this plot. It is important to note that no luminosity-size relation, such as the often used L $\propto$ $D^{\eta}$ (where $D$ is the diameter), follows the well-spread distribution of our measurements. In fact, Scoville et al. (2001), who found $\eta$ = 2, stated that this relation is biased by unresolved, superposed, or blended HII regions. Liu et al. (2013) found a value of $\eta$ = 3 but also made a careful analysis of the possible misinterpretation of the data, and concluded that the estimation of $\eta$ is strongly affected either by the algorithm used to define the regions or by observational biases. Here, the combination of the good spatial resolution, our emission region identification technique, and radius definition method reveal the complex nature of the luminosity-size distribution of NGC 628 HII regions. Multiple models for different ISM conditions are required to understand this distribution. From the comparison with the theoretical curves, we can see that the electron temperature has a very small effect on the HII region radius, whereas both the electron density and the volume filling factor can be held responsible for most of the scattering. Also, as the size of a bright HII region increases, the electron density and/or the volume filling factor decreases. From the general distribution of the low-luminosity regions and the models, we can conclude that they are either dense and small, with a higher electron density and volume filling factor, or (most probably) diffuse and large, with a smaller electron density and volume filling factor in agreement with the fact that low-mass HII regions are less dense and ionisation bounded. For bright regions, the size is globally smaller (between 30 and 100 pc) in comparison with the low-luminosity sample. The bright HII regions have either a high electron density and/or a high volume filling factor.

Figure  17 shows the size distribution function of the HII regions. We simply plotted the histogram of log(N($R_{HII}$)d$R_{HII}$) as a function of $R_{HII}$d$R_{HII}$, where N($R_{HII}$)d$R_{HII}$ is the number of regions with a radius in the range $R_{HII}$ to $R_{HII}$+d$R_{HII}$. The ability to recover the size of a region is limited by the spatial resolution of the data and the physical distance between regions. Therefore, the measured size of the smallest regions is only an upper limit. Although the physical meaning of this plot is not clear, it is interesting as it presents a relatively simple slope, $\rho$ = 0.016$\pm$0.001 (measured between a radius of 70 and 250 pc).

5.4 Luminosity Function of the HII Region Candidates

The total H$\alpha$ flux of an emission region is obtained by summing the flux of all the pixels up to the outer limit radius, and after subtraction of the DIG background. Using a distance of 9.006 Mpc for the galaxy (NED), the total flux is then converted into the H$\alpha$ luminosity. For some fainter regions, the central peak can be well detected but noise in the region’s boundary can produce a total luminosity below the initial detection threshold. Therefore, from the 4285 regions detected, 4158 have a log(LHα in erg s*-1*) above the detection limit of 35.65.

Considering that the majority of the emission regions selected are HII regions (this is well supported by the BPT diagrams discussed in Section 6), an HII region luminosity function was calculated using the sample of 4158 emission regions and is shown in Figure  18. It follows the relation N(L)dL = ALαdL, where N(L)dL is the number of regions with a luminosity in the range L to L+dL. The analysis was done by considering regions before and after the extinction correction. Bins of 0.05 in log(LHα in erg s*-1*) were selected due to the data quality, allowing us to show a smooth and well-sampled curve. The luminosity function slopes derived with the extinction corrected data, $\alpha_{cor}\,=\,-1.12\pm 0.03$, and without the correction, $\alpha_{obs}\,=\,-1.10\pm 0.04$, are identical. The median luminosity for the extinction corrected data is 36.95$\pm$0.02, and for the uncorrected data is 36.77$\pm$0.02. No break in the slopes is seen, as previously noticed at high luminosity, e.g. Bradley et al. (2006).

We looked for variations of the luminosity function with the galactocentric distance of the HII regions, separating them into two samples, one inside and one outside a radius $R_{G}$ of 6.5 kpc (in the case of the extinction corrected data only). This distance was first selected because it simply gives a nearly equal number of regions in both samples (2075 regions inside and 2083 regions outside 6.5 kpc). These two samples are also covering a different gas metallicity range, as already known from the metallicity gradient reported by several authors (e.g. Sánchez et al. 2011; Pilyugin et al. 2014; Berg et al. 2015) and as supported by the line ratios studied in Section 6. Figure  19 shows the resulting luminosity functions. We do not see a significant change in the slope of the function with respect to the radial position in the disc ($\alpha_{>\,6.5\,kpc}$ = $-1.11\pm 0.04$ and $\alpha_{<\,6.5\,kpc}$ = $-1.14\pm 0.05$). Regions in the central portion of the disc ($<$ 6.5 kpc) are in general slightly brighter with a median value of log Lmed = 37.01$\pm$0.02 compared to log Lmed = 36.86$\pm$0.02 in the external portion of the disc. Indeed, brighter regions (logL $>$ 37.5) and less faint regions (logL $<$ 37) are seen in the central portion of the disc ($R_{G}$ $<$ 6.5 kpc). This could results from a crowding effect in the central disc, where many small regions overlap and are measured as one bright region. To evaluate if the level of crowding is different in both subsamples, the distance between two peaks was evaluated for the closest, the second-closest, and the third-closest peak for all the 4285 peaks in the sample. We found a median value for the distance between peaks of 110.0 pc at $R_{G}$ $<$ 6.5 kpc and 115.5 pc at $R_{G}$ $>$ 6.5 kpc for the closest, 159.5 pc at $R_{G}$ $<$ 6.5 kpc and 168.3 pc at $R_{G}$ $>$ 6.5 kpc for the second closest, and 198.6 pc at $R_{G}$ $<$ 6.5 kpc and 211.0 pc at $R_{G}$ $>$ 6.5 kpc for the third closest peak. The error for all median values is $<$ 1 pc if evaluated by moving the galactocentric radius threshold by 500 pc. Although small, this difference supports the idea that a crowding effect could explain the small differences observed in the luminosity functions presented in Figure  19.

Figure  20 shows the luminosity function for two samples of HII region candidates separated according to their DIG background. A threshold for the DIG surface brightness of 2.05 $\times$10*-18* erg s*-1cm-2* over one pixel was used to separate the two samples. Although some very bright regions with a significant DIG background are not in the arms, the higher DIG background group contains regions mostly located in the core of the spiral arms (as previously found by Kreckel et al. 2016). As the bright regions have brighter DIG background, it suggests that the main ionisation source for the DIG is within the regions themselves, i.e. massive stars. When excluding spiral arm regions, we find that the DIG background level is correlated with the galactocentric radius, as higher DIG emission is measured in the regions closer to the galaxy centre (see Fig.  21). This suggests that the evolved star background, more important in the galaxy inner region, is also a source of ionisation for the inter-arm regions. The selected DIG threshold also ensures two samples with a similar size (2075 regions above and 2083 below the threshold). In this case, a significant difference is measured for both the median value and the slope of the luminosity functions (Fig.  20). The sample with a lower DIG background has lower values : log(L${}_{med})\,=\,36.75\pm$0.02, compared to 37.14$\pm$0.02 for the high DIG background regions, and a slope $\alpha_{low\,DIG}\,=\,-1.46\pm 0.06$, compared to $\alpha_{high\,DIG}\,=\,-1.11\pm 0.04$. The slope of the lower DIG background regions is similar to the slope found in M31 ($\alpha_{M\,31}\,=\,-1.52\pm 0.07$) by Azimlu et al. (2011). Therefore there is a clear relation with the DIG background and the average luminosity of the emission regions, here associating the DIG origin with the importance of the massive star population. But in the case of the lower DIG background regions, the change in the slope of the luminosity function and the correlation with the galactocentric radius suggest multiple origins for the DIG background.

6 LINE RATIOS OF THE EMISSION REGIONS

6.1 Line Ratios for Individual Spaxels

Line ratios were calculated using flux maps of the lines extracted with ORCS ($\S$  3.5). Line ratios of [OII]/[NII], [OIII]/[NII], [OIII]/[OII], [OII]/H$\beta$, ([OII]+[OIII])/H$\beta$, [NII]/H$\alpha$, [SII]/H$\alpha$, [SII]/[NII], [OIII]/H$\beta$, and [SII]$\lambda$6716/[SII]$\lambda$6731 were obtained after the subtraction of the stellar population contribution (§ 3.4) and the extinction correction (§ 5.1). Appendix B contains line ratio maps covering the whole galaxy as well as for four enlarged reference sections of the disc. A detection threshold of 3$\sigma$ for the individual lines was applied to produce these maps. The maps of the four enlarged reference sections show variations of the line ratios within emission regions themselves (the H$\alpha$ emission peaks are all indicated with a cross on the images), and variations of the DIG background. For example, in Figure  32, the [NII]/H$\alpha$ ratio can change significantly within a distance of 100 pc. As seen in Figure 33, the H$\alpha$ emission peak may be superimposed to low- or high-local values of the ratio and large fluctuations can be found in the DIG surrounding the regions. Particularly in the case of [SII]/H$\alpha$ (Fig.  35), there is a strong correlation between the position of the peaks and local minimum values. This behavior is expected as this ratio is known to increase when the DIG contribution is important (Hafner et al. 1999; Elwert & Dettmar 2005). But also, with a resolution of 35 pc, we could as well be seeing the ionisation structure within the HII regions (Levesque et al. 2010; Alarie, private communication).

Figure  22 shows the BPT diagram (as initially defined by Baldwin et al. 1981) of [OIII]/H$\beta$ vs [NII]/H$\alpha$ for all the spaxels measured over the whole galaxy disc. Two versions of the same diagram are shown, one colour coded according to the H$\alpha$ flux for each spaxel and the other according to the spaxel density. As expected, most spaxels lie in the HII region ionisation regime area of the BPT diagram and a few spaxels cross the transition zone and AGN limit. The general trend drawn by the bright H$\alpha$ spaxels (red curve in the Fig.  22) can be represented by the polynomial fit :

[TABLE]

where $x$ = log([OIII]/H$\beta$). Globally, the fit (going from high to low [OIII]/H$\beta$ ratio) follows an increase in the gas metallicity. We can also notice a turnover at log([OIII]/H$\beta$) = $-$0.5 and log([NII]/H$\alpha$) = $-$0.5, where the [NII]/H$\alpha$ ratio starts to decrease while the [OIII]/H$\beta$ continues to decrease. This turnover is expected when the metallicity reaches a value above 12+log[O/H] $\simeq$ 9 (Fig. 7 from Kewley & Dopita 2002) while considering that the number of ionizing photons and therefore the ionisation parameter ($q$ = $Q$H/$n_{H}$, the number of hydrogen ionizing photons over the hydrogen atom density) both decrease with increasing metallicity (number of ionizing photons as a function of the metallicity; Topping & Shull 2015, or ionisation parameter as a function of the metallicity; Pérez-Montero et al. 2014). A turnover may also be due to a variation in the relative abundances (e.g. N/O) in the galaxy. Therefore, the location of the polynomial fit and the possible turnover are expected to vary for galaxies with different metallicities, relative abundances, and other properties.

In Figure  22, the spaxels that are bright in H$\alpha$ are predominantly associated with massive HII regions, where it is most likely that the IMF is well represented (i.e. no stochastic sampling with the theoretical slope and mass limits). The dispersion observed along the fit in the bright H$\alpha$ spaxels could therefore be mainly attributed to different ages (as the ionizing spectrum gets less energetic with time with the death of the most massive stars). On the contrary, the dispersion of faint H$\alpha$ spaxels can be due to multiple causes: age, IMF under-sampling, DIG contamination, different sources of ionisation, SNR, etc.

In the right panel of Figure  22, the density envelop of the BPT diagram shows a concentration of spaxels around the ratio values [$-$0.5, $-$0.5], associated with some of the most intense peaks, and a more extended cloud of spaxels centered at [$-$0.8, 0.2]. Many of the spaxels in this clouds reach the transition and AGN zones, suggesting that young stellar clusters are not their main ionizing source. These spaxels are not associated with any particular structure (like the arms) in the galaxy, but, as shown in the following figures, are more related to the less bright regions.

Figure  23 shows the [OIII]/H$\beta$ vs [NII]/H$\alpha$ BPT diagram by taking into account the spaxel distance to its ionizing source (i.e. the distance between the spaxel and the emission peak within the zone of influence defining an emission region). Diagrams are displayed independently for very bright regions (i.e. with an integrated H$\alpha$ flux above 5$\times$10*-15erg  s-1cm-2*) and for moderately bright spaxels (i.e. with an integrated H$\alpha$ flux above 5$\times$10*-16erg  s-1cm-2*. In the case of very bright regions (Fig.  23 left panel), the further away they are from their emission peak, the greater is the dispersion around the best polynomial fit. In the case of the sample including moderately bright regions (Fig.  23 right panel), the relation with the distance to the emission peak is more ambiguous. In this case, for moderately bright regions, the effect of the distance is mixed with other parameters like the age, IMF under-sampling, DIG contribution, and other sources of ionisation, as seen in the previous paragraphs.

Figure  24 shows the [OIII]/H$\beta$ vs [NII]/H$\alpha$ BTP diagram as a function of the spaxel galactocentric radius $R_{G}$, for bright regions and faint regions. The general distribution of the colours in these plots reflects the metallicity gradient of the galaxy. The central spaxels (in blue) are more metal rich than those of the extended disc (in yellow-orange-red), as they have larger values of [NII]/H$\alpha$ and smaller values of [OIII]/H$\beta$ (as expected from the global metallicity indicators, e.g. Pettini & Pagel 2004). Nevertheless, the northern spiral arm, located at a distance larger than 15 kpc (red dots in the plots near log([OIII]/H$\beta$) $\simeq$ 0), behaves differently than this somewhat traditional linear-negative metallicity gradient. This could be due to a flattening or an inversion of the slope of the metallicity gradient at large galactocentric radius and/or to a change in the gradient of the relative abundance of N/O.

Figure  25 presents the [OIII]/H$\beta$ vs [SII]/H$\alpha$ BTP diagram. As shown in the left panel, the brightest spaxels are mostly found within the HII region ionizing regime. When considering only the brightest regions (the right panel), the spaxels that are closer to the emission peak display a smaller [SII]/H$\alpha$ ratio in general. Larger values for this ratio are a signature of DIG and/or of the ionisation structure, as mentioned already.

6.2 Line Ratios for Each Emission Region

For each region identified, as described in Section  4, global line ratios were evaluated by fitting the integrated flux that is encompassed within an aperture radius corresponding to the half-width $\sigma$ of the region pseudo-Voigt fitted intensity profile. All spaxels within this aperture have been corrected for the stellar population (§ 3.4) and the extinction (§ 5.1), but not for the DIG background. We adopted the half-width $\sigma$ radius instead of the region s outer limit in order to minimise effects related to the DIG background and to the crowding of the regions. All lines in the integrated spectrum of each region were fitted using ORCS as described in Section  3. Figure  26 shows the BPT diagrams of [OIII]/H$\beta$ vs [NII]/H$\alpha$, [OIII]/H$\beta$ vs [SII]/H$\alpha$, [OIII]/H$\beta$ vs [OII]/H$\beta$, [OII]/H$\beta$ vs [SII]/H$\alpha$, and [OII]/H$\beta$ vs [NII]/H$\alpha$ separating the regions according to the SNR. We defined the parameter SNRcross to consider a signal threshold which takes into account the different SNR of the lines involved in a diagram: first all the lines used in a diagram must have an SNR $>$ 3 and then, at least one of these lines must have an SNR greater than SNRcross. In the diagrams, different colours are used for the emission regions according to different values of SNRcross (as specified in Fig. 26). A larger SNRcross basically points to a brighter region. The bright HII regions of the CHAOS project (Berg et al. 2015) are also shown in our diagrams.

Again, most regions identified fall in the ionizing zone for the HII regions. The comparison with the CHAOS sample, where only the highest H$\alpha$ surface brightness knots have been selected (Berg et al. 2015), demonstrates SITELLE’s ability to reproduce line ratios previously observed for these regions, i.e. the CHAOS data are well superimposed to the SITELLE s regions with the highest SNRcross threshold selected. As these bright regions probably represent the most massive and youngest HII regions, our larger sample includes objects with a wider range of physical properties (masses, ages, stellar content, ISM properties, etc.) allowing us to see an unbiased distribution in the BPT diagrams.

The Appendix B also presents plots of line ratios as a function of the galactocentric radius. While the plots for [NII]/H$\alpha$ and [SII]/H$\alpha$ show decreasing ratios at large radius, the other ratios have a tendency to increase or stay flat at large radius. These behaviours indicate changes in the physical gas parameters, like the abundances, which will be discussed in a following paper. The comparison with the CHAOS bright HII regions sample (superposed in red on the SITELLE data in these plots) shows a larger dispersion for the SITELLE data since they contain a broader diversity of emission regions detected. For the bright regions, SITELLE is in good agreement with CHOAS for most line ratios, except for those involving the [SII] lines (Fig. B3 and B5). In these plots, the SITELLE’s ratios [SII]/H$\alpha$ and [SII]/[NII] are higher, but also the SNRcross is somewhat smaller than in Figure 26. The discrepancy with the CHAOS data may be explained by the fact that the [SII] lines are semi-strong lines, are more dependent on the electron density, and are more sensitive to the DIG (Hafner et al. 1999; Elwert & Dettmar 2005) and its background subtraction. Berg et al. (2015) mention a sky subtraction for the CHAOS data using a 2D background fit within the slit used to observe an HII region; this background may include a significant contribution from the DIG. In our case, we measured the line ratios using a limited aperture size (i.e. the half-width of the pseudo-Voigt fitted on the region intensity profile) to minimise the DIG without doing any subtraction. But it could also be due to a selection bias, as bright regions (with high ionisation parameter) are more likely to present low [SII]/H$\alpha$ ratios. Figure B19 presents the ratio [SII]$\lambda$6716/[SII]$\lambda$6731 as a function of the galactocentric radius where we are considering a higher SNR limit (from 5 to 10). For the brightest regions with SNR $>$ 10, the agreement with the CHAOS sample is then excellent suggesting indeed that the electron density in bright HII regions is often lower than the limit of sensitivity of the ratio (i.e.[SII]/[SII] $\simeq$ 1.4 when $n_{e}$ $>$ $10^{2}$ electrons/cm3). This does not exclude as well the fact that a larger SNR may be required when measuring semi-strong line ratios due to the propagation of the uncertainties of both lines.

7 Summary

In this first paper of a large program focusing on star formation and chemical evolution in disc galaxies, we have shown the high potential of SITELLE for the study of nearby extragalactic HII regions from observations on the well-known spiral NGC 628. Our results are as follows:

From our datacubes covering the whole disc of NGC 628, we identified 4285 HII region candidates. Ten strong and semi-strong emission lines have been sampled in about 4 billion spaxels using automated routines specifically created for SITELLE. We developed a new identification procedure for the regions, i.e. to find out all the emission peaks and define their zone of influence. Using both the spatial and spectral information, this new procedure optimizes the number of parameters extracted for each region (RA/DEC position, dust extinction, velocity, H$\alpha$ profile, DIG background, luminosity, size, morphological type, integrated line fluxes, etc.).

• 2)

By comparing the intensity profile of the regions, four morphological categories, from spherical to diffuse, were defined. These morphologies point out to different types of stellar aggregates (dense or sparse), different gas distributions (homogeneous to filamentary), and multiple sources of ionisation.

• 3)

A well-sampled luminosity function of the HII region candidates was obtained and we studied its variation as a function of the location of the regions along the galactic radius and the DIG background flux. Using the complete sample, a constant slope of $\alpha_{cor}\,=\,-1.12\pm 0.03$ (corrected for the dust extinction, using the Balmer lines, and for the galaxy disc stellar populations on the line of sight to each region) was found; no evidence of a break in the distribution was noticed. The slope remains the same while selecting subsamples of regions located at different galactic radii (which can be seen as different metallicity ranges). A significant change of the slope is seen when considering only regions with a low DIG background. We noticed that high DIG regions are close or belong to the spiral arms of the galaxy (in agreement with the findings of Kreckel et al. 2016). A correlation between the DIG background flux and the total H$\alpha$ luminosity of the bright regions was found, suggesting that ionizing photons primarily related to massive stars and escaping directly from these regions are the main source for the DIG. We also showed that the global DIG background (in fainter regions outside the spiral arms) increases slowly toward the galaxy centre. This behavior rather points to an additional source for the DIG which is related to evolved stars, as these are more numerous closer to the galaxy centre.

• 4)

The high spatial resolution of SITELLE allowed us to compare the integrated H$\alpha$ luminosity of the regions with their intensity profile width (i.e. $\sigma$ of a pseudo-Voigt fit). Diffused/extended high luminosity regions are not seen while low luminosity regions present a wide range of compactness. This is in agreement with the proposition that denser clouds are the progenitors of the most massive HII regions (Morrison MacLachlan 2015) and the fact that massive regions are density bounded rather than ionisation bounded. In less massive HII regions, the number of ionizing photons is subject to IMF fluctuations and the Strömgren radius of regions in a lower density medium is ionisation bounded. Therefore low luminosity HII regions present a wider range in size.

• 5)

The size of the regions was evaluated using the H$\alpha$ profile and a simple Strömgren sphere model. The size distribution function (i.e. number of regions versus their radius) shows a well-defined slope $\gamma$ = 1.81$\pm$0.02 for radii between 70 to 250 pc, as found in other studies.

• 6)

We looked at the distribution of the line ratios for the individual spaxels and for the integrated regions in different BPT diagrams. Bright spaxels show a distinct turnover in the [OIII]/H$\beta$ vs [NII]/H$\alpha$ diagram. A larger dispersion is seen for the fainter H$\alpha$ spaxels which can be due to their age, IMF under-sampling, DIG contamination, different sources of ionisation, SNR, etc. High-definition maps of the line ratios reveal the complex variations of the ionisation conditions within the emission regions themselves. A polynomial fit of each line ratios as a function of the galactocentric radius is now available for comparative studies.

• 7)

And finally, detailed information on the emission line regions of NGC 628 have been extracted from the SITELLE data. We have consolidated the information for all the 4285 regions identified in a catalogue.

Table  3 shows an example of the catalogue content for the HII region candidates in NGC 628. The objective of this new catalogue is to provide a complete database to study the vast parameter space covered by star-forming regions. Further work on NGC 628 will address HII region characteristics and nebular physics (abundances, mass, age, ionisation parameter, density) in a more detailed manner. In a following paper, we will use photoionisation models to investigate the ionisation conditions in the HII regions of NGC 628.

acknowledgments

This research is based on observations obtained at the Canada-France-Hawaii Telescope (CFHT) which is operated from the summit of Maunakea by the National Research Council of Canada, the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique of France, and the University of Hawaii. The observations at the Canada-France-Hawaii Telescope were performed with care and respect from the summit of Maunakea which is a significant cultural and historic site. The observations were obtained with SITELLE, a joint project between Université Laval, ABB-Bomem, Université de Montréal and the CFHT with funding support from the Canada Foundation for Innovation (CFI), the National Sciences and Engineering Research Council of Canada (NSERC), Fonds de Recherche du Québec - Nature et Technologies (FRQNT), and CFHT.

LRN, CR, and LD are grateful to the Fonds de recherche du Québec - Nature et Technologies (FRQNT) for individual and team financial support. LD thanks the Canada Research Chair program. CR and LD are grateful to the Natural Sciences and Engineering Research Council of Canada (NSERC). This research has made use of NASA’s Astrophysics Data System and of the VizieR catalogue access tool, CDS, Strasbourg, France.

Appendix A Flux and Velocity Maps

The flux map for all the emission lines seen in the three filters are shown in Figures 27 and 28. Emission lines from the same filter have been fitted simultaneously using the routine from ORCS $\S$  3.5) applied to all the spaxels individually. The subtraction of the galaxy stellar population spectrum on the line of sight for each spaxels was done ($\S$  3.4). The correction for the internal extinction has not been applied for these figures. For these figures, all spaxels have been kept regardless of the SNR measured for the lines.

The velocity map in Figure  29 was produced after a refinement of the wavelength calibration of the SN3 datacube using a sky line ($\S$  3.1). All the emission lines from this high-spectral resolution datacube (R $\simeq$ 1800; Tab. 1) have been used simultaneously to extract the velocity in each spaxel individually using the line fitting routine from ORCS ($\S$  3.5). The consideration of all the lines simultaneously, as well as the fact that the lines have a sinc profile, improve significantly the precision of the velocity measurements (Martin et al. 2016, 2018). For example, using all the lines in the SN3 filter for spaxels with SNRHα between 3 and 20 gives a velocity precision between 30 and 5 km s*-1*. All spaxels with SNR${}_{{\rm H}\alpha}>3$ have been kept in Figure  29.

Appendix B Line Ratios

Appendix B contains maps of the different line ratios, for individual spaxels covering the whole galaxy and for four enlarged reference sections in the disc (their exact location may be seen in Fig.  28). A detection threshold of 3$\sigma$ on the individual lines was applied to produce the maps. The stellar population contribution (§ 3.4) was subtracted and theextinction correction was done (§ 5.1) before the lines have been measured. The first set of figures (Fig.  30 and 31) presents the H$\alpha$/H$\beta$ emission line ratio along with the extinction map for the regions.

This appendix also presents plots of line ratios as a function of the galactocentric radius for the integrated line flux ratios in each emission region. The regions are encompassed within an aperture radius corresponding to the half-width $\sigma$ of the region pseudo-Voigt fitted intensity profile. All spaxels within this aperture have been corrected for the stellar population (§ 3.4) and the extinction (§ 5.1), but not for the DIG background. These plots are colour coded according to different SNR thresholds. The best polynomial fit obtained in the case of the highest SNRcross subsample is drawn over each plot (and its corresponding equation is written at the top of each plot). The HII regions from the CHAOS project (Berg et al. 2015) are superimposed to the SITELLE’s regions in these plots.