THE NuSTAR EXTRAGALACTIC SURVEYS: OVERVIEW AND CATALOG FROM THE COSMOS FIELD

We performed the NuSTAR data reduction using HEASoft v.6.15.1 and the NuSTAR Data Analysis Software NuSTARDAS v.1.3.1 with CALDB v.20131223 (Perri et al. 2013).²⁷ We processed level 1 products by running nupipeline, a NuSTARDAS module which includes all the necessary data reduction steps to obtain a calibrated and cleaned event file ready for scientific analysis. A further bad pixel file was also used in the filtering to flag outer edge pixels with high noise.

Light curves were produced using nuproducts in the low‐energy range (3.5–9.5 keV), with a 500 s interval binning. Twenty-one observations (fields 13, 20, 21, 55 , 56, 57, 58, 70, 83, 87, 88, 89, 91, 92, 96, 97, 110, 111, 116, 119, 120) were affected by abnormally high background radiation due to solar flares. Intervals with count rates >0.1 cnt s⁻¹ were cleaned by reprocessing the 21 observations, applying good time interval files created using standard HEASoft tools. The resulting loss of exposure time corresponds to 2%–9% of the exposure in these twenty-one observations, for a total of 21 ks. This is less than 1% of the total NuSTAR COSMOS exposure.

Exposure maps were created using nuexpomap, which computes the net exposure time for each sky pixel, for a given observation. In order to reduce the calculation time, we binned the maps using a bin size of 5 pixels. The exposure map accounts for bad and hot pixels, detector gaps, attitude variations and mast movements. Ideally, given the strong energy dependence of the effective area, exposure maps at each energy should be created and then summed together using weights based on the effective area at each energy. This is more complicated if the typical source spectrum is not constant with energy. This procedure is computationally expensive and it can be simplified by convolving the instrument response with a specific model for the incident spectrum. The energies at which the exposure maps were created were computed by weighting the effective area with a power law spectrum with Γ = 1.8 (see Section 4.1 for spectral choice motivation). The mean, spectrally weighted energies are 5.42 keV for the 3–8 keV band, 13.02 keV for the 8–24 keV band, and 9.88 keV for the 3–24 keV band.

The 121 vignetting corrected exposure maps were mosaicked using the HEASoft XIMAGE tool. The effective exposure time, corrected for vignetting, is plotted versus area in Figure 1. A small difference of <5% in exposure per area is seen between FPMA and FPMB, with FPMA being more sensitive.

As described in detail in Wik et al. (2014), the NuSTAR observatory has several independent background components which vary spatially across the field of view. The background spectrum can be decomposed into five components of fixed spectral shape (excluding instrumental line strengths) with individual normalizations that are position dependent. Below 20 keV, the background is dominated by stray light from unblocked sky emission leaking through the aperture stop. This component is by nature spatially non-uniform. Below 5 keV, there is a spatially uniform component across all detectors, due to emission from unresolved X-ray sources focused by the optics (denoted fCXB for focused Cosmic X-ray Background). A further low-energy background component is related to solar photons reflecting off the back of the aperture stop, producing a time-variable signal related to whether or not the instrument is in sunlight. Above 20 keV, the background is dominated by the detector emission lines produced by interactions between the spacecraft/detectors and the radiation environment in orbit, as well as several fluorescence lines.

We used nuskybgd (Wik et al. 2014) to produce accurate background maps for each observations taking into account all the components described above. Using this code, we extracted spectra (and response matrices) in four circular regions, covering each quadrant of the field of view (with radius of 28) and avoiding the gaps between the detectors for both FPMA and FPMB for each observation. The eight extracted spectra were jointly fitted with XSPEC v.12.8.1 (Arnaud 1996) to determine the normalization of all the above components in each observation. Because the fCXB component is more than 10 times fainter than the aperture background component, we first fit a fixed normalization to this component (using the nominal value from HEAO-1 measured normalization, Boldt 1987) and then we let it vary once the other components were constrained. Given the overall small number of counts in the background spectra (∼1000 counts), the Cash statistic (Cash 1979) was used for the spectral fitting. Background maps were then produced using the fitted normalizations.

To compare the generated background maps and the background value measured in the data images, we extracted the counts from the background maps to compare with the counts extracted in the same regions from the observed data for both FPMA and FPMB. We covered each field with 64 regions of 45'' radius. Figure 2 presents the normalized distribution of relative differences between background and observed data counts extracted in each region in each field (red: FPMA; blue: FPMB). Given the absence of bright NuSTAR sources in the COSMOS field, we find that removing detected sources when computing the background maps does not significantly change the overall background distribution. From Gaussian fitting of the distributions of the difference between source and background shown in Figure 2, we find centroids at $(\mathrm{Data}-\mathrm{Bgd})/\mathrm{Bgd}$ = 0.0023 and 0.0047 and standard deviations of 0.144 and 0.147 for FPMA and FPMB, respectively, showing a remarkable agreement between the generated background maps and the data.

Zoom In Zoom Out Reset image size

Any errors in the astrometric solution for the different exposures can introduce a loss of sensitivity due to the decrease of angular resolution when summing the observations into the mosaic. We therefore tested whether significant astrometric offsets are affecting our observations. The catalogs of detected Chandra sources in COSMOS (Elvis et al. 2009, F. Civano et al. 2015 submitted) was used as the reference for computing astrometric offsets. However, given the lack of multiple bright sources in individual observations required to perform accurate alignments, we used a stacking technique applied to FPMA and FPMB separately. At first, Chandra sources in each individual NuSTAR observation were stacked (by removing neighbors and keeping brighter sources) but the resulting signal in most fields had insufficient signal-to-noise ratio (S/N<1) to provide an accurate astrometric correction. We then stacked Chandra sources in contemporaneous and contiguous NuSTAR observations, i.e., during which the telescope did not move to observe another target between one COSMOS field and the next. This procedure assumes a stable alignment of the observatory. The astrometric offsets measured for stacked sources with high S/N (>3) were in the range 1''–7'', comparable to the NuSTAR pixel size and consistent with expected uncertainties (see also Section 4.3). Therefore we decided not to perform any astrometric corrections to our data. No significant offsets were found between FPMA and FPMB.

The 121 observations were merged using the HEASoft tool XSELECT into three mosaics: FPMA, FMPB, and the summed FPMA+B. Following Alexander et al. (2013), the FPMA, FMPB, and FPMA+B event mosaics were filtered in energy using the CIAO (Fruscione et al. 2006) tool dmcopy into three bands, 3–8, 8–24, and 3–24 keV. The high energy limit of 24 keV for the analysis has been imposed by the presence of relatively strong instrumental lines at 25–35 keV, whose parametrization is uncertain which can lead to spurious high residuals during the modeling of the background. Source detection at energies above 35 keV is possible but is beyond the scope of our paper, and will be the subject of a future work. In order to achieve the deepest sensitivity, we performed detection and analysis on the merged FPMA+B mosaic, as we are confident of the alignment of the two detectors (Section 3.4). The 3–8 and 8–24 keV band mosaics are shown in Figure 3 compared to the area covered by Chandra and XMM-Newton on the same field. The NuSTAR COSMOS survey is fully covered by both lower energy X-ray telescopes.

Zoom In Zoom Out Reset image size

To generate simulated maps, mock sources were assigned fluxes drawn randomly from the number counts estimated assuming the Treister et al. (2009) model in the 3–24, 3–8, and 8–24 keV bands, using an online number count Monte Carlo calculator.²⁸ Given that we are simulating the total population and not a sub sample of sources, we can assume that the number count model used in this work is consistent with other models available in the literature as Gilli et al. (2007), Akylas et al. (2012), and Ballantyne et al. (2011).

The minimum flux for the input catalog was 5 × 10⁻¹⁵ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ in the 3–24 keV band, which is a factor of ∼10 below the expected NuSTAR COSMOS limit. Hence, background fluctuations due to unresolved sources are included in the simulations. Fluxes for the 3–8 and 8–24 keV bands were computed using a power-law model with slope Γ = 1.8, the typical value for AGNs in this energy range (Burlon et al. 2011; Alexander et al. 2013), and Galactic column density N_H = 2.6 × 10²⁰ cm⁻² (Kalberla et al. 2005). The counts-to-flux conversion factors (CF) adopted here, computed using the response matrix and ancillary file available in the adopted CALDB, are CF = 4.59, 3.22 and 6.64 × 10⁻¹¹ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{\mathrm{counts}}^{-1}$ in the 3–24, 3–8 and 8–24 keV bands, respectively.

These sources were randomly added to a background map produced as described in Section 3.3, though without the fCXB component included. The point-spread function (PSF) used to add sources to the background map is taken from the NuSTAR PSF map available in the CALDB. It changes as a function if the off-axis angle (and azimuthal angle), and is the sum of the PSFs for all the observations covering any certain position. Then, to match the total number of counts in the observed NuSTAR mosaic, additional background counts (4%–6% of the total) were added by scaling a background map with only the fCXB component for which the normalization has been averaged over all the fields. These additional counts are due to the fact that we averaged the fCXB component normalization, to sources fainter than 5 × 10⁻¹⁵ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ which are not included in the simulation, and also to different spectral shapes of the true source population. A Poisson realization of the final map (sources plus background) was made. With the above procedure, we simulated a set of 400 mosaics in three bands (3–24, 3–8, and 8–24 keV) for both FPMA and FPMB, which were then summed to create simulated FPMA+B mosaics.

In order to fully exploit the large area and depth of the NuSTAR COSMOS survey, a dedicated analysis procedure for detection and source reliability was developed and applied to the simulated data set with the aim of validating it. The same procedure, described below, was then applied to the real data (see Section 5) as well as to the ECDFS data in Mullaney et al. (2015).

Following Mullaney et al. (2015), we used SExtractor (Bertin & Arnouts 1996) to obtain a large catalog of potential sources. The source detection was performed on false probability maps generated by convolving the data mosaic (either simulated or observed) and the corresponding background map (the mosaic where the normalization of the fCXB component has been averaged over all the fields) using a circular top-hat function with two smoothing radii (10'', 20'') to detect sources with different sizes, i.e., to take into account overlapping point sources across the mosaic. To convert the convolved maps into Poisson probability maps, we used the incomplete gamma function, igamma (available in IDL), so that ${P}_{\mathrm{random}}={igamma}(\mathrm{Sci},\mathrm{Bgd})$, where Sci and Bgd are the smoothed science and background mosaics. These probability maps give the likelihood that the signal at each position in the mosaic is due to random background fluctuations. We computed the logarithm of these maps and inverted them so that significant fluctuations are positive. These maps were then input to SExtractor using a detection significance of 10^−4.5, set to avoid the loss of any real but faint detection.

We performed aperture photometry at the positions obtained by SExtractor for each detected source. Total source counts were extracted from the data mosaic and background counts were extracted from the background mosaic in the 3–24, 3–8, and 8–24 keV bands, using a circular aperture of 20'' radius. With total and background counts, we computed the maximum likelihood (DET_ML) for each source (see Cappelluti et al. 2009; Puccetti et al. 2009; LaMassa et al. 2013 for similar approaches). The DET_ML is related to the Poisson probability that a source candidate is a random fluctuation of the background (P_random): DET_ML = −ln P_random. Sources with low values of DET_ML, and correspondingly high values of P_random, are likely to be background fluctuations. Detection and photometry were both computed in the three different bands separately.

The sources detected in each probability map were merged into a single list, and duplicate sources were removed using a matching radius of 30'', i.e., if there are two sources with a separation smaller than 30'' only the one with the higher DET_ML is kept in the catalog. Given the size of the PSF, a 30'' matching radius should not produce a large number of false matches (<few %). The mean number of sources detected in each of the probability maps (in each band) and the final number of sources after cleaning the list of multiple detections of the same source are reported in Table 1.

We then used the procedure described by Mullaney et al. (2015; their Section 2.3.2) to deblend the counts of sources that have been possibly contaminated by objects at separations of 90'' or lower. Deblended source and background counts were used to compute new DET_ML values for each source.

We thus obtained for each band a catalog of detected sources which was then matched to the list of input mock sources using a positional cross-correlation, with a maximum separation of 30''. We report the mean number of detected sources matched to an input catalog source for each band in the fourth line of Table 1. The fraction of matched sources with respect to detected sources is ∼66%, 72% and 60% in the 3–24, 3–8, and 8–24 keV bands, respectively.

In Figure 4, we plot the distribution of the separation between input source positions and detected source positions for the three bands. These distributions are flux dependent: the distribution of bright sources (green lines, ≥10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$) peaks at ∼5'', while the one of fainter (<10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$) sources peaks at 8''. In the 3–24 and 3–8 keV bands, ∼55% of the matches are within 10'' and ∼90% within 20''. These numbers are slightly lower in the 8–24 keV band, where ∼45% of the matches are within 10'' and ∼85% within 20'', because of the lower number of counts in this band, however the fraction of sources within 20'' is still very high. The 30'' matching radius was chosen a posteriori to avoid losing the tail of sources at large separations.

Zoom In Zoom Out Reset image size

The threshold for source detection must be set by balancing reliability and completeness. Reliability, which is an indicator of the number of spurious sources in the sample, is computed using the ratio between the number of detected sources matched to an input source and the total number of detected sources. Completeness is instead the ratio between the number of detected real sources and the number of input simulated sources. The analysis of the simulations allows us to choose a threshold in source significance, or DET_ML, to maximize the reliability of the sources in the sample, while simultaneously maximizing completeness. A lower DET_ML threshold gives higher completeness at the cost of lower reliability.

Figure 5 shows the cumulative distribution of reliability as function of DET_ML in three bands for the FPMA+B simulations. The horizontal lines represent 97% and 99% reliability (e.g., 3% and 1% spurious sources). For what follows, we use the significance level corresponding to 99% reliability: DET_ML(99%) = 15.27, 14.99, and 16.17 for the 3–24, 3–8, and 8–24 keV bands, respectively, corresponding to Poisson probabilities of log₁₀P_random = −6.63, −6.51, −7. Figure 6 shows the completeness in the three bands as a function of X-ray flux at the DET_ML threshold corresponding to 99% reliability. Table 2 gives the flux limits corresponding to four completeness fractions in the three bands. The mean numbers of detections above the 99% reliability DET_ML threshold found in the simulations are listed in Table 1 (last row). The fraction of detected sources above this threshold is between 35%–40% in all bands.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2.  Completeness as Function of Flux (for 99% Reliability Catalog)

Completeness	F(3–24 keV)	F(3–8 keV)	F(8–24 keV)
	($\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$)	($\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$)	($\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$)
90%	1.1 × 10−13	6.1 × 10−14	1.3 × 10−13
80%	1.0 × 10−13	5.3 × 10−14	1.1 × 10−13
50%	7.6 × 10−14	3.9 × 10−14	8.6 × 10−14
20%	5.9 × 10−14	2.9 × 10−14	6.4 × 10−14

Download table as:  ASCIITypeset image

The "sky-coverage" is the integral of the survey area covered down to a given flux limit. If at the chosen detection threshold the completeness is sufficiently high (with reliability also high), the number of detected sources should correspond to the number of input sources with DET_ML higher than the threshold value. In this case, the curves in Figure 6 represent normalized sky-coverages. The survey "sky-coverage" in the three bands is plotted in Figure 7.

Zoom In Zoom Out Reset image size

In Figure 4 (dashed lines), we plot the distribution of the separation between input and detected source positions restricted to the matches above the DET_ML threshold, for the three bands. The fraction of matches within 10'' and 20'' is significantly improved when considering only sources with DET_ML above the 99% reliability threshold, increasing for all the bands to 68% and 98%, respectively. Comparing to the distribution of the whole sample (solid line) it is possible to see that the tail at large separations is made by sources at low DET_ML values. By performing a two-dimensional Gaussian fitting of the distributions, we find a consistent width in all bands, and at both bright and faint fluxes of 66, which can therefore be associated with the positional uncertainty of the detections. This value is consistent with what is expected for the source positional uncertainties according to the minimum resolvable separation Rayleigh criterion, assuming instead of the first diffraction minimum of a point source, the size corresponding to 20% of the encircled energy fraction (∼10'').

To validate the aperture photometry performed and the count deblending applied, we computed the 3–24 keV (and 3–8 and 8–24 keV) band fluxes for all detected sources. We use the counts-to-flux CF reported in Section 4.1, computed using a power law model with Γ = 1.8 and Galactic column density. We then convert the fluxes from aperture (in 20'') to total assuming a factor, derived from the NuSTAR PSF, such that F_aperture/F_total = 0.32. We find that this value is approximately constant across the field of view (with only a few percent variation) because the size of the PSF core is constant. Therefore, this factor can be applied to all the sources to convert the aperture counts computed from the mosaic, i.e., using the counts from different positions on the detector, to total counts. The difference of the aperture correction between FPMA and FPMB is of the order of 4%, and would affect the flux estimates at the level of the statistical uncertainties. Moreover, this aperture correction factor is energy independent and can be applied to convert from aperture to total fluxes in the 3 bands used here.

In Figure 8, we present the input versus output fluxes in the 3–24 keV band for all sources above the detection threshold in all 400 simulations. The agreement at bright fluxes validates the aperture correction. The fluxes are within a factor of 1.5 of the input value down to the flux corresponding to almost 90% completeness (10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$). The spread of the distribution increases toward lower fluxes, becoming a factor of 2.5 wider at the flux corresponding to 50% completeness of the survey (7 × 10⁻¹⁴ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$). The spread at low fluxes, in particular at the flux limit, is expected, and is due to Eddington bias.

Zoom In Zoom Out Reset image size

Once we tested and optimized the source detection procedure on simulated data and defined all the matching radii on the basis of the highest rate of detected and significant sources, we applied the same procedure to the NuSTAR COSMOS data. Summarizing, for each band, we first ran SExtractor using the parameters defined in Section 4.2 and we obtained a catalog of detected sources merging all the outputs. We then extracted 20'' source and background counts at the detection position and corrected these counts for contamination of neighbors. Using the deblended counts, we computed DET_ML values for all sources. We define as detected all those sources found by SExtractor and detected above threshold those for which the DET_ML is above the threshold.

There are 81, 61, and 32 sources above the 99% reliability detection thresholds in the 3–24, 3–8, and 8–24 keV bands, respectively. These numbers agree, to within a few percent, with the numbers of detections expected from the simulations as presented in the last line of Table 1. We compared the number of detections to those predicted by the X-ray background synthesis models described in Ballantyne et al. (2011), which consider several luminosity functions (Ueda et al. 2003; La Franca et al. 2005; Aird et al. 2010) as well as the most recent one by Ueda et al. (2014). We combined the luminosity functions with the spectral model of Ballantyne (2014), and folded the results with the COSMOS sensitivity curves. The models under-predict the actual number of detected sources by 20%–30% in both the 3–8 keV and 8–24 keV bands, though it is still consistent within the uncertainties. We also compared with the number counts presented in Ajello et al. (2012) in the 15–55 keV band, properly converting the flux to the 8–24 keV band (assuming a Γ = 1.8 slope). Also in this case, the predicted counts (17–34 sources) is in agreement with the observed data. A more extended analysis of the NuSTAR number counts bands is described in F. A. Harrison et al. (2015, in preparation).

We then generated a master catalog by performing a simple positional matching (30'' matching radius) between the catalog of sources detected in the 3–24 keV band and in the 3–8 keV band, and then matched the resulting catalog (including matches but also unmatched sources in both bands) to the 8–24 keV band catalog. From the master catalog, we determine the number of significant sources (with DET_ML above the threshold) in the total sample. Table 3 lists the number of different combinations of sources above the threshold in at least one band. We include the combination of sources detected and above the threshold (labeled with capital F, S, and H) and detected but below the threshold (labeled with lower case f, s, and h). By summing all the possible combinations, the total number of sources above the threshold in at least one band is 91.

The counts for each detected source in a given band were computed by performing photometry in a 20'' radius aperture in the FPMA+B mosaic as well as in the background mosaic and then correcting the counts using the deblending code as per the simulations. If the source was either detected, but below the threshold, or undetected in a given band, we computed upper limits by extracting the counts (in the same 20'' radius aperture) at the position of the detected (but below threshold) source, or at the position of the significant detection in another band if the source was undetected. The method adopted for error determination is Gehrels (1986; 1σ equivalent errors are used). For non detections, 90% confidence upper limits were computed by following standard approach (see Narsky 2000).

Vignetting-corrected count rates for each source are obtained by dividing the best-fit counts derived from aperture photometry for each band by the net exposure time, weighted by the vignetting at the position of each source. Total fluxes were obtained by converting the count rates, assuming a power-law model with Γ = 1.8 and Galactic column density, using the conversion factor in each band and applying the aperture correction factor. The energy CF, and therefore the fluxes, are sensitive to the spectral shape: a flatter spectral slope (Γ = 1) would produce a difference of ∼20%, <5% and <15% in the 3–24, 3–8, and 8–24 keV fluxes, respectively. The energy conversion factor for the 3–24 keV band depends most strongly on the spectral shape because of the wider band considered. In Figure 9, the histograms of count rates and total fluxes in each band are presented (90% confidence upper limits have been included).

Zoom In Zoom Out Reset image size

The 91 NuSTAR detected sources were matched to both the Chandra (Civano et al. 2012) and XMM-COSMOS point-source catalogs (Brusa et al. 2010) to obtain lower energy counterparts and multiwavelength information. We use the nearest neighbor matching approach, as done in the simulations (Section 4.2), with a 30'' matching radius. Given that the NuSTAR survey overlaps also with the new Chandra COSMOS Legacy survey (F. Civano et al. 2015, in preparation), we also used the Chandra catalog for this project (S. Marchesi et al. 2015, in preparation). We applied a flux cut to both catalogs at 5 × 10⁻¹⁵ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ (2–10 keV), consistent with the limit in flux used in the simulations. At this flux, the number density of sources in the 2–10 keV band is 600 deg⁻², therefore the number of Chandra or XMM-Newton sources found by chance in the searching area around each NuSTAR source is <0.13. All the matches are therefore very likely to be real associations.

The cross-correlation returned 87 matches within 30''. The distribution of separations between the Chandra or XMM-Newton and the NuSTAR positions is presented in Figure 10. The fraction of matches within 10'' and 20'' are 56% and 97%, respectively, both consistent with the simulations (see Figure 4, dashed lines). Among the 87 sources with a Chandra and/or XMM-Newton counterpart, 14 have multiple matches within 30''. The distribution of separations considering the secondary counterpart (defining as secondary the source at larger separation) is shown in Figure 10 as a dashed line. In this case, the fraction of matches within 10'' drops to 52%, and to 90% for matches within 20''. Two of the NuSTAR sources with a multiple match show a significant iron Kα line in the NuSTAR spectrum at the energy as expected from the redshift of the primary Chandra and/or XMM-Newton counterpart, therefore these sources (NuSTAR J100142+0203.8 and J100259+0220.6, source ID 181 and ID 330 in the catalog) can be securely associated with their lower energy counterpart. Of the 12 remaining NuSTAR sources with multiple lower X-ray energy counterparts, only two (NuSTAR J095845+0149.0 and J095935+0241.3, source ID 134 and ID 257) have a primary and secondary counterpart with separations and fluxes in the 2–10 keV Chandra/XMM-Newton band within 3% one of the other; therefore both primary and secondary could be considered a possible association. For ten of the 14 sources with two possible counterparts, the separation between the NuSTAR position and the secondary candidate is 30% larger than the separation between the NuSTAR position and the primary. The flux of the secondary is also significantly fainter (50% or more) than the primary, making the primary association stronger. Hereafter, we consider the primary match to be the Chandra and/or XMM-Newton counterpart. We flag those NuSTAR sources with secondary counterparts, providing a supplementary catalog of matches.

Zoom In Zoom Out Reset image size

Of the 87 matched sources, 79 are associated with a Chandra and XMM-COSMOS source (41 to a C-COSMOS source and 38 to a Chandra COSMOS Legacy source), 7 with a Chandra only source (of which 4 are from the new Chandra COSMOS Legacy survey), and only 1 is matched with an XMM-COSMOS source outside the Chandra COSMOS Legacy survey area (see also Figure 3). Although the Chandra COSMOS Legacy survey overlaps with the XMM-COSMOS area, we consider here fluxes from the already published XMM-COSMOS catalog for the 38 sources detected in both. All 87 sources are detected in the 2–10 keV band in at least one of the three low energy surveys. When considering the Chandra and XMM-Newton sources in the central area of the NuSTAR survey that has uniform and deepest exposure, all those with 2–10 keV (Chandra or XMM-Newton) fluxes larger than 10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ are also detected by NuSTAR. At fluxes in the range 5 × 10⁻¹⁴−10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$, the fraction of Chandra and XMM-Newton sources detected by NuSTAR drops to 60% and becomes lower than 10% at fluxes below 2 × 10⁻¹⁴ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$. These fractions are in full agreement with the NuSTAR survey completeness presented in Table 2 and Figure 6. In Figure 11, we compare the 3–8 keV NuSTAR fluxes with the Chandra or XMM-Newton fluxes in the same band. We converted count rates from both the C-COSMOS and XMM-COSMOS catalogs in the 2–7 and 2–10 keV bands respectively into 3–8 keV fluxes, accounting for the slightly different energy range covered and the spectral model assumed here. In Figure 11, downward arrows are 90% confidence upper limits on the NuSTAR flux . Black contours represent the locus occupied by the simulated sample when comparing input and measured fluxes in the simulations in the 3–8 keV band (as shown in Figure 8 for the 3–24 keV band). The scatter around the 1:1 line is similar to that observed in the simulations, keeping in mind that flux variability could affect the real source distribution. The break observed at ∼10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ and the flattening of the distribution are both consistent with the simulation results and with Eddington bias affecting the sources at the flux limit. Similar behavior is found by Mullaney et al. (2015) in the ECDFS NuSTAR survey and also in the serendipitous source sample presented by Alexander et al. (2013).

Zoom In Zoom Out Reset image size

The redshift distribution of the 87 NuSTAR sources associated with a lower energy counterpart is presented in Figure 12. Spectroscopic redshifts and optical spectral classifications are available for 80 matched sources and photometric redshifts (and spectral energy distribution classifications) for 87 sources from either the XMM or C-COSMOS catalogs. According to the spectroscopic classifications, the sample is equally divided between broad line AGNs and narrow line AGNs. For the remaining seven sources without optical spectra, the spectral energy distribution fitting suggests six are best fitted by a narrow line AGN-like template and one is best fit by a broad line AGN-like template. Although we detect both broad line and narrow line AGNs up to z ∼ 2, the mean redshift of the broad line AGN is z ∼ 1, while the narrow line AGN and narrow line AGN-like sources peak at lower redshift (z ∼ 0.6). Narrow line AGNs and/or obscured sources are on average fainter (and less luminous) than broad line AGNs and/or unobscured ones, and therefore are generally detected at lower redshifts, explaining the different redshift distributions and implying that the volume sampled when surveying obscured sources is smaller than that sampled by unobscured sources.

Zoom In Zoom Out Reset image size

We refer to L. Zappacosta et al. (2015, in preparation) and A. Del Moro et al. (2015, in preparation) for a detailed X-ray spectral analysis of both bright and faint sources and for the absorption distribution in the sample. Here, to obtain a rough estimate of the obscuration level characterizing the COSMOS NuSTAR sources, we computed the hardness ratio, defined as HR = $\frac{H-S}{H+S}$, where H and S are the number of net counts in the 8–24 keV and 3–8 keV bands, respectively. We used the Bayesian Estimation of Hardness Ratios method (BEHR, Park et al. 2006) which is the most suitable tool to compute hardness ratios and uncertainties in the Poisson regime of low counts, whether the source is detected in both bands or not. The HR reported in the catalog is the mode value computed by BEHR. To compare the HR (in the 3–8 and 8–24 keV bands) computed above with X-ray spectral models, and to characterize the level of intrinsic obscuration, the redshift of each source needs to be taken into account. In Figure 13, the HR values are plotted for each source versus redshift. If the upper or lower value of the HR is at its maximum value (1 or −1, respectively), the HR values are considered to be lower or upper limits. The error computed with BEHR were estimated using the Gibbs sampler (a special case of the MCMC) to obtain information on the posterior distribution of the 3–8 and 8–24 keV counts and therefore on the HR (see Park et al. 2006 for more details). The errors and the upper and lower limits on the HR are derived from the MCMC draws. The limits do not necessarily correspond to a non detection in a given band, because BEHR computes HR directly using total counts and background counts, and relies on the combined statistics of both sub-bands.

Zoom In Zoom Out Reset image size

Even though spectral complexity is likely present in these sources (see Section 6 as an example), we compared the HRs with two sets of spectral models: a power law model with slope of Γ = 1.8 and column densities 10²², 10²³, 5 × 10²³, and 10²⁴ cm⁻² (dotted lines, from bottom to top) and the more complex MYTorus model (Murphy & Yaqoob 2009). For the latter, we assumed a uniform torus with opening angle with respect to the axis of the system fixed to 60 deg, corresponding to a covering factor of 0.5, and column densities of 10²⁴ cm⁻² and 3 × 10²⁴ cm⁻² (dashed lines, from bottom to top). Baloković et al. (2014) used the latter model in spectral analysis of particular heavily obscured AGNs observed with NuSTAR (NGC 1320), and found hardness ratios calculated from the best-fit models which are consistent with those considered here. As for comparison we also report in Figure 13 the hardness ratio evolution computed for the best fit spectra model of NGC 1320 (with 4 × 10²⁴ cm⁻², red solid line) from Baloković et al. (2014).

As shown in Figure 13, if we choose HR = −0.2, corresponding to the commonly adopted definition of obscured X-ray AGNs (10²² cm⁻²), to divide the sample in obscured and unobscured, we find that 50% of the NuSTAR sources are obscured. We caution that in the 3–24 keV NuSTAR passband, the HR for a modestly obscured AGN is very close to that expected for completely unobscured AGN (HR = −0.3) and so it is difficult to estimate the reliability and statistical uncertainty in this measure of the obscured fraction from NuSTAR HR alone. A more robust and accurate measure of the fraction including modestly obscured AGNs (down to 10²² cm⁻²) will require X-ray spectral analysis to measure N_H directly and comparison with lower-energy Chandra and XMM-Newton data, which will be addressed in future work (L. Zappacosta et al. 2015, in preparation). NuSTAR HR is more sensitive to higher obscuration (>10²⁴ cm⁻²), and can therefore be used to identify candidate CT sources. We refer to these sources as candidates as confirming their CT nature requires detailed analysis of their X-ray spectra to measure N_H directly or other independent means, which is again beyond the scope of this paper. We compute the fraction of candidate CT sources, defined as the number of sources with HR and redshift combination above the >10²⁴ cm⁻² line obtained using the MYTorus model (the magenta line in Figure 13), using the results from the BEHR MCMC analysis²⁹ and thus allowing for the uncertainties in the individual HR estimates. We also compute the fraction using the 10²⁴ cm⁻² line obtained with an obscured power law model (the red dashed line in Figure 13). The fraction of candidate CT AGN obtained is 13% ± 3% with the MYTorus model and 20% ± 3% with the obscured power-law model. These estimates are consistent with the fractions based simply on the best estimates (posterior mode) for HR (9% and 19%, respectively). We note that these values correspond to the observed fraction of CT candidates, combined over the entire luminosity and redshift range of our sample.

Two sources have lower limits on their HR: NuSTAR J100204+0238.5 (ID 557) is only detected above the threshold in the 8–24 keV band (and it is the solo source with only an 8–24 keV band detection in the whole sample) while source NuSTAR J100229+0249.0 (ID 249) is also detected above the threshold in the 3–24 keV band. Their HRs suggest high levels of obscuration, above 5 × 10²³ cm⁻². Source ID 249 is also part of a sample of candidate highly obscured AGNs from the XMM and C-COSMOS spectral analysis (Lanzuisi et al. 2015). A total of six sources detected by NuSTAR, highlighted in green in Figure 13 (NuSTAR ID 107, 249, 299, 129, 181, 216), are identified as candidate obscured AGN by the same Lanzuisi et al. (2015) spectral analysis. For all of them, the NuSTAR HR suggests columns densities exceeding 5 × 10²³ cm⁻², considering also the 1σ error bars.

Two sources (ID 330 and ID 557, labelled in the figure) have HR > 0.5 strongly suggesting obscuration exceeding 10²⁴ cm⁻² with both an obscured power law model and the MYTorus model. Both of these are candidate CT AGN. A more detailed analysis of ID 330 is presented in Section 6.

NuSTAR fluxes in the 3–24 keV band, where we have the highest number of detections, have been converted into 10–40 keV rest frame luminosities, assuming a power law model with Γ = 1.8 and a standard k-correction of (1 + z)^(Γ−2) to take into account the different bandpasses. The luminosities here are not corrected for absorption, although, the 10–40 keV band is not sensitive to obscuration up to columns of < few × 10²⁴ cm⁻². In Figure 14, the X-ray luminosity for the NuSTAR COSMOS sources is plotted versus redshift. Upper limits for eight sources not detected in the 3–24 keV band have been included as downward arrows. The COSMOS sample is compared here with the Swift-BAT 70 month all-sky survey sample (Baumgartner et al. 2013). The flux limit of the ECDFS survey (Mullaney et al. 2015) is also presented as a long-dashed line. The COSMOS NuSTAR survey reaches luminosities two orders of magnitude fainter than the Swift-BAT sample and extends to significantly higher redshift. Broad line and unobscured sources (blue squares) have on average higher luminosities, while the faint end of the luminosity distribution is dominated by narrow line or obscured AGNs. Given the large area covered by the COSMOS survey, we are able to sample also rare sources at very low redshift and faint luminosity, such as source ID 330, a spiral galaxy at z = 0.044 (see Section 6).

Zoom In Zoom Out Reset image size

To further study the nature of the NuSTAR sources, we compare the X-ray fluxes to the optical magnitudes of their counterparts. The X/O ratio (Maccacaro et al. 1988) is defined as $X/O=\mathrm{log}({f}_{{\rm{X}}}/{f}_{\mathrm{opt}})=\mathrm{log}({f}_{{\rm{X}}})+C+{m}_{\mathrm{opt}}/2.5$, where f_X is the X-ray flux in a given energy range, m_opt is the magnitude at the chosen optical band, and C is a constant which depends on the specific filter used in the optical observations. Usually, the r- or i-band flux is used (e.g., Brandt & Hasinger 2005). Originally, a soft X-ray flux was used for this relation, and the majority of luminous spectroscopically identified AGNs in the Einstein and ASCA X-ray satellite surveys were characterized by X/O = 0 ± 1 (Stocke et al. 1991; Schmidt et al. 1998; Akiyama et al. 2000; Lehmann et al. 2001). Harder X-ray surveys, performed with Chandra and XMM-Newton, found that a large number of X-ray detected sources have high (>10) X/O values, and later studies showed that high X/O is associated with large obscuration (Alexander et al. 2001; Hornschemeier et al. 2001; Fiore et al. 2003; Perola et al. 2004; Civano et al. 2005; Eckart et al. 2006; Brusa et al. 2007; Cocchia et al. 2007; Laird et al. 2009; Xue et al. 2011). High X/O sources are extreme in that their optical magnitude is faint due to a combination of high redshift and/or obscuration. At low X/O, the optical emission is dominated by the host galaxy. Given the correlation of X/O with redshift, sources with low X/O are also typically at low redshift. In Chandra and XMM-Newton surveys, low X/O sources have been dubbed optically dull or X-ray Bright Optically Normal Galaxies (XBONGs). Several studies have shown that these XBONGs could harbor highly obscured AGNs (Comastri et al. 2002; Civano et al. 2007), but so far no clear case has been reported.

Figure 15 presents the i-band magnitude plotted against the 3–8 keV (left) and 8–24 keV (right) fluxes for all NuSTAR detected sources. For both bands, the X/O = ±1 locus (yellow area) has been defined using C = 5.91, computed taking into account the width of the i-band filters in COSMOS (Subaru, CFHT, or for bright sources SDSS). The locus takes into account the spectral slope used to compute the X-ray fluxes (Γ = 1.8). The long dashed lines represent the region including 90% of the AGN population as derived in the 2–10 keV band in the C-COSMOS survey (Civano et al. 2012). The four sources not matched to a Chandra or XMM-Newton counterpart are plotted as upper limits with i < 27 (see Section 5.4). Even though this is the first time such a locus is presented above 10 keV, the agreement between detections and the AGN locus is remarkable. About 10% of the sources in both the 3–8 and 8–24 keV bands lie outside the locus, consistent with what was found in the Chandra 2–10 keV band (see, e.g., Civano et al. 2012). Flux upper limits are consistent with the locus moving to the left of the plot, and could increase the number of sources at high X/O. It is interesting that the two sources with HR > 0.5 (starred symbols in Figure 15) are located at extremes of the diagram, one with high X/O (ID 557) and z > 1, the other at very low redshift and very low X/O. Both are candidate highly obscured and perhaps CT AGNs (see Section 6).

Zoom In Zoom Out Reset image size

As mentioned in the Section 1, lower energy X-ray missions like Chandra and XMM-Newton are not sensitive to very obscured AGN at redshift z < 1.5. Therefore, any 8–24 keV detected NuSTAR source not associated with a low energy counterpart could harbor a CT AGN. Four sources (NuSTAR J100047+0139.2, J095820+0149.3, J095831+0150.1, and J095930+0250.1; ID 111, 135, 141, 245 in our catalog) are not matched with a Chandra or XMM-Newton point source within a 30'' radius. Given the 99% reliability cut applied to the sample, we expect ∼2 spurious sources in the sample out of 91 detected sources (summing the number of spurious sources expected in each band). The two spurious detections could be found among the sources with a lower energy counterpart or those without a counterpart. The probability of spurious associations is determined by the number density of sources above the 2–10 keV flux limit used for the matching, which is 600 deg⁻². Using this density, for each source we calculate the probability of finding a Chandra or XMM-Newton source by chance at a distance less than the observed separations between NuSTAR sources and their primary Chandra or XMM-Newton counterparts (see Figure 10). Over the full population, we expect a total of between 1 and 2 of the associations to be due to random chance.

Source NuSTAR J100047+0139.2 (ID 111) is coincident with an extended Chandra and XMM-Newton source centered at R.A. = 10:00:45.55 and decl. = +01:39:26.1 with redshift z = 0.220 and X-ray flux F_{0.5–2 keV} = 1.6 × 10⁻¹³ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ (Finoguenov et al. 2007; George et al. 2011; Kettula et al. 2013). Its strong detection in the 3–8 keV band (and no detection in the 8–24 keV band) with a flux of 4.5 × 10⁻¹⁴ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$ (in agreement with the flux in the same band measured in the XMM-Newton spectrum within the error bars) suggests that the NuSTAR emission also arises from a soft extended source, likely a galaxy cluster. Spectral analysis of this source, including both Chandra and XMM-Newton data, will be presented in D. zr. Wik et al. (2015, in preparation).

NuSTAR J095820+0149.3 (ID 135) is significantly detected in the 3–8 keV band and just below the DET_ML threshold in the 3–24 keV band. Source NuSTAR J095831+0150.1 (ID 141) is significantly detected in the 3–24 keV band and just below the threshold in the other two bands. The last source, NuSTAR J095930+0250.1 (ID 245), is significantly detected in the 3–24 keV band and just below the threshold in the 8–24 keV band.

Figure 16 presents the Hubble ACS F814W and Spitzer IRAC (3.6 μm) 30'' × 30'' cutouts around the NuSTAR positions for the three unmatched sources, excluding the extended one. Given their relatively bright X-ray fluxes (see red arrows in Figure 15), their optical infrared counterparts could be any of the objects labelled in Figure 16. By matching the NuSTAR position of sources ID 135, 141, 245 with the COSMOS photometric catalog (Ilbert et al. 2009), we find that within 30'' there are about 200 detected sources, of which about ∼30 have an i-band magnitude that is in the same range i = 16–25 as the 87 sources matched to a Chandra or XMM-Newton counterpart (see Figure 15).

Zoom In Zoom Out Reset image size

There are two extended optical sources, neither with a spectroscopic redshift, and a bright star, a few arcseconds away from the NuSTAR position of source ID 141. A source with a spectroscopic redshift of z = 0.836, identified as a galaxy is detected 15'' from source ID 141 position. For source ID 245, two bright spectroscopically identified sources are detected at separations of 11'' and 15'', respectively, one at z = 1.277 and the other at z = 0.358. Both of these are classified as galaxies. The one at z = 1.277 is also identified as an AGN using the infrared selection criterion of Donley et al. (2012). The catalog of IR selected AGNs of Donley et al. (2012) includes ∼1500 sources, detected with S/N > 3 in all four Spitzer-IRAC channels over 2 deg². Chandra stacking analysis of the individually non-detected AGN candidates leads to a hard X-ray signal indicative of heavily obscured to mildly CT obscuration. The number of expected IR selected AGNs in a circle of 11'' radius is 0.02. Therefore the proximity of the IR AGN to source ID 245 together with its redshift make it the likely counterpart of the NuSTAR source. The measured 12 μm luminosity of the IR AGN (∼10⁴⁵ erg s⁻¹), and the 2–10 keV luminosity as derived from the NuSTAR flux (∼6 × 10⁴⁴ erg s⁻¹) follows the relation between the observed infrared and intrinsic X-ray luminosities of AGNs (e.g., Fiore et al. 2009; Gandhi et al. 2009). The Chandra upper limit on the 2–10 keV luminosity is instead significantly fainter than the intrinsic value estimated from NuSTAR (<2 × 10⁴³ erg s⁻¹), implying very high obscuration.

For source ID 135, only one spectroscopically identified galaxy at z = 0.127 is detected ∼9'' away the NuSTAR position therefore a possible counterpart given the separation distribution as shown in Figure 10, and its optical magnitude of i_AB = 22.5, which places it within the AGN locus in Figure 15, top panel.

As a general conclusion, assuming these three unmatched sources are all narrow line AGNs with a mean redshift of z = 0.6, their X-ray luminosities in the 10–40 keV band would be 10⁴⁴ $\mathrm{erg}\;{{\rm{s}}}^{-1}$ or higher. Due to their relatively bright X-ray fluxes, and the significantly deeper flux limits in the 2–10 keV band of the Chandra COSMOS Legacy survey (∼3 × 10⁻¹⁵ $\mathrm{erg}\;{\mathrm{cm}}^{-2}\;{{\rm{s}}}^{-1}$), X-ray variability could explain their lack of low energy counterparts.

X-ray flux variability on timescales from hours to years by factors 10–100 is uncommon, although such events have been detected in the past (e.g., Ulrich et al. 1997; Uttley et al. 2002; McHardy 2013; Lanzuisi et al. 2014). Source ID 141 and 245 are only detected in the full NuSTAR band, so a direct comparison with the COSMOS 2–10 keV flux limit cannot be made. Source ID 135 is only detected in the 3–8 keV band with a flux a factor 10 brighter than the Chandra limit. This could be explained by variability.