unWISE Coadds: The Five-year Data Set

Before proceeding to coadd the five years of W1/W2 L1b images downloaded, we begin by computing a custom photometric zero-point for each exposure so as to place all exposures on the most uniform possible relative photometric calibration. Per-exposure zero-points are needed to specify the multiplicative rescaling of each L1b image during unWISE coaddition. In brief, we measure the zero-point variation in each band in day-long intervals using repeat observations of stars near the ecliptic poles, which are observed on a frequent basis thanks to the WISE scan strategy. We find that this procedure captures subtle zero-point variations with time that are not reflected in the per-exposure zero-point provided by the MAGZP keyword in each exposure's L1b header. Full details of our relative photometric calibration methodology can be found in Section 4 of Meisner et al. (2017b). Figure 1 shows the results of our photometric zero-point determination procedure as applied to the fourth year of NEOWISE-R exposures (i.e., those exposures newly incorporated into our coadds in this work).

Zoom In Zoom Out Reset image size

We ran the unWISE coaddition code without substantive modifications relative to the most recent algorithmic updates described in Meisner et al. (2017a, 2018b), but now employing five years of W1/W2 single-exposure imaging, which is a larger quantity of raw data than previously used for any full-sky unWISE coadd processing. Sections 5.2.1 and 5.2.2 provide a high-level overview of the full-sky coaddition results.

5.2.1. Five-year Full-depth unWISE Coadds

5.2.2. Five-year Time-resolved unWISE Coadds

Our most recent prior release of time-resolved unWISE coadds (Meisner et al. 2018c) incorporated four years of W1/W2 imaging, and typically offered 8 time-resolved coadd epochs per band at a given sky position spanning a 6.5 yr baseline. By adding in a fifth year of W1/W2 exposures, we generally obtain two additional time-resolved coadd epochs per band per sky location, while simultaneously extending the time baseline to 7.5 yr. In our five-year unWISE coadd data set, every coadd_id tile footprint has at least 9 time-resolved coadd epochs available per band. A median of 10 time-resolved coadd epochs per band is available over the full sky. The maximum number of time-resolved coadd epochs per band is 185, close to the ecliptic poles. Figure 2 shows the spatial distribution of the number of coadd epochs per band. This distribution of the number of coadd epochs matches expectations based on the time periods during which WISE was/was not observing, the WISE survey geometry, and our time-slicing algorithm. Over the entire sky, the total number of single-band epochal coadds in the five-year unWISE data set is 387,915.

Zoom In Zoom Out Reset image size

As discussed in Section 4, the unWISE coadds are optimized for forced photometry, in particular that performed to enable DESI's selection of luminous red galaxy and quasar targets as cosmological tracers (DESI Collaboration et al. 2016a; Dey et al. 2019). Because such forced photometry adopts source positions determined in the optical when measuring off of the full-depth unWISE coadds, it is critical that those coadds be accurately astrometrically calibrated, as any source mis-centering will lead to WISE flux underestimation. As a result, spatially coherent systematics in the full-depth unWISE coadd astrometry can imprint onto the DESI selection function, potentially biasing downstream power spectrum measurements. It is therefore very important to assess and optimize the quality of the full-depth unWISE coadd astrometry.

In previous work (Meisner et al. 2018b, 2018c), we have presented a detailed characterization of the unWISE time-resolved coadd astrometry, along with refined world coordinate system (WCS) parameters designed to optimize the time-resolved coadds for proper motion studies. However, until now we have not provided such an astrometric analysis for our full-depth unWISE coadds, which are constructed by simply propagating the L1b WCS without modification. Making use of high-fidelity astrometry from the recently released Gaia DR2 (Gaia Collaboration et al. 2018), we present an assessment of the five-year full-depth unWISE coadd astrometry and describe our corresponding creation of improved WCS solutions. To do so, we obtain centroids of sources in the five-year unWISE stacks from the full-sky unWISE Catalog (Schlafly et al. 2019), which provides completely independent astrometric measurements in each of the W1 and W2 bands. The unWISE Catalog astrometry is based on model fits, and for the purposes of this analysis is equivalent to simple flux-weighted centroiding, because the unWISE Catalog PSF models are constructed so as to always place the flux-weighted centroid at the origin.

6.1.1. Computing Per-coadd Offsets Relative to Gaia

We begin our astrometry evaluation process by computing per-coadd median positional offsets relative to Gaia for each full-depth unWISE image. Each full-depth coadd identified by its unique (coadd_id, band) pair is analyzed independently. We seek to compute median per-coadd offsets along the unWISE x and y pixel coordinate directions, and in R.A. and decl. For each unWISE coadd_id footprint, we gather the list of Gaia DR2 sources that fall within this footprint and have full astrometric solutions available. Without a full astrometric solution available, we would not be able to account for the motion of each calibrator source. The Gaia proper motion allows us to accurately propagate the Gaia coordinates to the mean epoch of each unWISE coadd before performing any unWISE versus Gaia comparisons. We use the unWISE -frames single-exposure metadata table corresponding to each coadd to compute the mean MJD of the contributing frames. We then propagate Gaia positions to the mean epoch of each full-depth unWISE coadd based on the Gaia proper motions. We ignore parallax because the full-depth unWISE coadds average together equal amounts of imaging on opposite sides of the parallactic ellipse as a result of the WISE survey strategy. The typical mean epoch of the five-year full-depth unWISE coadds is MJD ≈ 56970 (year ≈ 2014.85), whereas Gaia DR2 natively reports positions at epoch 2015.5.

For each full-depth unWISE coadd, we match the unWISE and Gaia sources using a radius of 2'', and compute offsets Δx = median(x_unwise − x_gaia), Δy = median(y_unwise − y_gaia) in terms of unWISE pixel coordinates. We also compute Δα = median((α_unwise − α_gaia) × cos(δ_gaia)) and Δδ = median(δ_unwise − δ_gaia). Figures 3, 4 shows a full-sky map of Δα (Δδ) in Galactic coordinates for W2, where the per-coadd values have been binned into N_side = 16 HEALPix maps (Górski et al. 2005) which are then smoothed with a 5° FWHM Gaussian kernel. The median number of unWISE-Gaia matches per coadd is 17,600 (13,500) in W1 (W2), varying between ∼6100 (4900) at very high Galactic latitude and a maximum of ∼196,000 (174,000) toward the inner Galaxy.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The large angular scale trends shown in Figures 3–4 mainly reflect Galactic rotation and the solar apex motion. This results from the fact that pre-hibernation L1b images were calibrated to 2MASS without a proper motion correction. As a result, our maps of Δα and Δδ look much like Figures 1⁸ and 2⁹ of the AllWISE explanatory supplement Section V.2.b.ii (Cutri 2013), albeit substantially reduced in amplitude because our full-depth coadds average together pre and post hibernation imaging, with post-reactivation L1b WCS determined in a way that does account for proper motion of the calibrator sources using UCAC4 (Zacharias et al. 2013).

To provide a better sense for the behavior of the measured per-tile unWISE-Gaia offsets, Figure 5 shows plots of their trends in W2 as a function of Galactic latitude. Δα ramps between ∼−70 mas and ∼+60 mas, while Δδ mainly shows an overall offset of ∼40 mas. Very similar trends are present in W1.

Zoom In Zoom Out Reset image size

6.1.2. Recalibrated Full-depth WCS

6.1.3. Validation of Recalibrated Full-depth WCS

To verify that improved astrometry results from our recalibrated WCS solutions, we analyze a sample of spectroscopically confirmed quasars drawn from the SDSS "DR14Q" quasar catalog (Pâris et al. 2018). These quasars serve as a useful validation sample because, for our purposes, they can safely be assumed stationary (negligible parallax and proper motion). We downselect to quasars that have astrometry available from Gaia DR2, retaining DR14Q quasars with a Gaia counterpart within 2''. ∼355,000 of ∼526,000 DR14Q quasars have such Gaia DR2 counterparts. Henceforward, for the DR14Q-Gaia sample of ∼355,000 quasars, we always adopt the Gaia DR2 positions as ground truth. We match the DR14Q-Gaia quasar sample against the unWISE Catalog using a 2'' radius, and compute per-object residuals in R.A. and decl. The R.A. residuals include a multiplication by cos(decl.) so that they are in true angular units rather than being simple R.A. coordinate differences. Note that the unWISE Catalog astrometry does not take into account the WCS recalibration described in Section 6.1.2. Figure 6 plots the unWISE Catalog W1 quasar astrometry residuals as a function of Galactic latitude. As anticipated, these trends are similar to those shown for the per-tile unWISE-Gaia offsets in Figure 5, although we should not expect exactly the same trends given that the quasars sample a footprint which covers only ∼1/3 of the sky. Figure 7 shows the corresponding W1 quasar astrometry residuals upon converting the unWISE Catalog pixel coordinate centroids to (R.A., decl.) with the recalibrated WCS solutions of Section 6.1.2. Figure 7 illustrates that systematic deviations of the quasar positions with respect to Gaia have been removed thanks to the recalibrated WCS solutions. Analogous plots of the W2 quasar astrometry residuals before and after recalibration show very similar results.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The Gaia-DR14Q quasars are not especially bright in WISE (median Vega magnitudes of W1 ≈ 16.3, W2 ≈ 15.2), and their Gaia astrometry may contain excess scatter because these objects have host galaxies which are likely often resolved by Gaia. Therefore, to gauge the bright-end scatter of the unWISE astrometry after recalibration, we assemble a sample of "low-motion" stars in Gaia DR2, which are meant to be brighter analogs of our Gaia-DR14Q quasar sample. Specifically, we select Gaia DR2 sources with full astrometric solutions, parallax < 3 mas, parallax_error < 1 mas, ∣pmra∣ < 1 mas yr⁻¹, pmra_error < 1 mas yr⁻¹, ∣pmdec∣ < 1 mas yr⁻¹, and pmdec_error < 1 mas yr⁻¹. We then cross-match these Gaia sources to the unWISE Catalog with a radius of 2''. We downselect to relatively WISE bright sources which can still be confidently assumed unsaturated, with unWISE Catalog mags between 10 and 11.6 (Vega) in the band under consideration. We further require $| {b}_{\mathrm{gal}}| \gt 20^\circ$ and unWISE Catalog fracflux > 0.999. This results in Gaia-unWISE comparison samples of ∼17,000 (18,000) sources in W1 (W2). For these samples the median offsets before recalibration are −9 mas (−8 mas) in dR.A. in W1 (W2), and 29 mas (34 mas) in dDecl. After recalibration, these median offsets are −1 mas (0 mas) in dR.A. in W1 (W2), and 0 mas (−2 mas) in dDecl. In W1, the bright-end scatter as measured by the robust rms is 58 mas (45 mas) in dR.A. before (after) recalibration and 52 mas (42 mas) in dDecl. In W2, the bright-end scatter is 52 mas (37 mas) in dR.A. before (after) recalibration and 47 mas (36 mas) in dDecl.

In summary, considering that the W1/W2 FWHM is ∼6'', we find that the bright-end scatter provided by our full-depth coadd astrometry is very good, even without employing recalibrated WCS solutions: ∼1/150 (1/110) of a FWHM with (without) recalibration. There are large angular scale systematics in the native full-depth unWISE coadd astrometry at the level of ∼150 mas peak-to-peak (∼0.025 FWHM), and these systematics can be entirely eliminated by using our recalibrated WCS solutions.

Because we expect the signature scientific application of our time-resolved unWISE coadds to be motion-based discovery of faint/cold brown dwarfs not detectable by Gaia (e.g., Kuchner et al. 2017), it is critical to understand and optimize the astrometric properties of these coadds. WISE brown dwarf motion searches are often limited by the need to scrutinize ∼50,000–1,000,000 finder charts, pinpointing the small number of real moving objects among a vastly larger candidate pool dominated by false positives (e.g., Kirkpatrick et al. 2014, 2016; Luhman 2014; Schneider et al. 2016). In such scenarios, both the candidate selection and vetting processes benefit to the extent that true astrophysical motion can be best distinguished from spurious apparent motion due to slight astrometric misalignments between WISE images at different epochs. Additionally, substellar object characterization can leverage high-fidelity WISE astrometry, for example in determining moving group membership to obtain age/mass constraints (e.g., Gagné et al. 2017).

In Meisner et al. (2018b) we presented an extensive astrometric evaluation of the time-resolved unWISE coadds and showed that a bright-end scatter of ∼50 mas per coordinate could be achieved with a coadd-level WCS recalibration following coaddition. That study was performed prior to the release of Gaia DR2, and hence we used the HSOY catalog (Altmann et al. 2017) to obtain calibrator sources with high-quality astrometry including proper motions. For our five-year set of time-resolved coadds, we employ essentially the same astrometric recalibration procedure as in Meisner et al. (2018b), but we now replace HSOY with Gaia DR2 as our calibrator catalog. Here we provide astrometric recalibrations based on Gaia DR2 for all time-resolved unWISE coadds, not just those coadds which are built from the most recently processed year of NEOWISE-R exposures. These recalibrated WCS solutions for the time-resolved coadds are derived using SCAMP (Bertin 2006).

In selecting our sample of Gaia astrometric calibrators, we require a full 5-parameter solution so that proper motion is available for propagation of coordinates to each time-resolved coadd's epoch. The number of Gaia DR2 SCAMP calibrators per time-resolved coadd is typically quite similar to the number previously available when using HSOY. The mean (median) factor by which the number of SCAMP calibrators per coadd increases for Gaia DR2 relative to HSOY is 1.13 (1.00). The larger mean number of SCAMP calibrators with Gaia DR2 arises primarily from crowded fields in the Galactic plane, where HSOY suffers from relying in part on ground-based catalogs. At high Galactic latitude, ($| {b}_{\mathrm{gal}}| \gt 30^\circ$), there are typically ∼2500 calibrator sources contributing to a given coadd's SCAMP solution.

For each of the 387,915 time-resolved coadds, we gather the set of overlapping Gaia DR2 calibrators and calculate the Gaia coordinates at the relevant epoch. We then compute a first order SCAMP solution for each coadd by comparing flux-weighted centroids measured from the time-resolved unWISE coadd versus Gaia calibrator positions. Only 966 SCAMP recalibration failures occurred, representing just ∼0.25% of time-resolved coadds. These rare failures generally correspond to time-resolved coadds that consist mostly of "empty" regions with zero integer frame coverage (see Figure 7 of Meisner et al. 2018b). We compared our new SCAMP recalibrations based on Gaia DR2 against those based on HSOY for tens of thousands of coadds where both are available, and the results are extremely similar. With recalibration based on Gaia DR2, our time-resolved coadds achieve an excellent bright scatter of ∼50 mas per coordinate, corresponding to ∼1/55 of the WISE pixel sidelength. The results of our time-resolved coadd WCS recalibrations are available as part of our five-year unWISE coadd data release (see Section 9 for details).

In this work we did not make use of the Gaia calibrator parallaxes when propagating Gaia coordinates to each time-resolved coadd's epoch, so our astrometric recalibrations do not provide absolute astrometry. Our Gaia DR2 calibrator sources typically have a parallax of ∼1–1.5 mas, so we did not foresee any clearcut means to conclusively validate that a parallax correction had in fact been successful. In the future we will endeavor to incorporate Gaia parallaxes into our unWISE astrometric recalibration analyses.

A first step toward creating the unWISE bitmask images is defining a sample of sources which are sufficiently bright to potentially cause artifacts that require flagging. We maintain separate bright source catalogs in W1 and W2. In Meisner et al. (2017a), we obtained these bright source samples via simple full-sky queries of the AllWISE catalog. Specifically, we selected all AllWISE sources with w1mpro < 9.5 (w2mpro < 8.3) in W1 (W2). We use Vega magnitudes and fluxes throughout this appendix unless otherwise noted. The equivalent thresholds in AB are 12.2 (11.64) in W1 (W2). Over the entire sky, this yields samples of 6365,819 (2221,020) sources in W1 (W2), with the density of such sources increasing dramatically near the Galactic plane. At $| {b}_{\mathrm{gal}}| \gt 18^\circ$ (high enough Galactic latitude to fall within DESI's footprint), this bright source sample has on average ∼29 (∼11) sources per deg² in W1 (W2).

Unfortunately, the AllWISE catalog lacks entries for some bright stars, and can report inaccurate (or "null") fluxes for some exceptionally bright objects that have far exceeded the WISE saturation threshold. To remedy this, the WISE team's ARTID artifact flagging module employed a custom "Bright Source List" that merged information from WISE itself, 2MASS, and even IRAS.¹⁴ Here, we take a similar approach to enhancing our AllWISE-based bright source list, but using only 2MASS K_S magnitudes in a simplistic manner.

To merge 2MASS K_S information into our AllWISE-based bright source list, we begin by selecting all 68,661 2MASS sources with K_S < 5. We then cross-match this 2MASS K_S bright sample with our AllWISE sample, using a 15'' radius. This radius is chosen to allow for matching of stars with proper motions of up to $| \mu | \approx 1\buildrel{\prime\prime}\over{.} 5\,{\mathrm{yr}}^{-1}$ given the ∼10 yr AllWISE-2MASS time baseline. This proper motion threshold seems reasonable given that only ∼275 sources with $| \mu | \gt 1\buildrel{\prime\prime}\over{.} 5\,{\mathrm{yr}}^{-1}$ are known, many of which (brown dwarfs, white dwarfs) are not sufficiently bright in W1 or W2 to require masking. For cases where an AllWISE bright source has a 2MASS K_S < 5 match within 15'', we replace our bright source catalog's w?mpro magnitude with the 2MASS K_S magnitude, if K_S < w?mpro. For 2MASS K_S < 5 sources with no AllWISE bright source match within 15'' or a match at >5'' separation, we instantiate a new object in our WISE bright source list at the 2MASS object's position and with its w?mpro value set to the 2MASS K_S magnitude. In W1 (W2), 4582 (877) bright sources have their WISE magnitudes replaced with 2MASS K_S, and an additional 2767 (2414) WISE bright source sample entries are instantiated based on 2MASS.

In the future we will investigate ways to more accurately predict W1 and W2 magnitudes of bright stars based on multi-band 2MASS JHK_S photometry, perhaps in combination with Gaia magnitudes and/or parallaxes. Also, given the availability of Gaia proper motions, improved 2MASS-AllWISE cross-matching could be enabled by propagating 2MASS positions to the AllWISE epoch. Lastly, unWISE coadds incorporate data spanning a considerable ∼8 yr time period, so that one could imagine augmenting our WISE bright source list with motion information to instantiate multiple entries for bright sources that are moving very rapidly.

Many of the unWISE mask bits are generated via a method that we refer to as "PSF model thresholding": on each tile's footprint in each band, we render a model image containing only objects from our bright source sample, then produce binary masks by flagging model image pixels brighter than a specified threshold as being "contaminated." Figure 9 illustrates an example of this methodology. The W1 and W2 PSF models themselves are shown in the left column of Figure 10.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The key inputs to this procedure are the positions and fluxes from our bright source catalog (Appendix A.1), and the PSF models used to render the bright star profiles. For this purpose, we employ the Meisner & Finkbeiner (2014) W1 and W2 PSF models. These PSF models natively exist in WISE detector coordinates, so we must perform several modifications before applying them to render models of unWISE coadds. The Meisner & Finkbeiner (2014) PSFs include terms that account for PSF variation as a function of position within the detector. Because we wish to render models of unWISE coadd images that sample many detector positions at each bright source's location, we begin by constructing a PSF in each band that averages over detector position. We then use the WCS of each coadd footprint to rotate the averaged PSF such that it has an orientation matching that of the detector on the sky. Typically there are two such rotations for each unWISE tile footprint, corresponding to northward and southward WISE scans—these two rotations of the PSF differ by 180° from one another. At very high ecliptic latitude, each bright star's location samples a substantial range of approach angles toward the ecliptic pole, meaning that the simple approximation of two discrete scan directions with PSF orientations separated by 180° no longer suffices. Currently, the unWISE mask bits constructed via PSF model thresholding do not account for this effect, which is only relevant over a very small fraction of the sky. With the PSF model averaged and rotated, for each bright source, we scale the PSF amplitude according to its magnitude listed in our bright source catalog and add it to the model image at the appropriate location. With the bright star model images in hand for each unWISE tile footprint, scan direction, and band we can proceed to apply various thresholding criteria, thereby flagging several types of artifacts. Full details of these thresholding steps are provided in subsequent sections.

One major limitation of the PSF model thresholding approach is that such flagging is limited in angular extent by the size of the PSF model. The Meisner & Finkbeiner (2014) W1/W2 PSFs extend quite far into the wings, with sidelengths of 149 (∼135 FWHM). Still, extremely bright stars can have diffraction spikes and "halos" reaching beyond the extent of our PSF models. Diffraction spikes become too long for the PSF models at parent magnitudes brighter than ∼4. At ∼0–0.5 mag, the circularly symmetric component of the PSF profile begins to extend beyond our PSF model. Therefore, we also produce geometric masks for diffraction spikes (Appendix A.14) and bright source circular halos (Appendix A.12), with radii that can extend far beyond the size of our PSF models.

Bits 0–3 are the original four mask bits from the Meisner et al. (2017a) unWISE bitmasks. They can be thought of as flagging the ∼0.6% (0.4%) of pixels most strongly affected by nearby bright source profiles at high Galactic latitude in W1 (W2). These masks were initially created to identify areas affected by difference imaging artifacts. Thus, these bits may provide sufficient bright star masking for unWISE-based difference imaging analyses, or other analyses in which the PSF modeling yields residuals comparably good to those from difference imaging.

Having rendered our bright source only model as described in Appendix A.2 for the relevant band, scan direction and unWISE tile footprint, generating mask bits 0–3 amounts to simply thresholding on the model pixel values. Based on examination of typical extragalactic sky regions, we adopt the same baseline threshold in W1 and W2: Γ₀ = 100 nanomaggies per pixel. We then scale the brightness threshold (Γ) toward larger values with increasing source density, in attempt to avoid masking an excessive fraction of area at low Galactic latitude. To implement this scaling of Γ with source density, we begin by constructing a map of the number of Gaia DR2 (Gaia Collaboration et al. 2018) sources per N_side = 32 HEALPix pixel. We interpolate over the raw Gaia sources counts in HEALPix pixels affected by small-scale source density spikes (e.g., globular clusters) and localized dust clouds. We then interpolate off of the resulting source density map at coordinates corresponding to the center of the unWISE tile footprint for which we are constructing a bitmask, and refer to the resulting value as n_src (units of Gaia sources per N_side = 32 HEALPix pixel). The threshold Γ in nanomaggies/pixel as a function of source density is then determined by:

$$\begin{eqnarray}\begin{array}{rcl}{\rm{\Gamma }} & = & {{\rm{\Gamma }}}_{0}+(\max ({n}_{\mathrm{src}},{10}^{5})-{10}^{5})\\ & & \times ({{\rm{\Gamma }}}_{\max }-{{\rm{\Gamma }}}_{0})/({n}_{\mathrm{src},\max }-{10}^{5}).\end{array}\end{eqnarray} \tag{ 1 }$$

So the threshold ramps linearly from its minimum value of Γ₀ to its peak value of Γ_max between n_src = 10⁵ and n_src = n_src,max. Γ remains at its baseline value so long as n_src < 10⁵, which is the case for ∼75% of the sky. n_src,max is 3.1 × 10⁶ Gaia sources per N_side = 32 HEALPix pixel. Γ_max is 10,000 nanomaggies/pixel (9250 nanomaggies/pixel) in W1 (W2).

Bit 0 (1) marks pixels which have values larger than Γ in the southward (northward) scan W1 model rendering, as does bit 2 (3) for W2. In constructing each of bits 0–3, we also dilate the binary map of pixels above threshold by a 3 × 3 square kernel. We refer to these bits as "core and wings" because they capture the cores of bright star profiles, extending somewhat out into the wings so as to also capture diffraction spikes and optical ghosts for sufficiently bright objects. The most noticeable scan direction dependence seen in bits 0–3 pertains to the W2 ghost location, appearing on opposite sides of the parent bright star in opposite scan directions.

Note that bits 0–3 mask regions containing the parent bright source centroids themselves. Also, bits 0–3 account properly for cases in which some portion of a parent bright source's PSF profile overlaps the tile under consideration despite the parent's centroid falling outside of the tile boundaries.

PSF modeling of bright sources will be compromised by attempting to fit saturated pixels near the profile core. To address this, unWISE mask bit 4 (5) marks "saturated" pixels in W1 (W2). Pixels flagged by bit 4 (5) are the subset of pixels flagged by bright source bits 0 OR 1 (2 OR 3) with model profile values which exceed thresholds of 85,000 (130,000) Vega nanomaggies in W1 (W2). Note that these per-band thresholds are fixed across the entire sky, with no spatial variation.

For the purposes of these unWISE mask bits, we have not intended to define the "saturation" threshold such that it truly corresponds to saturation of the WISE detectors. Instead, we are masking pixels which may be either fully saturated, potentially nonlinear, or be bright enough to experience nonlinearities to due interactions with outlier rejection during unWISE coaddition. At high Galactic latitude, the fraction of area masked by bit 4 (5) is 4.3 × 10⁻⁵ (3.0 × 10⁻⁵).

unWISE tile footprints are not mutually exclusive, and overlap by varying amounts depending on sky position. Typically, the overlap is ∼3' along each boundary. As a result, one may sometimes wish to determine the "best" unWISE tile to consult for a specific (R.A., decl.) sky location given the multiple tiles available to choose from. There are different possible ways to define "best." Here we define the best unWISE tile for a given (R.A., decl.) as the tile which maximizes the minimum distance of that (R.A., decl.) from any tile edge. For each pixel center's location within a given tile, one can compute whether this tile is the best tile for that pixel's (R.A., decl.), or instead whether that sky location would be better analyzed in some other tile. Pixels marked with the "center of pixel not primary" unWISE mask bit are sky locations that would be better analyzed in another tile. The region flagged by this bit is a border along the edges of the tile, with a typical width of ∼35 pixels.

A bright source with PSF wings that contaminate an unWISE tile's footprint despite the bright source's centroid itself falling outside of the tile boundary can be problematic for image modeling analyses. For instance, the unWISE Catalog pipeline treats each unWISE tile footprint independently, detecting and modeling only those sources with centroids inside of the tile boundaries. As a result, the catalog creation process would by default be predisposed to model flux from the wings of off-edge bright sources as sums of large numbers of fainter point sources. To avoid this outcome, unWISE mask bits 7 (W1) and 8 (W2) flag pixels affected by bright stars with centroids outside of the tile's footprint. The unWISE Catalog pipeline recognizes these bits and uses them to preemptively make deblending less aggressive in the affected regions.

Bits 7–8 are based on bits 0–3 described in Appendix A.3. For W1, if bit 0 and/or bit 1 is set and the parent bright star centroid is off the tile edge, then bit 7 is also set. For W2, if bit 2 and/or bit 3 is set and the parent bright star centroid is off the tile edge, then bit 8 is also set. Because the bit 7 (8) logic does not distinguish between bits 0–1 (2–3) in W1 (W2), information about scan direction is not retained in bits 7–8.

Certain image analyses may encounter problems in regions affected by resolved galaxies. For instance, the unWISE Catalog pipeline does not perform any galaxy model fitting, and by default deblends aggressively in an attempt to explain all flux above background as a sum of point sources. Flagging of resolved galaxies allows unWISE Catalog deblending to be made less aggressive in these regions; this avoids shredding resolved galaxies into large numbers of point sources. For our purposes, resolved galaxies are those with sizes ≳6'', which is the approximate W1/W2 PSF FWHM.

To flag resolved galaxies, we use an early version of the Legacy Survey Large Galaxy Atlas¹⁵ (LSLGA). This catalog contains ∼2.1 million galaxies with angular sizes of d₂₅ ≳ 7'', a size threshold coincidentally well-matched to the WISE FWHM. We visually inspected the largest (in terms of angular size) 300 LSLGA galaxies, manually removing a small number of SDSS (York et al. 2000) filter edge reflection artifacts and very low surface brightness dwarf galaxies, and also tweaking a few of the d₂₅ parameters by eye. Based on this slightly modified LSLGA catalog, we use unWISE mask bit 9 to flag elliptical regions about each resolved galaxy, with the ellipse size set by d₂₅ and the appropriate shape/orientation determined by the axis ratio and position angle. d₂₅ is taken to be the major axis of the elliptical mask, as this visually appeared to work well. We mask circular regions in cases where position angles and/or axis ratios are not available. We also floor the b/a axis ratio at 0.5 to avoid overly line-like masked regions. The median per-tile fraction of area flagged by the resolved galaxy bit is 0.08%.

This mask bit was inspired by (and largely copied from) the big_obj mask bit within the Schlegel et al. (1998, SFD) dust map data products. Its purpose is to flag regions affected by the LMC, SMC and M31. For the LMC and SMC, the unWISE big object mask bit is set for each pixel based on querying the SFD big_obj mask bit at that pixel center's (R.A., decl.) coordinates. For M31, we use an elliptical mask with a = 100', a/b = 2.82, and position angle of 35° east of north. This position angle was chosen to best cover M31's outskirts.

Bright source ghosts are well-captured by the Meisner & Finkbeiner (2014) PSF models (see Figure 10). As a result, we can employ our PSF model thresholding approach to specially flag regions affected by ghosts, making use of the thresholding criterion defined in Appendix A.3. When rendering the same bright source only model used to create bits 0–3, we use the ghost regions labeled in the right column of Figure 10 to keep track of which pixels fall within bright source ghosts. We then generate our ghost-specific mask bits (11–12 in W2 and 25–26 in W1) by flagging those pixels that are within ghost-affected regions and have model surface brightnesses larger than Γ_ghost = 0.15Γ, where Γ is defined in Equation (1). In other words, we apply a fainter surface brightness threshold for masking within regions affected by ghosts, and assign the results of this modified thresholding to ghost-specific mask bits.

Persistence artifacts referred to as "latents" within the WISE documentation are present in W1 and W2 imaging. When a bright source is imaged and saturates one or more WISE pixels, the sky location imaged by those same detector pixels in subsequent exposures will display a diffuse blob-like latent feature. Latents are more pronounced in W1 than W2, and therefore tend to appear as somewhat blue smudges (see Figure 11). Exceptionally bright stars can cause latents that remain noticeable for many exposures while decaying in amplitude over time. We refer to the number of exposures since imaging of the parent bright source as the "order" of a latent. Currently, our unWISE masking accounts for only first and second order latents i.e., those arising one and two exposures subsequent to imaging of the parent bright source. Because WISE scans at ∼07 per exposure, the sky position of a latent is significantly offset from its parent source along the scan direction. This makes unflagged latents particularly troublesome for rare object searches—latents are difficult to trace back to their parent bright sources by eye, and can also be misinterpreted as flux variable or moving objects because they appear at different sky locations during different WISE sky passes which have differing scan directions. Latent positions are deterministic, and could all be predicted exactly given the WCS and timestamp of every exposure in combination with a perfectly accurate bright source catalog.

Zoom In Zoom Out Reset image size

As a first step toward creating image-level unWISE latent masks, we generate a catalog of latents appearing in each unWISE tile footprint and each band. To make the computation of latent positions within each coadd efficient, we precompute per-band lookup tables containing the MJDs and full WCS parameters of all ∼25 million single-exposure W1/W2 images. For each unWISE tile footprint in each band, we use the unWISE -frames metadata table to determine the list of exposures contributing to the corresponding full-depth coadd. We then loop over these exposures, combining our WCS/MJD lookup table and bright source list (Appendix A.1) to determine the (R.A., decl.) positions of latents appearing in each contributing exposure. We aggregate these latent world coordinates on a per-coadd basis, and convert this list of (R.A., decl.) coordinates to pixel coordinates within the coadd under consideration. The result is one latent catalog per unWISE coadd_id per band. In addition to world and coadd pixel coordinates, each such catalog contains a variety of metadata that will enable us to subsequently render image-level latent bitmasks. These metadata include, for each latent, the parent bright source magnitude, the WISE scan direction, and the latent's order (1 for first latent, 2 for second latent). In practice, when computing each coadd's catalog of latents, we do not analyze every single contributing exposure. Instead, we sort the contributing exposures by MJD and consider only every sixth exposure. This is done to reduce the computational cost of generating the latent catalog. Because a WISE sky pass typically includes ≳12 exposures at each sky location, and the WISE scan direction at a given sky position varies by ≲05 over the time period of six exposures, analyzing only the latents from every sixth exposure should not cause any image-level under-flagging downstream.

Only a subset of the objects within our bright source list are sufficiently bright to create latents. Table 3 provides the magnitude thresholds that we employ for determining which bright sources are capable of causing first and second latents in each band. These thresholds are applied directly to the w?mpro values in our bright source list, and are taken to be constant over time and across the entire sky.

To determine the angular size of the region that should be masked around each single-exposure latent location, we define an effective parent magnitude (m_eff,l) that takes into account several factors: the parent bright source's magnitude, the absolute ecliptic latitude and the background level. The following terms contribute to the determination of m_eff,l:

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{\mathrm{bg},l}={f}_{\mathrm{bg},l}\cdot {\mathrm{log}}_{10}\left[\max ({intmed}/{{intmed}}_{0},1)\right],\end{eqnarray} \tag{ 2 }$$

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{\mathrm{arc},l}=2.5\cdot {\mathrm{log}}_{10}(\max (\delta \theta ,\delta {\theta }_{0})/\delta {\theta }_{0}),\delta {\theta }_{0}=1^\circ .\end{eqnarray} \tag{ 3 }$$

The effective parent magnitude used to compute the size of each single-exposure latent is then

$$\begin{eqnarray}&&{m}_{\mathrm{eff},l}=m+{{\rm{\Delta }}}_{\mathrm{bg},l}+{{\rm{\Delta }}}_{\mathrm{arc},l},\end{eqnarray} \tag{ 4 }$$

where m is the parent source's w?mpro value taken from our bright source list. intmed is the sky background level in WISE L1b DN. This sky background level is computed for each band on a per unWISE tile basis, by taking the median of the intmedian values in the unWISE -frames metadata table, restricted to frames that actually contributed to the coadd. The intmed value used in Equation (2) is that of the unWISE tile for which the bitmask image is being generated. intmed₀ is a fiducial background level appropriate at high Galactic latitude and low ecliptic latitude. In Equation (2), we adopt intmed₀ values of 25 DN (60 DN) in W1 (W2). The f_bg,l prefactor is set to 2.38 (3.78) in W1 (W2), and Δ_bg,l is capped at 4.3 mag. We use Δ_bg,l to account for the fact that an increased background level tends to reduce the region over which a latent's profile is non-negligible.

Δ_arc,l accounts for reduced latent surface brightness in the unWISE coadds due to increased arcing of latent imprints as $| \beta |$ becomes larger (see bottom row of Figure 11). Because of the continuous variation of the WISE scan direction at high ecliptic latitude, the coadd-level imprint of a bright source's multiple single-exposure latents is spread across an area which is $\sim \delta \theta /\delta {\theta }_{0}$ larger than it would have been in the ecliptic plane, with δθ given by

$$\begin{eqnarray}&&\delta \theta =0\buildrel{\circ}\over{.} 78/\cos (\beta ).\end{eqnarray} \tag{ 5 }$$

Equation (5) can be thought of as the range of ecliptic longitude spanned by a WISE exposure (078 on a side), as a function of ecliptic latitude. This approximately corresponds to the range of ecliptic pole approach angles sampled as a function of β. δθ₀ is set to 1° because this is the approximate ecliptic pole approach angle spread within a given scan parity (northward or southward) for regions near the ecliptic plane. We cap δθ at 180°, at which point latents will have been spread out into complete rings surrounding their parent stars.¹⁶ Note that β in Equation (5) refers to the location of the parent bright source.

To flag pixels affected by latents in a given unWISE coadd footprint, we begin by using nearest neighbor interpolation to mark each pixel corresponding to a single-exposure latent listed in the relevant latent catalog and having m_eff,l ≤ m_thresh, where m_thresh is the magnitude threshold for the type of latent being considered (see Table 3). We then apply a binary dilation about each of these marked pixels using a circular kernel, with a radius that depends on the difference (${m}_{\mathrm{eff},l}-{m}_{\mathrm{thresh}}$). Figure 12 shows the dilation radius in pixels as a function of (${m}_{\mathrm{eff},l}-{m}_{\mathrm{thresh}}$). We determined the shape of this function by assuming that for a parent source of magnitude m_thresh, only the centermost pixel of the profile is sufficiently saturated to yield a non-negligible latent imprint. Then, using the Meisner & Finkbeiner (2014) PSF model profile, we determined the radius of saturation that would result from making the total parent flux larger by a factor of ${10}^{-({m}_{\mathrm{eff},l}-{m}_{\mathrm{thresh}})/2.5}$.

Zoom In Zoom Out Reset image size

As listed in Table 2, there are eight unWISE mask bits for latents, one for each possible combination of WISE band (W1 or W2), scan direction (north or south), and latent order (first or second).

Many unWISE mask bits associated with bright sources flag the location of the parent bright object itself, a feature which may at times be considered undesirable. To address this situation, mask bit 21 (22) marks a 3 × 3 pixel box surrounding the centroid of each object in the W1 (W2) bright source sample described in Appendix A.1.

To offer a set of unWISE mask bits similar to the AllWISE CC flags, we include bits marking circular "halos" around bright sources. Bit 23 (24) flags circular halos around bright sources in W1 (W2). Our halo radius formula, adapted from that of AllWISE, depends on the parent source's brightness, the sky background level, and absolute ecliptic latitude. The sky background level is relevant because the halos of bright stars become negligible at smaller radii when the sky background is higher, for example in the Galactic plane. The WISE coverage increases toward the ecliptic poles, so that for fixed parent star brightness and sky background level, the halo will tend to be appreciable relative to typical pixel noise out to a larger radius at higher absolute ecliptic latitude.

To construct our sample of halo parent sources, we downselect the bright star sample of Appendix A.1 to sources with w?mpro < 8 in the relevant band.¹⁷ We then use the following formulae to determine the halo radius, in arcseconds:

$$\begin{eqnarray}&&{r}_{h}=B\cdot {10}^{a\cdot {m}_{\mathrm{eff},h}+b},\end{eqnarray} \tag{ 6 }$$

$$\begin{eqnarray}&&B={c}_{\mathrm{bg}}\cdot {\mathrm{log}}_{10}({intmed})+{d}_{\mathrm{bg}},{B}_{\min }\leqslant B\leqslant {B}_{\max },\end{eqnarray} \tag{ 7 }$$

$$\begin{eqnarray}&&{m}_{\mathrm{eff},h}=m-2.5\cdot {\mathrm{log}}_{10}\left[\sqrt{1/\max [\cos (\beta ),0.2]}\right].\end{eqnarray} \tag{ 8 }$$

The parameters of these equations are provided in Table 4. m is the w1mpro (w2mpro) value from our bright source list for W1 (W2). For each halo parent source, the intmed background level of the nearest tile center is adopted (see Appendix A.10 for a description of how per-tile intmed values are calculated). m_eff,h is an effective magnitude appropriate for use in the halo radius computation given the values of the other parameters. m_eff,h is brighter than m by an offset which tracks the signal-to-noise increase of a fixed flux source as its ecliptic latitude, and therefore WISE coverage, increases.¹⁸ The halo radius parameters were tuned based on inspection of five-year full-depth unWISE coadds—different halo radii may be preferable when deeper or shallower unWISE coadds are being analyzed. Our halo radius functional forms and parameters are substantially rooted in those used by ARTID. Figure 13 shows our halo radii as a function of parent magnitude at β = 0 and with fiducial high Galactic latitude sky background levels.

Zoom In Zoom Out Reset image size

At high Galactic latitude ($| {b}_{\mathrm{gal}}| \gt 18^\circ$), the typical fraction of area flagged by the unWISE circular halo masks is 1.2% (0.7%) in W1 (W2). We caution that circular masks are likely better implemented at the catalog level than at the image level, although image-level halo bitmasks do allow end users to conveniently dilate the halos as desired. Halos due to parent bright sources with centroids falling outside of a tile's boundaries but nevertheless partially impinging within that tile's footprint are always taken into account.

Given the spatial extent of the PSF models we employ, we can use our PSF thresholding approach to make "PSF-based" diffraction spike masks out to a radius of up to ∼105. The methodology for doing so is analogous to that described in Appendix A.9 for generating ghost-specific mask bits. We track which pixels in our bright source only models are within diffraction spikes, based on the PSF regions labeled in the center column of Figure 10. Within these diffraction spike regions, we apply a reduced brightness threshold of Γ_spike = 0.05Γ, thereby generating mask bits specific to the diffraction spikes (bit 27 for W1 and 28 for W2). Because the WISE diffraction spikes are quite nearly symmetric under 180° rotation, the northward and southward scan direction diffraction spike masks are virtually identical. We therefore do not report PSF-based diffraction spikes using separate mask bits for each scan direction. Instead, each band's PSF-based diffraction spike bit is the OR of PSF-based diffraction spike flagging computed in the two scan directions.

Diffraction spikes are captured to some extent in bits 0–3, and to a large extent in bits 27–28. However, for sufficiently bright stars, the diffraction spikes can reach beyond the boundaries of the Meisner & Finkbeiner (2014) W1/W2 PSF models, which are limited to ∼15' on a side. To handle such situations, we have implemented geometric diffraction spike mask bits for w?mpro ≲ 6 sources, extending up to ∼1° from the parent's centroid.

Our geometric diffraction spike masks consist of straight lines emanating outward from the bright source centroid at angles of 45°, 135°, 225°, and 315° from ecliptic north. In a manner similar to that used for determining circular halo size (Appendix A.12), we consider multiple factors when calculating the geometric diffraction spike length: the parent source brightness, the sky background level, and the absolute ecliptic latitude. These factors are taken into account by computing an effective magnitude (m_eff,sp) for each object in our bright source sample. The following are terms that contribute to our computation of m_eff,sp:

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{\mathrm{bg},\mathrm{sp}}=2.5\cdot {\mathrm{log}}_{10}\left[\max ({intmed},{{intmed}}_{0})/{{intmed}}_{0}\right],\end{eqnarray} \tag{ 9 }$$

$$\begin{eqnarray}&&{{\rm{\Delta }}}_{\mathrm{cov},\mathrm{sp}}=-2.5\cdot {\mathrm{log}}_{10}\left[\sqrt{1/\cos (\min (\beta ,80^\circ ))}\right].\end{eqnarray} \tag{ 10 }$$

The effective magnitude is then computed as

$$\begin{eqnarray}&&{m}_{\mathrm{eff},\mathrm{sp}}=m+{{\rm{\Delta }}}_{\mathrm{bg},\mathrm{sp}}+{{\rm{\Delta }}}_{\mathrm{cov},\mathrm{sp}}+{{\rm{\Delta }}}_{{fl},\mathrm{sp}}\end{eqnarray} \tag{ 11 }$$

Where m is the w1mpro (w2mpro) value from our bright source list in W1 (W2). We only generate geometric diffraction spike masks for bright sources with m_eff,sp < 6. Δ_bg,sp is a penalty that increases a bright source's magnitude (makes it considered to be effectively fainter) in regions of relatively high background. intmed is the sky background level in WISE L1b DN, defined the same way as in our latent effective magnitude computation. When computing geometric diffraction spike radii, we floor intmed at a fiducial extragalactic, low ecliptic latitude value intmed₀, taken to be 25 DN (60 DN) in W1 (W2). This forces Δ_bg,sp to be a strictly non-negative correction. Δ_cov,sp acts to decrease a bright source's effective magnitude (make it be considered brighter) with increasing absolute ecliptic latitude. At higher $| \beta |$, WISE provides larger frame coverage and hence reduced background noise in the unWISE coadds. The functional form of Δ_bg,sp results in a correction that tracks the decrease in background pixel noise with increasing $| \beta |$. In computing Δ_bg,sp, we cap this ecliptic latitude correction at its $| \beta | =80^\circ$ value, to avoid applying excessively large corrections very nearby the ecliptic poles.

At high ecliptic latitude, each unWISE coadd averages together frames with a significant spread in approach angles toward the ecliptic pole. Diffraction spikes therefore begin to take on a "flared" appearance at high $| \beta |$, and are ultimately washed out into nearly disk-like patterns in the immediate vicinity of the ecliptic poles. This results in decreased diffraction spike surface brightness at high ecliptic latitude, which we account for with the Δ_fl,sp term in Equation (11). The value of Δ_fl,sp is given by the right-hand side of Equation (3). In the case of geometric diffraction spikes, δθ (given by Equation (5)) is capped at 90° because the flaring of diffraction spikes will form a complete "disk" around bright stars once the spread in ecliptic pole approach angles reaches this value (due to the fact that the single-exposure WISE diffraction spikes emanate outward at azimuthal angles spaced at 90° intervals).

Analytic functions of m_eff,sp determine the geometric diffraction spike radius for each bright source. These functional forms and parameters are substantially rooted in those used by ARTID. The geometric diffraction spike radius is given by

$$\begin{eqnarray}&&{r}_{\mathrm{sp}}={B}_{L}\cdot {10}^{{a}_{L}\cdot \max ({m}_{\mathrm{eff},\mathrm{sp}},-2)+{b}_{L}}\cdot T({m}_{\mathrm{eff},\mathrm{sp}}),\end{eqnarray} \tag{ 12 }$$

$$\begin{eqnarray}&&T({m}_{\mathrm{eff},\mathrm{sp}})=1-(6-\max ({m}_{\mathrm{eff},\mathrm{sp}},-2))/16.\end{eqnarray} \tag{ 13 }$$

r_sp has units of arcseconds and T(m_eff,sp) is a tapering function that modulates the ARTID-like prescription in Equation (12). Given that Equation (13) floors m_eff,sp at −2 and is only applied for m_eff,sp < 6, we always have 0.5 ≤ T(m_eff,sp) < 1. T(m_eff,sp) ramps linearly with effective magnitude in the range −2 ≤ m_eff,sp < 6. The parameters of Equation (12) are listed in Table 5, and the geometric diffraction spike radius as a function of effective magnitude is shown in Figure 14.

Zoom In Zoom Out Reset image size

We initially mark geometric diffraction spikes as narrow, ∼1 pixel wide lines of the appropriate orientation and length. We then dilate these narrow versions of the geometric diffraction spike masks by an amount which depends on the parent bright source effective magnitude, dilating more aggressively for brighter parent sources. Specifically, we dilate with a 5 × 5 pixel square kernel for m_eff,sp > 0, a 9 × 9 pixel square kernel for −2 < m_eff,sp ≤ 0, and a 13 × 13 pixel square kernel for m_eff,sp ≤ −2.

Interior to each unWISE tile footprint, geometric diffraction spikes due to parent bright sources with centroids falling outside of the tile boundaries are taken into account. Although we reduce the geometric diffraction spike radii to account for flaring at high $| \beta |$, this flaring is not reflected in the morphology of the regions masked. We hope to implement this functionality in a future release of the unWISE bitmasks; this limitation of the present bitmasks is only relevant over a small fraction of the sky near the ecliptic poles.

For some applications, the full content of the native unWISE bitmasks (Table 2) may not be necessary. In particular, the splitting of artifacts into multiple bits based on scan direction is only relevant for time domain applications. Our unWISE bitmasks are already being employed by catalogs which do not require scan direction dependent masking information—the DESI pre-imaging surveys (Dey et al. 2019) and the unWISE Catalog (Schlafly et al. 2019). For the purposes of these catalogs, we have developed a scheme to collapse the native 31-bit unWISE mask information into a smaller number of summary bits. The definitions of our 8 summary bits in terms of the native unWISE mask bits from Table 2 are provided in Table 6.

In the preceding sections, we have suggested a number of detailed unWISE bitmask improvements that we may implement in the future. In this section we propose a few additional upgrades that could be also be performed.

• 1.  
  Homogenizing procedures for different mask bits. More work could be done to homogenize the separate procedures used to create different unWISE mask bits. For example, our PSF thresholding scales back masking in crowded regions based on a map of Gaia source density, whereas the circular halo and latent mask bits do so based on sky background level.
• 2.  
  Overflagging in the Galactic plane. Despite our efforts to limit masking in high source density regions, we still flag an excessively large fraction of area in some parts of the Galactic plane, mostly toward the Galactic center. Considering tiles that have centers within 1° of the Galactic plane, those with 120° < l_gal < 240° typically have ∼10%–15% of area masked. This fraction rises sharply toward the Galactic center, reaching ∼50% at (l_gal, b_gal) = (±45°, 0). At the Galactic center, the fraction of area masked peaks at ∼80%–90%. Users should indeed be cautious of essentially all unWISE measurements made near the Galactic center given the extreme crowding that results from the ∼6'' FWHM W1/W2 PSF. Nevertheless, further tuning of the ways in which unWISE masking is scaled based on source density and sky background level could be attempted.
• 3.  
  More finely pixelized bitmasks. The unWISE bitmasks are pixelized at the native WISE pixel scale of 275 pixel⁻¹. Each such pixel covers an on-sky area equivalent to >100 pixels of optical data from the Mosaic3 (Dey et al. 2016) and Dark Energy Camera (Flaugher et al. 2015) instruments used for DESI pre-imaging. Thus, our current pixel size may limit the utility of applying our masks to WISE forced photometry based on detections in much higher resolution optical imaging. To address this, we could generate the unWISE bitmasks using a smaller pixel size.
• 4.  
  W3 and W4. In the future we could include information about W3 and W4 in our unWISE bitmask products.
• 5.  
  Nebulosity. Nebulosity due to Galactic dust can be problematic for source detection and modeling analyses. The unWISE Catalog processing uses a neural network to classify whether or not a sky region is affected by nebulosity (Schlafly et al. 2018), and reports this information via its flags_info image product. This nebulosity flagging could potentially be added into our unWISE coadd bitmask images.