unWISE: UNBLURRED COADDS OF THE WISE IMAGING

The Wide-field Infrared Survey Explorer (WISE) mission is described in detail by Wright et al. (2010); I only briefly summarize it here. The instrument measures four infrared bands centered on 3.4 μm (W1), 4.6 μm (W2), 12 μm (W3), and 22 μm (W4), using HgCdTe arrays for the shorter W1 and W2 bands, and Si:As for W3 and W4. Survey operations began 2010 January 14. Coverage of the full sky was completed 2010 July 17. Observations continued until 2010 August 6, when the solid hydrogen cryogen ran out and the longest band, W4, became unusable. "Three-band" observations using the remaining W1, W2, and W3 channels continued until 2010 September 29, and the "NEOWISE post-cryo" program (Mainzer et al. 2011) continued observations in the shortest W1 and W2 bands until 2011 February 1. The WISE telescope is 40 cm and provides a field of view 47'× 47' which is imaged simultaneously via dichroics on the four detectors. Exposure time was 7.7 s (W1 and W2) and 8.8 s (W3 and W4). The sky was scanned from ecliptic pole to pole, with roughly 90% overlap from scan to scan. Over 99% of the sky has 11 or more exposures in W3 and W4, and 23 or more exposures in W1 and W2. Median coverage is 33 exposures in W1 and W2, 24 in W3, and 16 in W4. Coverage increases toward the ecliptic poles, with over 3000 exposures in each band at the poles; see Figure 1.

Zoom In Zoom Out Reset image size

The pixel scale is roughly 275 per pixel, and the delivered images are 1016 by 1016 pixels. The W4 images were binned 2 by 2 on board. The four-band primary mission collected ∼1.5 million four-band image sets. The three-band phase collected ∼400,000 three-band image sets, and the NEOWISE post-cryo phase collected ∼900,000 two-band image sets.

The AllWISE Data Release¹ was announced 2013 November 13, and includes data taken during all three mission phases. The data products in the AllWISE Release include a source catalog of nearly 750 million sources, a database of photometry in the individual frames at each source position, and "Atlas Images": coadded matched-filtered images. The previous All-Sky Data Release² also included calibration individual exposures ("level 1b" images), which I use here. The algorithms are described in considerable depth in the Explanatory Supplement documents, and the Atlas Image coadds are described by Masci & Fowler (2009). I briefly review relevant parts here.

The Atlas Images in the AllWISE Release are intentionally convolved by a model of the point-spread function (PSF), producing what are sometimes called "detection maps" because they are close to optimal for detecting isolated point sources. This operation, however, blurs the images, making them less than ideal for some purposes. In this work, I produce coadds that do not include this blurring step; see Figure 2 for an example. As a result, my coadds maintain the full resolution of the original images and are appropriate for making photometric measurements. While in principle I typically advocate making photometric measurements by simultaneously forward-modeling the individual exposures, the WISE instrument is remarkably stable and has a nearly isotropic PSF, so creating a coadd results in little loss of information about the non-variable sky.

Zoom In Zoom Out Reset image size

The WISE Atlas Images (and my coadds) start with the calibrated "level 1b" individual exposures. These have been calibrated astrometrically and photometrically, flat-fielded, and had static artifacts (bad pixels, etc.) flagged. The first step is selecting the frames to be included in the coadd. The criteria used (which I also follow) are to omit W3 and W4 frames within 2000 s of an annealing event (the W3 and W4 Si:As arrays were periodically heated—"annealed"—to remove latent images); frames with a bad frame-level quality-assurance score from the WISE image-processing pipeline; and frames taken during some on-orbit experiments. For each frame, the photometric calibration is applied to convert the images to common units. Next, background level matching is performed by matching the median of a new image to that of the existing coadd, in the pixels where they overlap. Frames suspected of being affected by moon glow are subjected to tests on their pixel distributions; a rapidly varying background due to moon glow increases the pixel rms significantly above that of an image not affected by the moon. Downstream, the coadding process uses median statistics, so is quite robust to outliers; some moon-contaminated frames can be tolerated, but moon-contaminated frames are rejected if they form too large a fraction of the frames. Next, all frames to be included in the coadd are projected and interpolated onto the coadd frame in order to detect outliers. A "pixel overlap-area weighting method," or top-hat interpolation kernel, is used in this step. The list of interpolated pixel values for each coadd pixel is stored so that robust statistics (the median and approximate median-of-absolute differences) can be computed pixel-wise. A set of heuristics is used to flag outlier pixels using the computed robust statistics. An additional set of heuristics expands masked regions so that regions containing several masked pixels become completely masked. Finally, after outlier masking, the final interpolation is performed. This uses an approximation to "pixel response function" (PRF) weighting, equivalent to first smoothing the image by its PSF and then resampling that image. In practice, a sub-pixelized PRF representation and nearest-neighbor interpolation are used.

The primary difference between my coadds and those of the AllWISE Release is in the interpolation kernels we use. This belies a difference in philosophy as well. The WISE Explanatory Supplement documents describe pixels as though they are "little boxes." This picture leads to confusion and heuristics with no theoretical basis in sampling or information theory. Unfortunately, this "pixels are little boxes" view has persisted for some time in astronomical image processing, perhaps most strikingly in the form of the "drizzle" process predominantly used for Hubble Space Telescope (HST) images (Fruchter & Hook 2002).³ In contrast, I think of detector pixels as taking (noisy) delta-function samples of an underlying smooth, continuous function: the pixel-convolved PSF convolved with the scene. Once the pixels have been read out, it is not necessary or helpful to think of them as "little boxes"; they are simply samples of an underlying two-dimensional function. As I show below, the WISE pixels well-sample the PSF and thus this underlying function, so the Nyquist–Shannon sampling theorem applies; a discrete set of samples contains sufficient information to reconstruct the original signal.

Briefly, the Nyquist–Shannon sampling theorem, in one dimension, states that if a continuous function contains no frequency components above a band-limit frequency B, then a discrete set of samples spaced 1/2B apart completely characterize the function; the original function can be reconstructed exactly through sinc interpolation of the samples. A good review article is Jerri (1977). In the context of astronomical images, an image satisfies the sampling criterion—it is "well-sampled"—if the underlying function contains no spatial frequency components (Fourier modes) shorter than 2 pixels. In that case, a discrete set of image pixels (samples) are sufficient to completely characterize the underlying function. As mentioned above, the underlying function is the astronomical scene convolved by the pixel-convolved PSF; since the scene can contain (nearly) point sources, the well-sampled condition depends on the pixel size relative to the pixel-convolved PSF. Given a well-sampled image, we can resample it to a different pixel grid, without loss of information—as long as the target pixel grid also well samples the image—through sinc interpolation. The sinc function is the Fourier transform of the ideal low-pass filter, leaving all frequencies below the band limit untouched and removing all frequencies above it. The sinc interpolation kernel has infinite extent, so in practice we taper the sinc kernel to give it finite extent. In this work, I use the popular Lanczos kernel as a tapered approximation to the sinc kernel.

Figure 3 shows the loss of information resulting from using the top-hat kernel compared to the Lanczos kernel.

Zoom In Zoom Out Reset image size

Figure 4 shows the WISE PSF models from Meisner et al. (2014) and their Fourier transforms. These models are constructed by combining a large number of bright stars (to get the wings) and moderately bright stars (to get the core, which is saturated in bright stars). The PSF profiles roughly follow the envelope of an Airy profile scaled to be critically sampled, suggesting that the WISE images are approximately critically sampled. The FWHMs quoted in the WISE Explanatory Supplement are 61 and 56 for the major and minor axes in W1, 6.8 and 6.1 for W2, 7.4 and 6.1 for W3, and 12.0 and 11.7 for W4; for comparison, the rule-of-thumb of ∼2.1 pixels across the FWHM (which is derived from the width of the critically sampled Airy profile) would be 57. The Fourier transform of the PSF models shows a small amount of power at the band limit, comparable to that of a Gaussian whose FWHM equals that of the critically sampled Airy. In practice, this indicates that we can treat the WISE images as well-sampled; if there is indeed power outside the band limit, we will incur a small aliasing error when we resample the images.

Zoom In Zoom Out Reset image size

In the case of the WISE Atlas Images, in addition to resampling to the coadd pixel grid, there is also an extra convolution by the PSF. In order to avoid the issues with top-hat resampling (as shown in Figure 3), the Atlas Image pixels were chosen to have twice the resolution of the input images. In contrast, I do not perform any additional convolutions on the input images, and I produce coadds at the original nominal pixel scale. The resulting PSF widths are roughly 85 FWHM for the AllWISE Release Atlas Images in W1, W2, and W3 bands, versus roughly 6'' for my coadds; the W4 PSF is roughly twice as wide.

The overall procedure follows that of the AllWISE Release, though in general each step is simpler. I have used the high-quality data products provided by the WISE pipeline as far as possible. I use the same R.A., decl. centers and axis-aligned orientations as the Atlas Image tiles. My coadds are 2048 × 2048 pixels at the nominal native pixel scale, 275 per pixel, rather than the 4095 × 4095 images at 1375 per pixel chosen in the AllWISE Release. Since the WISE images are well sampled, there is nothing to be gained in my formulation by producing coadds at double the resolution of the inputs. For the sake of consistency, my W4 coadds are on the same 275 per pixel grid as the other bands, even though the input images have only half the resolution.

I produce a coadd in two phases. In the first, I accumulate an initial mean and variance estimate for the coadd pixels and use this to mask outlier pixels in the input frames. I mask and patch these pixels, then produce the final coadd.

I select input frames from all three mission phases that touch the R.A., decl. box of the target coadd tile. I remove any frames with a WISE pipeline-assigned frame-level quality score of zero. For W3 and W4, I remove frames within 2000 s of an anneal. For W4, I remove scans in the range 03752a–03761b inclusive (as do the WISE team), in which the W4 bias level was changed for testing purposes.

I drop frames that appear to be contaminated by scattered moonlight, following the procedure used by the WISE team. After all other cuts have been made, I examine frames that are not within the moon mask (i.e., unlikely to be strongly affected by moon glow). I use the robust estimate of the pixel standard deviation measured in the WISE pipeline (the median minus 16th percentile), computing the median and median absolute difference of this quantity in the uncontaminated frames. This roughly characterizes the "typical" background variations in the coadd tile. I reject any frame within the moon mask whose pixel standard deviation is greater than the median plus 5 median absolute differences times a Gaussian correction factor of 1.4826.

I then read each frame in turn, along with its uncertainty image and masked-pixels image. I scale each intensity image and its uncertainty image, using the WISE pipeline zero points, to have a zero point of 22.5. That is, a source with integrated flux of unity in these images (and thus the coadds) corresponds to a Vega magnitude of 22.5. I estimate an image-wide standard deviation, σ_i, based on the median of unmasked pixels in the uncertainty image. The uncertainty images (produced by the WISE pipeline) include read noise and pixel-wise Poisson noise from sources, but I want a single image-wide scalar to avoid biases; the median yields roughly the read noise plus sky background noise. I use the inverse variance as the weight w_i for image i; weight w_i = σ_i⁻².

I find that many of the W3 and W4 images have spatially varying background levels that appear to be instrumental rather than astrophysical. I apply a simple approximate spatial median filter on the input frames in order to remove some of this background variation. A full median filter is too expensive so I compute the median in subimages of side length 101 pixels, then interpolate between the subimage medians using a piecewise parabolic kernel. This will of course also remove astrophysical spatially varying backgrounds, which leads to severe artifacts in and around extended structures such as large galaxies and nearby star-forming regions. This is discussed below and illustrated in Figures 12 and 13.

In order to resample the images, I patch masked pixels, in a crude manner. For each masked pixel, I take the average of its four neighboring pixels, ignoring masked neighbors. I iterate this process until no pixels remain to be patched. The resulting images have reasonable values in all masked pixels. As a result, I do not need to handle masked pixels specially in the interpolation step; if I did not patch the images, then a single masked pixel would invalidate an image region the size of the support of the interpolation kernel (dozens of pixels). Patching the images allows me to perform the interpolation with a small approximation error. I propagate the masks (via nearest-neighbor resampling) to avoid including pixels that are dominated by patched values in the coadd.

Next, I estimate and subtract a (scalar) background level from each frame. I have taken a somewhat different approach than the WISE team, based on finding the mode of the pixel value distribution. Details are given in Section 2.1 below. For W3 and W4, which have already had a median-filtered background subtracted, this additional background subtraction has little effect.

Next, I resample the input frames onto the output coadd frame. As mentioned above, the sinc interpolation kernel can be used to perform this interpolation without loss of information. The sinc has infinite extent, so I use the third-order Lanczos kernel as an approximation (see Section 2.1 below). For each input frame, I compute a number of interpolated values: the image itself (Lanczos-interpolated), I_i, a binary image indicating which coadd pixels are touched by input pixels (nearest-neighbor interpolated), M_i, and a binary image indicating which input pixels were unmasked, or "good" (also nearest-neighbor interpolated), G_i; see Figure 5. From these, I accumulate three images:

where V_w is the weighted variance image, K_w the weighted coadd, and W the sum of (masked) weights. Recall that I_i and M_i are images, while w_i is a scalar (the image-wide inverse variance described above). The first-round coadd C₁ and the per-pixel sample standard deviation S₁, are

see Figure 6.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Using the per-pixel sample variance computed in the first-round coadd, I detect pixels that are significant outliers. In order to avoid outliers themselves biasing the estimated sample mean and variance, I have to remove each frame's contribution from the sample mean (coadd) and sample variance,

Furthermore, since the sample variance estimate given a small number of samples can be noisy, I apply a prior on the variance, $P_i^2$, to get a regularized variance estimate $\hat{S}_i^2$:

where e is a multiplicative flux error (set to 3%; see Section 3), and ν is the effective weight of the prior (set to 5); this is similar to a Gamma prior. Using these estimates of the mean and standard deviation, I compute the per-pixel chi value,

and mask any pixel i with |χ_i| > 5 as an outlier. See Figure 7. This readily removes artifacts such as cosmic rays and satellite trails, and will also remove transient or fast-moving sources. I mask the four-connected neighbors of masked pixels as well, resample these masks back to the original image frame via nearest-neighbor interpolation, and record them for future use. I discard any frame that has more than 1% of its pixels masked. I then patch the frame's resampled image and accumulate it into the second-round coadds.

Zoom In Zoom Out Reset image size

During the second round, I patch outlier pixels, and add the outliers to the pixel masks:

(where ∧ denotes binary AND) and then proceed to accumulate "unmasked" and "masked" values:

and at the end, the coadd products are

where the difference between the "masked" C_m and "unmasked" C_u is that the masked versions ignore masked pixels in the input images (including pixels flagged as outliers), while the unmasked versions use the "patched" pixel values for masked pixels rather than ignoring them altogether. The per-pixel sample standard deviations S_m and S_u are scaled so that they are approximate error estimates for the coadd intensity pixels.

Finally, I re-estimate and subtract a scalar background level from the final coadds. This might seem to be unnecessary, since I subtracted the background from the input frames going into the coadd. But the coadds have considerably higher signal-to-noise; a pixel with one sigma of source flux in the individual frames will appear to be a blank pixel affected only by the background, but in the coadd it will have several sigma of source flux, and will no longer be considered part of the "blank" sky that defines the background level. In this way, the coadd's pixel distribution preferentially shifts pixels with positive flux out of the "blank sky" peak into the "source" tail, meaning that the peak of the remaining "blank sky" pixel distribution is shifted negative.

The output products include, for each coadd tile, (1) co-added intensity images that exclude masked outliers and masked input pixels, C_m, as well as ones that replace these masked pixels by patching, C_u; (2) pixel-wise inverse-variance maps W_m and W_u; (3) per-pixel sample standard deviation maps S_m and S_u; (4) number-of-frames coverage maps N_m and N_u; (5) FITS tables summarizing the frames considered for inclusion in the coadd; and (6) mask images for the input-frame pixels flagged as outliers.

2.1. Implementation Details

An example coadd image is shown in Figure 8.

Zoom In Zoom Out Reset image size

Figure 9 shows the distribution of pixel values with respect to their error estimates. I find that the AllWISE Release coadd uncertainty values appear to underestimate the variation in pixel values, at least at full-tile scales: the distribution of pixel values in full Atlas Image tiles is considerably broader than expected for a flat background. This could be due to spatial variations in the background, but since the same variations are not seen in my coadds, I expect this could be due to variations in the AllWISE Release background estimates, rather than true variations in the background levels. The method involves matching the background of an incoming frame with the existing coadd; thus "the final background level in a co-add is dictated by the first frame that was (re)projected."⁵ An example is shown in Figure 10.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 11 examines the behavior of my per-pixel sample standard deviation measurements. This can be seen as a measurement of the uncertainty of the coadd based on the data values encountered during coadding. Interestingly, there are considerable differences between this measurement of uncertainty and the "error propagation" approach used in the AllWISE Release; this could be explained by multiplicative systematics of about 3%. I suspect the All-Sky Release Atlas Images are also affected by this systematic scatter.

Zoom In Zoom Out Reset image size

Figure 12 demonstrates the effect of applying a spatial median filter to estimate and remove a varying background from the input frames before coadding them. As expected, the median filter removes large-scale structure from the coadds, but also depresses the background levels around bright sources. A more extreme case of this is shown in Figure 13.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

As outlined above, I perform two rounds of coadding, where the first round is used to reduce the impact of artifacts. In my implementation, I perform the resampling once, keeping the resampled input frames in memory until the second round. The memory required thus increases with the number of pixels touching the coadd tile. In practice, this has not proven to be problematic: the tile containing the north ecliptic pole (2709p666) has nearly 5000 W1 exposures at its peak and is touched by approximately 10¹⁰ input-frame pixels in 20,000 exposures. Producing this coadd tile required peak memory of about 100 Gb and took about 12 hr on a single CPU core. The tile shown in Figure 8, 1384p454, is more typical: it is covered by 432 exposures in W1, incorporates a total of ∼3 × 10⁸ input pixels, and took peak memory of 2.3 Gb and 630 s on a single CPU core.

The method presented here includes a number of heuristics to handle the imperfections in real data. These heuristics include noise estimation for image weighting; masked-pixel patching; background estimation; and masking of transients and defects in the imaging.

I weight the images in these coadds based on an estimate of an image-wide (scalar) noise estimate. I did not want to use pixel-wise error estimates (which incorporate a Poisson noise term from the source as well as the sky) since these lead to biased estimates (pixels with fewer counts get a smaller error estimate and hence more weight). Using an error estimate based only on the baseline background and read-noise should result in weights that optimize the detectability of faint (background-dominated) isolated point sources. I chose a simple heuristic for estimating this scalar error estimate: I take the median of the unmasked uncertainty pixels from the WISE single-frame pipeline. Since these uncertainty pixels include a component due to Poisson noise from the signal as well as read noise, the median value is slightly biased high. I could instead attempt to estimate the mode of the uncertainty distribution in the same way I estimate the background level.

Since I resample the input frames using the Lanczos kernel, masked pixels in the input images would lead to large holes (as large as the support of the kernel) in the resampled images. Since several percent of the pixels in typical WISE images are masked, this would have resulted in unacceptably large losses, so I chose to "patch" pixels. My crude approach is simply to average neighboring valid pixels. Other approaches for patching defects include the linear predictive filter used by the Sloan Digital Sky Survey (Stoughton et al. 2002). Patching pixels can be seen as redistributing the weights of the interpolation kernel from the masked pixel to its neighbors. I mask the resampled pixel closest to the patched pixel, so it is not included in the "masked" coadds, C_m. Nearby resampled pixels also include some contribution from the patched pixel; I chose to include these since their effect should be rather small. A more principled (but much more expensive) approach to patching images would be to fit a forward model to a small image patch and take the model prediction (with the model fit to unmasked pixels) as the patched pixel value. This would take into account the PSF, nearby sources, and many nearby pixels, resulting in predictions that were consistent with the modeled astrophysical scene.

For the W1 and W2 bands, I assume flat backgrounds and estimate a single scalar background level. For W3 and W4, I found it helpful to fit a varying background, and I chose to use a spatial median filter for this purpose. If I had a better model for the backgrounds (from astrophysical and instrumental sources) and could fit or marginalize over this model, that would come closer to producing ideal coadds. A median filter should generally remove trends on spatial scales of the order of the filter size and preserve features smaller than the filter, but will introduce biases near bright sources and crowded regions; see Figure 13 for an example. It would perhaps be better to avoid including in the median estimate any pixels that appear to be "contaminated" by sources (for example, by skipping pixels that are significantly different than their neighbors), though this is just a further heuristic rather than a real solution to the problem of estimating backgrounds. I am currently exploring possible approaches to reducing the highly depressed background level around bright structures (as shown in Figure 13). One option would be to also median-smooth on a larger spatial scale and add those filtered values back in, though this would likely still result in artifacts around structures in a range of scales. Another option would be to find the median or some other statistic of the stack of median-smoothed values, adding this back to the coadd. This way, extended structures that are consistent in many individual exposures would appear in the coadd. Yet another approach would be to match the large-scale variations to those found in other data sets, as done by Meisner et al. (2014); indeed, their maps, with point sources and other small-scale structures removed, could be ideal for repairing the large-scale structures damaged by my median filtering.

The AllWISE Atlas Images were created using a "background-matching" approach for removing frame-to-frame variations in the background level. In this approach, as a frame is being accumulated into the coadd, its median value is matched to that of the existing coadd, in the region where they overlap. This approach tends to produce backgrounds that are smooth on large scales (see Figure 13), even in the presence of bright objects, since it uses no explicit estimation of the background level. However, it has the disadvantages that it is order-dependent—the first few frames accumulated into the coadd have disproportionately large influence on the final background—and also depends on the extent of the coadd tile: shifting the tile boundaries results in different backgrounds.

I attempt to mask transients and defects in the individual exposures by building an initial coadd to estimate the mean and variance of the input frame pixels contributing to each coadd pixel. This approach of storing only summary statistics is in contrast to the approach used in the AllWISE Release of keeping in memory the full list of input frame pixels contributing to each coadd pixel, and then performing operations on those pixel lists. Using robust statistical estimates, the WISE team's approach should allow nearly 50% of the input frames to be contaminated without affecting the results. This power presumably comes at significant computational cost and complexity. My simpler and cheaper approach of accumulating only the mean and variance (and removing an input frame's contribution to these estimates before evaluating that frame) should be expected to work well at flagging a single moderate outlier in each pixel, or multiple extreme outliers if the number of contributing frames is large. If there is more than a single outlier, the sample variance estimate will be inflated by the outlier, and moderate outliers will no longer deviate from the mean enough to be flagged. Despite these expected shortcomings, I have found that my simple approach seems to work well in practice. Possible improvements would be to detect pixels whose sample variance is larger than would be expected from a Poisson noise model, and do more careful outlier detection there.

In attempting to remove artifacts such as satellite trails and cosmic rays, we also remove fast-moving sources and clip sources that vary in brightness. One would expect variable sources to induce a larger sample variances than a non-varying source of the same mean magnitude, and this might be a simple way of flagging such sources. In general, though, these coadds should be used for the non-transient sky only.

One might ask why I felt it necessary to write a new code to produce these coadds when publicly available codes such as SWarp (Bertin et al. 2002) already exist. The core coaddition (image warping and interpolation) could indeed have been implemented using SWarp, but the detection and masking of transients and defects (which I do by creating first-round coadds) would have required first running SWarp to do the interpolation and first-round coadds, then detecting and masking transients in the "temporary" files written by SWarp, then running the second-round coadds on these altered files. I felt it would be simpler to reimplement the warping and interpolation functionality.

The WISE team should be commended on their timely, high-quality data releases. One shortcoming I have attempted to address in this work is that the Atlas Images are blurred by convolution by the PSF. I present a new set of coadded images that preserve the resolution of the original exposures. I use a theoretically justified resampling method with some approximations. In the process of validating my results, I found that my pixel-wise sample variances, and likely the AllWISE Release Atlas Images as well, include a ∼3% systematic uncertainty in addition to the expected statistical uncertainty. These images, and the code used to produce them, are publicly available at http://unwise.me.

It is a pleasure to thank David W. Hogg (NYU), David Schlegel (LBL), Tod Lauer (NOAO), and Aaron Meisner (CfA), and the anonymous referee for helpful comments on the manuscript; Rachel Mandelbaum (CMU) for enlightening discussions; David Schlegel and Stephen Bailey (LBL) for data-management assistance; Michael Blanton (NYU) for the spatial median-filtering code; and Roc Cutri (IPAC) and the rest of the WISE team for very helpful responses to numerous queries at the WISE Help Desk.

This publication makes use of data products from the Wide-field Infrared Survey Explorer, which is a joint project of the University of California, Los Angeles, and the Jet Propulsion Laboratory/California Institute of Technology, and NEOWISE, which is a project of the Jet Propulsion Laboratory/California Institute of Technology. WISE and NEOWISE are funded by the National Aeronautics and Space Administration.

This research used resources of the National Energy Research Scientific Computing Center (NERSC), which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.

This research has made use of NASA's Astrophysics Data System.