Radio Weak BL Lac objects in the Fermi era

TL;DR

This paper investigates the existence of radio weak BL Lac objects (RWBLs) in the Fermi era, proposing a new method using WISE mid-IR colors to distinguish them from similar objects and confirming two candidates as genuine RWBLs.

Contribution

It introduces a novel IR color-color diagnostic to identify RWBLs and confirms two new RWBLs associated with Fermi sources using multifrequency data.

Findings

Identified two new genuine RWBLs using IR color diagnostics.

Demonstrated the effectiveness of WISE colors in classifying BL Lac objects.

Confirmed spectral behavior consistent with BL Lac population for the candidates.

Abstract

The existence of "radio weak BL Lac objects" (RWBLs) has been an open question, still unsolved, since the discovery that quasars could be radio-quiet or radio-loud. Recently several groups identified RWBL candidates, mostly found while searching for low energy counterparts of the unidentified/unassociated gamma-ray sources listed in the Fermi catalogs. Confirming RWBLs is a challenging task since they could be confused with white dwarfs (WDs) or weak emission line quasars (WELQs) when there are not sufficient data to precisely draw their broad band spectral energy distribution and their classification is mainly based on a featureless optical spectra. Motivated by the recent discovery that Fermi BL Lacs appear to have very peculiar mid-IR emission, we show that it is possible to distinguish between WDs, WELQs and BL Lacs using the [3.4]-[4.6]-[12]$\mu$m color-color plot built using the WISE magnitudes when the optical spectrum is available. On the basis of this analysis, we identify WISE J064459.38+603131 and WISE J141046.00+740511.2 as the first two genuine RWBLs, both potentially associated with Fermi sources. Finally, to strengthen our identification of these objects as true RWBLs, we present multifrequency observations for these two candidates to show that their spectral behavior is indeed consistent with those of the BL Lac population.

Figures

Full text

††slugcomment: pre-proof version

1. Introduction

BL Lac objects constitute the most elusive class of radio-loud active galactic nuclei (AGNs; Urry & Padovani, 1995). Hosted in elliptical galaxies (see e.g., Falomo et al., 2014), their emission extends from radio frequencies up to the highest TeV energies and is characterized by rapid and large-amplitude flux variability, flat radio spectra, even below $\sim$1GHz (see e.g., Massaro et al. 2013a, ; Massaro et al. 2013b, ), peculiar infrared (IR) colors (Massaro et al. 2011a, ; D’Abrusco et al., 2012) and significant and variable optical polarization (i.e., $P\sim 10\%$) (see e.g., Andruchow et al., 2005; Smith et al., 2007, and references therein).

BL Lacs show featureless optical spectra, sometimes having weak emission lines of rest-frame equivalent width ($EW_{r}$) lower than 5Å or absorption lines due to their host galaxies, both superimposed to a blue continuum (see e.g., Laurent-Muehleisen et al., 1999). BL Lac objects are one of the two subclasses of the radio-loud AGNs known as blazars (see e.g., Massaro et al., 2009, for a recent review). The other subclass is constituted by the flat spectrum radio quasars (FSRQs), characterized by properties similar to those of BL Lacs but having a quasar-like optical spectra (i.e., with broad emission lines; see e.g., Véron-Cetty & Véron, 2000) and radio spectral indices111Spectral indices, $\alpha$, are defined by flux density, S${}_{\nu}\propto\nu^{-\alpha}$ lower than 0.5 between $\sim$300 MHz and a few GHz (see Massaro et al. 2014a, , and references therein).

Over the entire electromagnetic spectrum both BL Lacs and FSRQs feature a typical double-bumped spectral energy distribution (SEDs) with the first component peaking between the IR and the X-rays and the second one in the gamma-rays (see e.g., Padovani & Giommi, 1995; Paggi et al., 2009; Abdo et al. 2010b, ). Mid-IR colors clearly helps to determine the energy peak location of the low energy component (Massaro et al. 2011a, ). BL Lac SEDs and their broad-band features were described in terms of emission arising from particles accelerated in relativistic jets directly pointing towards the observer (Blandford & Rees, 1978).

BL Lacs are traditionally discovered using radio, X-ray surveys and/or their combination (see e.g., Stickel et al., 1991; Sambruna et al., 1996; Giommi et al., 2005), since in the optical energy range the lack of spectral features combined with their flux variability makes a color-color selection extremely challenging (Jannuzi et al., 1993; Londish et al., 2002; Collinge et al., 2005; Plotkin et al., 2008). However, recent observations to search for optical polarization in BL Lacs selected using deep optical survey as the Sloan Digital Sky Survey (see e.g., SDSS DR12 Alam et al., 2015, and references therein) and the Two-Degree Field Quasar Redshift Survey (2QZ Smith et al., 2005, and references therein) validated this strategy of investigating proper motion and broad-band color properties to minimize contamination by featureless Galactic objects while searching for quasi-featureless spectra (Smith et al., 2007).

At the current moment, the most efficient survey for discovering new BL Lacs is the Fermi gamma-ray one (Abdo et al. 2010b, ; Ackermann et al. 2011a, ; Nolan et al., 2012), in particular when combining gamma-ray observations with mid-IR colors analyses (Massaro et al. 2012a, ; Massaro et al. 2012b, ) and follow up optical spectroscopic observations (Paggi et al., 2013; Shaw et al., 2013; Massaro et al. 2015a, ; Massaro et al. 2015b, ). This occurs because BL Lacs are the largest known population emitting at MeV-GeV energies, constituting more than 30% of the associated sources in the latest release of the Fermi catalog (The Fermi - Large Area Telescope Third Source Catalog, 3FGL Acero et al., 2015). Moreover, recent optical spectroscopic campaigns aimed to confirm the nature of blazar candidates of uncertain type listed in the 3FGL showed that a significant fraction of those observed to date are also BL Lacs (see e.g. Álvarez-Crespo et al. 2016a, ; Álvarez-Crespo et al. 2016b, ; Álvarez-Crespo et al. 2016c, ).

However, nowadays, all the known Fermi BL Lacs have radio counterparts and radio surveys are therefore crucial to associate $\gamma$-ray objects with their counterparts. Thus, one of the major open issues for this rare AGN class is proving the existence of radio-weak BL Lac objects (RWBLs). We indicate them as radio-weak sources instead of using the term “radio-quiet”, generally adopted for the AGN unification scenario. This designation is more useful because radio follow up observations could reveal some radio emission. However, the current lack of counterparts in the major radio surveys available to date place them in the tail of the distribution of the IR-to-radio flux ratios for the whole BL Lac population.

At the current moment only a handful of sources have been claimed to be RWBLs on the basis of their broad band properties as featureless optical spectra combined with X-ray emission and and their lying within the positional uncertainty regions of unidentified/unassociated gamma-ray sources (UGSs; Paggi et al., 2013; Ricci et al., 2015; Landoni et al., 2016).

Their existence will have a crucial impact (i) on the estimates of their source counts and luminosity functions (Ajello et al., 2012, 2014), a task already complicated by the challenge of measuring their redshifts (Falomo et al., 2014), (ii) on the gamma-ray association methods all currently based on radio surveys (see e.g, Ackermann et al., 2015, and references therein), if proved that RWBLs could be associated to Fermi sources, as well as (iii) on the unification scenario of AGNs (Urry & Padovani, 1995), to name a few implications. An accurate accounting is crucial to compute correctly their contribution to the extragalactic gamma-ray background for which, together with FSRQs (Ajello et al., 2012) and radio galaxies (Inoue, 2011; Di Mauro et al., 2014; Massaro & Ajello 2011b, ) they provide the most significant ones. Additional implications are related to jet physics. Radio emission in extragalactic radio-loud AGNs permits us to investigate the “large-scale” (i.e. from pc to Mpc size) evolution of relativistic jets and their interactions with the environments (see e.g., Liuzzo et al., 2013; Nulsen et al., 2013; Fabian, 2012, and references therein). The existence of RWBLs could be crucial to address different phases of jet formation in cases when relativistic jets are not able to form extended radio structures (Blandford & Königl, 1979; Begelman et al., 1984), as it seems from low radio frequency (i.e., below $\sim$1 GHz) spectral studies (see Massaro et al. 2013b, ; Nori et al., 2014; Giroletti et al., 2016, and references therein).

There are two major problems to solve before claiming the existence of RWBLs, even when featureless optical spectra are available. They can be confused with or misinterpreted as hot white dwarfs (WDs; see e.g. Kepler et al., 2015) and/or weak emission line quasars (WELQs; Diamond-Stanic et al., 2009; Plotkin et al., 2010, and references therein) when optical spectroscopic data are available. Optical polarimetric observations could shed light on their nature. However, unfortunately, optical polarization in BL Lacs can be variable and it is not always observed thus not representing the most definitive stamp to confirm their nature (see e.g., Londish et al., 2004).

In the present paper we explore the use mid-IR colors as a diagnostic tool to disentangle between BL Lacs, WDs and WELQs. This new tool could potentially open a new window to search for RWBLs. This investigation is motivated by the recent discovery that Fermi BL Lacs show peculiar mid-IR colors (see also, Massaro & D’Abrusco, 2016, and references therein), result made possible by observations performed with the WISE all-sky survey (Wright et al., 2010), covering a region well isolated from that of other Galactic and extragalactic objects. The main explanation underlying this distinction is the BL Lac non-thermal IR emission (see e.g., D’Abrusco et al., 2012, 2013).

This manuscript is organized as follows. In § 2 we present all the samples used to perform our mid-IR color-color analysis. In § 3 we compare the mid-IR colors of BL Lacs with those of WELQs and WDs, while in § 4 we search for RWBLs. § 5 is dedicated to the analysis of the IR variability for the best candidates selected on the basis of the mid-IR colors ad the optical spectra. Then in § 6 we also investigate the SWIFT observations available for these candidates while in § 7 we discuss on the SEDs. Summary and conclusions are given in § 8

Unless stated otherwise, for our numerical results, we use cgs units. WISE magnitudes used here are in the Vega system and are not corrected for the Galactic extinction. As shown in our previous analyses, such corrections mostly affect the magnitude at 3.4$\mu$, for sources lying at low Galactic latitudes and it ranges only between 2% and 5% (see e.g. D’Abrusco et al., 2014), thus being negligible. In addition this allows us to compare plots in the present analysis with those of previous investigations (Massaro et al. 2013c, ).

2. Sample selection

To perform our analysis we built three main samples. The first two composed of WDs and WELQs, respectively; both used in comparison with the third one of Fermi blazars. We searched for mid-IR counterparts of all samples in the latest release of the WISE catalog222http://wise2.ipac.caltech.edu/docs/release/allwise/. All our final samples of WDs, WELQs and gamma-ray blazars only include sources having a WISE counterpart within 3″.3 that is detected at 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m so allowing us to build at least one of the mid-IR color-color plots used in our previous blazar analyses (D’Abrusco et al., 2014; Massaro & D’Abrusco, 2016). The choice of 3″.3 association radius is the same adopted in our previous investigations (D’Abrusco et al., 2013) and it was computed on the basis of Montecarlo procedure corresponding to a probability of having spurious associations of the order of $\sim$1% (see also D’Abrusco et al., 2014, for additional details).

2.1. BL Lac objects and Flat Spectrum Radio Quasars

The first sample used in our analysis lists all the blazars, both BL Lacs and FSRQs, included in the Roma-BZCAT333http://www.asdc.asi.it/bzcat/ (Massaro et al. 2015c, ) and having a counterpart in the latest release of the Fermi catalog (Acero et al., 2015). This sample has been recently used to perform the analysis on the IR–gamma-ray connection (Massaro & D’Abrusco, 2016). According to the nomenclature proposed in the Roma-BZCAT hereinafter we will refer to BL Lac objects as BZBs and to the FSRQs as BZQs. We choose this sample because it has the largest fraction of BZBs with a mid-IR counterpart in WISE corresponding to 603 out of 610 (i.e., $\sim$99%) within our association radius; this fraction is similar to that occurring for the BZQs with 419 out of 426 (i.e., $\sim$98%).

2.2. White Dwarfs

Hot WDs could be misclassified as BZBs because they could show almost featureless optical spectra and X-ray radiation (see e.g., Heise, 1985; Bilḱová et al., 2010). Thus in the present analysis we want to show that the mid-IR colors of the WISE counterparts of WDs is different from that of BZBs. We expect that for WDs mid-IR emission arises from the presence of a companion brown dwarf or the presence of a dusty disk444http://www.stsci.edu/$\sim$debes/wired.html, both expected to have thermal IR colors (Debes et al., 2011) while those of BZBs are dominated by the non-thermal synchrotron radiation of particles accelerated in their jets.

To achieve this goal we considered Sloan Digital Sky Survey DR7 White Dwarf Catalog (SDSSDR7WD) (Kleiman et al., 2013). This new catalog of spectroscopically confirmed WDs from the Sloan Digital Sky Survey (SDSS) Data Release 7 (DR7) listing more than 19000 stars. Crossmatching the SDSSDR7WD catalog with the latest release of the WISE catalog we found 4125 have a unique mid-IR potential counterpart within 3″.3. However, only for 255 out of 4125 sources their mid-IR counterparts were detected by WISE in all three of its short-wavelength bands, and so we only do a comparison with this subset.

We assume that mid-IR correspondences found for the SDSSDR7WD are the real counterparts of the WDs. However in the worst scenario if none of the WDs listed in the SDSSDR7WD has mid-IR emission this will automatically confirm that their mid-IR spectral behavior is different from that of Fermi BZBs known to date. We also computed the chance probability of having spurious associations within 3″.3 when crossmatching the AllWISE source catalog with the SDSSDR7WD; the result is less than 4% (see e.g., D’Abrusco et al., 2013, for more details on the method adopted for the pchance calculation).

To verify the reliability of our choice for the WD sample we also used the latest release of the McCook & Sion White Dwarf Catalog555http://www.astronomy.villanova.edu/WDCatalog/index.html available to date (see McCook & Sion, 1999, for details). This catalog has been also used to search for dusty disks around WDs (Hoard et al., 2013) using WISE observations. However in this sample the number of sources detected at 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m for which is possible to built the mid-IR colors is only 91. The chance probability of spurious associations between the McCook and Sion WD catalog and the AllWISE survey is lower than 1%.

2.3. Weak emission line quasars

In the last two decades an intriguing and also elusive subclass of radio-quiet quasars was discovered: the weak emission line quasars (WELQs), featuring exceptionally weak or completely missing broad emission lines in the ultraviolet (UV) rest-frame energy range (see e.g, Collinge et al., 2005; Shemmer et al., 2006; Diamond-Stanic et al., 2009).

The original criterion to define WELQs with respect to other members of the same class (i.e., normal quasars) is based on the $EW_{r}$ distribution of the Ly$\alpha$ $\lambda$1216+N V $\lambda$1240 blend for a sample of quasars in the Sloan Digital Sky Survey at redshift $z>$3, being lower than 15$\AA$ as occurs for the $>3\sigma$ weak tail of the $EW_{r}[Ly\alpha+NV]$ distribution. Using this definition, a first sample of $\sim$70 WELQs was defined (Diamond-Stanic et al., 2009). Since at low $z$ it is not always possible to apply the original criterion on the $EW_{r}[Ly\alpha+NV]$ because these emission lines lie outside the optical energy range covered by the SDSS spectrograph an equivalent criterion was defined selecting WELQs on the basis of the $EW$ of Mg II $\lambda$2800, C III] $\lambda$1909, and/ or C IV $\lambda$1549 (see e.g., Meusinger & Balafkan, 2014).

When WELQs were discovered the possibility of being classified BZBs was explored and discarded (see e.g. Plotkin et al., 2015, for a recent review), however a comparison of their mid-IR properties with those of confirmed BZBs has not yet been completely performed. Here we selected a sample of WELQs to show that their mid-IR colors are generally similar to those of normal quasars but not to BZBs. Additional analyses based on the SEDs and/or broadband spectral indices (i.e., radio-to-optical $\alpha_{ro}$ and optical-to-Xray $\alpha_{ox}$) were also used to distinguish between WELQs and BL Lacs (Plotkin et al., 2010; Lane et al., 2011; Wu et al., 2012) when multifrequency observations are available. Thus in the following sections we also investigate the broad band SEDs of RWBL candidates selected on the basis of the combination of the mid-IR colors and the optical spectra.

We considered two catalogs of WELQs, spectroscopically selected from the SDSS. The first catalog lists 73 WELQs in the redshift range between 3.03 and 5.09 (Diamond-Stanic et al., 2009, hereinafter DS9 subsample) while the second lists 46 sources with 0.614$<z<$3.351 (Meusinger & Balafkan, 2014, hereinafter M14 subsample). We excluded 6 sources out of these 119 WELQs, namely: SDSS J074451.36+292006.0 and SDSS J092145.37+233548.1 from the M14 subsample and SDSS J141318.86+450522.9 from DS9 subsample because they are all listed in the Roma-BZCAT and SDSS J104831.29+211552.2, SDSS J115326.70+361726.3, SDSS J144204.03+132916.0 again from the M14 subsample since they show a flat radio spectrum that is widely interpreted as jet signature so suggesting that their quasi-featureless optical continuum marks non-thermal synchrotron emission.

Only 90 out of the remaining 113 sources arising from the combination of the DS9 and the M14 subsamples have a counterpart in the AllWISE catalog. These 90 sources thus constitute our WELQ reference sample to carry out the comparison with BZB mid-IR colors. We estimated the chance probability of spurious associations between the SDSSDR7WD catalog and the WISE survey as being lower than 1%.

3. BZBs, WDs and WELQs at mid-IR wavelengths

In Fig. 1 we show the mid-IR colors built with the in the 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m magnitudes of the BZBs and the BZQs in our Fermi blazar sample in comparison with generic infrared sources selected at high Galactic latitudes ($|b|>$50°), to highlight their peculiar spectral behavior. Since we know that BZQs generally lie in the region occupied by “normal” quasars, at the top end of the so called “WISE Gamma-ray Strip” (Massaro et al. 2011a, ), to highlight the mid-IR color region where BZBs reside, we will only consider RWBL candidates those WISE sources lying in the dashed polygon shown in Fig. 1 where the contamination of the BZB region by BZQ is lower than 5%, lacking a radio counterpart. It is worth noting that the selected region also includes TeV BL Lacs.

Fig. 2 shows the comparison between WISE sources associated with WDs and our blazar sample. It is evident that only 4 out of the 255 sources in the SDSSDR7WD catalog have mid-IR colors consistent with the main BZB region of the “WISE Gamma-ray Strip”, while the remaining $\sim$99% of the WDs lie well separated from this region. Three out of four sources are extragalactic, while WISE J122859.87+104032.8 has a clear WD optical spectrum, making it the only WD contaminant of the BZB region in the mid-IR color-color diagram. Moreover WISE J122859.87+104032.8 has a negative value of the $u-r$ a situation that never occurs for the known BZBs (Massaro et al. 2012c, ; Massaro et al. 2014b, ). WISE J172633.50+530300.3 also shows a featureless optical spectrum with a peculiar continuum shape; this could be due to the presence of a nearby source clearly visible in the SDSS image but not in the WISE one. This implies a relatively small contamination of the mid-IR color selection.

Both WISE J001736.91+145101.9 and WISE J121348.83+642520.0 are BL Lacs with a radio counterpart in addition to WISE J121348.83+642520.0 listed in the latest release of the Roma-BZCAT catalog.

To complete our comparison between WDs and BZBs we also show the mid-IR color-color plot with the 91 sources selected out of the McCook & Sion WD catalog in Fig. 2. There are also in this case only two objects that could be confused as BZBs on the basis of their mid-IR colors, namely WISE J085506.13+063904.1 and WISE J115642.21+130603.9, both lying in the footprint of the SDSS. Thus in both cases we checked the optical images using the SDSS Navigate tool666http://skyserver.sdss.org/dr12/en/tools/chart/navi.aspx. We found that about 6″6̇ from the location of WISE J085506.13+063904.1 there is a clear WD, namely SDSS J085506.62+063904.7, for which an optical spectrum is available being most likely the star listed in the McCook & Sion WD catalog and lacking a WISE counterpart. This mis-association is due to the positional uncertainty of the sources in the McCook & Sion WD catalog. Thus we conclude that WISE J085506.13+063904.1 is a simple contaminant unrelated to the WD selection. Then for the latter, WISE J115642.21+130603.9 the lack of its optical spectrum combined with the uncertainties of the SDSS photometric observations does not allow us to verify its nature. However, this object is classified as galaxy thus unlikely to be a WD that can contaminate our RWBL selection.

The comparison between mid-IR colors of BZBs and WDs prove the usefulness of the mid-IR color-diagram to distinguish the two classes while searching for RWBLs, once optical spectra are available.

Fig. 3 shows the mid-IR comparison between the WELQs and the BZBs. Here too is clear the neat separation between these two classes. It is also worth noting that sources in the DS9 subsample appear to lie in a different region of the color parameter space with respect to the M14 subsample. The reason underlying this effect is due to their different redshift distributions. We note that there is also a marginal difference between the mid-IR colors of the WELQs in the M14 sample and that of the BZQs but it is less useful to search for RWBLs hidden among this sample.

In this case we highlight the fact that two WELQs lie below the threshold marking the BZB region of the color-color plot and additional three objects are close to it. These two sources are WISE J001514.87-103043.4 and WISE J150427.69+543902.2; both showing MgII$\lambda$2798 emission line of rest-frame equivalent width $EW_{r}$ consistent with the 5Å threshold, usually adopted to classify BL Lacs, within 3$\sigma$, thus being “border line cases”. The remaining three sources, belonging to the M14 subsample (see Fig. 4, where they are shown in magenta), all have a radio counterpart at 1.4 GHz in the Faint Images of the Radio Sky at Twenty centimeter (FIRST; Helfand et al., 2015) catalog, but the lack of additional radio data did not allow us to verify the flatness of their radio spectra, although they are definitively not RWBLs. The location of all five sources with respect to that of the BZBs in the mid-IR color-color plot is shown in Fig. 4.

Finally, we conclude that both WDs and WELQs have distinct mid-IR colors with respect to the BL Lac population and [3.4]-[4.6]-[12]$\mu$m color-color diagram can be used as diagnostic tool to confirm the source nature having featureless optical spectra. Sources that appear to be contaminants of the mid-IR selection procedure are ruled out when WISE data are combined with optical spectroscopic observations.

4. Radio Weak BL Lacs

In the literature there are five sources tentative classified as RWBLs.

2QZ J215454.3-305654 (Londish et al., 2002) for which several multifrequency follow up observations have been performed to date was the first. However its nature is still uncertain since it could be a RWBL or a “lineless” quasar (Londish et al., 2004). 2. 2.

Paggi et al. (2014), carried out optical spectroscopic observations of blazar candidates, potential UGS counterparts and discovered WISE J064459.38+603131. This mid-IR sources lies within the positional uncertainty region of 2FGL J0644.6+6034.7 (3FGL J0644.6+6035). This object shows a single optical emission line of EW consistent with the 5Å threshold, arbitrarily set in the literature for the BL Lac class definition, similar to the two cases discussed at the end of § 3. 3. 3.

In 2015 WISE J173052.85-035247.2 was discovered during the optical spectroscopic follow up observations of UGSs (Ricci et al., 2015). This source could be the potential counterpart of 2FGL J1730.6-0353 (a.k.a. 3FGL J1730.6-0357). 4. 4.

More recently Landoni et al. (2016) claimed that also SDSS J004054.65-0915268 could be a tentative RWBL. 5. 5.

Finally, Marchesini et al. (2016) discovered the last known candidate of RWBL, WISE J141046.00+740511.2. This source was also found while carrying out spectroscopic observations of UGS listed in the First Fermi-LAT Catalog of Sources above 10 GeV (1FHL; Ackermann et al., 2013) selected on the basis of the SWIFT X-ray observations (Landi et al., 2015). It has been indicated as potential counterpart of 1FHL J1410.4+7408 (a.k.a. 3FGL J1410.9+7406).

In Fig. 4 we show the positions of the five RWBLs overlay to the “WISE Gamma-ray Strip” subregion of BZBs. SDSS J004054.65-0915268 (a.k.a WISE J004054.63-091526.8) has an upper limit at 12$\mu$m and being a high $z$ source it seems to be more consistent with the location of the DS9 subsample of WELQs rather than the BZB population. A similar situation occurs for 2QZ J215454.3-305654 (a.k.a. J215454.32-305653.2) that only a follow up observations at 12$\mu$m could potentially reveal if the source is consistent with the mid-IR colors of the BZBs. WISE J064459.38+603131.7, being selected on the basis of its mid-IR colors is definitively consistent with the BZB population on the 3.4$\mu$m - 4.6$\mu$m - 12$\mu$m diagram, while WISE J173052.85-035247.2 appears to have such colors more similar to the WELQs listed in the M14 subsample.

The remaining case of WISE J141046.00+740511.2 is the most intriguing one. This object is potentially associated to a 1FHL object and its potential WISE counterpart lies in the region of the “WISE Gamma-ray Strip” where all known TeV BL Lacs are located (Massaro et al. 2013d, ; Arsioli et al., 2015). The lack of any radio counterpart in both the NRAO VLA Sky Survey (NVSS; Condon et al., 1998) and the FIRST catalog (Helfand et al., 2015) that cover the field of WISE J141046.00+740511.2 makes this source the most promising candidate to be the first, genuine, RWBL. Nevertheless is it worth emphasizing that the distance between this mid-IR source and the gamma-ray position of 3FGL J1410.9+7406 is $\sim$90″well below the average of the associated Fermi blazars listed in the Roma-BZCAT corresponding to $\sim$130″. At the angular separation of $\sim$90″, the association probabilities of counterparts listed in the latest releases of the Fermi catalogs is typically grater than 99.5% (Ackermann et al., 2015).

Applying this mid-IR color-based analysis, in combination with optical spectra, we conclude that WISE J064459.38+603131 and WISE J141046.00+740511.2, lying within the positional uncertainty regions of two distinct UGSs listed in the 3FGL, are the first two RWBLs detected so far. Both of them can also be associable with 1FHL sources thus strengthening our interpretation, since the high energy sky above 10 GeV is mostly dominated by BZBs. In the following we will also show that additional multifrequency observations can strengthen our results.

5. mid-IR variability

BL Lac objects are known to be among the most variable extragalactic sources at all frequencies, showing variations either in spectral shape or in flux/luminosities. This peculiar behavior could in principle mine the mid-IR color selection proposed here and in our previous works. However this situation does not seem to occur.

BL Lacs appear to be variable in the mid-IR band (see e.g., D’Abrusco et al., 2012) although their color variations still reside within the boundaries of the “WISE Gamma-ray Strip” (see Buson et al., 2016, for more details). Here we show an example of two BZBs: 5BZBJ1542+6129 and 5BZBJ1959+6508, the latest also detected at TeV energies (Aharonian et al., 2003), for which more than 50 WISE observations are available, when both sources are detected at least in the first three mid-IR filters at 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m. In Fig. 5 it is clear how both BZBs during their single-epoch WISE observations move in the mid-IR color-color plot in the same region occupied by the BL Lac objects. In addition their colors are also consistent with the boundaries chosen for our RWBL selection procedure.

According to the WISE catalog one of our best RWBL candidate: the WISE J141046.00+740511.2 is flagged as variable in the mid-IR. This is in agreement with the fact that such emission it is unlikely due to dust emission, as occurs in normal quasars or in white dwarfs (Debes et al., 2011). Fig. 6 shows the WISE magnitude at 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m as function of time for two different epochs; mid-IR variability, not only between the two observing periods, but also in each single set of observations is found. We also present the lightcurve for the second RWBL candidate WISE J064459.38+603131 that clearly shows infrared variability being detected more times at 12$\mu$m during the WISE single-epoch observations with respect to WISE J141046.00+740511.2.

In Fig. 7 we also show their mid-IR color variability as function of time. In the same figure we also present the single epoch observations when both sources are detected in the first three WISE filters overlaid to the “WISE Gamma-ray Strip”. As for the other two BZBs mentioned above these variations are clearly consistent with the region chosen to select RWBLs for WISE J141046.00+740511.2 while for WISE J064459.38+603131 there is partial agreement. It is worth mentioning that the location of integrated mid-IR colors for both WISE J064459.38+603131 and WISE J141046.00+740511.2 could appear different from that where they exhibit variability in the mid-IR color diagram. The underlying reason is that there are many WISE single-epoch observations when sources (see Fig. 6) are not detected that contribute to the color change for the integrated values. It is also worth noting that for the two RWBL candidates the uncertainties on the WISE magnitudes are greater than for the two bright BZBs selected above as example.

On the basis of this investigation of mid-IR variability and spectral behavior we conclude that the evidence supports the interpretation of these sources as RWBLs.

6. SWIFT observations

Several analyses were carried out to study the nature of WELQs in the X-rays (see e.g., Shemmer et al., 2009; Wu et al., 2012). These observations were also used to build their broad band spectral indices and/or their SEDs (see e.g., Plotkin et al., 2010) to verify their potential similarities with BL Lac objects. Since the two sources identified here as RWBLs were discovered during a campaign aimed to search for UGS counterparts and all the Fermi UGSs are in a SWIFT observational program (Stroh et al., 2013)777http://www.swift.psu.edu/unassociated/, we reduced the SWIFT data available to date to build the RWBL SEDs.

WISE J064459.38+603131 was observed by SWIFT both with the X-Ray Telescope (XRT), for $\sim$ 3 ksec, and with the Ultraviolet/Optical Telescope (UVOT) instruments on board but for the latter the only ultraviolet filter used was M2. On the other hand, for WISE J141046.00+740511.2, there are optical data in the U band and all the UV filters available, together with about 11 ksec of X-ray data.

Here adopted the standard XRT and UVOT data reduction procedure, extensively used in our previous analyses (see e.g. Massaro et al. 2008a, ; Maselli et al., 2013, 2016, and references therein) while only the very basic details are described below. XRT data were processed with the XRTDAS software package (v.3.0.0) developed at the ASI Science Q9 Data Center (ASDC) and distributed within the HEASOFT package (v.6.18) by the NASA High Energy Astrophysics Archive Research Center (HEASARC). All the XRT observations were carried out in the photon counting (PC) readout mode. Event files were calibrated and cleaned applying standard filtering criteria with the XRTPIPELINE task and using the latest calibration files available in the Swift CALDB. Events in the energy range 0.3–10 keV with grades 0-12 were selected in the analysis according to the standard guidelines. Exposure maps were also created using XRTPIPELINE. The SWIFT observations available for both sources have variable exposures, we merged cleaned event files obtained with the above procedure using XSELECT task and considering only SWIFT observations with telescope aim point falling in a circular region of 12 arcmin radius centered in the median of the individual aim points, to get a uniform exposure. The corresponding merged exposure maps were then generated by summing the exposure maps of the individual observations with XIMAGE package. Source detection of in the XRT images was performed using the detection algorithm DETECT.

We collected 13 counts for WISE J064459.38+603131 in the merged XRT event file while there are 36 counts for WISE J141046.00+740511.2, making both sources well detected above the 5$\sigma$ threshold. However the poor number of counts available, corresponding to a count rate of 0.004 cts/s and 0.003 cts/s, respectively, did not permit us to carry out a detailed spectral analysis. Thus assuming an absorbed power-law model with a photon index of 2.2 and with the Galactic hydrogen column density (Kalberla et al., 2005) and estimated the unabsorbed X-ray flux of in the 0.5-10 keV for WISE J064459.38+603131 and of for WISE J141046.00+740511.2 using WebPIMMS888https://heasarc.gsfc.nasa.gov/cgi-bin/Tools/w3pimms/w3pimms.pl.

Details on the UVOT reduction procedure used here are described in Massaro et al. (2008b) and Maselli et al. (2008). We performed the photometric analysis using a standard task UVOTSOURCE available in the HEASOFT package (v.6.18). Counts were extracted from a 5 arcseconds radius circular region in the V, B, and U filters as well as for the UV filters (UVW1,UVM2, and UVW2). The count rate was corrected for coincidence loss and the background subtraction was performed by estimating its level in an offset annular region centered on the source position with inner radius of 7.5 arcseconds and outer radius of 15 arcseconds. Since estimate of magnitude errors is complex due to possible instrumental systematics and calibration, here we only showed a typical uncertainty of 5% for all filters.

The values of the observed AB magnitude for WISE J064459.38+603131 in the UV band at 1.36$\times$1015Hz is 19.69 mag (M2 UVOT filter) while for WISE J141046.00+740511.2 the AB magnitudes measured in the UVOT filters are: 19.23 mag (U), 19.85 mag (W1), 20.24 mag (M2) and 20.31 mag (W2), respectively.

In Fig. 8 we show the XRT event file in the 0.3-10 keV energy range and one of the UVOT images for both the RWBL selected. We overlaid to the SWIFT images the positional uncertainty region, at 95% level of confidence, from the 2FGL and the 3FGL catalogs for the nearby UGSs. WISE J064459.38+603131 lies within the error ellipse of the 2FGL, as occurred when it was selected (Massaro et al. 2012a, ) but not within the one of the 3FGL, while WISE J141046.00+740511.2 that is associable to the 3FGL source, it was not found a potential counterpart of the 2FGL UGS.

SWIFT observations, both in the UV and in the X-ray energy ranges are used in the following to build the SEDs for both the RWBLs selected.

7. Spectral Energy Distributions

For the two RWBL selected we also built the SEDs to check if they are similar or not to those of WELQs (Lane et al., 2011). Analysis of the SEDs is equivalent to comparing broad band spectral indices that could give an additional indication of the differences between BL Lacs and WELQs (see e.g Wu et al., 2012).

In the case of WISE J064459.38+603131, the source has an upper limit in the NVSS catalog of 0.13 mJy/beam while WISE J141046.00+740511.2 1.69 mJy/beam. The former is detected in all the WISE bands and in the 2MASS. From the SWIFT observations we also got an estimate of the UV flux density in the M2 filter and the X-ray flux. The latter source still has a WISE counterpart detected at all mid-IR frequencies but it does not have a 2MASS correspondence. SWIFT observations provided us flux density estimates in the optical (U) and in all the UV filters (W1, M2, W2) together with the X-ray one.

We corrected the all the infrared, optical and ultraviolet flux densities for the Galactic absorption following the Draine (2003) and the reddening corrections of Schlafly & Finkbeiner (2011), with the only exception for the WISE magnitudes at 12$\mu$m and 22$\mu$m.

In Fig. 9 we show the SEDs for both WISE J064459.38+603131 and WISE J141046.00+740511.2. For the former SED the infrared emission appear to be compatible with a log-parabolic shape (Massaro et al., 2004, 2006) with the only exception of the UVOT flux density in the M2 filter. This could be due to the presence of the host galaxy and/or to variability since UV and IR observations were not taken simultaneously and both of them are integrated over different epochs. The latter source: WISE J141046.00+740511.2 has also an SED well described by a log-parabolic shape but in this case the SWIFT optical and UV data seems to mark the presence of an host galaxy. In both panels of Fig. 9 we also report the unabsorbed X-ray flux computed from the source counts of the XRT event file. Since this X-ray flux is an estimate we neglect second order corrections to convert the X-ray flux into a $\nu F_{\nu}$ SED value.

The agreement between the multifrequency data and the broad band log-parabolic SED support the interpretation of WISE J064459.38+603131 and WISE J141046.00+740511.2 as RWBLs.

8. Summary and Conclusions

• •

The mid-IR colors of BL Lacs appear to be clearly different from those of WDs and WELQs. This strongly indicate that the WISE color-color diagram built with the mid-IR magnitudes at 3.4$\mu$m, 4.6$\mu$m and 12$\mu$m together with optical spectroscopic observations are a powerful diagnostic tool to distinguish these three source classes. We only found few contaminants within the WD and the WELQ catalogs used in our analysis. However, for all these cases, the combination of mid-IR and optical spectroscopic data allowed us to exclude them as BL Lacs.

• •

Applying this mid-IR color-based analysis we conclude that WISE J064459.38+603131 and WISE J141046.00+740511.2, lying within the positional uncertainty regions of two distinct UGSs, are the first two RWBLs detected so far. Both of them are also associable with 1FHL sources thus strengthening our interpretation, since the high energy sky above 10 GeV is mostly dominated by BZBs. Both sources appear to be variable in the mid-IR as expect for BZBs and they also show a broad-band SEDs, built thanks to the WISE and the SWIFT observations, consistent with a log-parabolic shape.

• •

In addition to them 2QZ J215454.3-305654 could also potentially be a RWBL, but lacking of the mid-IR detection at 12$\mu$m it is not possible to confirm its nature using the mid-IR color-color plot. In addition the two contaminants: WISE J001514.87-103043.4 and WISE J150427.69+543902.2, found in the WELQ sample used here, deserve a deeper investigation to clarify their nature.

We highlight that the strategy used to search for potential counterparts of UGSs via X-ray follow up observations has been crucial in both cases to discover WISE J064459.38+603131 and WISE J141046.00+740511.2, thus emphasizing that combining X-ray and mid-IR observations could be crucial for future searches of low energy counterparts of gamma-ray sources. However, we also stress that such discovery has been done thanks to optical spectroscopic follow up campaigns that finally revealed the source nature. Since both these two sources are not only associable with 3FGL objects but also with Fermi sources detected above 10 GeV their existence challenges current gamma-ray association methods mostly based on radio surveys that do not allow to catch these potential counterparts. The existence of RWBLs will also strong impacts on jet physics since they could allow us to study sources whose emission is dominated by relativistic jets that do not tend to form extended radio structures as their parent population of radio galaxies.

Finally, we remark that deeper radio follow up observations than those available to date could potentially reveal their low energy emission at GHz frequencies. However the ratio between the infrared and the radio emission of so called RWBLs, being below the flux limit of current radio surveys is unexpectedly unusual. In Fig. 10 we show the distribution of such ratio to highlight how the two selected RWBLs lie in the tail of the BL Lac population being below the 0.5 value.

We thank the anonymous referee for useful comments that led to improvements in the paper. F.M. gratefully acknowledges the financial support of the Programma Giovani Ricercatori – Rita Levi Montalcini – Rientro dei Cervelli (2012) awarded by the Italian Ministry of Education, Universities and Research (MIUR). H.A.S. acknowledges partial support from NASA Grants NNX15AE56G and NNX14AJ61G. This research has made use of data obtained from the high-energy Astrophysics Science Archive Research Center (HEASARC) provided by NASA’s Goddard Space Flight Center. The NASA/IPAC Extragalactic Database (NED) operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. Part of this work is based on the NVSS (NRAO VLA Sky Survey): The National Radio Astronomy Observatory is operated by Associated Universities, Inc., under contract with the National Science Foundation and on the VLA Low-frequency Sky Survey (VLSS). This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, funded by the National Aeronautics and Space Administration. TOPCAT999http://www.star.bris.ac.uk/$\sim$mbt/topcat/ (Taylor, 2005) for the preparation and manipulation of the tabular data and the images.