© M. C. i Bernadich et al. 2022

1. Introduction

Ultraluminous X-ray sources (ULXs) are extragalactic point-like X-ray sources whose luminosity exceeds that of the Eddington limit of an accreting stellar-mass black hole (StBH, M_BH < 10 M_⊙, L_Edd ≈ 10³⁹ erg s⁻¹), and that are not the central source of their host galaxy. These objects have been the subject of active research since the advent of X-ray astronomy (see Kaaret et al. 2017 for a review) and for a time, due to the scaling of the Eddington limit with mass, it was believed that they hosted accreting BHs of masses between those of StBHs and supermassive BHs. Today they are more commonly associated with super-Eddington accretion onto common compact objects such as StBHs and neutron stars (NSs), with a few exceptions where BHs of higher mass need to be invoked (Bachetti 2016).

The current understanding comes from three main lines of evidence. Firstly, predictions of the formation rate for intermediate-mass BHs (IMBHs, 100 M_⊙ < M_BH < 10⁶ M_⊙) with donor companions fall short of explaining the observed abundances of ULXs (King 2004; Madhusudhan et al. 2006). Secondly, ULXs are most commonly associated with star-forming galaxies (Swartz et al. 2011) or even young elliptical galaxies with recent star formation events rather than old elliptical ones with no star formation (Kim & Fabbiano 2004, 2010). This is consistent with ULXs belonging to a high-luminosity extrapolation of the X-ray binary population whose abundances correlate with the star formation rate (SFR) and the stellar mass (M_*) of their host galaxies (Gilfanov 2004; Grimm et al. 2003; Mineo et al. 2012). And finally, but perhaps most importantly, the direct observation of regular pulsations in some ULXs unmistakably points towards accretion onto a NS in binary systems. Famous examples of neutron star powered ULXs (NS-ULXs) are M82 X-2 (Bachetti et al. 2014), NGC 7793 P13 (Fürst et al. 2016; Israel et al. 2016), NGC 5907 ULX1 (Israel et al. 2017), and NGC 300 ULX1 (Carpano et al. 2018). Since some of these objects present L_X > 10⁴⁰ erg s⁻¹, they are clear cases of super-Eddington accretion. In line with these observations, population synthesis models even suggest that NS-ULXs constitute a significant fraction of the ULX population (Wiktorowicz et al. 2017).

The prospect of super-Eddington accretion in ULXs has sparked a great deal of interest for the modelling of ULX systems. Many works have shown that sufficiently extreme accretion rates to power ULXs can be achieved during certain phases of regular HMXB and LMXB evolution. These rates can actually explain most of the ULX population (King 2002; Rappaport et al. 2005; Wiktorowicz et al. 2015; Pavlovskii et al. 2017). Even so, many aspects of the physics of super-Eddington accretion itself are poorly understood, especially in cases where the material falls onto strongly magnetised neutron stars, and very rich theoretical research is being performed to model this phenomenon (see Kaaret et al. 2017 for a comprehensive discussion of this topic, which otherwise lies outside of the scope of this Introduction).

Nonetheless, ULXs are still good sources to look for IMBHs as a handful of especially bright candidates still remain. ESO 243-49 HLX-1 has been detected with luminosities of L_X ≈ 10⁴² erg s⁻¹ and with spectral state transitions typical of sub-Eddington Galactic BH binaries (GBHBs, Remillard & McClintock 2006), being consistent with a mass of ∼ 10³ M_⊙ (Servillat et al. 2011). M82 X-1 shows a similar behaviour (Kong et al. 2007), with the addition of quasi-periodic oscillations (QPOs, Remillard & McClintock 2006) also characteristic of GBHBs that grant it a mass estimate of 400 M_⊙ (Pasham et al. 2014). And just as in GBHBs, powerful radio jets are also present in the IMBH candidates ESO 243-49 HLX-1 (Webb et al. 2012; Cseh et al. 2014) and NGC 2276-3c (Mezcua et al. 2015).

The relevance of ULXs as proxies of recent star formation, in accretion physics, and in IMBH searches has been so far revealed. However, ULX studies are limited by their small numbers and large distances. For this reason, a great deal of effort has been put into building compilations of ULX candidates from X-ray surveys. Early works relied on ROSAT data (e.g. Roberts & Warwick 2000; Liu & Bregman 2005; Liu et al. 2006; Colbert & Ptak 2002), yielding up to ∼100 candidates. More recent works have been based on XMM-Newton and Chandra data. Catalogues based on XMM-Newton have been highly successful, yielding up to ∼300 ULX candidates (Walton et al. 2011; Earnshaw et al. 2019), but they are somewhat hampered by the low angular resolution of the observatory; works based on the Chandra telescope have provided better identification of individual sources thanks to its crisp angular resolution (e.g. Swartz et al. 2004, 2011; Wang et al. 2016; Kovlakas et al. 2020).

Of particular interest to us are the works of Earnshaw et al. (2019) and Kovlakas et al. (2020) since they both constructed ULX samples from the largest and most recent XMM-Newton and Chandra samples available at the time. Earnshaw et al. (2019) is, in fact, our main predecessor and inspiration regarding methodology. They identify 384 ULX candidates within the fourth XMM-Newton data release¹ (3XMM-DR4, Rosen et al. 2016), using the Third Reference Catalogue of Bright Galaxies (RC3, de Vaucouleurs et al. 1991) and the Catalogue of Neighbouring Galaxies (Karachentsev et al. 2004) as references for the host galaxies. Their study focuses mostly on the spectral properties of ULXs, and they find that ULXs tend to be somewhat harder in late-type galaxies, and that their hardness ratios in general are like those of X-ray binaries below the Eddington threshold. Kovlakas et al. (2020) find 629 ULX candidates within the Chandra Source Catalog 2.0² (CSC2) using the Heraklion Extragalactic CATaloguE (HECATE, Kovlakas et al. 2021) as the reference for host galaxies, which comes with valuable information such as the SFR and M_*. Their study focuses heavily on the scaling properties of the ULX content of galaxies with the SFR and M_* parameters, and finds that SFR is the determinant parameter for late-type galaxies, while M_* dominates over the abundance in late-type galaxies, in accordance with recent population synthesis models (Wiktorowicz et al. 2017).

In this work, we used the ninth XMM-Newton data release³ (4XMM-DR9, Webb et al. 2020), its stacked version (4XMM-DR9s, Traulsen et al. 2020), and the HECATE catalogue to build the largest catalogue of extragalactic non-nuclear point-like X-ray sources based on XMM-Newton data. Then, we used the catalogue to preliminarily explore the distribution of ULX spectral and variability properties, and the dependence of ULX abundances galaxies with the SFR and M_* on XMM-Newton’s sky.

This paper is organised as follows: in Sect. 2 we explain the specifics of selection in archival data, in Sect. 3 we describe our methodology for identifying ULX candidates and other interloping objects, in Sect. 4 we present our main results, and in Sect. 5 we discuss some of the limitations and implications of our work and compare it with previous works.

2. Data samples

2.1. X-ray samples

We used the Fourth XMM-Newton Serendipitous Source Catalog, Ninth Data Release, (4XMM-DR9, Webb et al. 2020) as the basis for our ULX catalogue. It lists 810 795 individual detections of 550 124 unique sources discovered across 11 204 observations, and it represents an increase of 177 396 sources with respect to 3XMM-DR4, the resource used in the previous ULX study based on XMM-Newton data (Earnshaw et al. 2019). The columns we used during the construction are the following: DETID and SRCID, the identification number of detections and sources; SCRA, SCDEC, and SCPOSERR, the source J2000.0 sky coordinates and their positional uncertainty; SCEXTENT, the measured extension of the source; SCDETML, the detection likelihood of sources taken from the highest value of DETML among their detections; EP8FLUX and EP8FLUXERR, the measured X-ray flux of detections and its uncertainty in the 0.2–12 keV band derived from the European Photon Imaging Camera⁴ (EPIC) photon counts; SCEP8FLUX and SCEP8FLUXERR, the averaged X-ray flux of sources and its uncertainty in the 0.2–12 keV band derived from the EPIC photon counts; SCHRi, where i runs from 1 to 4, the source hardness ratios from the source count rates in the respective energy bands; and SCSUMFLAG, the summary quality flag of a unique source, which runs from 0 to 5 (0 means all detections of the source are trustworthy; 5 means at least one detection is most likely spurious), a value given by the highest value of SUMFLAG among all detections of the source.

Since the XMM-Newton FWHM resolution is of 6″ (ESA: XMM-Newton SOC 2019), we included only sources with SC_EXTENT < 6″ to work with point-like sources alone. We also considered only sources with SC_DET_ML > 8. As these parameters are based on the worst result for all the detections of a source, it can happen that an otherwise point-like source is counted as extended if in one of its detections it is measured as such. We nevertheless excluded them from the catalogue because other source parameters that consist of averaged detection parameters are most likely unreliable. This left 452 602 sources available for the catalogue.

Additionally, we used the stacked version of 4XMM-DR9. This is the second XMM-Newton serendipitous source catalogue built from overlapping observations (4XMM-DR9s, Traulsen et al. 2020), which compiles variability information for 218 283 sources detected in 6604 overlapping XMM-Newton observations, an improvement of 146 332 sources with respect to the first version (Traulsen et al. 2019). The strength of 4XMM-DR9s is that, using the overlapping observations, extra source variability parameters that are not included in the bare 4XMM-DR9 are computed. We used the following columns from 4XMM-DR9: SRCID, the source identification number; RA, DEC, and RADECERR, the source J2000.0 sky coordinates and their uncertainty values; NCONTRIB, the number of times a source contributes to the computation of variability, usually equal to the number of times it has been observed; VARPROB, the probability of a source not showing inter-observation variability, computed from the reduced χ² of EPIC flux variability; and FRATIO, the F_max/F_min ratio of the highest to the lowest observed fluxes in a single source.

2.2. Galaxy sample

We used the Heraklion Extragalactic CATaloguE (HECATE, Kovlakas et al. 2021) compilation of 204 733 galaxies within 200 Mpc. Built from the HyperLEDA catalogue (Makarov et al. 2014) and other databases such as NED, SDSS, and 2MASS, it is a much more complete galaxy compilation than the Third Reference Catalogue of Galaxies (RC3, de Vaucouleurs et al. 1991), traditionally used during previous ULX studies (Swartz et al. 2011; Wang et al. 2016; Earnshaw et al. 2019). A more detailed description of its contents is available in Kovlakas et al. (2021), but here we mention the columns we used: PGC, the Principal Galaxy Catalogue identification number, originating from the Principal Catalogue of Galaxies (Paturel et al. 1989) and still used in HyperLEDA; RA and DEC, the J2000.0 coordinates of the galactic centre; R1 and R2, the minor and major D25 isophotal radii; PA, the north-to-east position angle of the major axis; D and DERR, the galaxy distance and its uncertainty; T, the Hubble type (T_H) of a galaxy; logSFRHEC, the decimal logarithm of the SFR estimate; and logMHEC, the decimal logarithm of M_*.

The SFR values in HECATE are based on infrared calibrations, which tend to overestimate the SFR in early-type galaxies (Kovlakas et al. 2021). Therefore, SFR values were only inspected for galaxies with T_H ≥ 0.

2.3. The interloper samples

Our X-ray catalogue is inevitably populated by a fraction of non-ULX contaminating objects such as background active galactic nuclei (AGNs) or foreground stars. Therefore, we built a filtering pipeline (see Sect. 3.1) to identify interlopers in external catalogues and databases of already-known objects. These consist of the second Gaia data release (Gaia DR2, Gaia Collaboration 2018), the Tycho-2 catalogue of bright stars (Tycho2, Høg et al. 2000), the 14th Sloan Digital Sky Survey data release (SDSS-DR14, Blanton et al. 2017), the 13th edition of the Véron-Cetty & Véron catalogue of QSOs and AGNs (VéronQSO, Véron-Cetty & Véron 2010), the SIMBAD database (Wegner et al. 2000), the Panoramic Survey Telescope and Rapid Response System database⁵ (Pan-STARRS1, Flewelling et al. 2020) from the Pan-STARRS1 surveys (Chambers et al. 2016), and the NASA/IPAC Extragalactic Database⁶ (NED, Helou et al. 1991).

All the cross-matches were performed using the TOPCAT⁷ and STILTS⁸ software tools to manipulate the tables in all of the following steps (Taylor 2005).

3. Methods

3.1. Automatic filtering process

Correlation of sources with galaxies. We selected all 4XMM-DR9 sources that overlapped with the isophotal ellipses of HECATE within their positional uncertainty. Catalogue entries matched to galaxies with available measurements of R1, R2, and PA were labelled with MATCH_FLAG = 0. When PA was unknown, detections of sources matched with the minor isophotal circle were labelled as MATCH_FLAG = 1, while detections of sources within the annulus drawn by R1 and R2 were given MATCH_FLAG = 2. When both PA and R2 were unlisted, sources matched to the circle drawn by R1 were labelled as MATCH_FLAG = 3. Additionally, we stored in nGalaxies the number of galaxies a source has been matched with. Detections of sources with n_Galaxies > 1 are listed more than once and treated independently for each galaxy, being identifiable by the unique DETPGCID and SRCPGCID detection and source-galaxy identity numbers. In this way a comprehensive catalogue of 50 446 entries was built, of which 49 816 have n_Galaxies = 1 and 48 206 present MATCH_FLAG = 0.

Identification of central sources. At low redshifts AGNs typically present luminosities of L_X > 10⁴¹ erg s⁻¹ (Brandt & Alexander 2015). However, some have been observed to overlap with the ULX luminosity regime (Ho et al. 2001; Eracleous et al. 2002; Ghosh et al. 2008; Zhang et al. 2009). Therefore, we identified all sources whose position overlapped within three times their uncertainty with a circle of radius 3″ around the centre of their host galaxy. A total of 3658 entries in our catalogue were thus labelled as central sources. We also checked for sources with more than one potential host galaxy (n_Galaxies > 1) that had been flagged as central-source candidates in one of their iterations, and flagged them as central-source candidates in all of their other galaxy associations to indicate their potentially nuclear nature. This highlighted 137 further entries, all of which are composed of optically coincident galaxies, whether related or not. Figure 1 illustrates how essential it is to identify these central sources before classifying them as ULXs, as otherwise they would easily slip through and contaminate our sample.

Fig. 1. X-ray XMM-Newton (top) and optical Pan-STARRS1 (bottom) images of galaxy NGC 3631. Non-nuclear sources are indicated with red crosses, while the green cross is the central source. The Pan-STARRS1 image covers a slightly smaller sky area than the XMM-Newton one.

Identification of foreground stellar objects. X-ray emission from stars is ubiquitous across the Hertzsprung-Russell diagram, and there are manifold reasons why it can lead to the misclassification of Galactic objects as extragalactic ULXs, including active coronae, hot stellar winds, or binary accretion. We performed a cross-match of unflagged objects with Gaia DR2 making use of the query tool provided by CDS, Strasbourg, and then with Tycho2 to identify potential contaminating foreground stars. In both cases a positional overlap within three times the positional uncertainty was required on both sides of the cross-match. Assuming that all matched objects were stars, a very stringent constraint of log(F_X/F_V) <  − 2.2 was used to decide whether the matched object could explain the X-ray luminosity of the XMM-Newton source, where F_X is the source X-ray flux provided by SCEP8FLUX and F_V its flux in optical light. This limit is established in consistency with the work of Freund et al. (2018), where log(F_X/F_BOL) = − 2.2 is considered the maximum luminosity ratio for early-type stars, the ones presenting the largest share of X-ray luminosity. However, because the bolometric flux F_BOL was unavailable in most existing catalogues, we approximated it as the optical flux of the star, F_V. We used the formula from Maccacaro et al. (1988) to write

$$\begin{array}{c}\hfill log\left(\frac{F_{\text{X}}}{F_{\text{V}}}\right)=log{F}_{\text{X}}+\frac{m_{\text{V}}}{2.5}+5.37<-2.2,\end{array}$$(1)

where m_V is the optical magnitude of the studied object. In cases where m_V is not listed, we used the G-band magnitude m_G instead, which typically has the same value as the V-band magnitude for bright enough sources. From Gaia DR2, we used the listed G-band magnitude values Gmag, and for Tycho2 we used the VT-magnitude measurements VTmag. This resulted in 2257 entries flagged as stars from Gaia DR2 and 96 entries from Tycho2.

Optical counterparts of ULXs are typically very faint, presenting m_V ≳ 21 (Kaaret et al. 2017), which leads to log(F_X/F_V)≳0 for a ULX candidate with L_X = 10³⁹ erg s⁻¹ at a distance of 20 Mpc. Therefore, we were safe from accidentally disregarding genuine ULX candidates as stellar objects. Objects with 0 > log(F_X/F_V) >  − 2.2 were considered in an upcoming section of the pipeline.

Identification of background QSOs. From the 2–10 keV log N − log S distribution above a flux of 10⁻¹⁴ erg cm⁻² s⁻¹ presented in Mateos et al. (2008) we expected a background source density of 300 deg⁻² in the field. With an accumulated galaxy sky area outside of the local group (D > 1 Mpc) of 2.25 deg², this implied that around 670 background contaminants lie in the line of sight of galaxies where we expected to find most of our ULX candidates. This clearly motivated the need to identify such objects in currently available catalogues. Therefore, all remaining unflagged objects found to overlap within three times their positional uncertainties during the cross-match with the SDSS-DR14 and VéronQSO were flagged as QSOs. The first query was performed with SDSS-DR14 and highlighted 140 entries, while the second highlighted 135 more.

Identification of miscellaneous objects inSIMBAD. We cross-matched our remaining unflagged objects with the SIMBAD database with the query tool provided by CDS, Strasbourg. In this step of the pipeline sources were treated differently depending on the nature of the matched objects, indicated as main_type in the SIMBAD catalogue⁹. Every object overlapping within three times their positional uncertainty and with either main_type = Star or with their maintype containing an asterisk (*), and whose optical magnitude Vmag holds equation (1) was flagged as a stellar contaminant. Objects with main_type = AGN, AGN_Candidate, QSO, QSO_Candidate, or SN were also flagged, regardless of their optical magnitude. Supernovae in particular were highlighted as their emission is dominated from the expanding envelope. This step highlighted 5745 entries in the catalogue.

Identification of miscellaneous optical objects. As previously stated, ULX optical counterparts typically present m_V ≳ 21 (Kaaret et al. 2017), which leads to log(F_X/F_V)≳0 for a ULX candidate with L_X = 10³⁹ erg s⁻¹ at a distance of 20 Mpc. Objects with larger optical luminosity are usually AGNs or foreground stars. Sources with log(F_X/F_VOL) <  − 2.2 were already labelled as stars in a previous section of the pipeline. In this section we looked for optically bright objects of extragalactic origin. We performed a cross-match with the query tool provided by CDS, Strasbourg, with the Pan-STARRS1 catalogue to find all objects that overlapped within a radius of three times their positional uncertainties and that held

$$\begin{array}{c}\hfill log\left(\frac{F_{\text{X}}}{F_{\text{V}}}\right)=log{F}_{\text{X}}+\frac{m_{\text{V}}}{2.5}+5.37<0\end{array}$$(2)

in order to highlight all sources brighter in the optical than in the X-rays, using the G-band magnitude gmag as an estimate for m_V. Up to 3151 entries were flagged as Pan-STARRS1 extragalactic objects.

3.2. Identification of ULX candidates

To identify ULX candidates, we computed individual detection and average source luminosities from the X-ray fluxes as listed in 4XMM-DR9 (Sect. 2.1) and from the distances as listed in HECATE (Sect. 2.2), along with the corresponding error propagation. These are listed as Luminosity and LuminosistyErr for detections and SCLuminosity and SCLuminosistyErr for sources. In our catalogue we consider any X-ray source to be within the ULX luminosity regime if this first condition is met:

(1) It has at least one detection with luminosity above the Eddington limit within the uncertainty (Luminosity + LuminosityErr > 10³⁹ erg s⁻¹) in at least one of the potential host galaxies.

This check identified 5943 entries of diverse nature, corresponding to 3280 sources in 2729 galaxies; 2205 of these objects were flagged as central, while only 856 objects in 552 galaxies did not have counterparts. Furthermore, some of them present large uncertainties in their data. Therefore, we imposed additional quality conditions for a source to qualify as a ULX candidate:

(2) It has no identified counterpart.

(3) It has SC_Luminosity > SC_LuminosityErr, as otherwise it would indicate that its source parameters are unreliable.

(4) It has SC_SUM_FLAG < 1, which means that we only considered sources for which none of the individual detections was flagged as probably spurious.

(5) It has a single potential host galaxy (n_Galaxies = 1), as otherwise their distances and luminosities are unreliable.

This originally left 730 ULX candidates in 490 galaxies. Additionally, we created a sub-category of bright ULX candidates, which follows that Luminosity + LuminosityErr > 5 × 10⁴⁰ erg s⁻¹, into which only 130 in 122 galaxies were included. We also found seven sources with more than one potential host galaxy (n_Galaxies > 1), but whose distances differed by only a small fraction and thus did not affect their ULX status. These sources were included to the ULX count, one of them being a bright candidate.

Sources that were below the ULX luminosity regime in all of their detections but still hold conditions 2–5 are simply referred to as quality sources. Many of these sources and low-luminosity interlopers identified in Gaia DR2, Pan-STARRS1, and SIMBAD are found in nearby galaxies, as shown by Fig. 2.

Fig. 2. TOPCAT astrometric map of M51 showing all 4XMM-DR9 sources within its isophotal ellipse. The interlopers are composed of Pan-STARRS1 and SIMBAD objects, and all of them but one are below the ULX luminosity regime if assumed to be at the same distance as the others.

3.3. Manual inspection of candidates and contaminants

Our procedure revealed 50 ULX candidates with L_X > 10⁴¹ erg s⁻¹, well within the AGN luminosity regime (Brandt & Alexander 2015). However, we had to consider the possibility of some of them being contaminants in disguise that survived the filtering pipeline. With the aim of identifying potential counterparts, we inspected the available optical and X-ray images from the Pan-STARRS1 and XMM-Newton surveys for 112 bright ULX candidates. Here we discuss some of the relevant results.

Source SRCID = 200559905010001 (Fig. 3), with detections of L_X = (4.3 ± 1.3)×10⁴² erg s⁻¹ and L_X = (8.2 ± 2.3)×10⁴² erg s⁻¹, is most likely a blend of the two AGNs in the interacting galaxy pair NGC 5256. Source SRCID = 206701401010002 presents ⟨L_X⟩=(4.5 ± 0.9)×10⁴¹ erg s⁻¹ and is found in the interacting galaxy pair II Zw. This pair has been thoroughly studied in the X-rays (Inami et al. 2010; Iwasawa et al. 2011), and the XMM-Newton source is most likely a blend of sources A, C, and D from Goldader et al. (1997). Source SRCID = 201389514010019 (Fig. 4) clearly points towards a background galaxy, even though it is listed as belonging to galaxy NGC 2528. According to NED, this galaxy is at z = 0.13. An additional eight sources were revealed as spurious detections around an oversaturated region on the X-ray image. All of the mentioned sources were flagged accordingly.

Fig. 3. Pan-STARRS1 optical image of the location of source SRCID = 200559905010001. The X-ray source (not shown here) is a blend of the two cores of the interacting pair NGC 5256.

Fig. 4. Pan-STARRS1 optical image of the location of source SRCID = 201389514010019 in galaxy NGC 2528 (large spiral on the right). The source lies clearly on top of a background galaxy at z = 0.13.

Unfortunately, 19 of the sources did not have Pan-STARRS1 images available, so we inspected their positions on NED to look for possible counterparts. This is particularly relevant for objects with an ecliptic latitude below −30°, as this area is not entirely covered by the Pan-STARRS1 survey. Fifteen additional sources showed Pan-STARRS1 counterparts of unclear nature, so we double-checked them with NED. The results were the following.

Sources SRCID = 200936502010001 and SRCID = 200029702010002 were found to be the central sources of their host galaxies, NGC 5128 and UGC 01841, due to a discrepancy in the value for the central coordinates with HECATE. Source SRCID = 201241101010001 was identified as the QSO MR 0205. Finally, four remaining sources were also identified with distant background galaxies.

There is another side to the manual inspection of sources. As ULXs are more commonly associated with star-forming regions (Kaaret et al. 2017), some genuine ULXs are located within or next to optically bright H II regions of spiral and irregular galaxies. This caused confusion in the filtering step involving Pan-STARRS1, and otherwise good candidates were flagged as interlopers. Therefore, we performed manual inspection for 128 Pan-STARRS1 associations that held SC_SUM_FLAG ≤ 1, SC_Luminosity +SC_LuminosityErr ≥ 10³⁹ erg s⁻¹ and SC_Luminosity > SC_LuminosityErr with the intention of recovering genuine ULX candidates.

A total of 72 sources were matched to bright nebulous features in spiral galaxies (see Fig. 5 for an example). A NED-based inspection of ten further sources also suggested an H II region or having only X-ray counterparts, and therefore they were re-flagged as clean. Source SRCID = 201109302010027 in particular was not matched to any potential contaminant in NED, but it showed two possible X-ray counterparts. It was kept as an interloper due to its possibly blended nature.

Fig. 5. Pan-STARRS1 optical image of the location of source SRCID = 205562803010007 (green marker) in galaxy NGC 2903. The red marker indicates the position of the Pan-STARRS1 counterpart identified by the pipeline. According to NED, the optical counterpart constitutes an H II region.

Another set of 11 sources were identified with infrared sources, also in NED, while four more sources where found to coincide with background galaxies. These sources were kept in the group of interlopers.

Finally, we searched for possible extragalactic objects in Gaia DR2. Shu et al. (2019) identify 3 175 537 AGN candidates from 641 266 363 Gaia DR2 sources, a fraction of 0.5%, which suggested further possible undetected extragalactic interlopers in our catalogue. Therefore, we performed a final cross-match of the remaining unmatched sources with Gaia DR2, this time applying condition (2) to decide whether to flag the object as a possible interloper. Twenty-three ULX were matched with a Gaia DR2 source. We inspected all of them manually in the NED database and found that 14 of them did not have a known counterpart, so they were re-flagged as clean. In contrast, seven were found to have infrared counterparts, SRCID = 203049401010029 was matched to a UV source, and SRCID = 207843705010011 was found to be the central source of NGC 7632, and therefore they were accordingly confirmed as interlopers.

After the manual inspections, the final version of the catalogue yielded 1452 detections of 779 ULX candidates in 517 galaxies, and 163 detections of 94 bright ULX candidates of quality in 94 galaxies.

3.4. Luminosity and complete sub-samples

As low-luminosity sources become increasingly more difficult to detect with increasing distance, our catalogue presents a bias towards brighter sources at large distances, hampering any potential study of the ULX population properties. Following the methodology of Earnshaw et al. (2019), we created the luminosity sub-samples to ensure the completeness of the set of sources above a luminosity threshold L_min inside of a radius D_max.

We established a minimum flux for a source to be detected F_min from which we compute the maximum distances D_max at which sources with luminosity L_min can be seen. Then we selected all objects with D < D_max in each case to make sure that every source with L > L_min is accounted for, mitigating this way the bias towards brighter sources in the L > L_min regime. From these objects we then selected all quality sources and ULX candidates whose host galaxy’s isophotal dimensions are known (MATCH_FLAG = 0) and that lie at least 25° away from the Galactic plane. This last condition ensures minimal photoelectric absorption from the Milky Way.

The most controversial aspect of this method is the choice of F_min since there is no hard detection threshold for the XMM-Newton observatory as it depends on exposure time and background intensity. We observed in Fig. 6 that the distribution of all detection fluxes in 4XMM-DR9 peaks early at F_X ≈ 10⁻¹⁴ erg cm⁻² s⁻¹, with 77% of the detections at higher fluxes, so to work around the problem we adopted this value as the baseline for XMM-Newton’s sensitivity.

Fig. 6. Distribution of detection fluxes within 4XMM-DR9, regardless of their nature. Most detections (77%) have been detected with fluxes higher than FX ≈ 10−14 erg cm−2 s−1.

With this in hand, we built the luminosity sub-samples 38ss, 39ss, and 5x40ss, with corresponding L_min = 10³⁸, 10³⁹, and 5 × 10⁴⁰ erg s⁻¹ and D_max ≈ 9, 29, and 204 Mpc. We then selected all ULX candidates with n_Galaxies = 1 belonging to 39ss to build the complete ULX sub-sample (cULXss), and we did the same with the bright ULX candidates and the 5x40ss, building the complete bright ULX sub-sample (cbULXss).

The effectiveness of this method is illustrated in Fig. 7, where the bias towards brighter sources has been significantly downplayed by excluding sources with D > D_max from the sub-sample. As shown in Table 1, these two sub-samples contain a reduced amount of 292 and 69 sources. The sub-samples are also helpful in an additional way. Due to the clustering nature of ULXs, objects at D > 10 Mpc are prone to blending due to the limited angular resolution XMM-Newton. The cut at D_max = 29 Mpc also helps in the mitigation of this bias.

Fig. 7. Distance-luminosity dispersion of all quality sources within 38ss, 39ss, and 5x40ss, and within cULXss and cbULXss, with an imposed sensitivity of 10−14 erg cm−2 s−1. The distance cuts correspond to Dmax ≈ 9, 29 and 204 Mpc.

4. Results

4.1. Contents of the catalogue

Tables 1 and 2 summarise the content of our final catalogue. The catalogue contains 50 446 entries corresponding to 23 262 sources, of which 23 120 are associated with a single galaxy (n_Galaxies = 1). X-ray sources with Gaia DR2 associations form the largest group of contaminants, containing 41.8% of the total number of sources. This is due to the 8437 matches with objects of high optical–to–X-ray flux ratio, 8270 of which are located in the Magellanic Clouds and M31, indicating a possible stellar nature. Nonetheless, none of them fall within the ULX X-ray luminosity regime so no further regard is paid to them. Central objects follow, constituting 9.7% of the catalogue. They are followed by 8.4% of objects directly identified as AGNs or QSOs in SDSS_DR14, VeronQSO and SIMBAD. A further 7.1% have been identified as stars in Gaia DR2, Tycho2, or SIMBAD. Only 4.9% of the objects have been selected as potential interlopers from the Pan-STARRS1 survey, and a residual fraction is constituted by supernovae, infrared sources, background galaxies, and a single ultraviolet source.

A total of 3274 objects (14.1%) have been detected at least once within the ULX regime (L_X + ΔL_X > 10³⁹ erg s⁻¹). Of these objects, 2208 consist of central sources, only 779 (3.3%) qualify as ULX candidates according to the thresholds established in Sects. 3.2 and 3.3, and 287 are interlopers of various kind. Of the ULX candidates, 516 have been detected at least once within the ULX regime with high certainty (L_X − ΔL_X > 10³⁹ erg s⁻¹), and 666 with high likelihood (L_X > 10³⁹ erg s⁻¹); 761 have at least one detection well above the NS Eddington limit (L_X − ΔL_X > 1.8 × 10³⁸ erg s⁻¹). A sub-set of 94 objects also qualify as bright ULX candidates (L_X + ΔL_X > 5 × 10⁴⁰ erg s⁻¹), and all of them are the only candidate in their host galaxy. Among some of the candidates we find well-known objects such as the NS-ULXs NGC 7793 P13 and NGC 5907 ULX-1; the IMBH candidate NGC 2276-3c; and M51-ULS-1, host of the extragalactic exoplanet candidate M51-ULS-1b (Di Stefano et al. 2021). Other known objects such as the IMBH candidates M82 X-1 and ESO 243-49 HLX-1 also appear in the catalogue, but with SC_SUM_FLAG > 1.

Only sources and candidates with n_Galaxies = 1 are considered during the remainder of this section due to their more reliable luminosities and host galaxy associations.

4.2. ULX distribution and abundances

In Fig. 8 it can be seen that all ULX candidates are found at D > 1 Mpc, while bright candidates are found at D ≳ 20 Mpc. Figure 8 also illustrates how at D ≳ 100 Mpc the bright candidates dominate almost completely due to the bias in favour of brighter sources, as discussed in Sect. 3.4 and possibly source blending. This is consistent with the information in Table 1, which shows that median distances of ULX and bright ULX candidates are 28.4 Mpc and 127 Mpc. Sources of quality cluster in the nearest galaxies for the same reason. They constitute a population of low-luminosity X-ray sources that dominate the serendipitous X-ray sky in the local group of galaxies. For instance, M31, M33, and the Magellanic Clouds alone already concentrate 17% of all sources of quality. The X-ray content of the local group is not directly relevant to the ULX population as it is mostly composed of typical X-ray binaries and supernovae remnants and it has already been thoroughly studied (e.g. Sturm et al. 2013; Pietsch 2008). Therefore, for the remainder of this paper we focus most of our attention to sources hosted by galaxies at D > 1 Mpc.

Fig. 8. Distance distribution of all sources of quality, ULX candidates, and bright ULX candidates in our catalogue.

In Table 3 it can be seen that objects from the cULXss are most commonly associated with early spiral galaxies (0  ≤  T_H < 5) and late spiral galaxies (5 ≤ T_H < 9). However, while 68.1% of the cULXss candidates are found in spiral galaxies, only 50.7% of the objects in the cbULXss dataset follows the same trend. In Fig. 9 it can be seen how the overwhelming majority of cULXss objects concentrate in spiral galaxies in the range of 0 ≤ T_H < 10, with a spike in elliptical galaxies (T_H < −3). For the cbULXss the distribution is similar, but the spike in elliptical galaxies is more prominent. In both samples lenticular galaxies (−3 ≤ T_H < 0) and irregular galaxies (T_H ≥ 9) are disfavoured. Figure 10 also shows that while ULX candidates of quality are more frequently found in late-type galaxies (T_H ≥ 0) than in early-type galaxies (T_H < 0), this distinction becomes less clear for the bright ULX candidates.

Fig. 9. Hubble type distributions of cULXss and cbULXss. In both cases there is a spike in ULX abundances in the earliest galaxies.

Fig. 10. Luminosity distribution of the sources of quality and ULX candidates in our catalogue, divided between late-type galaxies (LTG) and early-type galaxies (ETG). A distance cut of D > 1 Mpc has been imposed to avoid contamination from the Local Group, heavily dominated by low-luminosity sources.

In addition to looking at how ULX candidates distribute themselves among galaxy morphological types, it is also interesting to investigate the fraction of galaxies that contain any ULX and general ULX frequencies. For this, we built the set of all galaxies in the field of view that would be included in the cULXss or cbULXss if hosting a ULX was not a requirement, with the extra condition of having R1 > 3″ to exclude galaxies that can only contain central sources. In essence, this consists of the baseline of galaxies that would belong to the complete sub-samples if they hosted at least one ULX candidate, including the ones that actually do. This is shown in Table 3 as the row of potential galaxies, and it is used implicitly in Table 4 to compute the mean SFR and M_* values. In Table 3 it is seen that 21% of the potential cULXss galaxies contain at least one ULX candidate. For elliptical galaxies this gets reduced to 12%, for lenticular galaxies it is 21%, for early spirals it grows to 49%, and for late spirals to 37%. The general ULX rate is ∼0.4 ULX/galaxy, being the highest at ∼0.98 ULX/galaxy in early spiral galaxies. Figures 11 and 12 show these results in more detail for the different galaxy morphological sub-types, with ULX frequencies peaking at ∼1.8 ULX/galaxy for T_H ∼ 3 galaxies, of which ∼70% host ULX candidates. These results are mostly in agreement with earlier works, which either find a larger fraction of ULXs in late-type galaxies (Swartz et al. 2011) or find higher abundances in late-type galaxies (Earnshaw et al. 2019; Kovlakas et al. 2020). However, Kovlakas et al. (2020) find that the galaxies with a higher chance of hosting ULXs are T_H ∼ 5 spirals instead.

Fig. 11. Distribution of ULX frequencies in cULXss and cbULXss across the morphological sub-types of host galaxies. Frequencies from cbULXss are multiplied by 100 to aid visibility.

Fig. 12. Distribution of fraction of host galaxies in cULXss and cbULXss across their morphological sub-types. Values from cbULXss are shown in % to aid visibility.

To find an explanation for these distributions, we need to explore the properties of host galaxies, and more particularly the relationship between ULX abundances, host SFR, and M_*. In spiral galaxies with high specific-SFR (sSFR = SFR/M_*) values and that are dominated by young stellar populations, the ULX population is typically associated with HMXB evolution timescales of τ ∼ 100 Myr (Wiktorowicz et al. 2017), the rates of which are expected to scale with the SFR. In Table 4 we show that spiral galaxies are typically those with the highest SFRs, as given in HECATE. Early spirals have higher SFR absolute values than late spirals (∼1.99 vs. ∼0.80 M_⊙ yr⁻¹), but as expected they also have lower sSFR values. Even though ULXs are more abundant in early spirals, the ULX frequency per SFR is higher in late spirals (∼0.49 vs. ∼1.01/M_⊙ yr⁻¹). The same holds for the number of ULXs per unit mass, and the whole picture is exaggerated for irregular galaxies. Figure 13 also shows the distribution of individual ULX/SFR galaxy rates for early and late spiral galaxies in the cULXss that contain available SFR information. The distributions peak at ∼0.65 and ∼0.55 ULX/M_⊙ yr⁻¹ each. Wiktorowicz et al. (2017) predict the existence of 400 ULXs in a galaxy with solar metallicity, M_* = 6 × 10¹⁰ M_⊙, and SFR = 600 M_⊙ yr⁻¹ for a period of 100 Myr, which implies 0.67 ULX/M_⊙ yr⁻¹, staying well within the range of values provided by our sample.

Fig. 13. Histogram of the number of ULX per star formation rate for early and late spiral galaxies from cULXss. The large values on the right are due to galaxies with low SFR that host a single ULX, and the most extreme cases have been cut out of the picture.

For elliptical and lenticular galaxies, different considerations need to be applied. As expressed by Wiktorowicz et al. (2017), there is a ULX subpopulation constituted by LMXBs that reach Roche-lobe overflow at τ ∼ 1 Gyr. This evolutionary path explains the presence of ULXs in galaxies with very low star-forming activity. In this case ULX frequencies are not expected to depend significantly on their current SFR values.

On the other hand, objects in the cbULXss provide a slightly different picture. All of them are located only beyond the D ≳ 20 Mpc mark and have extremely low abundances, as shown in Table 3. Additionally, all of them are the only bright ULX candidate in their host galaxy. From the Wiktorowicz et al. (2017) predictions for the rate of objects with L_X > 10⁴⁰ erg s⁻¹, we estimate a typical number of 0.023 bright ULX/M_⊙ yr⁻¹ in star-forming galaxies. However, we only measure 0.0037 bright ULX/M_⊙ yr⁻¹ and 0.0055 bright ULX/M_⊙ yr⁻¹ for early and late spiral galaxies, which is an order of magnitude lower. This can be partially explained by considering that our luminosity cut of L_X > 5 × 10⁴⁰ erg s⁻¹ is slightly higher, and that the ULX luminosity distribution in general decays exponentially in most studies. Nonetheless, a total of only 36 objects are involved in these estimations. A larger sample will need to be used in the future to extract more reliable conclusions.

4.3. Spectral properties of the catalogue

The X-ray spectral properties of the identified sources can provide additional information regarding their nature (e.g. Earnshaw et al. 2019; Walton et al. 2011). Therefore, we proceed with the inspection of the spectral properties of the objects within our catalogue. We use the hardness ratios, defined as

$$\begin{array}{c}\hfill {\mathit{HR}}_i=\frac{R_{i+1}-R_i}{R_{i+1}+R_i},\end{array}$$(3)

where R_i and R_i + 1 are the count rates in the consecutive energy bands i and i + 1. With this quantity, spectral hardness can be compared between different sources, with larger values indicating harder sources. The 4XMM-DR9 computes the source HR_i values from the count rates in the 0.2–0.5, 0.5–1.0, 1.0–2.0, 2.0–4.5, and 4.5–12.0 keV energy bands, from the EPIC photon counts, HR₁ and HR₂ being the most useful thanks to larger photon counts in lower energy bands.

From Table 5 it is apparent how the contaminant populations have distinct average spectral properties from each other. Supernovae, AGNs, and QSOs tend to be the objects with the hardest spectra. Comparing the interloper hardness ratios with the ULX hardness ratios shown in Table 6 also shows that, as expected, the typical spectra of the ULX population is distinct from those of stellar objects, supernovae, infrared sources, and nuclear sources, while bearing the closest spectral resemblance with the AGN and QSO populations. These properties are made more apparent by looking at Fig. 14, where it can be seen that AGNs and QSOs cluster around slightly harder spectra than the cULXss and cbULXss samples, and that there is even less overlap of the ULX population with the supernovae and the stellar populations. Figure 15 shows that the ULX population has the greatest resemblance to the population of X-ray objects just below the ULX luminosity threshold (10³⁸ < L_X < 10³⁹ erg s⁻¹).

Fig. 14. Hardness ratio dispersion of cULXss and cbULXss against the abundance contours of some populations of interlopers. AGNs and QSOs include both confirmed and candidate objects. Only hardness ratios with uncertainties lower than 0.2 are accepted for the sake of reliability.

Fig. 15. Hardness ratio dispersion of cULXss and cbULXss against the abundance contours of lower luminosity objects present in 38ss. Only hardness ratios with uncertainties smaller than 0.2 are accepted for the sake of reliability.

There are also noteworthy spectral differences within the ULX sample itself. Table 6 shows a division between ULX spectral properties in relation to their host galaxies for the cULXss. ULXs in lenticular and elliptical galaxies present the softest spectra in all cases for faint and bright ULXs. The hardest spectra are found in late spiral galaxies. This notion is clearer in Fig. 16 where the hardness ratio distributions of ULXs in late-type and early-type galaxies are seen to cluster around slightly different values. However, Fig. 17 illustrates that the data comes with a large dispersion, despite the increasing trends in the values of HR₁ and HR₂ with the morphological sub-class of each galaxy up until the irregular galaxies. This information points towards photoelectric absorption being higher in spiral galaxies than in elliptical and irregular ones, in agreement with what has been reported in previous XMM-Newton-based works (Walton et al. 2011; Earnshaw et al. 2019), but also towards unaccounted factors being more determinant for the individual properties of ULXs.

Fig. 16. Hardness ratio vs. luminosity dispersion of objects in cULXss (squares) and cbULXss (diamonds), distinguishing between objects hosted in late-type galaxies (TH ≥ 0) and early-type galaxies (TH < 0). Density contours have been drawn to help the eye. Only hardness ratios with uncertainties lower than 0.2 are accepted for the sake of reliability.

Fig. 17. HR1 and HR2 across morphological sub-types of host galaxies from cULXss sources. The error bars show the 1σ dispersion of the values, and the data points have been slightly moved to the left (HR1) and to the right (HR2) for ease of visibility.

Looking at Table 6, a division between the cULXss and cbULXss is seen. On average, low-luminosity ULXs present harder spectra than the bright ones. It is unlikely that the host galaxies of bright ULX candidates contain less neutral gas in general, and perhaps their generally softer spectra is an indication of bright candidates having a different physical nature than the rest of the candidates. Their emission could also result from the blending of lower-luminosity sources, as suggested in Sect. 5.1.3, perhaps with the diffuse emission of the stellar population in the host galaxies. Within the cbULXss itself, it also appears that candidates in late-type galaxies cluster around slightly harder spectra than those hosted by early-type galaxies, as shown in Fig. 16.

4.4. ULX variability

The differences between ULX candidates also extend to their variability properties. We find 147 cULXss counterparts within three times their positional uncertainty in the 4XMM-DR9 objects with NCONTRIB  ≥  2. This number goes down to 31 for the cbULXss sample. For all these sources, we have at our disposal variability information drawn from their multiple overlapping observations (Traulsen et al. 2019).

We pay attention to VARPROB, which is the probability of a source not showing variability between observations, in essence a long-term variability measure. For our analysis we put the ULX candidates into five bins of variability probability. We consider objects with VAR_PROB > 1σ, where 1σ = 0.6827, to be most likely constant, and we collect them in the 1σ constant bin. Thereafter, we collect the objects with 1 − 1σ > VAR_PROB > 1 − 2σ, 1 − 2σ > VAR_PROB > 1 − 3σ, and VAR_PROB < 1 − 3σ into the 1σ variable, 2σ variable, and 3σ variable bins, being 2σ = 0.9545 and 3σ = 0.9973, while the remainder are left in the uncertain variability bin. We also register the FRATIO values of the matched objects, which corresponds to the ratio between the highest and the lowest detected flux from a source, and compute the average value not only for every population time, but also for every VARPROB bin.

The results are shown in Table 7, where it can be seen that ULX candidates from the cULXss hosted by late spiral galaxies are those with the highest chances of showing variability between different observations as they present the lowest VARPROB average value, with 61.2% of them belonging to the 3σ variable bin, and only 4.1% being in the 1σ constant bin. In contrast, ULX candidates in lenticular galaxies present the highest VARPROB average value, falling to 20% of them in the 1σ constant bin. However, the ULX population with the lowest fraction of candidates in the 3σ variable bin is that of candidates hosted by elliptical galaxies. ULX candidates hosted by early spiral galaxies present intermediate variability properties.

Regarding the measured variation of fluxes between different observations, Table 7 also shows that for ULX candidates hosted by late spiral galaxies the FRATIO values span two orders of magnitude on average, for the entire set and for those in the 3σ variable bin. ULX candidates in lenticular galaxies show the lowest FRATIO values, being followed by those in elliptical galaxies. In these cases the inter-observation luminosity variability does not go beyond one order of magnitude. ULXs in early spiral galaxies repeat their role as a bridge between late spiral-hosted and lenticular-hosted candidates. Not surprisingly, candidates in the 1σ constant bin present average values of FRATIO of the order of 1 regardless of their host galaxy, and while some exceptions are seen, in general FRATIO increases with decreasing VARPROB.

Some trends are similar for the cbULXss. As shown in the lower half of Table 7, bright candidates hosted by late spiral galaxies are those with the highest fraction of ULX candidates in the 3σ variable bin, and they once again have the highest FRATIO average values. However, they are typically less variable than their cULXss counterparts, as seen from the higher VARPROB and lower FRATIO average values for the entire set of candidates regardless of the morphological type of the host galaxy. The differences between early spiral-hosted, elliptical-hosted, and lenticular-hosted bright ULXs are also blurrier, most likely due to the small size of the cbULXss.

It should be noted that despite these general trends, the dispersion in the VARPROB and FRATIO values is too large to find a reliable dependence on other quantities. It is tempting to propose a relationship between variability and the star formation activity in galaxies, as coincidentally late spiral galaxies also present the highest sSFR (excluding irregular galaxies from the analysis), while both elliptical and lenticular galaxies present the lowest sSFR values. Despite these hints of a correlation between variability and galaxy morphology, a further investigation of the variability properties in the context of the sSFR of the galaxies does not show a clear correlation.

Sutton et al. (2013) make a thorough study of the variability of several ULXs, stating that an increase in X-ray brightness usually goes hand in hand with a softening of the source due to accretion winds narrowing their funnel shape and obscuring the inner regions of the disc. We do not find such relationship for many cases either. For instance, NGC 5204 X-1, one of their notable examples, which appears in our cULXss with ⟨L_X⟩=(6.0 ± 0.6)×10³⁹ erg s⁻¹, presents the opposite behaviour within XMM-Newton observations, ranging from HR₁ ≈ 0.33 and HR₂ ≈ −0.02 at L_X ≈ 5 × 10³⁹ erg s⁻¹ to HR₁ ≈ 0.42 and HR₂ ≈ 0.02 at L_X ≈ 8 × 10³⁹ erg s⁻¹. Nonetheless, the methodology and sample data in both of our studies diverge greatly, and in particular our methodology is rather limited when making statements about individual sources.

5. Discussion

5.1. Limitations of our study

5.1.1. Limitations of the catalogue

With 779 identified candidates our catalogue provides the largest ULX sample built to date, only followed by the 629 candidates provided by the Chandra-based catalogue from Kovlakas et al. (2020), the 470 candidates from Walton et al. (2011), and the 394 candidates from Earnshaw et al. (2019). Of our 779 candidates, 94 are bright candidates with at least one detection with L_X > 5 × 10⁴⁰ erg s⁻¹. However, to avoid biases towards bright sources and blending, we have only considered the 292 candidates from the cULXss and 62 from the cbULXss for the population study. While this ensures the veracity of our results for the whole cULXss, the small size of the cbULXss hampers our study severely. In addition, as suggested in Sect. 5.2.2, the bright population is more prone to contamination by unidentified interlopers and source blending, hampering the study even further.

Aside from the manual inspection of sources during the filtering process, our study is also strictly limited to collective properties of our samples. Both the interloper and the ULX populations are heterogeneous samples at the astrophysical level, and while some trends are identified, individual objects can differ greatly. For instance, we identify in Sect. 4.3 that ULX candidates in late-type galaxies tend to present harder spectra than those hosted by early-type galaxies. However, it is postulated that ULXs present hard or soft states depending on the viewing angle and the accretion rate (Sutton et al. 2013), and therefore it is most likely that the spectrum of an individual ULX is determined by these factors rather than the morphological type of the host galaxy.

From a population study such as the one presented in this work we cannot draw any conclusions on the nature of individual sources. However, we can make connections between their collective properties (spectra, variability, frequency) and the nature of the stellar populations in their host galaxy. For example we see a trend for ULXs in late-type spirals to exhibit harder colours, indicative of more significant photoelectric absorption. However, we do not see any trend of variability indicators with the sSFR. The variability and spectral trends can be better addressed with systematic observations of individual objects (e.g. Sutton et al. 2013; Kaaret et al. 2017).

5.1.2. Selection biases

The main limitation of 4XMM-DR9 is that it has been constructed from observations of galaxies or clusters that were already of interest to astronomers. As such, the characteristics of our ULX sample may deviate from those of the true ULX population. In Fig. 18 it can be seen how galaxies included in 4XMM-DR9 present a slight bias towards earlier types with respect to the full content of HECATE. This bias only gets amplified further along in the construction of our catalogue.

Fig. 18. Hubble type distributions of galaxies with listed values in HECATE, within 4XMM-DR9 XMM-Newton observations, and within our catalogue. Portions have been normalized to sum to 1 in all of the three cases due to the large differences in size between the sets.

Choosing to use 4XMM-DR9 as the only X-ray catalogue of reference also limits the size of the X-ray sample available to us as it covers only 2.7% of the sky area and 8.6% of all the galaxies present in HECATE. Fortunately, other recent works such as Kovlakas et al. (2020) already take into account non-overlapping parts of the sky, even though they suffer the same bias as our own individually.

The most straightforward way to tackle these issues is by studying the ULX content of blind all-sky surveys. A good opportunity for this is presented by the still ongoing eROSITA All Sky Survey (eRASS, Merloni et al. 2012; Predehl et al. 2021). By observing the entirety of the sky instead of selected patches, and accounting for the average eROSITA sensitivity at the end of the eRASS, we assess that around ∼1000 new ULX candidates will be discovered in HECATE galaxies of at least 20″ in size, eROSITA’s FWHM angular resolution.

Further biases in our catalogue stem from our own design choices, particularly from focusing our studies on source parameters. Occasionally, quality flags and some measures of source parameters pick the worst possible value among all the values assigned to individual detections. This means that if a legitimate source has one detection with SUM_FLAG > 1, this bad value percolates into SCSUMFLAG and renders the source ineligible as a ULX candidate under our criteria. Ignoring the SCSUMFLAG value, the total number of ULX candidates of quality increases to 847. This choice was made because a bad detection often leads to unreliable source parameters, but the true extent of this effect would need to be checked individually for every source. Nonetheless, all of these sources are still available to the user of the catalogue.

Another bias we introduced was the selection of point-like sources with high detection likelihood, which excluded 1632 sources within HECATE’s isophotal ellipses, in contrast to the 23 262 accepted ones in the final catalogue. However, we think that the exclusion of extended sources is justified as ULXs are point-like sources.

5.1.3. XMM-Newton angular resolution

The angular resolution of XMM-Newton is a limiting factor due to the clustered nature of star-forming regions, where ULXs tend to be found (Anastasopoulou et al. 2016; Kaaret et al. 2017). H II regions, both Galactic and extragalactic, present sizes of the order of ∼10 pc, at most ∼100 pc (Tremblin et al. 2014; Hunt & Hirashita 2009, plus references therein). This implies that we should expect contamination from neighbouring X-ray sources in the ULX candidate parameters at D ≳ 10 Mpc, or even blending of sources into spurious ULX candidates. The cULXss deals with this issue by establishing the cut at D_max ≈ 29 Mpc, so we expect little contamination to our ULX population study. However, it is very telling how at D_max ≳ 100 Mpc the bright ULX candidates completely overtake over the normal ULX candidates. It was explicitly revealed in Sect. 3.3 that some of the manually investigated objects suffer from this problem.

Other catalogues built with higher resolution telescopes such as Chandra are able to partially alleviate this problem, but the issue reappears at D ≳ 40 Mpc as stated in Kovlakas et al. (2020). Individual follow up of interesting sources is once again essential to untangle the nature of these objects.

5.2. Discussion of contents

5.2.1. Candidates at the ULX luminosity threshold

As expressed in Sect. 3.2, our criterion for the selection of ULX candidates consists in whether a source has at least one detection whose luminosity is above this established threshold within its 1σ uncertainty. This rather lax condition allows for typically low-luminosity sources to be considered ULX candidates. However, it manages to take into account that some variable objects may behave like ULXs with intermittence and also like NS-ULXs close to the Eddington limit. While in this work we have used the common definition for the ULX luminosity threshold in the literature (L_X > 10³⁹ erg s⁻¹, Kaaret et al. 2017), we note that it is but an approximation, as the Eddington limit lies at 1.8 × 10³⁸ erg s⁻¹ for a 1.4 M_⊙ NS, and at 1.3 × 10³⁹ erg s⁻¹ for a 10 M_⊙ StBH.

As explained in Sect. 4.1, only 516 of our ULX candidates have at least one detection with luminosity above the ULX threshold with more than 1σ significance (L_X − ΔL_X > 10³⁹ erg s⁻¹) and 150 others are also likely to be ULXs within a 1σ uncertainty (L_X > 10³⁹ erg s⁻¹). This leaves the catalogue with 106 objects that qualify as ULXs only within the uncertainty of their brightest detection (L_X + ΔL_X > 10³⁹ erg s⁻¹, L_X < 10³⁹ erg s⁻¹), and are therefore unlikely to be ULXs according to the common definition. However, from this sub-set, a vast majority of 95 have been observed with luminosity well above that of the Eddington limit for a 1.4 M_⊙ accreting NS (L_X − ΔL_X > 1.8 × 10³⁹ erg s⁻¹), and therefore still deserve a place in this catalogue in consideration of a more comprehensive ULX luminosity threshold.

5.2.2. The most luminous ULXs

Our catalogue contains 94 ULX candidates of quality with at least one detection with luminosity L_X > 5 × 10⁴⁰ erg s⁻¹. This sample extends to 109 sources in total if we ignore the SCSUMFLAG column. We showed in Sects. 4.3 and 4.4 that, while the cbULXss follows the same trends as the cULXss, in the sense that candidates hosted by late-type galaxies tend to present harder spectra and higher variability, it is also true that the bright candidates generally have softer spectra and more subdued variability than their cULXss cousins.

A tempting way to explain this difference is to invoke different physical origins of the sources. Hyperluminous X-ray sources (HLXs) are ULXs with L_X > 10⁴¹ erg s⁻¹, and although the community has recently shifted towards explaining them as typical X-ray binaries at the high end of the luminosity distribution due to very high mass transfer rates (Bachetti 2016; Wiktorowicz et al. 2017), they are still the best objects where to look for high mass BHs or even IMBHs (Kaaret et al. 2017). This search is further motivated by the first direct confirmation ever of an IMBH in the gravitational wave event GW190521 (Abbott et al. 2020a,b), which is proof that black holes of these masses do exist and are to be considered a potential explanation for the brightest ULXs. A handful of ULXs such as ESO 243-49 HLX-1 (Kong et al. 2007; Servillat et al. 2011; Pasham et al. 2014) or NGC 2276-3c (Mezcua et al. 2015) constitute promising candidates of being such kind of objects, three of which are included in our catalogue.

In our catalogue, ESO 243-49 HLX-1 corresponds to source SRCID = 202045402010003, and has an average luminosity of ⟨L_X⟩=(1.7 ± 0.4)×10⁴¹ erg s⁻¹ and a maximum luminosity of L_X = (7 ± 2)×10⁴¹ erg s⁻¹. Unfortunately, one of its six detections has SUM_FLAG = 3, so it was not included in our list of ULX candidates of quality. M82 X-1 has also been recovered as source SRCID = 201122902010001, but this source is actually a blend of M82 X-1 and the NS-ULX M82 X-2 (Bachetti et al. 2014), as they are only 0.52″ apart and are therefore indistinguishable by XMM-Newton. Furthermore, this source also has one detection with SUM_FLAG = 3. NGC 2276-3c (SRCID = 200223402010005) is the only one of these three that appears as a bright ULX of quality. However, it is also constituted by three blended sources in XMM-Newton, one of them being a genuine IMBH candidate (Mezcua et al. 2015).

Nonetheless, it is known that not all bright ULX candidates need to explained with the presence of accreting massive BHs. For example, another object classified as a bright ULX candidate in our catalogue is NGC 5907 ULX-1, a NS-ULX that has been observed with a luminosity of L_X > 10⁴¹ erg s⁻¹ (Israel et al. 2017).

This implies that there must be other explanations to their different properties. Many times the reason may not even be astrophysical. As seen in Fig. 8, almost all of them are located at D ≳ 20 Mpc, which makes them prone to be the result of blending of several low-luminosity sources. The fact that they are always found alone in their host galaxies is also a strong indication of source blending. This is reinforced by the explicit detection of confused sources among manually inspected objects in Sect. 3.3. It is also possible that there is significant contribution of unidentified nuclear sources to our cbULXss. Table 8 clearly shows that indeed a much larger share of background interlopers was identified for the bright ULX candidates than in the general sample, implying that the unidentified contaminants may also constitute a larger fraction.

Nonetheless, from the aforementioned objects, only NGC 2276-3c has ⟨L_X⟩> 10⁴¹ erg s⁻¹. In fact, only 25 belong to the HLX class within their uncertainty, 33 if SCSUMFLAG is ignored. From the ones with SC_SUM_FLAG ≤ 1, only three match with the brightest objects from Earnshaw et al. (2019), the study closest to ours in terms of sampling and methodology. NGC 4077 and IC 4596 are included in Earnshaw et al. (2019) as potential IMBH candidates. IC 4320 also appears in Earnshaw et al. (2019) and Walton et al. (2011), and has already been investigated by Sutton et al. (2012). The remaining 22 do not seem to appear in previous literature as far as we are aware. If extended to other objects qualified as contaminants, we also find a match for the HLXs in IC 4252 and UGC 6697. The former is considered a central source in our catalogue, while the latter holds n_Galaxies = 2, SC_SUM_FLAG = 2 and has been identified with a Pan-STARRS1 counterpart, and therefore they do not meet the requirements to be considered sources of quality. Others also included in Earnshaw et al. (2019), such as UGC 1934 and NGC 2276 belong to our bright ULX candidates, but do not qualify as HLXs.

5.2.3. Effectiveness of the filtering pipeline

Table 5 shows that a total number of 1943 sources are matched with AGNs, QSOs, or background galaxies, and 1652 with stellar objects. Additionally, in Table 8 we show that 74 and 27 objects in each of these groups would have been classified as ULX candidates if it were not for their association. In total, they would constitute 11.5% of the ULX candidates of quality. It is therefore natural to ask whether the exclusion of these objects is accurate or, to the contrary, sufficient at all.

In the case of stellar objects, we applied extra conditions in addition to positional coincidences based on the X-ray flux and optical magnitude of the sources. This way, an X-ray source was only classified as a stellar contaminant if the stellar counterpart was bright enough to explain the X-ray emission or, at least, to contaminate it significantly. In any case, most of the objects classified as stellar contaminants concentrate in nearby galaxies devoid of ULX candidates. Nonetheless, the best way to attest the efficacy of the pipeline is by comparing the spectral values of ULXs to those of stellar contaminants, shown in Tables 6 and 5, and seeing that the set of cULXss and the cbULXss have distinguishable spectral properties with respect to the bulk of stellar objects.

On the other hand, the direct exclusion of any X-ray source coincident with a known QSO, AGN, or SDSS-DR14 source may raise more concerns. To determine whether this method has overshot or has been insufficient, we can compare the number of expected background sources to the identified ones. In Sect. 3.1 it is mentioned that around 670 background contaminants are expected in the set of galaxies at D > 1 Mpc from the log N − log S presented in Mateos et al. (2008), which is equal to ∼12% of the sources at that distance if unflagged. In our catalogue, 538 sources at D > 1 Mpc (∼10%) were actually classified as possible QSOs or AGNs. Therefore, it can be concluded that our filtering pipeline does a good job of identifying most of the background contaminants, but also that ∼2% of our sources may still be unidentified background sources, including ∼15 ULX candidates.

5.3. Comparison with other works

5.3.1. Previous serendipitous XMM-Newton ULX catalogues

Earnshaw et al. (2019) followed much of the methodology presented in Walton et al. (2011). As Walton et al. (2011) built their sample from 2XMM-DR1¹⁰ and Earnshaw et al. (2019) from 3XMM-DR4, their work can be seen as an update on the ULX content of the XMM-Newton survey. Likewise, our work can also be seen as an update on the ULX content of XMM-Newton with respect to Earnshaw et al. (2019). Therefore, we expect to recover most of the ULX candidates presented their work and Walton et al. (2011).

A source positional cross-match recovers 1008 sources from the total of the 1314 listed in Earnshaw et al. (2019). As Earnshaw et al. (2019) do not preserve sources classified as interlopers, 80% of the matches correspond to clean sources on our side. The remaining sources consist of 11 central sources, 30 stars, 1 QSO, 16 SIMBAD objects, 21 objects from Gaia DR2, 32 extragalactic objects from Pan-STARRS1, and 1 background galaxy manually identified in NED.

The missing 306 sources are explained mostly by the HECATE updates to the isophotal radii of galaxies with respect to RC3 (de Vaucouleurs et al. 1991; Corwin et al. 1994) used in Earnshaw et al. (2019). Figure 19 shows how most of the newly updated galaxy dimensions often differ slightly from the older values, both towards larger and smaller values. This leads both to the loss and new inclusion of some peripheral sources if the new sizes are smaller or larger, respectively. This phenomenon is clearly illustrated in Fig. 20, and shows how smaller isophotal ellipses lead to the exclusion of seven sources from Earnshaw et al. (2019). For the same reason, we only recover 271 galaxies as sometimes all the sources present in one galaxy are lost. As a rough estimate for the loss of sources, we take 154 galaxies present in both Earnshaw et al. (2019) and our catalogue and that are smaller in HECATE than RC3 and compute that, on average, 26% of the area sky area is lost. On the other hand, for 112 galaxies that are larger, n_Galaxies = 1, and excluding a source that we found to be matched in different galaxies for both catalogues, the mean sky area increment is 39%. This gives an intuition into the 23% of missing sources from Earnshaw et al. (2019), but it also indicates that some extra sources have been included. Nonetheless, the newer dimension values are considered more accurate due to the increase in photometric data (Makarov et al. 2014).

Fig. 19. Dispersion of the HECATE to RC3 semi-major and semi-minor axes ratios for galaxies appearing in this work and in Earnshaw et al. (2019).

Fig. 20. TOPCAT astrometric maps of sources found in galaxies NGC 4649 (large ellipses at the centre) and NGC 4647 (smaller ellipses at the upper right). The isophotal ellipses from RC3 and HECATE are also shown.

In addition to this, the discrepancies can be further explained by the improvement of detection algorithms in the 4XMM editions of the serendipitous XMM-Newton catalogues with respect to the 3XMM editions. As such, a fraction of the spurious detections present in 3XMM-DR4 are expected to be properly flagged or not included in 4XMM-DR9.

Regarding the ULX content, 260 sources classified by us as ULX candidates in our catalogue have a counterpart in Earnshaw et al. (2019). This implies the recovery of 68% of their 384 candidates. This can be explained by the more strict filtering pipeline used in our analysis in addition to the reasons given above. We also recover 337 of the 470 objects in Walton et al. (2011), 72% of the total.

Our findings regarding the spectral and abundance properties of the cULXss are a confirmation of the findings from Walton et al. (2011) and Earnshaw et al. (2019). Both of them find an overabundance of ULXs in spiral galaxies and a clear tendency for candidates hosted by late-type galaxies to present slightly harder spectra. In addition, Earnshaw et al. (2019) also find that the spectral properties of the ULX population resemble mostly that of the X-ray population in the range of 10³⁸ < L_X < 10³⁹ erg s⁻¹, while the AGN population presents a slightly different distribution. Finally, we recover five of the HLX candidates in their catalogues.

5.3.2. Previous serendipitous Chandra ULX catalogues

As of the writing of this paper, Kovlakas et al. (2020) is the most recent and largest Chandra-based ULX catalogue. Their work goes in parallel to our own, presenting 629 ULX candidates in 309 galaxies out of 23 043 sources in 2218 galaxies. Therefore, our catalogues complement each other to a great extent. Most remarkably, we use the same reference list of galaxies, HECATE, which leads to a very interesting comparison of results.

If we restrict our comparison to mutually shared galaxies at D > 1 Mpc, their numbers decrease to 11 359 sources in 849 galaxies. Kovlakas et al. (2020) focused on galaxies with D < 40 Mpc, and selected as ULX candidates all non-nuclear objects in AGN galaxies or all objects in non-AGN galaxies that have L_X > 10³⁹ erg s⁻¹ and present negative pileup and unreliability flags in Chandra. Using these conditions, 341 of these sources consist of ULX candidates in their catalogue. These 849 galaxies contain 3130 sources in our catalogue, out of which 457 are ULX candidates on our side. Only 2091 sources have a direct counterpart on both sides, including 301 ULX candidates on our side and 144 from Kovlakas et al. (2020). Finally, only 74 of the ULX candidates coincide on both sides.

Two remarkable discrepancies are noticed. Firstly, the higher density of sources in Kovlakas et al. (2020) in comparison to our own catalogue, and secondly, the lower yield of ULXs in Kovlakas et al. (2020) for the same sky area. The first can be easily explained by Chandra’s sharper resolution, which allows the detection of many fainter sources that are harder to resolve with XMM-Newton. The second point can be explained by the difference in our ULX filtering criteria. We choose as ULX candidates those objects that, in addition to having good data and being clean of interlopers, find themselves within the ULX luminosity regime in at least one of their detections, including the corresponding uncertainties. This implies the inclusion of objects at the ULX luminosity threshold into our ULX candidate set, as discussed in Sect. 5.2.1. By contrast, Kovlakas et al. (2020) use only the average source luminosity as their criterion. Using their criterion, our set of ULX candidates in the same area would be reduced from 457 to 347, closer to their 341 candidates. Additionally, our data are prone to source blending due to XMM-Newton’s poorer resolution, leading to slight luminosity overestimates.

There is an additional, more subtle, factor that contributes to the lack of coinciding candidates. As illustrated in Fig. 21, the estimated luminosity of matched sources are in good agreement between the two catalogues, but with a typical dispersion of around one order of magnitude. This alone explains why there are ULX candidates on our side that are not accounted for as such in Kovlakas et al. (2020) and vice versa, despite having a positional match on the sky. The source of the spread can be both intrinsic variability in the sources and uncertainties arising during the flux measurement. Nonetheless, these discrepancies add to the mutual completeness of the two catalogues. Even in overlapping observations, the two catalogues should be considered to have a complete view of the total ULX content.

Fig. 21. Dispersion of source luminosities for coinciding D > 1 Mpc sources in our catalogue (horizontal axis) and Kovlakas et al. (2020) (vertical axis). ULX candidates in our catalogue are highlighted in blue.

In addition to quantity, we perform similar qualitative work using the same parameters as provided by HECATE. Kovlakas et al. (2020) find an overabundance of ULX in spiral galaxies as well as a spike in elliptical galaxies, in agreement with our results. The most interesting common result is that the number of ULX candidates per SFR in early spiral galaxies is the closest to the one predicted in Wiktorowicz et al. (2017), of 0.67 ULX/M_⊙ yr⁻¹. In our case, we estimate a peak value of 0.65 ULX/M_⊙ yr⁻¹ for early spiral galaxies, while Kovlakas et al. (2020) present a value of 0.51 ULX/M_⊙ yr⁻¹. We also agree with Kovlakas et al. (2020) in that a ULX population in early-type galaxies exists despite their low SFR values, indicating that it is composed of an older population of LMXB. However, regarding the spiral galaxies we see that in the XMM-Newton sample, galaxies with T_H ∼ 3 present the highest ULX frequencies and fraction of hosts, while Kovlakas et al. (2020) find it at galaxies with T_H ∼ 5 (Sc galaxies in the original paper). This difference will most likely get blurred as ULX catalogue sizes grow in the future.

5.4. Future prospects

Our catalogue consists of the largest ULX collection built to date, and it is based on the most recent XMM-Newton serendipitous catalogue available at the time. Nonetheless, the available catalogues keep improving, both in size and in quality. For instance, in late 2020, the tenth data release of the XMM-Newton serendipitous catalogue, 4XMM-DR10¹¹, was made available to the public, adding 25 034 new sources from 443 new observations, an increase of almost 5% with respect to 4XMM-DR9. As the available samples keep increasing, so will the samples of ULX candidates available to the community. This will allow for more solid statements regarding the spectral and variability properties of the bright ULX population, and to better identify it as a sub-group of the ULX population.

A particularly promising resource that will soon be available is that of the eRASS survey (Merloni et al. 2012), which has been active for more than a year as of the writing of this paper (Predehl et al. 2021). Based on the area observed by eROSITA on the western side of the galactic hemisphere as of April 27, 2020, the preliminary study of Bernadich (2020) finds a total number of 132 ULX candidates, of which only 30 are matched to our catalogue within three times their positional uncertainty. This is to be expected as sources from the eRASS in general are spread all across the sky, and the observations overlap little with the XMM-Newton field. Taking into account the area covered by eROSITA and its average sensitivity at the end of the survey by 2024, the discovery of around 1000 new ULX candidates is expected.

X-ray catalogues are not the only ones of relevance to the search for ULXs. Information from other wavelength ranges helps to investigate the nature of the observed objects. A great deal of work is being done on that side too. For instance, also in late 2020, the Gaia Early Data Release 3 was made available to early access, providing positional and apparent brightness information for 10⁸ new sources (Gaia Collaboration 2021), with a full release planned by 2022¹². Given our methodology, catalogues like this will be of great help in identifying further interlopers in future works.

5.5. Conclusions

We have built an expanded ULX catalogue from the latest available XMM-Newton serendipitous source catalogue, 4XMM-DR9, and the HECATE list of galaxies. A total number of 23 262 point-like sources of high detection likelihood have been included within the isophotal ellipses of HECATE galaxies, of which 3274 have at least one detection within the ULX luminosity regime (L_X > 10³⁹ erg s⁻¹). However, most of these sources consist of contaminating interlopers, 2208 of them being the nuclear source of their host galaxy. We have built a filtering pipeline to identify such possible interlopers from other available catalogues and databases, such as Pan-STARRS1, Tycho2, SDSS-DR14, VéronQSO, SIMBAD, Pan-STARRS1, and NED. In the end, 779 objects in 617 galaxies qualify as ULX candidates, 94 of which with at least one detection with luminosity over the L_X > 5 × 10⁴⁰ erg s⁻¹ threshold. Around 30 of these objects qualify as HLXs, with ⟨L_X⟩> 10⁴¹ erg s⁻¹, many of which require individual follow-ups to untangle their physical nature.

We show that ULX candidates from our complete sub-sample are preferably found in late-type galaxies, in agreement with the notion presented in the recent ULX census (Kovlakas et al. 2020; Earnshaw et al. 2019) and reviews (Kaaret et al. 2017, and references therein). In particular, early and late spiral galaxies present the highest ULX frequencies, with 49% of early spiral galaxies and 37% late spiral galaxies hosting at least one ULX. These galaxies are also the ones with the highest SFR per unit of stellar mass according to the HECATE values, confirming the correlation between ULX and star formation. We compute 0.49 ULX/M_⊙ yr⁻¹ for early spirals and up to 1.01 ULX/M_⊙ yr⁻¹ for late spirals. These numbers range across the values reported in previous literature (Wiktorowicz et al. 2017; Kovlakas et al. 2020).

From the hardness ratio measurements provided by 4XMM-DR9, we have shown that our ULX sample has spectral properties that are distinct from those of the interloper population. At the same time we show that ULX candidates hosted by late-type galaxies tend to have harder spectra, most likely due to higher photoelectric absorption. These results are overly similar to those already presented by previous XMM-Newton-based ULX studies, such as Walton et al. (2011) and Earnshaw et al. (2019). The same trend is seen in the brighter candidates, but they have softer spectra in general.

From the stacked version 4XMM-DR9s, built from overlapping XMM-Newton observations, we also see that ULX candidates hosted by late spiral galaxies are the ones with the highest probability of showing inter-observational variability, with the highest amplitude modulations, while candidates hosted by elliptical and lenticular galaxies fall on the other end. As far as we are aware, this is the first study of this kind. Other variability studies such as Sutton et al. (2013) focus on a more thorough study of individual objects rather than on the whole sample, and therefore the results are not fully comparable. To fully understand the variability properties of a population as heterogeneous as ULXs, individual follow-ups are essential. Again, the bright candidates distinguish themselves by showing milder variability in general, even though similar trends surface among them.

We are not able to make solid statements regarding the spectral and variability properties of bright candidates due to the small size of the complete bright sub-sample. The expected growth of available X-ray samples in the future will mitigate this problem. Of particular promise is the all-sky blind eRASS survey being performed by the eROSITA observatory (Predehl et al. 2021).

Acknowledgments

This work was supported by the German DLR under project 50OX1901. Konstantinos Kovlakas and Andreas Zezas acknowledge support from the European Research Council under the European Union’s Seventh Framework Programme (FP/2007-2013)/ERC Grant Agreement n. 617001, and the European Union’s Horizon 2020 research and innovation programme under the Marie Skłodowska-Curie RISE action, Grant Agreement n. 691164 (ASTROSTAT). The authors would also like to thank all peers who helped in the construction of the catalogue by offering their own tools and experience, and the anonymous referee who helped to improve this paper. This research made use of the cross-match service provided by CDS, Strasbourg.

References

All Tables

All Figures

	Fig. 1. X-ray XMM-Newton (top) and optical Pan-STARRS1 (bottom) images of galaxy NGC 3631. Non-nuclear sources are indicated with red crosses, while the green cross is the central source. The Pan-STARRS1 image covers a slightly smaller sky area than the XMM-Newton one.
In the text

	Fig. 2. TOPCAT astrometric map of M51 showing all 4XMM-DR9 sources within its isophotal ellipse. The interlopers are composed of Pan-STARRS1 and SIMBAD objects, and all of them but one are below the ULX luminosity regime if assumed to be at the same distance as the others.
In the text

	Fig. 3. Pan-STARRS1 optical image of the location of source SRCID = 200559905010001. The X-ray source (not shown here) is a blend of the two cores of the interacting pair NGC 5256.
In the text

	Fig. 4. Pan-STARRS1 optical image of the location of source SRCID = 201389514010019 in galaxy NGC 2528 (large spiral on the right). The source lies clearly on top of a background galaxy at z = 0.13.
In the text

	Fig. 5. Pan-STARRS1 optical image of the location of source SRCID = 205562803010007 (green marker) in galaxy NGC 2903. The red marker indicates the position of the Pan-STARRS1 counterpart identified by the pipeline. According to NED, the optical counterpart constitutes an H II region.
In the text

	Fig. 6. Distribution of detection fluxes within 4XMM-DR9, regardless of their nature. Most detections (77%) have been detected with fluxes higher than FX ≈ 10−14 erg cm−2 s−1.
In the text

	Fig. 7. Distance-luminosity dispersion of all quality sources within 38ss, 39ss, and 5x40ss, and within cULXss and cbULXss, with an imposed sensitivity of 10−14 erg cm−2 s−1. The distance cuts correspond to Dmax ≈ 9, 29 and 204 Mpc.
In the text

	Fig. 8. Distance distribution of all sources of quality, ULX candidates, and bright ULX candidates in our catalogue.
In the text

	Fig. 9. Hubble type distributions of cULXss and cbULXss. In both cases there is a spike in ULX abundances in the earliest galaxies.
In the text

	Fig. 10. Luminosity distribution of the sources of quality and ULX candidates in our catalogue, divided between late-type galaxies (LTG) and early-type galaxies (ETG). A distance cut of D > 1 Mpc has been imposed to avoid contamination from the Local Group, heavily dominated by low-luminosity sources.
In the text

	Fig. 11. Distribution of ULX frequencies in cULXss and cbULXss across the morphological sub-types of host galaxies. Frequencies from cbULXss are multiplied by 100 to aid visibility.
In the text

	Fig. 12. Distribution of fraction of host galaxies in cULXss and cbULXss across their morphological sub-types. Values from cbULXss are shown in % to aid visibility.
In the text

	Fig. 13. Histogram of the number of ULX per star formation rate for early and late spiral galaxies from cULXss. The large values on the right are due to galaxies with low SFR that host a single ULX, and the most extreme cases have been cut out of the picture.
In the text

	Fig. 14. Hardness ratio dispersion of cULXss and cbULXss against the abundance contours of some populations of interlopers. AGNs and QSOs include both confirmed and candidate objects. Only hardness ratios with uncertainties lower than 0.2 are accepted for the sake of reliability.
In the text

	Fig. 15. Hardness ratio dispersion of cULXss and cbULXss against the abundance contours of lower luminosity objects present in 38ss. Only hardness ratios with uncertainties smaller than 0.2 are accepted for the sake of reliability.
In the text

	Fig. 16. Hardness ratio vs. luminosity dispersion of objects in cULXss (squares) and cbULXss (diamonds), distinguishing between objects hosted in late-type galaxies (TH ≥ 0) and early-type galaxies (TH < 0). Density contours have been drawn to help the eye. Only hardness ratios with uncertainties lower than 0.2 are accepted for the sake of reliability.
In the text

	Fig. 17. HR1 and HR2 across morphological sub-types of host galaxies from cULXss sources. The error bars show the 1σ dispersion of the values, and the data points have been slightly moved to the left (HR1) and to the right (HR2) for ease of visibility.
In the text

	Fig. 18. Hubble type distributions of galaxies with listed values in HECATE, within 4XMM-DR9 XMM-Newton observations, and within our catalogue. Portions have been normalized to sum to 1 in all of the three cases due to the large differences in size between the sets.
In the text

	Fig. 19. Dispersion of the HECATE to RC3 semi-major and semi-minor axes ratios for galaxies appearing in this work and in Earnshaw et al. (2019).
In the text

	Fig. 20. TOPCAT astrometric maps of sources found in galaxies NGC 4649 (large ellipses at the centre) and NGC 4647 (smaller ellipses at the upper right). The isophotal ellipses from RC3 and HECATE are also shown.
In the text

	Fig. 21. Dispersion of source luminosities for coinciding D > 1 Mpc sources in our catalogue (horizontal axis) and Kovlakas et al. (2020) (vertical axis). ULX candidates in our catalogue are highlighted in blue.
In the text