RX J0806.4–4123 has the second-longest period among the XTINSs, P = 11.37 s. Its spin-down power is erg s −1 , its inferred dipole magnetic field is B = 2.5 × 10 13 G, and its total X-ray luminosity at a distance of 250 pc is L X = 2.6 × 10 31 erg s −1 (Haberl et al. 2004 ; Haberl 2007 ; Posselt et al. 2007 ; Kaplan & van Kerkwijk 2009 ). Previously, we used the Very Large Telescope (VLT) in the near-infrared (NIR) to carry out a direct-imaging search for substellar companions around isolated XTINSs (Posselt et al. 2009 ), and we noticed a slight H-band flux enhancement at the location of RX J0806.4–4123. Motivated by an additional Herschel 160 μm detection very close to the neutron star position (Posselt et al. 2014 ), we followed up on the NIR source with the Gemini telescope. We indeed found indications of the VLT source in another NIR band. The two independent faint flux enhancements together, H VLT ≈ 23.3 ± 0.5 (1σ) and J Gemini ≈ 24.8 ± 0.5 (1σ), resulted in a combined NIR detection (H and J bands) significance of 3.1σ (Posselt & Luhman 2016 ). From the UV-optical spectral slope (measured by Kaplan et al. 2011 ), one would expect the XTINS to have a NIR magnitude on the order of 28, much fainter than the NIR detection. The NIR flux enhancements can be interpreted as a substellar companion, a fallback disk, or even a very unusual pulsar wind nebula (PWN). In this paper, we report the results of recently obtained Hubble Space Telescope (HST) NIR imaging follow-up observations that were carried out to confirm the previous weak NIR detection.

The X-ray pulsar RX J0806.4–4123 belongs to a group of seven nearby (<1 kpc), X-ray thermal isolated neutron stars (XTINSs, also dubbed the "Magnificent Seven") whose defining properties are as follows: a lack of detected non-thermal X-ray emission, radio-quietness, relatively large pulse periods (3–11 s), and inferred dipole magnetic fields on the order of ~10 13 G. The observational properties of the XTINSs (e.g., Kaplan & van Kerkwijk 2009 ; Haberl 2013 ) place this neutron star population between those of the rotation-powered pulsars and the magnetars. Since their X-ray luminosities exceed their spin-down energies by a factor of at least 10, explaining the XTINS properties requires an additional energy source other than the energy loss of a rotating dipole. Currently, it is commonly assumed that the most promising model for the heating of magnetars and XTINSs is the decay of their large magnetic fields. Viganò et al. ( 2013 ) presented simulations of a magnetothermal evolution model capable of unifying the diverse neutron star populations. In this model, magnetars and XTINSs are hot because they are heated by the decay of their large magnetic fields. In an alternative model, the observational properties of the XTINSs and magnetars are explained by formation of and influence by supernova fallback disks (e.g., Ertan et al. 2014 ). However, the fallback disk model alone cannot explain the giant outbursts of magnetars. Acceptance of this model is also hampered by the lack of fallback disk detections. So far, there is only one example of a possible (passive) fallback disk around the 3.9 kpc distant anomalous X-ray pulsar 4U 0142+61, detected with Spitzer IRAC by Wang et al. ( 2006 ). Ertan et al. ( 2007 ) reported that an active irradiated disk model can also reproduce the observed optical and infrared fluxes of 4U 0142+61.

There is no significant K s source at the position of RX J0806.4–4123, though there may be a very faint (less than 1σ) flux enhancement (see Figure 3 ). Formally, we measured the flux in an aperture with a radius of 7 pix and obtained a magnitude 22.4 ± 0.4 (1σ), resulting in an upper 3σ limit of K s < 21.2 (corresponding to <2.2 μJy) at the position of the neutron star. Such an upper limit is consistent with ~3σ detections in the surrounding area (radius r ~ 30''), where such faint sources are found to have magnitudes of 20.8 to 21.4, depending on the varying background emission.

We selected 16 stars (magnitudes K s = 11.92–14.62) from the 2MASS point-source catalog (Cutri et al. 2003 ) with quality flag AAA and no obvious nearby neighbors in the combined image to establish a photometry scale. We applied the daophot and photcal IRAF packages and used an aperture radius of 7 pix for the 2MASS stars (the FWHM is about 3 pix). We determined a zero-point of 25.50 ± 0.02. With this calibration, all but the three brightest sources (brighter than 12.35) have measured flux values within 2.2σ of their 2MASS magnitudes. The reason for the deviation of the brightest reference source is likely a slightly nonlinear behavior of the detector in some of the longer exposures due to imperfect screening. Based on the other 13 calibration stars, reliable K s magnitudes are 13.0 and fainter.

For the data reduction, we used the f2 package within the Gemini IRAF Package (version 1.12) with IRAF (version 2.16). We applied the usual corrections (dark current subtraction, bad pixel mask correction, flat-fielding). After constructing sky frames with nisky , the flat-divided sky frames were subtracted from the respective flat-divided science images. Taking into account the different exposure times, we combined the data from individual nights, as well as the data from all nights. Using Gaia data release 1 (Gaia Collaboration et al. 2016a , 2016b ) and the Graphical Astronomy and Image Analysis Tool (GAIA) of the JAC Starlink Project (Draper et al. 2016 ), we calibrated the astrometry of the reduced images to an absolute astrometric accuracy of 0 13 (1σ).

RX J0806.4–4123 was observed in the K s band with Gemini South equipped with FLAMINGOS-2 (Eikenberry et al. 2012 ). The observations were carried out from 2017 January 14 to 2017 February 6 in Fast Turnaround observing mode. The program ID was GS-2016B-FT-23. We obtained 492, 12, and 226 dithered exposures with integration times of 10 s, 12 s, and 15 s, respectively, resulting in a total exposure time of 2.3 hr. The field of view of the Hawaii-2 array has a diameter of 6 1 and a pixel scale of 0 18 per pixel at the f/16 telescope focal surface. The airmass of the individual exposure ranged from 1.02 to 1.32, the average airmass was 1.08

Our results led us to conclude that the NIR source at the location of RX J0806.4–4123 is definitely an extended one. The maximum extension of the target source is reached at a position angle of ~20° east of north. Along the line of maximum extension, the angular size of the full major axis is at least 0 8. At r = 0 4 from the center, i.e., at a distance equal to the semimajor axis, the flux has decreased to 23% of the peak flux value (for comparison, a similar flux decrease for a reference point source, e.g., s4 in Figure 6 , results in an "extension" of 0 3). In the perpendicular direction, the full minor axis was measured to be ≈0 5. If one assumes that there is indeed an additional point source, the flux of that point source represents ~20% of the total F160W flux according to our PSF subtraction. In order to estimate a firm upper limit of the contribution of a possible point source to the total flux, we slightly varied the position of the subtracted point source. Using 18 positions within r = 0 15 of the maximum NIR flux, we estimated a standard deviation of 2.5% in flux. In addition, we manually increased the subtracted point-source flux from the fit value to larger ones until the residual "hole" was inconsistent with background fluctuations of the source-free region. This criterion was met when the point source constitutes 42% of the total flux. Therefore, the estimated conservative upper limit on the contribution of a point source is 50%.

Figure 2. Comparison of images before (only sky subtraction) and after PSF (and sky) subtraction. We show the target (left) and a similarly faint source (right). The images employ the same limits, logarithmic scale, and color map. R maxmin is defined as the ratio of the maximum−minimum in the PSF-subtracted image (i.e., the residual range; right subpanels) to the maximum−minimum in the image without PSF subtraction (left subpanels) and was obtained from image statistics in the inner 0 7 × 0 7 regions. For the target, R maxmin is about a factor of 4 larger than it is for s4, underlining that there is substantial "residual flux" after PSF subtraction in case of the target. For other reference source examples in the same observation, see Figure 6 .

To test the hypothesis that the source is extended, we subtracted the point-spread function (PSF) from the emission. In a drizzled image, pixels are correlated due to the drizzling process, and different PSFs are expected for individual sources due to the different sampling and signal-to-noise properties of each object (Anderson 2014 ). Therefore, we apply the "bundle" approach by Anderson ( 2014 ), which uses the flat-fielded FLT images derived from each of the 6 exposures. Accurate (4 times supersampled) "effective" PSFs were constructed for the FLT domain and are available at the Space Telescope Science Institute's WFC3 webpages 6 (Anderson 2016 ). We chose 12 apparent isolated point sources as reference sources. Most of them are similarly faint as the target, but we also included a few moderately bright objects not afflicted by diffraction spikes. Using a pixel raster of 41 × 41, we derived bundles for the target as well as for the 12 reference sources. We created stack images of each source (i) after sky subtraction but without PSF subtraction, and (ii) after PSF and sky subtraction. Comparing these "before" and "after" images allowed us to assess the goodness of the PSF subtraction. In general, the FLT -based PSF subtraction works well. In the resulting images, residuals are visually negligible for objects as faint as the target. For brighter objects, the relative residuals (with respect to unsubtracted image) are negligible as well. We refer to the Appendix for the study of several reference sources. There are clear extended-emission residuals in the case of the target; see Figure 2, in which we also show the PSF subtraction result for a similarly faint source for comparison.

To measure the source flux we used aperture measurements with a sky annulus for the background and noise estimate. We used the daophot and photcal IRAF packages. In order to account for systematic errors and an apparent small extension, we measured values for four different aperture radii (r = 0 4, 0 5, 0 6, 0 7), six different sky annuli (inner radii from r = 0 96 to 1 6, with different widths ranging from r = 0 6 to 1 5, where the latter is limited by a nearby diffraction spike). We also checked 10 slightly different aperture center positions (maximal difference ~0 3). The largest systematic effect on the measured magnitude came from increasing the aperture radius from 0 4 to 0 6 (up to ~0.3 brighter magnitude), indicating that the source may have an extension of up to r = 1 2. The effect of changing the aperture center position reached ~0.1 in magnitude. The magnitude effect of different sky annuli was much smaller (~0.03). Considering these different measurements, we obtained the following median values and 1σ uncertainties (including systematic errors): a Vega magnitude of 23.7 ± 0.2, and an ST magnitude of 27.1 ± 0.2, corresponding to a flux density of 0.40 ± 0.06 μJy.

RX J0806.4–4123 was observed in X-rays in 2002, 2005, and 2015 with Chandra. From these observations the proper motion was found to be rather small, <50 mas yr −1 , (Haberl et al. 2004 ; Motch et al. 2009 ; B. Posselt et al. 2018, in preparation) RX J0806.4–4123 was also observed and detected with HST in the UV (ACS/SBC F140LP; STmag = 23.61 ± 0.11) and in the optical (ACS/WFC F475W, STmag = 27.92 ± 0.22) in 2009 and 2010 (Kaplan et al. 2011 ). The neutron star's position in X-rays, in particular in the optical F475W band, are consistent with an unambiguous HST detection in the F160W band; see Figure 1 .

We obtained NIR images (program GO-14745) of RX J0806.4–4123 with the HST Wide Field Camera 3 (WFC3; Kimble et al. 2008 ) on 2016 September 11 (MJD 57642.98), employing the F160W filter with a pivot wavelength of 1.537 μm. The total exposure time was 5417 s. We used the WIDE-6 POS-TARGs dither pattern (optimizes the subsampling of the pixels; HST Instrument Science Report 2016–14 5 ) for our 6 individual exposures in MULTIACCUM mode (linear SPARS timing sequence). We stacked and processed our images within PyRAF using AstroDrizzle (version 2.1.8) of the DrizzlePAC software (Fruchter et al. 2010 ; Gonzaga et al. 2012 ). Since our target source is very faint, we used the inverse variance map (IVM) weighting scheme for the final combination of the data (Gonzaga et al. 2012 ). We experimented with the drizzle parameters for a balance of spatial resolution and sampling noise. We chose a final pixel scale of 0 0642 pixel −1 and a pixel fraction of 0.7. We registered our astrometry to the Gaia data release 1 (Gaia Collaboration et al. 2016a , 2016b ); employing the common 117 reference sources in the field of view, the rms of our astrometry fit is 0 019.

The NIR source at the location of RX J0806.4–4123 is unambiguously confirmed with our HST F160W observations. The measured F160W flux, 0.40 ± 0.06 μJy, is clearly in excess of the expected neutron star flux, considering the extrapolation of its UV-optical flux (see Figure 4). In addition, the object has a size of at least 0 8 (Section 2) and an inhomogeneous brightness distribution, consisting of a brighter "core" and elongated "halo" emission (Figure 2). In order to estimate the likelihood of an unassociated background source, we measure a source density of 0.07 sources per square arcsec in the F160W image, which includes all sources with STmag ≤ 27.4 (conservative limit since the target has STmag = 27.1 ± 0.2). The separation of the flux peak of the extended emission from the neutron star position in the 2010 HST F475W image is 0 15. Recent Chandra observations indicate a proper motion of RX J0806.4–4123 of ≈20 mas yr−1 (B. Posselt et al. 2018, in preparation; the previous limit was μ < 86 mas yr−1 (2σ), Motch et al. 2009). Considering a (conservative) total position uncertainty of 0 4 (including proper motion 2σ < 50 mas yr−1, F475W position error <0 1), and N = 1 trials, we estimate a by-chance probability of 3.7%.7 Therefore, we regard an unassociated background source as unlikely.

Figure 4. HST F160W flux measurement and the 3σ upper limit in the K S band from the Gemini telescope are marked in blue; the other limits and measurements are from Posselt & Luhman (2016), Posselt et al. (2014), and Kaplan et al. (2011). The shown measurement error bars indicate 3σ uncertainties. The red dashed line shows the predicted flux of the fallback disk model discussed in Section 3.2. This is the same model that was found to fit the X-ray to optical spectrum by Ertan et al. (2017; their Figure 5). The plotted disk model was calculated for where is the inclination of the observer with respect to the disk midplane normal. The NIR flux would be the same for if about half of the flux comes from scattered light from the outer disk, as outlined in Section 3.2. The black dotted and dashed–dotted lines indicate different power-law spectral slopes. Download figure: Standard image High-resolution image Export PowerPoint slide

Even though the by-chance probability for an unassociated background source is low, a background galaxy might explain the detected NIR emission. The optical and UV data obtained with HST (Figure 2 by Kaplan et al. 2011) do not show extended emission, and the fluxes are consistent with those measured for the other six of the Magnificent Seven, i.e., it is reasonable to assume that there is no significant flux contribution from a galaxy to the F475W and L140P flux measurements of RX J0806.4−4123. The colors of galaxies indicate that any putative galaxy counterpart must be either very red (e.g., Daddi et al. 2000; Franx et al. 2003; van Dokkum et al. 2006; Guo et al. 2012), or, considering the lack of F475W flux, be at redshift z 3.5 (e.g., Giavalisco 2002). The same result follows from simulations using the mock galaxy catalog by Williams et al. (2018). These simulations also show that the probability of detecting galaxies with consistently red observed colors (either star-forming or non-star-forming) in the field is low—only ~0.2 per square arcmin are expected to have F160W fluxes within the measured flux uncertainties while also being sufficiently faint in F475W. The analysis of actual, sufficiently deep observations, results in similar numbers. Analysis of the data on the Hubble Ultra Deep Field (HUDF; Rafelski et al. 2015) yield 0.39+/−0.17 galaxies per square arcmin if red and faint sources similar to our target are considered. The uncertainty results from the large Poisson error in the very small 12.8 arcmin2 field of view.8 Furthermore, we note that high-redshift galaxies are in general more compact, with smaller sizes than low-redshift galaxies (van der Wel et al. 2014; Shibuya et al. 2015). The detected ~1'' NIR extension would indicate an unusually large high-redshift galaxy. In the mock galaxy simulations, for example, <5% of galaxies with consistent observed colors are large enough to also be consistent with the detected NIR extension. Thus, we conclude that a high-z galaxy is an unlikely counterpart of the extended NIR emission at the position of RX J0806.4−4123.

The detected extended emission could also be due to coincidence with emission of the interstellar medium (ISM). Scattered stellar radiation, hot dust, or line emission from the ISM gas are potential origins. We discuss these possibilities in more detail in Section 3.1.

If the NIR emission is associated with RX J0806.4−4123, then L F160W ~ 6 × 1027 erg s−1, assuming a distance of 250 pc. Using the X-ray flux and spin-down power of the NS, the NIR efficiency is , and the luminosity ratio . The F160W flux value, 0.40 ± 0.06 μJy, and a size of 0 5 × 0 8, translate into an average brightness temperature9 T B = 282 K. If our source were a blackbody with this temperature, then one would expect a very high IR flux, e.g., Jy which is clearly above the current flux limits by Spitzer. However, thermal emission with a temperature, T = 950 K and a filling factor10 of 7.3 × 10−11 would be consistent with the F160W and J-band flux measurements, as well as the upper limits obtained for K S and 4.5 μm. It follows from these estimates that thermal emission is, in principle, consistent with the (N)IR measurements if the emitting area is a very small fraction of the apparent area of the detected extended source. Thermal emission could come from warm gas or X-ray irradiation and subsequent cooling of solid bodies such as dust or a substellar companion. At a distance of 250 pc, an extension of ~0 8 corresponds to 3 × 1015 cm or 200 au. This rules out a substellar companion as the sole origin of the NIR excess. However, the existence of a warm substellar companion cannot be entirely ruled out, since the conservative upper limit on the contribution of a point source to the total flux is 50%. Such a flux limit translates into a substellar companion mass limit of <10 Jupiter masses (using a bolometric correction of 3.5 (Golimowski et al. 2004) according to the calculations presented by Burrows et al. (2001), assuming a distance of 250 pc and an age of 10 Myr). Warm gas and heated dust, either in the ISM or cirumstellar matter, are discussed in Sections 3.1 and 3.2.

Alternatively, the F160W measurement can also be explained by non-thermal emission. Assuming a power-law F ν ∝ να for non-thermal emission, the current spectral constraints imply either a rather steep slope α < −2.5, or an exponential cutoff between the F160W and F475W bands. Charged particles producing curvature and/or synchrotron emission would be the likely source. With respect to an association to a neutron star, the most obvious origins are magnetospheric emission (for the possible flux contribution from a point source) or a pulsar wind and the shocks it produces. We discuss such a PWN in Section 3.3.

3.1. Emission from the ISM?

In principle, an ISM density enhancement along the line of sight could scatter NIR emission from ambient light and produce the observed extended NIR emission. The so-called cloudshine is a large-scale (typically, several arcminutes) phenomenon that is caused by ambient star light around large and dense molecular clouds (Lehtinen & Mattila 1996; Foster & Goodman 2006). There are no other extended NIR features and no indication of a molecular cloud (e.g., Dame et al. 2001; Dobashi et al. 2005) in the field of RX J0806.4–4123, which one would expect if the NIR emission were due to cloudshine. Instead of the interstellar radiation field and a large, dense molecular cloud, one could postulate that the NIR emission could originate from local point-source emission that is scattered by a small (i.e., unaccounted for) ISM enhancement at a slightly larger distance than the one of RX J0806.4–4123.11 However, in this case the expected thermal NIR flux from the XTINS surface alone is too low to produce the large observed NIR excess by scattering. An additional NIR emitter would be required, as well as an unlikely "wall" in the ISM. If the NIR excess flux is due to genuine emission from the ISM, it could originate either from ISM gas or dust in the diffuse ISM or in an ISM density enhancement. There are a few emission lines of ISM gas in the H and the other NIR-bands. Strongly forbidden rotational-vibrational emission lines of H 2 are observed, for example, in reflection nebulae or shocked regions in star-forming regions and around AGB stars (e.g., Shull & Beckwith 1982). These H 2 emission lines are due to collisional excitation or UV fluorescence, and in accordance with their detection in dense environments they require high densities and significant molecular fractions (Black & van Dishoeck 1987; Hollenbach & McKee 1989). The low N H value of RX J0806.4–4123 implies not only a low density, but also a very low molecular fraction, , because inefficient self-shielding from cosmic rays leads to dissociation of H 2 (Savage et al. 1977; Burgh et al. 2007). Therefore, the NIR excess flux around the neutron star is unlikely to be produced by H 2 emission lines. Interstellar dust emission peaks in the far-infrared. It can also emit NIR photons if it is hot enough. The required high temperatures (1000–2000 K) in thermal equilibrium could only be reached very close to the XTINS (unresolved by HST) by typical (size ~0.1 μm) dust grains. A non-equilibrium process, the so-called stochastic heating of very small grains 10 Å or large polycyclic aromatic hydrogen (PAH) molecules by UV photons, was found to produce the large-scale NIR emission of reflection nebulae (e.g., Sellgren 1984; Draine & Li 2001). There are two reasons why very small grains/PAH molecules of the ISM cannot be responsible for the extended emission around RX J0806.4–4123. First, there is not a sufficient amount of very small grains in the ISM surrounding the XTINS to produce the observed NIR flux. Using the constraints on the XTINS spectrum in UV (Kaplan et al. 2011), and our flux measurement at 1.6 μm, one can estimate the required number density, n vsg , of very small grains (~10 Å), which are stochastically heated by ~10 eV UV photons (Sellgren 1984). The value, n vsg ~ 10−3 cm−3, is at least four orders of magnitude larger than what one can expect for the typical ISM with a hydrogen number density n H ~ 1 cm−3, dust-to-gas ratio of 6%, and typical grain size and composition distribution (all from, e.g., Zubko et al. 2004; we used their BARE-GR-S dust model, which matches observations of the ISM in the solar neighborhood). The lack of known "NIR halos" around other strong UV sources such as white dwarfs supports our estimate that the diffuse ISM does not provide enough small grains for substantially enlarged NIR flux close to these objects. The second reason why very small ISM grains are an unlikely origin for the observed NIR emission is the strong X-ray emission of the XTINS. Voit (1991) estimated that X-ray photons efficiently evaporate grains of sizes 10 Å. Overall, we conclude that the observed NIR emission around RX J0806.4–4123 is unlikely to come from either the gas or dust of the photon-heated ISM. We discuss the interaction of potential pulsar wind particles with the ISM in Section 3.3.

3.2. Emission from Circumstellar Material?

If the NIR emission comes from the immediate vicinity of RX J0806.4–4123, it could be due to a surrounding disk or torus structure, similar to what has been inferred for the magnetar 4U 0142+61 from the Spitzer detection by Wang et al. (2006) or for white dwarfs from NIR detections (e.g., Melis et al. 2010). The stochastic heating of very small grains, which was discussed in the previous section for the ISM, is a possibility for circumstellar matter, too. However, close to a neutron star, one expects quick destruction (evaporation by X-rays) and removal (by Poynting–Robertson drag; e.g., Figure 10 by Posselt et al. 2014) of such small grains. Hence, regular replenishment, e.g., by collisions of asteroids or large grains, would be required to keep such a "debris disk" bright in the NIR. Considering an active fallback disk model, Ertan et al. (2017, 2014) showed that such a disk is able to produce the X-ray properties and the optical excess of the XTINSs. Here, we use their model to check whether such a fallback disk can reproduce the measured NIR excess. The model employs dissipative heating, which is most relevant in the active inner part of the disk, for RX J0806.4–4123, from the inner disk radius, to 5 × 1011 cm. The model also considers X-ray irradiation, which is most relevant for the passive outer part of the disk. For details on the model and the model parameters for RX J0806.4–4123, we refer to Ertan et al. (2017, 2014). The outer radius of the passive disk is not well known and depends on the initial disk conditions and the X-ray luminosity in the early phase of evolution. However, the local irradiation temperature is very small for large radii—for RX J0806.4–4123, it decreases to <30 K for r > 1013 cm. For studies of the flux at infrared and shorter wavelengths, it is sufficient to use an outer radius cm. The extension of the NIR source at the distance of RX J0806.4–4123 is r ~ 1015 cm (or ~100 au). This raises the question of whether there is any possibility of the disk model explaining not only the observed NIR flux but also the observed source extension. Drawing on the findings in protoplanetary disk studies, favorable properties of the disk and viewing geometry can perhaps provide such an explanation. For protoplanetary disks, it was shown that in the case of a flared disk NIR emission can be scattered and produce extended emission (e.g., Mulders et al. 2013; Anthonioz et al. 2015). First, we check whether the disk model can reproduce the observed spectral energy distribution. Then we outline a potential explanation for the observed source extension. Using and cm, we obtain the model spectrum given in Figure 4. The disk model curve in Figure 4 is within the error range of the data, illustrating that the fallback disk model can reproduce the observed NIR flux. Within this disk model, the work done by the magnetic torque heats the inner rim of the disk, which then produces the optical and UV fluxes. In contrast, the NIR emission of the model spectrum is produced from the surface of the irradiated inner disk within a few times 109 cm of the XTINS. One implication is that a completely passive disk could also produce the NIR emission, but not the optical-UV emission. One possible way to explain the observed source extension is based on the assumption of disk flaring and scattering off the outer passive disk surface. The thickness of a disk can be described with the pressure scale height, , where is the sound speed and is the local Keplerian angular velocity of the disk. In flared disks the flaring is described by the aspect ratio, h/r. The outer radii of protoplanetary disks around T Tauri stars reach several hundred au and the h/r values can approach unity (see, e.g., Chiang & Goldreich 1997; Armitage 2010; Dullemond & Monnier 2010 for reviews). Chiang & Goldreich (1997) (CG97) calculated self-consistently the aspect ratio and the temperature profiles of these disks. We follow the result of CG97 for our estimates, but adapt the model parameters (luminosity, mass) to the properties of RX J0806.4–4123. CG97 found for an optically thick inner disk region. Outside this region there is a radially isothermal optically thin region along which the aspect ratio increases even faster with increasing r as . The disk photosphere, defined as the disk surface layer corresponding to the optical depth of unity to the irradiating flux, is estimated by CG97 to be H 4h. In the luminosity regime of RX J0806.4–4123, erg s−1; using (CG97) with M = 1.4, M ⊙ , we estimate the aspect ratio of the disk for RX J0806.4–4123 as , where is the radial distance in au. For , we find H/r 0.1. This could be taken as a lower limit, since the outermost part of the passive disk could be optically thin to its own emission, along which H increases faster with increasing r, as r1/2. In CG97, the radially isothermal disk region has a temperature ~20 K. For RXJ0806.4–4123, such a temperature is reached at r ~ 8 × 1012 cm. Thus, a large fraction of the outer passive disk of RXJ0806.4–4123 could be radially isothermal with H/r ∝ r1/2. In this case, the result is H/r ≈ 0.3 at 100 au. Hence, depending on the extension of the radially isothermal region, H/r could be between 0.1 and 0.3 at the outermost disk. For our estimates below we take H/r = 0.2 at . The NIR luminosity illuminating and being scattered off the outer disk surface, , is a small fraction of the total NIR luminosity emitted from the inner disk, . As seen from the inner disk, the outer passive disk subtends a solid angle from the midplane of the disk (θ = 90°) to the surface of the outer disk ( ). For H/r 0.2 at r 100 au, this angle is . Assuming an albedo of β 0.9, we estimate that up to 4 % of the total NIR emission from the inner disk could illuminate and be scattered by the outer disk surface, . The observed ratio of the NIR flux scattered off the outer disk to the NIR flux from the inner disk (r ~ 109 cm) depends on the details of the viewing and the scattering geometry, the actual h/r profile, temperature, and the extension of the currently passive disk. For instance, the scattering geometry matters because the scattered emission is unlikely to be purely isotropic. The angular distribution of photons scattered by a dust grain depends on the grain composition and grain size (Draine 2003; Mulders et al. 2013). Large grains in particular show very anisotropic scattering, where one scattering produces a conical beam with a small solid angle Ω sc . The NIR light emitted from the inner disk is estimated to be scattered from a thin layer of the disk surface. In order to estimate the contribution of the scattered component to the observed NIR flux, we use a simplified picture and represent the total scattered emission by two narrow fan beams opening from the two surfaces of the disk. In the case of single scatterings, the solid angle of one fan beam can be estimated by integration over the azimuthal angle as . Inside this solid angle, the ratio of the observed scattered flux to the observed inner disk flux, , is related to the ratio of the luminosities as where i is the inclination angle between the line of sight and the normal to the midplane of the disk. If we assume, for example, , and (i.e., i 75°), we obtain and , which are consistent with the observation. For this scenario, half of the H-band emission is produced at the (unresolved) inner disk, while the remaining half is the contribution from the scattered light. In this example, the inclination of the observer and the inclination of the scattering disk surface are very similar, and the observer remains inside the beaming solid angle of the near-side disk surface. Instead of the simple fan-beam geometry for the scattered NIR radiation, the NIR light could actually be scattered in all directions. Our explanation remains valid if the scattered flux density is still concentrated near the disk surface. If the disk is observed at small inclination angles (close to face-on view), the contribution of the scattered light to the observed flux decreases depending on the efficiency of the "beaming," while the flux from the inner disk increases. In other words, the observed manifestations of any large disks around XTINSs could vary significantly, depending on the optical properties of the disks, the inclination , and the observed wavelength. Equation (1) shows that for an observer inside the solid angle of the beamed scattering, , the observed could be significantly greater than . In our simplistic model of an assumed scattering anisotropy, the observer sees not the entire disk but only a fraction of the disk surface close to the line of sight, which, at a given outer radius, reduces the size of the observed extended emission in comparison to the actual size of the disk. However, the observed size increases with increasing contribution from multiple scatterings. A given light ray can go through multiple scatterings, both local and at different radii, which can significantly change its azimuthal direction. The details of the extended emission depend on the actual geometry and the scattering properties of the disk and require detailed modeling to quantify the "beaming" of the NIR emission and its observable azimuthal distribution. We also note that we extrapolated the results from protoplanetary disks to our case of a potential neutron star disk. The composition of such a disk is very likely quite different from that of a protoplanetary disk. For example, one could expect a higher metallicity and possibly a higher dust-to-gas ratio in supernova fallback material. In addition to thermal continuum, such a disk could also show line emission. Modeling of different disk geometries, densities, and compositions is beyond the scope of this paper, in particular because we currently have a good measurement in only one photometric band of the (N)IR. As an alternative to a disk, scattering of the emission from the inner disk by the dust grains in an additional outer cold dusty belt or torus structure around RX J0806.4–4123 could also explain the extended NIR emission. Our previous Herschel 160 μm detection for the region of RX J0806.4–4123 resulted in an estimate of dust grain locations at cm to cm for a reasonable dust temperature range of 10–22 K and in dependence of the grain size, a μm , in μm (Posselt et al. 2014). We had speculated about the existence of a belt structure to explain the large ratio of the luminosity at 160 μm to the XTINS' X-ray luminosity (for details, see Posselt et al. 2014). The association of the Herschel 160 μm detection with RX J0806.4–4123 is not as firm as the Hubble NIR detection because of the worse spatial resolution of Herschel. Nevertheless, it is intriguing that the spatial dimensions of the potential torus are of about the same order as the NIR extension. A dusty torus alone cannot explain the NIR emission because the scattered NIR emission from RX J0806.4–4123 would be expected to be much fainter. The implied two-component structure is not uncommon for debris disks around main sequence stars (e.g., Su et al. 2013) and may also be a possibility around a neutron star. Overall, a disk can be the source of the detected NIR flux, while the observed source extension can be due to scattering by cold dust grains located at larger separations from the neutron star, either in a disk with favorable properties or in a surrounding dusty torus.

3.3. A Pulsar Wind Nebula?

Extended emission around a pulsar can also be a PWN. Most PWNe are known from high-energy and radio observations (e.g., Reynolds et al. 2017). About a dozen have been detected in H α (Brownsberger & Romani 2014), and two have recently been discovered in the ultraviolet (e.g., Rangelov et al. 2017). A few high-energy PWNe were also detected in the (near-)infrared, e.g., the PWNe of Vela (Shibanov et al. 2003) and PSR B0540-69 (Mignani et al. 2012), but there is only one serendipitous discovery of an infrared-only PWN with Spitzer at wavelengths longer than 8 μm (Wang et al. 2013). The NIR emission around RX J0806.4–4123 could come either from line emission in the shocked ambient medium or from synchrotron emission of the shocked pulsar wind. In the case of a shocked ambient medium, possible origins are emission lines of H 2 and [Fe ii]. These shock-excited lines are frequently observed in the dense environments of star-forming regions (e.g., Herbig-Haro objects), evolved stars with eruptive mass loss, and supernova remnants (e.g., Reipurth & Bally 2001; Lee et al. 2014). Such a dense environment, however, is unlikely for RX J0806.4–4123 based on the X-ray absorption in the XTINS spectrum and 3D models of the ISM (Section 3.1). Hence, in the case of the PWN scenario, the emission is probably due to synchrotron emission from the shocked pulsar wind. At low , particles cannot be accelerated to sufficiently high energies to produce X-ray emission, but emission at lower energies could be possible. According to Reynolds et al. (2017), the maximum synchrotron photon energy E max of the accelerated electrons in the shocked pulsar wind can be estimated as where ζ ~ 1 is a numerical factor, is the spin-down power in units of 1030 erg s−1, B −6 is the magnetic field at the location of the shocked pulsar wind in units of μG, and σ is the magnetization, which is defined as the ratio of the Poynting flux to the particle enthalpy flux. Expected PWN magnetic field values are of the order of 10–100 μG depending on distance, while σ is unknown but usually assumed to be 1 or even σ 1 (Kennel & Coroniti 1984; Reynolds et al. 2017). The spin-down power of RX J0806.4–4123 is (Kaplan & van Kerkwijk 2009). Thus, E max ~ 1 eV (corresponding to λ ≈ 1.2 μm) is possible for σ 1. The X-ray efficiencies, , of observed PWNe vary greatly, reaching values between 10−6 and 10−2 (Reynolds et al. 2017). The similarly defined NIR efficiency of RX J0806.4–4123, η F160W = 4 × 10−3, is within the PWN efficiency range seen at other wavelengths. For a pulsar moving supersonically with a total velocity v PSR through the ambient medium with density ρ a , the stand-off radius R s of the bow shock apex can be described by , where is the fraction of spin-down power carried away by the wind, and the factor f takes into account a possible anisotropy of the pulsar wind (f = 1 for an isotropic wind). Using an observed source extension radius r ~ 0 4, which corresponds to R s ~ 1.5 × 1015 cm (Section 2), ξ w f 1, and (m p is the mass of a proton), we obtain where v ⊥ is the transverse velocity component. The proper motion of RX J0806.4–4123 is unusually low, mas yr−1 (2σ), with recent Chandra observations indicating ≈20 mas yr−1 (B. Posselt et al. 2018, in preparation). At a distance of 250 pc, these values correspond to v ⊥ < 100 km s−1 (2σ) and 24 km s−1, respectively. This does not significantly exceed the speed of sound in the typical diffuse ISM (c s ~ 10 km s−1 for T ~ 8000 K). The constraint on the velocity and number density in Equation (3) depends on the actual R s , but it is uncomfortably low, particularly considering that an anisotropic wind and a typically assumed fraction of the spin-down power emitted with the wind (ξ w ~ 0.5) will lower the value further. R s could, however, be smaller, since it is difficult to establish the exact location of the pulsar in the extended NIR emission. Hence, a PWN can currently not be ruled out as a possible origin of this emission. If RX J0806.4–4123 were a usual rotation-powered pulsar, the PWN would be an obvious explanation. Such an explanation is, however, somewhat surprising for a member of the Magnificent Seven, which show only thermal X-ray emission and no radio or γ-ray emission. A part of the detected NIR emission (the point-source contribution that can be at maximum 50%) could in principle come from the magnetosphere of RX J0806.4–4123. However, the optical-UV spectrum of the neutron star, F ν ∝ να, with α = 1.62 ± 0.20 (1σ) (Kaplan et al. 2011), is consistent with, and thus likely dominated by, a Rayleigh–Jeans tail of thermal emission (for which α = 2). If there is a PWN, then a lack of detectable non-thermal emission from the neutron star at higher energies than NIR is unprecedented. It is noteworthy that the XTINS spectra of other members of the Magnificent Seven, in particular RBS 1774,12 show much stronger deviation from a Rayleigh–Jeans tail, α = 0.53 ± 0.08 (1σ) (Kaplan et al. 2011). However, as discussed by Kaplan et al. (2011), a magnetospheric origin of the deviations from Rayleigh–Jeans is not the only possible explanation, in particular considering the high optical luminosities of these XTINSs. The current uncertainties in the models of strongly magnetized neutron star atmospheres at long wavelengths could be, for example, another explanation for the unusual spectra of the Magnificent Seven. If the extended NIR emission around RX J0806.4–4123 is indeed a PWN, it could indicate an interesting new avenue to study the properties of pulsar winds without very energetic particles originating from the (inner) magnetosphere. It is interesting that in comparison to typical spectral slopes of X-ray-detected PWNe, α X ≈ −0.5 (Kargaltsev & Pavlov 2008), the slope of the NIR PWN, α PWN −2.5, is rather steep according to Figure 4. Such a steep slope could be explained if the NIR photons correspond to the highest electron energies in the PWN.

3.4. Is RX J0806.4–4123 Special?