Two new catalogs of blazar candidates in the WISE infrared sky

TL;DR

This paper introduces two comprehensive catalogs of radio-loud blazar candidates selected via WISE infrared colors, aiming to improve identification of gamma-ray blazars and enhance understanding of their infrared and radio properties.

Contribution

The paper presents two new catalogs of blazar candidates with improved selection techniques and cross-matching, expanding the known sample and aiding gamma-ray source identification.

Findings

The WIBRaLS2 catalog contains 9541 sources with confirmed blazar-like properties.

The KDEBLLACS catalog includes 5579 candidate BL Lacs based on infrared color similarity.

Approximately 90% of previously identified candidates were confirmed as blazars.

Abstract

We present two catalogs of radio-loud candidate blazars whose WISE mid-infrared colors are selected to be consistent with the colors of confirmed gamma-ray emitting blazars. The first catalog is the improved and expanded release of the WIBRaLS catalog presented by D'Abrusco et al. (2014): it includes sources detected in all four WISE filters, spatially cross-matched with radio source in one of three radio surveys and radio-loud based on their q22 spectral parameter. WIBRaLS2 includes 9541 sources classified as BL Lacs, FSRQs or mixed candidates based on their WISE colors. The second catalog, called KDEBLLACS, based on a new selection technique, contains 5579 candidate BL Lacs extracted from the population of WISE sources detected in the first three WISE passbands ([3.4], [4.6] and [12]) only, whose mid-infrared colors are similar to those of confirmed, gamma-ray BL Lacs. KDBLLACS members area also required to have a radio counterpart and be radio-loud based on the parameter q12, defined similarly to q22 used for the WIBRaLS2. We describe the properties of these catalogs and compare them with the largest samples of confirmed and candidate blazars in the literature. We crossmatch the two new catalogs with the most recent catalogs of gamma-ray sources detected by Fermi LAT instrument. Since spectroscopic observations of candidate blazars from the first WIBRaLS catalog within the uncertainty regions of gamma-ray unassociated sources confirmed that ~90% of these candidates are blazars, we anticipate that these new catalogs will play again an important role in the identification of the gamma-ray sky.

Figures

Tables

Table 1. Table 2: Members of the WIBRaLS catalog, split according to their WISE  spectral types and classes, and provenance of the radio counterpart.

	BZB-like	Mixed	BZQ-like	Total
Class A	55	5	61	121
Class B	273	124	317	714
Class C	1086	579	1746	3411
Class D	2330	0	2965	5295
NVSS	3024	557	4093	7664
FIRST	36	7	21	64
SUMSS	694	144	975	1813
Total	3744	708	5089	9541

Table 2. Table 3: Sample of rows of the catalog of WIBRaLS sources.

AllWISE nameaafootnotemark:	R.A.bbfootnotemark:	Dec.ccfootnotemark:	W1-W2ddfootnotemark:	W2-W3eefootnotemark:	W3-W4fffootnotemark:	$s_{\mathrm{BZB}}$ggfootnotemark:	$s_{\mathrm{MIX}}$hhfootnotemark:	$s_{\mathrm{BZQ}}$iifootnotemark:	Classjjfootnotemark:	Typekkfootnotemark:	Radio counterpartllfootnotemark:	Sradiommfootnotemark:	$q_{22}$nnfootnotemark:
J000011.09-433316.4	0.0462203	-43.5545611	1.31	2.85	2.67	0.0	0.0	0.39	D	BZQ	SUMSSJ000011.2-433317	70.9	-1.41
J000020.40-322101.2	0.0850076	-32.3503443	1.35	3.26	2.47	0.0	0.0	0.52	C	BZQ	NVSSJ000020-322059	520.9	-1.94
J000029.07-163620.2	0.1211645	-16.6056221	0.39	2.01	2.45	0.53	0.0	0.0	C	BZB	NVSSJ000029-163621	89.1	-1.24
J000047.05+312028.2	0.1960452	31.3411703	0.82	2.48	2.4	0.32	0.0	0.0	D	BZB	NVSSJ000047+312027	46.4	-1.43
J000056.54-402206.4	0.235613	-40.368453	0.71	2.69	2.6	0.28	0.0	0.0	D	BZB	SUMSSJ000056.7-402208	76.0	-1.36
J000101.04+240842.5	0.2543718	24.1451458	1.23	3.06	2.13	0.0	0.0	0.34	D	BZQ	NVSSJ000101+240842	46.6	-1.26
J000105.29-155107.2	0.2720486	-15.8520035	1.18	3.42	2.14	0.0	0.0	0.33	D	BZQ	NVSSJ000105-155106	347.5	-2.01
J000108.11-373857.1	0.2838199	-37.6492076	1.32	2.72	2.29	0.0	0.0	0.39	D	BZQ	SUMSSJ000108.0-373901	23.1	-1.12
J000118.01-074626.9	0.3250683	-7.7741395	0.95	2.58	2.19	0.64	0.13	0.0	C	BZB	NVSSJ000118-074626	208.4	-1.23
J000131.63+165413.8	0.3818138	16.9038342	1.21	2.92	2.45	0.0	0.0	0.59	C	BZQ	NVSSJ000131+165416	63.1	-1.08
J000132.22+135258.4	0.3842501	13.8829081	1.07	3.0	2.48	0.0	0.0	0.23	D	BZQ	NVSSJ000132+135258	74.2	-1.59
J000132.34+240230.3	0.3847546	24.0417769	0.85	2.14	1.84	0.81	0.0	0.0	B	BZB	NVSSJ000132+240231	359.2	-1.28
J000132.74-415525.2	0.3864515	-41.9236926	0.68	2.17	2.09	0.48	0.0	0.0	D	BZB	SUMSSJ000133.1-415526	11.6	-0.84
J000132.83+145607.9	0.3868025	14.9355453	0.92	2.85	2.77	0.0	0.1	0.31	D	BZQ	NVSSJ000132+145609	314.6	-1.69
J000137.07+431543.9	0.4044828	43.2622002	0.72	2.61	2.32	0.76	0.0	0.0	B	BZB	NVSSJ000137+431544	61.8	-0.9

Table 3. Table 4: Break-down of the number of WISE  sources identified as BL Lacs candidate at each step of the selection described in Section  4 , split by provenance of the radio counterparts.

	NVSS	FIRST	SUMSS	Total
WISE selection	10166	5532	2128	17826
duplicates removal	10166	2099	1680	13945
close sources removal	100084	2049	1671	13804
$q_{12}$ selection	5310	327	305	5942
$\Vert b\Vert$ selection	4947	327	305	5579

Table 4. Table 5: Sample of rows of the catalog of KDEBLLACS.

AllWISE nameaafootnotemark:	R.A.bbfootnotemark:	Dec.ccfootnotemark:	W1-W2ddfootnotemark:	W2-W3eefootnotemark:	Radio namefffootnotemark:	S20cmggfootnotemark:	q12iifootnotemark:
J000007.63+420725.5	0.0318093	42.1237527	0.35	1.92	NVSSJ000007+420722	18.2	-1.5
J000010.29-363405.2	0.042887	-36.5681267	0.53	2.33	NVSSJ000010-363407	6.2	-1.3
J000056.22-082742.0	0.2342813	-8.4616809	0.43	1.94	NVSSJ000056-082747	14.4	-1.6
J000116.37+293534.5	0.3182368	29.5929424	0.61	2.46	NVSSJ000116+293534	3.5	-1.02
J000126.44+733042.6	0.3601711	73.5118347	0.77	2.22	NVSSJ000126+733042	23.6	-1.63
J000137.86-103727.3	0.4077672	-10.6242584	0.40	2.07	NVSSJ000137-103727	10.2	-1.41
J000147.28+455015.2	0.4470018	45.8375759	0.78	1.91	NVSSJ000147+455016	4.2	-1.22
J000236.06-081532.4	0.6502775	-8.2590058	0.71	2.11	NVSSJ000236-081533	28.3	-1.61
J000302.99-105638.1	0.7624893	-10.9439389	0.67	2.36	NVSSJ000302-105637	18.2	-1.61
J000311.94-070144.3	0.7997588	-7.0289838	0.61	1.98	FIRST J000311.9-070144	4.17	-1.01

Equations

$$q_{24}=\log (S_{24\mu m}/S_{1.4\mathrm{GHz}})$$

$$q_{22}=\log (S_{22\mu m}/S_{\mathrm{radio}})$$

$$q_{12}=\log (S_{12\mu m}/S_{\mathrm{radio}})$$

Full text

BL Lacertae objects: general - catalogs - galaxies: active - radiation mechanisms: non-thermal

1 Introduction

Blazars represent one of the most extreme classes of active galactic nuclei (AGNs). These radio-loud sources are characterized by flat radio spectra even at low radio frequencies (i.e., below $\sim$1 GHz) (Massaro et al., 2013a, b; Nori et al., 2014; Giroletti et al., 2016), superluminal motions (see e.g., Vermeulen, & Cohen, 1994; Lister, & Homan, 2005; Lister et al., 2009, and references therein), peculiar infrared colors (Massaro et al., 2011; D’Abrusco et al., 2014, hereinafter Paper I), high optical polarization (see e.g., Agudo et al., 2014; Pavlidou et al., 2014; Angelakis et al., 2016; Hovatta et al., 2016) and, not least, rapid and irregular variability (see e.g., Homan et al., 2002) at all frequencies with uncorrelated amplitudes and different time scales ranging between minutes to weeks (Homan et al., 2002). Blazars can reach bolometric luminosities up to $10^{49}$ erg s*-1* during $\gamma$-ray flaring states (Orienti et al., 2014).

The emission from blazars, according to the unification scenario of radio-loud AGNs, arises from a relativistic jet, closely aligned along the line of sight, that, in some cases, can outshine the host galaxy and the other AGN emission components (see e.g., Blandford, & Rees, 1978; Urry, & Padovani, 1995).

A close connection between the $\gamma$-ray and the infrared (IR) properties of blazars has recently emerged from the investigation of the all-sky surveys carried out with Wide-Field Infrared Survey Explorer (WISE) (Wright et al., 2010) and Fermi satellites (D’Abrusco et al., 2012; Massaro et al., 2016). The observed correlations across more than 10 decades of frequency has been proven to be an extremely powerful tool to identify new $\gamma$-ray emitting blazars among Fermi sources with unknown or uncertain lower-energy counterparts (Massaro et al., 2012a; D’Abrusco et al., 2013).

With a sky density of $\sim\!0.1$ sources per square degree, blazars dominate the $\gamma$-ray sky in the MeV-TeV energy range and represent $\sim\!40\%$ of the sources in the Third Fermi -Large Area Telescope (LAT) Source Catalog (3FGL, Acero et al., 2015; Massaro et al., 2015). Blazars are mainly divided into two sub-classes based on their optical spectra: (i) BL Lac objects characterized by featureless optical spectra or showing emission and/or absorption lines of equivalent widths $EW<$ 5Å  (Stickel et al., 1991; Falomo et al., 2014), and (ii) Flat Spectrum Radio Quasars (FSRQs) with a typical quasar-like optical spectrum. Hereinafter, we adopt the nomenclature used in the Roma-BZCAT catalog (Massaro et al., 2009, 2015), where BL Lacs and FSRQs are referred to as BZBs and BZQs, respectively.

The Spectral Energy Distributions (SEDs) of blazars are characterized by two broad components, a low-energy one peaking between the IR and X-ray bands and a high energy one whose emission ranges between the X-rays and the $\gamma$-rays. The former component is interpreted as synchrotron radiation from relativistic particles accelerated in a jet while the component at higher energies, according to leptonic models, is due to inverse Compton emission with seed photons that can have different origins (for more details see e.g., Böttcher, 2007, 2012), while hadronic models invoke the synchrotron emission by protons or secondary particles produced in proton-photon interactions (Dermer, & Schlickeiser, 1993; Mücke, & Protheroe, 2001; Murase et al., 2012, and references therein).

The BL Lac population is also divided, according to their SEDs, in “Low-frequency peaked BL Lac objects” (LBLs) when the peak of the first component lie in the IR-to-optical energy range and “High-frequency peaked BL Lac objects” (HBLs) when the synchrotron peak falls in the UV-to-X-rays energy range (Padovani, & Giommi, 1996). Another, more recent classification distinguishes blazars as low-synchrotron peaked (LSP), intermediate-synchrotron peaked (ISP) or high-synchrotron peaked (HSP) based on the peak frequency $\nu_{S,\mathrm{peak}}$ of the synchrotron component of their SEDs (Abdo et al., 2010; Ackermann et al., 2015). Here we will adopt the LBLs/HBLs sub-classification for the BZBs as it does not strictly depend on the exact location of the peak frequency.

The discovery of the peculiar WISE IR colors of blazars has been used to search for blazar-like sources within the positional uncertainty regions of the unidentified/unassociated $\gamma$-ray sources (UGSs) that could be their potential counterparts (see e.g., Massaro et al., 2015). Several procedures based on WISE data have been developed to investigate the nature of the UGSs listed in all Fermi source catalogs as well as to verify the nature of blazars candidates of uncertain type (BCUs; Massaro et al., 2012b; Cowperthwaite et al., 2013; Álvarez Crespo et al., 2016a). However, all these methods require spectroscopic confirmation of the natures of the blazar candidates selected.

One of the largest sources of candidate blazars used to identify unassociated $\gamma$-ray sources observed by the Fermi LAT has been the catalog of WISE Blazar-Like Radio-Loud Sources (WIBRaLS; Paper I). This catalog contains WISE sources detected in all four WISE bands, whose mid-IR colors are similar to those of confirmed Fermi blazars. WIBRaLS sources were also required to a) have a radio counterpart from one of three major surveys, namely the National Radio Astronomy Observatories Very Large Array VLA Sky Survey (NVSS, Condon et al., 1998) the VLA Faint Images of the Radio Sky at Twenty-cm Survey (FIRST, White et al., 1997; Helfand et al., 2015) and the Sydney University Molonglo Sky Survey Source Catalog (SUMSS, Mauch et al., 2003), and b) to be radio-loud, that is to have an observed ratio between radio and 22 $\mu$m mid-infrared flux densities $>3$ (Paper I).

Since the publication of the WIBRaLS catalog, extensive optical spectroscopic campaigns whose goal is to verify the nature of WIBRaLS1 blazar candidates that can be spatially associated to $\gamma$-ray sources observed by Fermi have been carried out (see e.g., Massaro et al., 2014; Paggi et al., 2014; La Mura et al., 2015; Landoni et al., 2015; Massaro et al., 2015; Ricci et al., 2015; Ajello et al., 2017a; Álvarez Crespo et al., 2016a, b, c; Peña-Herazo et al., 2017; Paiano et al., 2017a, b; Landoni et al., 2018; Marchesi et al., 2018; Paiano et al., 2019; Marchesini et al., 2019). The total number of WIBRaLS candidate blazars which have been spectroscopically followed up and reported upon in either one of the papers listed above is, to date, 159 split in 126 candidate BZBs, 16 candidate BZQs and 17 candidates with no spectral classification available (Mixed). The analysis of the optical spectra confirmed that $\sim\!93\%$ of candidate BZBs have featureless optical spectra typical of BL Lacs and $\sim\!52\%$ of the candidate BZQs show FSRQ spectra. Only $\sim\!3\%$, $\sim\!12\%$ and $\sim\!11\%$ of the spectra of observed candidate BZBs, BZQs or Mixed blazars cannot be classified as belonging to blazars, yielding a weighted average efficiency of the selection $\sim\!95\%$. Checks of spectra already published in the literature for 28 additional candidate blazars from WIBRaLS1 confirmed a BL Lacs or FSRQs nature for 27 of them.

Since the publication of the first release of the WIBRaLS catalog (Paper I), new versions of some of the main datasets used to define its selection method have been released:

• •

The ROMA BZCat, which contains the list of bona fide, spectroscopically confirmed blazars used to define the locus occupied by $\gamma$-ray emitting blazars in the WISE color space, has reached its 5th release (Massaro et al., 2009). This release contains $\sim$ 3600 sources, vs $\sim$ 3050 sources in the version used by D’Abrusco et al. (2014).

• •

The LAT 3-year Point Source Catalog (3FGL) (Acero et al., 2015) of $\gamma$-ray sources detected by Fermi  containing $\sim$ 3000 sources, has also become available. This recent update to the catalog of Fermi LAT sources is more than twice as large as the 2-year 2FGL (Nolan et al., 2012) catalog ($\sim$ 2250 members) used to extract the first version of the WIBRaLS catalog.

These two new datasets yield jointly a larger sample of confirmed $\gamma$-ray blazars that can be used to characterize more accurately their WISE mid-IR properties and hence to identify more effectively blazars that may be detected in the $\gamma$-ray energy range. In this paper we describe a new release of the WIBRaLS catalog that, by taking advantage of these most recent data available, maximizes the legacy value of WISE observations for the investigation of blazars.

The completeness of the WIBRaLS catalog is a function of the blazar spectral class, and decreases significantly for BL Lacs, which often are not detected in the $[22]$ $\mu$m WISE band. BZBs, in particular HBLs, have lower detection rate in the $[22]$ WISE band (Paper I), since their emission in the mid-infrared at 22 $\mu$m may be lower than the limiting sensitivity in the fourth WISE band.

For this reason, in this paper we will also present a new, complementary catalog of candidate BZBs selected with a novel technique that employs the WISE colors obtained from the first three WISE bands only and has been applied to AllWISE sources not detected in W4. We focus only on the BZB spectral class for two main reasons: i) the region of the IR color-color space occupied by BZB is less contaminated by spurious IR sources than that of BZQs, ii) during our campaign of spectroscopic follow-up of candidate BL Lacs selected with WISE colors, we found that a large fraction of UGSs and BCUs are classified as candidate BZBs, providing an indication that they are the most elusive counterparts of Fermi sources (70.5 % and 65.4 %, respectively, see Massaro, & D’Abrusco, 2016).

The paper is organized as follows: in Section 2 we provide a brief introduction to the WISE and radio catalogs used to extract both samples of candidate blazars. Sections 3 and 4 describe the selection methods of the new WIBRaLS catalog and the catalog of BL Lac candidates selected using two WISE colors only, respectively. The comparison between the two catalogs presented in this paper and the literature is described in Section 5. Finally, our conclusions are summarized in Section 6.

We use cgs units unless otherwise stated and spectral indices, $\alpha$, are defined by flux density $S_{\nu}\propto\nu^{-\alpha}$ indicating as flat spectra those with $\alpha\!<\!0.5$. WISE magnitudes used here are in the Vega system and are not corrected for the Galactic extinction. As shown in our previous analyses (D’Abrusco et al., 2013, 2014), such correction affects only the magnitude at 3.4 $\mu$m for sources lying at low Galactic latitudes, and it ranges between 2% and 5% of the magnitude, thus not affecting significantly our results.

2 Data

2.1 Infrared data

All the candidate blazars discussed in this paper were sources extracted from images produced by the WISE space telescope (Wright et al., 2010). WISE has observed the whole sky from 2009 to 2011 in the four bands W1, W2, W3 and W4 centered on 3.4, 4.6, 12 and 22 $\mu$m, respectively.

The AllWISE source catalog has superseded the previously available WISE All-Sky catalog; it was produced by combining WISE single exposures images from the first two years of the mission with the post-cryogenic phase data (Mainzer et al., 2011), that observed the sky in three (W1, W2 and W3) and then only two (W1 and W2) bands. The final result is a two-fold improvement in the depth-of-coverage in the first two bands thanks to the additional observations that have increased the sensitivity of the stacked images, and improved photometric accuracy in all filters thanks to updated background measurements. The AllWISE sensitivities are 54, 71, 730 and 5000 $\mu$Jy in the W1, W2, W3 and W4 passbands respectively, with angular resolutions of 6″.1, 6″.4, 6″.5 and 12″. Due to the WISE survey strategy, the limiting sensitivity of the catalog of sources extracted from the WISE images is not uniform on the sky (compare with Figure 8 at the Explanatory Supplement to the AllWISE Data Release Products111http://wise2.ipac.caltech.edu/docs/release/allwise/expsup/sec4_2.html).

The AllWISE Source Catalog contains astrometry and photometry in the IR for 747,634,026 objects; only 25,882,082 of these sources ($\sim$3.5%) were seen in all the four bands, increasing up to 99,118,890 sources ($\sim$13.3%) detected in the first three bands W1, W2 and W3 only. These two catalogs represent the parent samples for the updated WIBRaLS (Section 3) and the new catalog of BL Lacs candidates (Section 4), respectively.

2.2 Radio data

We have crossmatched the catalog of AllWISE sources according to the procedure described in Section 3.2 with the NVSS, the FIRST and the SUMSS radio surveys. A brief description for each of these three radio surveys is given below.

• •

NVSS (Condon et al., 1998) is a 1.4 GHz continuum survey performed using the Very Large Array (VLA) and covering the entire sky north of $-40$ deg declination, i.e. 82% of the celestial sphere with a beam size of 45″ FWHM. The result is a catalog of over 1.8 million discrete sources brighter than S $\sim$ 2.5 mJy in the entire survey.

• •

FIRST (White et al., 1997; Helfand et al., 2015) is a project designed to produce the radio equivalent of the Palomar Observatory Sky Survey over 10,000 square degrees of the North and South Galactic Caps using the NRAO Very Large Array (VLA). The beam size varies between 5″.4 FWHM for circular beam and 6″.8 along the major axis for elliptical shape, as a function of the declination of the observation. The survey area has been chosen to coincide with that of the Sloan Digital Sky Survey (SDSS, see e.g., Gunn et al., 1998) and at the mv $\sim$23 limit of SDSS, $\sim$40 % of the optical counterparts to FIRST sources are detected. At the 1 mJy source detection threshold, there are $\sim$90 sources per square degree.

• •

SUMSS (Mauch et al., 2003) is a radio imaging survey of the southern sky carried out with the Molonglo Observatory Synthesis Telescope (MOST) operating at 843 MHz, with a beam size $\approx$ 43″ FWHM. The catalog covers approximately 3500 deg2 with declination $\delta<$ $-$30∘, about 43% of the total survey area. The survey has a limiting peak brightness of 6 mJy/beam at declinations $\delta\leq$ $-$50∘, and 10 mJy/beam at $\delta>$ $-$50∘. The SUMSS is therefore similar in sensitivity and resolution to the NVSS, with $\sim\!$7000 sources found in the overlap region.

3 The second release of the WIBRaLs catalog

In this paper, a new release of the catalog of WISE Blazar-like RAdio-Loud Sources (WIBRaLS) has been built up by following an improved version of the procedure described by Paper I and including new samples of confirmed blazars associated with Fermi $\gamma$-ray sources. Schematically, the steps we followed to extract the WIBRaLS catalog are the following:

We select WISE sources detected in all four WISE bands, with IR colors similar to those of associated Fermi blazars (Section 3.1). We call these sources “WISE blazar-like sources”. 2. 2.

The WISE blazar-like sources are positionally cross-matched with sources extracted from either one of the radio surveys NVSS, SUMSS and FIRST (Section 3.2) and only those with a radio counterpart are retained. 3. 3.

Among the WISE blazar-like sources with a radio counterpart, only radio-loud sources are selected as members of the WIBRaLS catalog (Section 3.3).

In the following, we provide descriptions of each step above.

3.1 WISE color selection of blazar-like sources

We define the blazar-like locus using the sample of confirmed Fermi blazars, based on the 5th version of the ROMA-BZCat catalog (Massaro et al., 2015) and the 3FGL catalog (Acero et al., 2015), associated with AllWISE counterparts detected in all four WISE filters. This locus thus includes newly identified 4-band sources and extends the locus presented in Paper I.

The BZCat v.5.0 contains 3561 bona fide blazars with spectroscopic confirmation. Optical spectra are also used to classify members of the BZCat as BZBs or BZQs, according to the total equivalent width of all emission and absorption features, while blazars whose properties are intermediate between BZB and BZQ are tagged as Uncertain (BZU) and other blazar-like objects whose emission is mostly contaminated by light from the host galaxy are labeled as BZGs.

All BZCat sources classified as BZU and BZG were discarded, leaving 1151 BZBs and 1909 BZQs. In order to spatially crossmatch the positions of these BZCat sources with AllWISE sources detected in the four filters, we adopted the maximum radial distance of 3″.3 obtained by conservatively combining a nominal uncertainty of 1″on the radio positions of BZCat sources with positional uncertainty in the WISE W4 passband (see D’Abrusco et al., 2013, for details). This sample was then spatially crossmatched with $\gamma$-ray sources included in the 3FGL catalog (Acero et al., 2015), the latest release of sources detected by the LAT instrument on board NASA’s Fermi spacecraft, based on the first 48 months of survey data. The crossmatch has been performed by taking into account the position angles and the lengths of the semi-major and semi-minor axes of the 95% confidence region of each Fermi source, available in 3FGL, and the positional uncertainty of the AllWISE counterparts of the BZCat sources. The final sample used to define the WISE three-dimensional locus consists of 901 confirmed $\gamma$-ray emitting blazars, split in 497 BZBs and 404 BZQs, and is more than twice as large as the analogous list in Paper I, which contained 447 sources.

Following Paper I, we define a model of the locus in the three-dimensional Principal Components (PCs) space generated by the three independent WISE colors W1-W2 ($[3.4]$-$[4.6]$), W2-W3 ($[4.6]$-$[12]$) and W3-W4 ($[12]$-$[22]$). The locus is modeled as a set of three coaxial cylinders whose axes lie on the direction of the first PC (PC1). The two extremal cylinders are populated by blazars of similar spectral classes, namely BZB and BZQ respectively, while the intermediate cylinder contains significant ($\geq\!25\%$) fractions of both BZBs and BZQs. The length along the PC1 axis of the two extremal cylinders is defined so that they contain at least 75% of blazars classified as candidate BZBs and BZQs, respectively, while excluding the locus sources with PC1 value smaller than the 1%-st percentile and larger than the 99%-th percentile of the distribution of locus PC1 coordinates. The mixed cylinder contains fractions of BZBs and BZQs each smaller than 75%. The radii of the three cylinders, which lie in the plane generated by PC2 and PC3 and orthogonal to PC1, are defined to contain 95% of the sources whose PC1 coordinates fall in each of the three cylinders. Sources located within the locus are selected on the base of the value of their “score”, a quantitative measure of the distance of a generic source to the locus model, defined as follows.

The WISE colors of a WISE source and the associated uncertainties are projected into the PC three-dimensional space, where they define an “uncertainty ellipsoid”. The position and orientation of this ellipsoid - relative to the model of each cylinder separately - are used to calculate a numeric value defined between 0 and 1. This normalized distance is then weighted by the volume of the error ellipsoid, so that two WISE sources placed in the same position relative to the locus model but with different uncertainties are assigned different scores. A more detailed description of the score and its properties can be found in Paper I.

The score is also used to classify candidate blazars compatible with the model of the locus according to their reliability. All WISE sources with non-zero score for either cylinder are ranked according to the decreasing compatibility with the locus model in the classes A, B, C and D, defined by the 90%-th, 60%-th, 20%-th and 5%-th percentiles of the score distribution for the sources of the locus sample. These classes are defined to facilitate quick prioritization of candidate blazars as targets of follow-up observation across sources of different spectral types (BZB-like, BZQ-like or Mixed) and regardless of the specific score distributions for each spectral type. Nonetheless, classes do not replace scores as a quantitative indicator of the degree of compatibility of each candidate blazar with the locus model in the WISE color space; all scientific analysis on the final list of WIBRaLS should be performed using scores. The lower 5% threshold corresponds to the fraction of locus sources that are located outside the locus model by definition. The values of the score thresholds used to define the three classes are reported in Table 1.

To all sources with score larger than the 5-th percentile threshold, for any one of the three cylinders, are assigned the corresponding type (BZB, BZQ or Mixed) and selected as WISE blazar-like sources. The Mixed type does indicate a specific spectral class, since the Mixed cylinder contains comparable fractions of both BZBs and BZQs.

The total number of WISE blazar-like sources selected is 526,681, split in 156,506 BZB candidates ($\sim\!30\%$), 348,805 ($\sim\!60\%$) BZQ candidates and the remaining 21,370 ($\sim\!4\%$) sources located in the Mixed region. The WISE blazar-like sources can also be split in 7,807 Class A sources ($\sim\!1\%$), 27,986 Class B sources ($\sim\!5\%$), 149,052 Class C sources ($\sim\!28\%$) and the remaining 341,836 ($\sim\!65\%$) belonging to the class D.

3.2 Radio counterparts

Following Paper I, we determined the optimal radii for the association of the WISE blazar-like sources with their potential radio counterpart in the NVSS, SUMSS and FIRST surveys.

We adopted a modified version of the procedure illustrated by Best et al. (2005) and Donoso et al. (2009). They computed the optimal radius for the spatial crossmatch of NVSS and FIRST radio detections with optical sources in SDSS by setting a threshold on the fraction of spurious associations (i.e. the contamination) obtained for different values of maximum crossmatch radius. In this paper, similarly at what done in Paper I, the optimal association radius is indeed set as the radial distance that provides a given fixed efficiency of the selection $e_{\mathrm{thr}}\!=\!99\%$, corresponding to a contamination $c_{\mathrm{thr}}\!=\!1\%$, where $c(\vartheta)\!=\!100\%\!-\!e(\vartheta)$. The efficiency or purity of the selection is defined as the number of sources around real radio positions $n_{\mathrm{real}}(\vartheta)$, reduced by the number of sources around mock radio positions $n_{\mathrm{mock}}(\vartheta)$ and divided by the number of “real” matches.

We estimated $n_{\mathrm{real}}(\vartheta)$ by counting the number of WISE sources detected in all four bands within circular regions of radius $\vartheta$ between 0″ and 60″  centered on a sample of $5\!\cdot\!10^{4}$ sources randomly extracted from each of the three radio surveys. The value of $n_{\mathrm{mock}}(\vartheta)$ was calculated by averaging over one hundred mock realizations of the coordinates of each real radio source, generated by moving the real position in a random direction and by a random radial distance in the $[60,120]$″ range.

The maximum cross-match radii evaluated for the NVSS and SUMSS surveys are $\vartheta_{\mathrm{NVSS}}\!=\!10\arcsec.4$, $\vartheta_{\mathrm{SUMSS}}\!=\!7\arcsec.4$ respectively, identical to those determined for the 1st release of the WIBRaLS catalog. This approach does not work for the FIRST survey because of the very high density of FIRST sources. For this reason, we assume as optimal search radius for FIRST the value previously determined $\vartheta_{\mathrm{FIRST}}\!=\!3\arcsec.4$ (Paper I), that was obtained by combining conservative estimates of the positional uncertainties of FIRST and AllWISE sources.

The number of WISE blazar-like sources associated with at least one radio counterpart in either of the three surveys within the maximum radial distances discussed above is 32630, split in 18903 with a NVSS counterpart, 1040 with a FIRST counterpart and 3323 with a SUMSS counterpart. In order to exclude WISE sources associated to distinct radio sources associated to the emission of lobes of the same radio galaxy, we searched and removed all duplicate radio sources whose positions would fall within 6′from each other. This radius matches the typical positional uncertainty of $\gamma$-ray sources detected by LAT in the 3FGL (Acero et al., 2015). In these cases, the AllWISE source with the largest score among all the duplicates was retained. After this step, the number of WISE blazar-like sources with a radio counterpart in the NVSS, FIRST and SUMSS catalogs is 18693, 10227 and 3287 respectively, for a total of 32207 sources. This sample includes 5547 sources with a counterpart from both the NVSS and FIRST surveys and 709 sources with one counterpart from both NVSS and one SUMSS. The number of unique WISE blazar-like sources associated to a radio counterpart (i.e., the previous sample after removing sources with radio counterpart listed in two catalogs) is 25951, with $\sim\!72\%$ (18693) has a NVSS counterpart, $\sim\!18\%$ (4687) has a FIRST counterpart and the remaining $\sim\!10\%$ (2578) are associated to a FIRST source. Moreover, 7101 sources ($\sim\!27\%$) of this list are classified as candidate BZBs, 2186 ($\sim\!8\%$) as Mixed candidates and the remaining 1664 ($\sim\!64\%$) as candidate BZQs 222The selection of candidate BZQs using solely their WISE photometry can be contaminated by normal quasars as they share the same region of the WISE colors space (Wright et al., 2010). Also, the existence of radio counterparts to WISE -selected BZQ candidates alone, lacking radio spectral characterization, does not remove the degeneracy (Stern et al., 2005, 2012), with 1366 class A ($\sim\!5\%$), 3302 class B ($\sim\!13\%$), 9265 class C ($\sim\!36\%$) and the remaining 12018 ($\sim\!47\%$) as class D sources.

3.3 Radio-loudness selection

Blazars are radio-loud AGNs but not every radio-loud AGN is a blazar. Blazars are generally hosted in elliptical galaxies, whose emission at radio frequencies is dominated by synchrotron emission from particles accelerated in the AGN relativistic jet pointed along the line of sight. In order to distinguish blazars from other radio sources, when lacking radio spectral information, we adopted the approach described in what follows. Padovani et al. (2011) as well as Bonzini et al. (2013) suggested that AGN-powered radio-sources can be identified using $q_{24}$, a modified definition of the so-called $q$ parameter (Helou et al., 1985), i.e. the logarithm of the ratio of far-IR to radio flux density, to overcome the dearth of accurate flux density measurements at far-IR frequency. The $q_{24}$ parameter is defined as:

[TABLE]

where $S_{24\mu m}$ is the observed flux density at 24 $\mu$m and $S_{1.4GHz}$ is the flux density measured at 1.4 GHz. Following D’Abrusco et al. (2014), we adopted a similar criterion to select “radio-loud sources” among the WISE blazar-like sources with a radio-counterpart by using the parameter $q_{22}$ defined as:

[TABLE]

where the flux density in the 24 $\mu$m band of the Multi-band Imaging Photometer used on Spitzer (MIPS) is replaced by the WISE W4 ($[22]$ $\mu$m) band. This approach is possible thanks to the similarities of the two wavebands (see Paper I for additional details). For the $S_{\mathrm{radio}}$ in Equation 2 we used the radio flux density at 1.4 GHz for sources with a NVSS or FIRST radio counterparts. Given the lack of flux measurement at 1.4 GHz for SUMSS sources, for those sources we used flux densities at 843 MHz.

Due to the remarkable flatness of radio spectra of blazars (see e.g., Healey et al., 2007; Massaro et al., 2013d), replacing the flux density at 1.4 GHz with the same quantity at 843 MHz produces a small effect of the value of the $q_{22}$ parameter, whose size was estimated using the sample of WISE blazar-like sources with a radio counterpart detected in both NVSS and SUMSS (508 sources). The distribution of the difference between the values of $q_{22}(1.4~{}\mathrm{GHz})$ and $q_{22}(843~{}\mathrm{MHz})$ is fairly constant across the interval of $q_{22}(1.4~{}\mathrm{GHz})$ covered by our sample, with $\Delta q_{22}\!\approx\!-0.07$, confirming a fairly flat radio spectrum for the sources in this sample. Based on this finding, for the WISE blazar-like sources associated with SUMSS counterparts only, we have used the corrected $q_{22}$ value defined as $q^{{}^{\prime}}_{22}(843~{}\mathrm{MHz})\!=\!q_{22}(843~{}\mathrm{MHz})\!+\!\Delta q_{22}$.

Both the $q_{24}$ and $q_{22}$ parameters show a dependence on the redshift of the source, as their values decrease for larger redshifts (Bonzini et al., 2013; D’Abrusco et al., 2014). Figure 1 shows the distribution of $q_{22}$ values for the confirmed $\gamma$-ray emitting blazars in the locus sample with reliable redshift measurements, as a function of the redshift and color-coded according to their spectral classification from the BZCat. Locus sources (mostly BZBs) with uncertain or unknown redshifts were not used, reducing the number of locus members to 563 (split in 159 BZBs and 404 BZQs). In Figure 1, the same trend reported by Bonzini et al. (2013) and Paper I is observed. If redshift estimates were available for all WISE blazar-like sources associated to a radio counterpart, we could have determined different $q_{22}$ thresholds in different redshift bins. Since redshifts are not available, we used fixed values of $q_{22}$.

In this paper, we improve over the approach used in Paper I by determining different $q_{22}$ thresholds for WISE blazar-like sources with radio counterpart classified as candidate BZQs, BZBs and Mixed. The thresholds for candidate BZBs and BZQs, calculated as the 95% of the $q_{22}$ distribution of all blazars in the locus sample classified as BZB and BZQ, are $q_{22}^{\mathrm{BZB}}\!\leq\!-0.61$ and $q_{22}^{\mathrm{BZB}}\!\leq\!-0.97$, respectively. The threshold for Mixed sources is set to the mean value of the thresholds for candidate BZQs and BZBs, $q_{22}^{(\mathrm{Mixed})}\!\leq\!-0.79$. These thresholds are shown as horizontal lines in Figure 1.

We evaluated the effect of the unknown underlying redshift distribution of the WISE blazar-like sources with radio counterpart on our selection based on fixed $q_{22}$ thresholds. We followed the same strategy used in Paper I. We computed the $q_{22}$ for all sources in the locus sample after varying their observed redshifts over an equally spaced grid covering the $[0,4]$ range with bins of 0.05 width, which includes the interval $[0,\sim\!1.25]$ covered by the observed redshifts. We assumed power-law spectral energy distribution with slope constrained by the observed flux densities at 22 $\mu$m and at 1.4 GHz. The fractions of candidate BZQs and BZBs that satisfy the $q_{22}$ conditions based on the observed redshift distribution are $\sim\!91\%$ and $\sim\!93\%$, respectively. These fractions are slightly lower than the $\sim\!94\%$ fraction of sources recovered with the $q_{22}\!\leq\!-0.5$ threshold used in Paper I.

The total number of WISE blazar-like sources with a radio counterpart that satisfies the radio-loudness criteria based on the $q_{22}$ parameter and, thus, belong to the second release of the WIBRaLS catalog is 9541. The break-down of the catalog according to WISE spectral type, class and provenance of the radio counterpart is given in Table 2, while the basic parameters for a subset of sources in the catalog are displayed in Table 3. The number of sources of the locus sample of confirmed $\gamma$-ray emitting blazars used to determine the WIBRaLS selection, that are found in the final WIBRaLS2 is 666, $\sim 26\%$ less than the original size of the locus sample (901 sources). The exclusion of these 235 locus sources is the result of a) the definition of the WISE model (Section 3.1), that excludes 2% of the sources based on their location along the PC1 axis, and 5% of each spectral class because the radii of the three cylinders are defined to contain the 95% of associated sources, and b) the definition of the $q_{22}$ thresholds (Section 3.3), which remove 5% of the remaining locus sources for each spectral class, by definition.

4 Catalog of WISE $\gamma$-ray BL Lac candidates

As mentioned briefly in Section 2.1, only $\approx 2.5\%$ of the AllWISE sources are detected at $[22]$ $\mu$m. For this reason, the WIBRaLS catalog, that requires its members to be detected in all four WISE bands, will certainly be incomplete for BZBs whose SEDs peak in the X-rays and are not bright enough to be detected in the W4 WISE band. This occurrence is schematically displayed in Figure 2, where the typical shapes of the SEDs of HBLs and LBLs are plotted together with the WISE sensitivity limits for the four WISE filters.

We designed a novel method that relies only on the colors obtained with the W1, W2 and W3 WISE bands to produce a list of candidate BL Lacs, extracted from the AllWISE survey. The steps adopted for our procedure can be summarized as follows:

We select all AllWISE sources detected only in the first three WISE bands, having IR colors similar to those of known Fermi BZBs belonging to the ROMA-BZCat v5.0. 2. 2.

From this sample, we further select only those sources with a radio counterpart in either one of the NVSS, SUMSS and FIRST radio surveys. 3. 3.

We identify as candidate BZBs sources whose values of the radio-to-infrared flux density ratio are compatible with those of confirmed Fermi BL Lacs.

Details on these steps are given in the following sections.

4.1 WISE color selection of BZB-like sources

The training set used to identify the WISE mid-IR colors of BL Lac objects was built selecting all the Fermi sources belonging to the 3FGL and associated with BZBs in the latest release of the Roma-BZCAT catalog (Massaro et al., 2015). We only considered sources whose AllWISE counterparts are not detected in the W4 band. The total number of unique sources in the training set is 93. Among these sources, 34 are associated with a NVSS radio counterpart only, 3 are associated with a FIRST source only and the remaining 56 have both a NVSS and FIRST counterpart.

We selected WISE sources with IR colors similar to those of our training set by adopting the same procedure used in previous analyses (see e.g. Massaro et al., 2011; Paggi et al., 2013; Massaro et al., 2013a), based on the Kernel Density Estimation (KDE). KDE is a non-parametric procedure that estimates the Probability Density Function (PDF) of a multivariate distribution with no assumption on the properties of the parent population. The KDE depends on only one parameter, i.e., the bandwidth of the kernel of the density estimator, which is qualitatively analogous to the window size for one-dimensional running average.

We applied the KDE to the 2-dimensional distribution of training set sources in the WISE W2-W3 vs W1-W2 color-color plane to determine its PDF. Then, we selected infrared sources in the AllWISE catalog not detected in the W4 passband and detected in the other three filters, whose colors are located within the isodensity contour enclosing 90% of the BL Lac training set. Only sources whose colors uncertainties ellipses are entirely contained within the 90% isodensity contour were retained. Figure 3 shows the distribution of the training set in the WISE color plane, together with the iso-density contours determined by the KDE method, with the 90% contour displayed as a thick black line. The projections on this color-color diagram of the WIBRaLS2 catalog is also shown for reference.

The total number of WISE sources extracted applying this method is $\sim$14406, corresponding to $\sim\!0.01\%$ of the parent sample of AllWISE sources not detected only in W4.

4.2 Radio counterparts

We further select possible BL Lacs candidates by searching for radio counterparts of the sample of sources selected based on their WISE colors. This procedure obviously misses “radio weak BL Lacs” but to date they are extremely rare and their associations with Fermi sources has not yet been verified (Massaro et al., 2017; Bruni et al., 2018).

We crossmatched the color-selected WISE sources with the radio surveys NVSS, SUMSS and FIRST using the same maximum association radii established as described in Section 3.2. We found a total of 17826 sources with at least one radio counterpart: 5532 are associated with a FIRST source within 3″.4, 10166 with a NVSS source within 10″and the remaining 2128 with a SUMSS source with 7″.4 (see Table 4). After removing from our list sources with SUMSS or FIRST counterparts that are also associated with a NVSS source, there remain 2099 WISE sources with a unique FIRST counterpart, 1680 WISE sources with a unique SUMSS counterpart and the 10166 associated with a NVSS radio source, for a total of 13945 sources. Similarly at what was done for the WIBRaLS catalog in Section 4.2, in order to avoid possible contamination from radio lobes originating from radio-galaxies, we also checked for sources that had from another radio source closer than 6′, and, for FIRST radio counterparts only, sources with side lobe probability $\geq$ 0.05. These constraints reduced the number of sources to 2049 sources with a unique FIRST counterpart, 1671 sources associated to a unique SUMSS source and 10084 w ith a radio counterpart in NVSS, for a total of 13804 candidates.

4.3 Infrared-to-radio ratio selection

The last step of the procedure to select WISE BL Lac candidates is based on the characterization of the distribution of their infrared-to-radio ratios $q_{12}$, similar to the $q_{22}$ parameter used to select WIBRaLS sources 3.3. We define the parameter $q_{12}$ as:

[TABLE]

i.e., the logarithm of the ratio between the WISE flux density measured in the W3 passband and the radio flux density. Flux density at 20 cm has been used to calculate $q_{12}$ of sources with a unique counterpart in the NVSS and FIRST surveys. Given the flatness of the radio spectrum of BL Lacs (Healey et al., 2007; Massaro et al., 2013a) and our estimate of the $\Delta q_{22}\!\approx\!-0.07$ for the WISE -selected sources with both SUMSS and NVSS radio counterparts (see Section 3.3), we assume a flat spectrum ($\alpha_{\mathrm{radio}}\!=\!0$) and use the flux density at 36 cm to calculate the $q_{12}$ parameter for sources with SUMSS counterparts. 333We estimate that for radio spectra with $\alpha_{\mathrm{radio}}\!=\!\pm 0.1$, the effect on the values of $q_{12}(843~{}\mathrm{MHz})$ calculated at 36 cm for SUMSS counterparts would be a factor ranging from $\approx\!0.95$ to $\approx\!1.05$. In the case of a much larger spectral index $\alpha_{\mathrm{radio}}\!=\!0.5$, the corrective factor would still be $\approx\!0.77$.

We consider BL Lac candidates those sources whose $q_{12}$ values are consistent with the observed $q_{12}$ distribution of the training set sources. The $q_{12}$ values for sources with a FIRST or NVSS radio counterpart (59 and 90 respectively) have been used separately to determine the intervals of acceptable $q_{12}$ values for both surveys, as in both surveys flux densities are measured at the same wavelength (20 cm).

Figure 4 shows the distributions of the $q_{12}$ values of the BL Lacs training set sources (light blue histograms), control samples including random WISE sources with radio counterparts (solid gray histograms) and of the candidate BL Lacs selected according to their WISE colors as described in Section 4.1 (dashed gray histograms). FIRST (upper panel), NVSS (mid panel) and SUMSS (lower panel) radio counterparts are shown separately; the light blue histogram is missing from the SUMSS panel as no source in the training set has a radio counterpart from the SUMSS survey. It is also interesting to notice that the peaks of the distributions of $q_{12}$ values for the control sample and the WISE -selected candidates with SUMSS radio associations (lower panel in Figure 4) are significantly shifted towards smaller values of $q_{12}$ relative to the NVSS and FIRST distributions. This difference is due to relative shallowness of the SUMSS survey, whose sensitivity allows the detection of sources with minimum density flux at 36 cm $S_{\mathrm{radio}}^{\mathrm{SUMSS}}(36~{}\mathrm{cm})\!=\!5$ mJy, larger than the $\sim\!2$ mJy and $\sim\!0.2$ mJy for NVSS and FIRST, respectively.

We defined the lower and upper limits of $q_{12}$ values used to select the candidate BL Lacs as the 10-th and 90-th percentiles of the $q_{12}$ distribution of the training set (red lines in Figure 4). The $q_{12}$ intervals used to select the candidates are between -1.85 and -1 for NVSS training set sources and between -1.64 and -1.01 for sources with FIRST counterparts. Given the lack of training set BL Lacs with SUMSS counterpart, we have conservatively defined the upper and lower thresholds for the $q_{12}$ selection of SUMSS BL Lac candidates as the highest and lowest values of the lower and upper thresholds on the $q_{12}$ distribution of training set sources with a FIRST or NVSS counterpart, respectively. As a result, the interval of allowed $q_{12}$ values for SUMSS counterparts ranges between -1.64 and -1. After applying the $q_{12}$ selections to the sample of 13804 WISE-radio selected sources, we obtain 5941 candidates, split in 327 associated to a FIRST source, 5310 with a NVSS counterpart and the remaining 305 associated with a SUMSS source.

Then, as final step, we discarded all sources in the catalog located at Galactic latitudes $|b|\!<\!10\deg$, since $\gamma$-ray sources at low galactic latitudes suffer from a higher detection threshold due to a higher Galactic diffuse emission background (Ackermann et al., 2015). This constraint reduces the total number of sources in the catalog of KDE-selected BL Lac candidates to 5579, due to the exclusion of 363 sources with a NVSS radio counterpart. In what follows, the sources in the catalog of KDE-selected candidate BL Lacs will be called KDEBLLACS. The break-down of the number of BL Lacs candidates selected at different stages of the procedure as a function of the different radio surveys is displayed in Table 4, while the properties available for a sample of sources in the final KDEBLLACS catalog is shown in Table 4

5 Discussion

The full characterization of the two catalogs presented in this paper would require optical, spectroscopic follow-up observations to confirm the nature of the candidates, and will be discussed in a future paper. In this Section, we examine the global spatial and WISE photometric properties of the WIBRaLS2 and KDEBLLACS catalogs, and discuss how they compare with the most recent catalogs of Fermi $\gamma$-ray sources.

The Aitoff projections of the sky positions in galactic coordinates of the WIBRaLS2 and KDEBLLACS catalogs are shown in Figure 5. The coverage of the WIBRaLS2 catalog is mostly uniform across the sky (left panel in Figure 5), with the exception of a region along the galactic plane where radio sources are not available from any of the three surveys used. Two regions of higher density can be observed where the FIRST and SUMSS surveys overlap with the NVSS coverage, respectively north and south of the Galactic plane. In the right panel, the sky distribution of the KDEBLLACS catalog features prominently a lower density of sources associated to SUMSS radio counterparts and the lack of sources due to the Galactic latitude selection described in Section 4.3.

It is interesting to compare the regions of the WISE color space occupied by the WIBRaLS2 and KDEBLLACS catalogs. The location of the sources belonging to the above catalogs in the WISE color space is displayed in Figure 6. The two samples occupy partially overlapping but distinct regions of the three dimensional WISE color space. Since KDEBLLACS are not detected in the W4 passband by definition, their positions along the W3-W4 axis cannot be established, but the three-dimensional volume of the space potentially occupied by these sources can be determined by using the upper limits on their W4 brightness (i.e., lower limits on their W4 magnitude) available in the AllWISE catalog.

The two catalogs of blazar candidates described in this paper are complementary because their members can differ in spectral properties (see discussion in Section 4) and in brightness. Figure 7, that displays the histograms of the magnitude values for the three WISE filters W1, W2 and W3 for the WIBRaLS2 and KDEBLLACS samples together, confirms that the KDE-based method selects sources that are increasingly fainter than those in the WIBRaLS catalog as the wavelength increase going from filter W1 to W3.

5.1 Comparison with literature

One of the main goals of the production of catalogs of candidate blazars, like the two discussed in this paper, is the discovery of multi-wavelength counterparts of unassociated $\gamma$-ray sources (UGSs) observed by Fermi. While the thorough analysis that is necessary to reliably associate low energy candidate blazars to UGS (see, for example D’Abrusco et al., 2013; Massaro et al., 2015; Schinzel et al., 2015) is beyond the goals of this paper, as a consistency check and in order to assess the potential of the two catalogs here discussed to improve the characterization of the currently known $\gamma$-ray sources, we compared the WIBRaLS2 and KDEBLLACS catalogs with the most recent catalogs of sources detected by Fermi LAT.

We spatially crossmatched the WIBRaLS2 catalog with the 3FGL catalog (Acero et al., 2015) using the 98% elliptical uncertainty regions for 3FGL sources and a fixed positional uncertainty of 1″on the WIBRaLS2 radio coordinates. We found a total of 1049 matches, including the 666 sources in the locus sample selected as WIBRaLS2 members (see Section 3.3). Among the remaining 373 sources, 49 are unassociated; the 320 sources that are associated or identified with a known multi-wavelength counterpart in 3FGL are all classified, according to the CLASS1 parameter in 3FGL (Acero et al., 2015), as BL Lacs, FSRQs or BZU (Blazar candidate of Uncertain type), except for 3 rdg (radio-galaxies) and 3 ssrq (soft-spectrum radio quasars). The crossmatch between the KDEBLLACS catalog and 3FGL, using the same positional uncertainties, returned 186 matches (all distinct from the WIBRaLS2 crossmatches), with 57 unassociated sources. $\sim 95\%$ of the remaining 131 3FGL sources crossmatched with a KDEBLLACS source that are associated or identified, are classified as BZU or BL Lacs according to the 3FGL CLASS parameter. The total number of unique crossmatched sources from either the WIBRaLS2 or KDEBLLACS catalogs increases to 1757 (1404 and 353, respectively) when the preliminary LAT 8-year Point Source List (FL8Y)444https://fermi.gsfc.nasa.gov/ssc/data/access/lat/fl8y/ (which contains a total of 5523 $\gamma$-ray sources) is used, for a total of 152 unassociated sources. Out of the 1605 associated or identified sources, $\sim 99\%$ are classified as either BL Lacs, FSRQs or BZU.

Following the same approach described above, we also crossmatched the WIBRaLS2 and KDEBLLACS catalogs with the Third Catalog of Hard Fermi-LAT Sources (3FHL; Ajello et al., 2017a). The spatial crossmatch returns a total of 807 distinct matches, split in 647 from WIBRaLS2 and 160 from KDEBLLACS, with 33 unassociated 3FHL sources and the remaining 774 composed by BL Lacs (520, $\sim 68\%$ of the total), and blazars of uncertain type (BZU, $\sim 16\%$ of the total). The $83\%$ of the WIBRaLS2 members crossmatched with 3FHL sources are classified as candidate BZBs or MIXED candidates, based on their WISE colors.

It is also useful to compare the two catalogs of candidate blazars presented in this paper with the largest catalog of candidate HSPs available in the literature to date, namely the second WISE High Synchrotron Peaked blazar (2WHSP) catalog (Chang et al., 2017), with 1691 entries. The 2WHSP is an expansion of the 1WHSP catalog (Arsioli et al., 2015), and contains HSPs candidates drawn from the AllWISE catalog that can be associated to radio and X-ray counterparts. The 2WSHPs candidates are further selected by requiring that their radio-to-IR and IR-to-X-ray broad-band spectral slopes are consistent with those of known, confirmed HSPs sources (Chang et al., 2017), and that the peak frequency of the synchrotron emission component of their SEDs $\nu_{\mathrm{peak}}$ is $\!>\!10^{15}$ Hz (Chang et al., 2017). The 2WHSP contains 1691 unique candidates or confirmed HSPs: 460 of these sources are associated to AllWISE sources detected in four bands and 717 sources are associated to sources detected in the first three (W1, W2 and W3) filters. The remaining sources cannot be compared to the catalogs discussed in this paper because not detected in both the W3 and W4 passbands. We determined that 248 of the 460 2WSHPs sources detected in all AllWISE filters are also in the new WIBRaLS catalog, while 267 of the 717 sources not detected in W3 have been selected as KDEBLLACS. The main cause of difference between our catalogs and the 2WHSP is the maximum spatial radius used to crossmatch AllWISE sources with radio counterparts from one of the three radio surveys FIRST, NVSS and SUMSS. Chang et al. (2017) report that radio counterparts were selected within 0.3′ and 0.1′ for NVSS/SUMSS and FIRST catalogs, respectively. These radii are significantly larger than the radii used in this catalog, discussed in Sections 3.2 and 4.2 (10″.4, 3″.4 and 7″.4 for NVSS, FIRST and SUMSS respectively). As a consequence, only 311 AllWISE sources with a unique radio counterpart contained in the sample used to extract the KDEBLLACS can also be found in the 2WSHP catalog. Another possible source of difference between the KDEBLLACS and the 2WHSP catalogs is the extent of the region in the WISE color-color diagram where 2WHSP sources are located, which is significantly larger than the area where the BL Lacs candidates are selected from (see Figure 3).

5.2 Comparison with the WIBRaLS1 catalog

The second release of the WIBRaLS catalog, WIBRaLS2, described in this paper, contains 5025 candidate blazars also found in the first WIBRaLS (WIBRALS1 hereinafter) catalog (Paper I). Of the total 7885 members of WIBRALS1, 2830 sources are not included in WIBRaLS2 because their $q_{22}$ values (for $\sim 95\%$ of them) are larger than the new, more stringent thresholds adopted for BZQs in WIBRaLS2 (see Section 3.3). The increase in the size of the training set of confirmed $\gamma$-ray blazars used to define the WIBRaLS locus has led to an increase of the volume in WISE color space that has produced 4711 sources in the WIBRaLS2 not included in WIBRaLS1.

6 Summary and conclusions

In this paper, we present two new catalogs of blazar candidates, selected on the basis of their WISE mid-IR colors, their association to a radio counterpart and their radio-loudness. These new catalogs contain a combined total of 15196 candidate blazars, including both candidate FSRQs and BL Lacs, distributed on $\sim 90\%$ of the sky.

The second, enhanced release of the WIBRaLS catalog supersedes the original WIBRaLS1 (Paper I) catalog, which has been extensively employed to associate Fermi unidentified sources with their low energy counterparts. WIBRaLS2 contains candidate blazars drawn from AllWISE sources detected in all four WISE passbands with colors similar to those of spectroscopically confirmed, $\gamma$-ray emitting blazars, that are associated to radio counterparts and identified as radio-loud. Spectral classification as candidate BZBs, BZQs or Mixed candidate blazars, following the ROMA-BZCat (Massaro et al., 2015) terminology, and derived on the WISE color properties is also provided. WIBRaLS2 contains 9541 candidate blazars, a $\sim\!25\%$ increase over the first version of WIBRaLS.

The KDEBLLACS catalog complements the WIBRaLS2 by identifying BL Lacs candidates that, because of their typical low infrared-to-radio ratios and/or WISE brightness (cp. HBLs), have gone undetected in the W4 WISE filter and, as a consequence, cannot be considered for selection in WIBRaLS2. The KDEBLLACS members are required to be detected in the three WISE filters W1, W2 and W3 and are selected based on their positions in the W2-W3 vs W1-W2 WISE color-color diagrams using the KDE technique. They are also associated to a radio counterpart and identified as radio-loud according to the $q_{12}$ radio-to-MIR spectral parameter. We select 5579 sources in the KDEBLLACS.

Previous samples of candidate blazars selected on the basis of the IR colors, eventually combined with radio and/or multifrequency observations (see e.g. Massaro et al., 2013e; Arsioli et al., 2015; Maselli et al., 2015; Massaro et al., 2016) have been used by the Fermi-LAT collaboration for the preparation of the 3FGL, the Third Catalog of Active Galactic Nuclei (3LAC; Ackermann et al., 2015) and the Second and the Third Catalog of Hard Fermi-LAT Sources (2FHL and 3FHL; Ackermann et al., 2016; Ajello et al., 2017a, respectively).

This work will contribute to a more comprehensive understanding of the diverse and fascinating $\gamma$-ray sky as observed by Fermi. In particular, the community will benefit from the WIBRaLS2 and KDEBLLACS catalogs of candidate blazars presented in this paper and the subsequent programs of follow-up spectroscopic observations needed to confirm their nature and, possibly, determine their redshifts by using them to: (i) improve our knowledge of the luminosity function of BL Lacs (see e.g., Ajello et al., 2014); (ii) select potential targets for the Cherenkov Telescope Array (CTA) as shown by Massaro et al. (2013c); Arsioli et al. (2015) ; (iii) obtain more stringent limits on the dark matter annihilation in sub-halos (see e.g., Zechlin, & Horns, 2012; Berlin et al., 2014); (iv) search for counterparts of new flaring $\gamma$-ray sources (see e.g., Bernieri et al., 2013) and of high energy neutrino emission (see e.g., IceCube Collaboration et al., 2018); (v) test new $\gamma$-ray detection algorithms (see e.g., Campana et al., 2015, 2016, 2017); (vi) and, finally, perform population studies of the remaining UGSs (see e.g., Acero et al., 2013).