In this study we examine 51 archival NICER observations and six archival NuSTAR observations of the neutron star ultracompact Xray binary UCXB 4U 0614091 which span over 5 yr. The source displays persistent reflection features so we use a reflection model designed for UCXBs with overabundant carbon and oxygen XILLVERCO to study how various components of the system vary over time. The flux of this source is known to vary quasiperiodically on a timescale of a few days so we study how the various model components change as the overall flux varies. The flux of most components scales linearly with the overall flux while the power law representing coronal emission is anticorrelated as expected. This is consistent with previous studies of the source. We also find that during observations of the highsoft state the disk emissivity profile as a function of radius becomes steeper. We interpret this as the corona receding to be closer to the compact object during these states at which point the assumed power law illumination of XILLVERCO may be inadequate to describe the illumination of the disk.
['X-ray binary stars', 'Neutron stars', 'High energy astrophysics', 'Accretion', '1811', '1108', '739', '14', 'Astrophysics - High Energy Astrophysical Phenomena']
A View of the Longterm Spectral Behavior of Ultracompact XRay Binary 4U 0614091
{'A View of the Long-Term Spectral Behavior of Ultra Compact X-Ray Binary 4U 0614+091': "D. L. MOUTARD, 1, 2 R. M. LUDLAM, 1 E. M. CACKETT, 1 J. A. GARC'IA, 3 J. M. MILLER, 2 D. R. WILKINS, 4  \n1  \nDepartment of Physics & Astronomy, Wayne State University, 666 West Hancock Street, Detroit, MI 48201, USA 2 Department of Astronomy, University of Michigan, 1085 S. University, Ann Arbor, MI 48109, USA 3 Cahill Center for Astronomy and Astrophysics, California Institute of Technology, 1200 E. California Blvd, MC 290-17, Pasadena, CA, 91125, USA 4 Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, 452 Lomita Mall,Stanford, CA 94305-4085, USA", 'ABSTRACT': 'In this study, we examine 51 archival NICER observations and 6 archival NuSTAR observations of the neutron star (NS) ultra-compact X-ray binary (UCXB) 4U 0614 + 091, which span over 5 years. The source displays persistent reflection features, so we use a reflection model designed for UCXBs, with overabundant carbon and oxygen (XILLVERCO) to study how various components of the system vary over time. The flux of this source is known to vary quasi-periodically on a timescale of a few days, so we study how the various model components change as the overall flux varies. The flux of most components scales linearly with the overall flux, while the power law, representing coronal emission, is anti-correlated as expected. This is consistent with previous studies of the source. We also find that during observations of the high-soft state, the disk emissivity profile as a function of radius becomes steeper. We interpret this as the corona receding to be closer to the compact object during these states, at which point the assumed power law illumination of XILLVERCO may be inadequate to describe the illumination of the disk.', '1. INTRODUCTION': "In low-mass X-ray binary (LMXB) systems, a neutron star (NS) or black hole (BH) typically interacts with a main sequence star ≲ 1 M ⊙ . This interaction takes place via Rochelobe overflow, where the stellar companion fills its Roche lobe, until the material begins to fall towards the NS or BH, and eventually circularizes into an accretion disk. Typical LMXBs have orbital periods on the timescales of hours or days, but some display a much shorter orbital period (Avakyan et al. 2023). This subclass of LMXB is called an ultra-compact X-ray binary (UCXB) and is generally comprised of a NS or BH, and a H-poor degenerate companion like a white dwarf (WD) or He star. The orbital period required to be considered a UCXB is typically < 80 minutes, which would not be possible for larger main sequence or Hrich stars without directly colliding with the compact accretor (Savonije et al. 1986; Bahramian & Degenaar 2023). Because UCXB companion stars are chemically distinct from a main sequence star, the abundances in the accretion disk will also differ quite greatly from those LMXBs which are not ultra-compact. These systems typically have negligible amounts of hydrogen and overabundance of oxygen (Nelemans et al. 2004, 2006). Aside from being useful avenues of studying accretion physics and probing the nature of compact objects in a unique way, these sources will also produce  \nCorresponding author: D. L. Moutard  \nmoutard@umich.edu  \ngravitational waves on the order of mHz. This will be detectable by future multi-messenger missions such as LISA (Chen et al. 2020), so understanding and classifying these systems well in advance is useful.  \nLMXB systems are typically defined as having a few primary components. These are a disk which produces thermal X-rays (Shakura & Sunyaev 1973), a corona close to the compact object which produces hard, non-thermal X-rays via the Compton scattering of disk photons, and in the case of NS systems, thermal photons from the NS itself (Syunyaev et al. 1991). These contribute to the overall continuum of the X-ray spectrum, but coronal photons will also illuminate the disk and be reprocessed by the material therein. This leaves a signature as the photons are reprocessed, often in the form a relativistically broadened iron K-shell fluorescence emission around 6-7 keV(Fabian et al. 1989). In the case of UCXBs, the overabundance of oxygen means that this feature may be screened as interactions occur more readily with more abundant elements in the disk, and instead a more dominant O VIII feature appears in the softest X-rays, at around 0.67 keV (Koliopanos et al. 2013, 2021). These additional features are considered part of the reflection spectrum, and are believed to come from the innermost region of the disk (Tanaka et al. 1995; Bhattacharyya & Strohmayer 2007). Based on the broadening of these lines we can determine several useful quantities such as the location of the inner disk with respect to the compact object and the inclination of the system.  \nOne such UCXB is 4U 0614 + 091. This source was discovered in 1978 by the Uhuru mission and considered to be  \na binary system (Swank et al. 1978). Later observations of type-I X-ray bursts determined that the compact object was a neutron star (Brandt et al. 1992), followed by later spectral studies verifying the presence of elements that make the companion likely to be WD (Juett et al. 2001; Nelemans et al. 2004). The orbital period was measured by Shahbaz et al. (2008) using optical data to be approximately 50 minutes. This was done by directly observing modulations in the optical light curve caused by the uneven heating of the donor star.  \nMoutard et al. (2023) used reflection modeling to study 4U 0614 + 091 as the flux changed. The reflection model used therein is XILLVERCO, which is a modification of the XILLVER reflection table (Garc'ıa et al. 2013). The key difference is that XILLVERCO has abundances designed for CO WD, with reduced amounts of hydrogen and helium (Madej et al. 2014). The source appears to have a quasi-periodic flux, varying on the timescale of ∼ a week. The study concluded that the flux of most of the spectral components are correlated with the overall flux, but found that in the lowest flux state, the non-thermal flux was at its highest , and the disk appeared to truncate to about twice the innermost stable circular orbit ( R ISCO ). R ISCO is the closest a particle can stably orbit a compact object before falling onto it, and is equal to 6 R g ( R g = GM/c 2 , with 6 R g = 12.4 km for a non-rotating 1 . 4 M ⊙ NS). A parallel is drawn to BH LMXB systems, where typically at lower luminosities, the disk truncates (Done et al. 2007). The trend is less clear in NS systems, as the complicated magnetic fields affect the location of the inner disk, but since both NSs and BHs are compact accretors, this comparison is not unwarranted.  \nIn this paper, we extend the study from Moutard et al. (2023), supplementing the 4 NICER observations with additional archival observations. Our goal is to study the long term behaviors of the system as the overall flux varies, and to test further analogies between NS and BH accretion systems. While the previous study benefited from simultaneous NICER and NuSTAR observations, the observations only cover roughly 1 period of the quasi-periodic flux oscillations. Because of this, it is difficult to draw conclusions about the system more generally. By extending the study we are able to draw more robust conclusions about model parameters and their relation to the variation in flux over time. This will also serve as a test of the limits of XILLVERCO. The model currently assumes the bulk of the illumination is from the nonthermal corona, but it is well understood that in higher luminosity states, the continuum shifts to be dominated by softer, thermal X-rays (Lin et al. 2007). We choose NICER because of the large number of archival observations. NICER also provides excellent coverage in the band containing the O VIII feature, which is key to understanding UCXBs. We also supplement this with 6 archival NuSTAR observations. NuSTAR can be useful for constraining the overall continuum shape, as the band extends into the hard X-rays. Although there are many observations of 4U 0614 + 091 using Swift, the source is quite bright. Therefore, the observations are taken in windowed mode to account for additional photon pileup. These  \npileup effects can introduce additional uncertainties to the emission spectrum of a source (Miller et al. 2010; Ng et al. 2010), so we opt to use NICER and NuSTAR, which do not suffer from such pileup effects. These instruments also allow us to directly compare to our previously published literature. Section 2 of this paper discusses the details of the archival observations and the selection criteria, Section 3 provides details on the modeling and reports the relevant results of fitting spectra to these models, and Section 4 discusses the physical interpretations of these results. The final section, Section 5, summarizes this work.", '2. OBSERVATIONS AND DATA REDUCTION': "To study 4U 0614 + 091 we primarily use archival NICER observations because the large number available makes it possible to extend this study over the course of several years. We limited archival NICER observations to those with > 5 ks before filtering, leaving us with 51 total observations. This exposure cut is chosen to guarantee all observations (regardless of state) have above 10 6 total counts, ensuring a comparable data quality to even the lowest flux observations from Moutard et al. (2023). We then calibrate all of these observations using NICERDAS 2023-08-22 v011a and CALDB 20221001. This is done by first using NICERL2 for geomagnetic prefiltering, then followed by NIMAKETIME with COR SAX > 4 . 0 to manage particle overshoots and KP < 5 to filter out high particle background. We then construct good time intervals (GTIs) using NIEXTRACT-EVENTS. To construct spectra, backgrounds, and responses we use NICERL3SPECT with the 3C50 background model. We bin the spectra using the optimal binning scheme with a required count rate of 30 counts per bin to ensure χ 2 statistics are valid (Kaastra &Bleeker 2016). We construct light curves using NICERL3LC, which normalizes the light curves to the number of active detectors. These light curves are generated using several Xray bands, including the total band (0.5-6.8 keV), the super soft band (0.5-1.1 keV), the soft band (1.1-2.0 keV), the intermediate band (2.0-3.8 keV), and the hard band (3.8-6.8 keV). We compile these light curves into hardness-intensity and color-color diagrams (Figure 1), which demonstrate that the hard color of this system is largely unchanged over the course of observations. One observation (OBSID = 1050020115) contained a type-I X-ray burst, which was filtered out using a custom good time interval (GTI). The observation date of all 51 NICER observations are shown overlaid on a long term MAXI light curve in Figure 2, with insets that display 20 day intervals.  \nWe also look at 6 NuSTAR observations of 4U 0614 + 091. These are generated by first defining source regions in a 100 '' circle centered on the source, and another 100 '' circle centered off the source to define a background region using ds9. The reduction of these observations is done using NUSTARDAS v2.1.2 and CALDB 20230816. More information on all observations can be seen in Table A.1 in Appendix A. We plot the the NuSTAR light curves using 1 second time bins, and detect no bursts in any observation, so no further filtering was done.  \nFigure 1. The hardness intensity diagram for the soft color (1.1-2.0 keV/0.5-1.1 keV) and the hard color (3.8-6.8 keV/2.0-3.8 keV), as well as the color-color diagram for all 51 observations. The observations are binned into 150s segments for clarity. We see that the hard color stays relatively constant throughout all observations, but the soft color varies more distinctly over time.  \n<!-- image -->", '3. RESULTS AND ANALYSIS': 'In this section we describe the specifics of the models used in this work to study the reflection spectrum. We begin by first modeling only the continuum components, then proceed to add a reflection component to the model.', '3.1. NICER Continuum Modeling': 'Before using reflection models, we first test a continuum model similar to the one used in Moutard et al. (2023). The models are constructed and fit in XSPEC VER. 12.13.1 (Arnaud 1996). The first component is TBABS, a multiplicative model component used to measure absorption due to neutral hydrogen in the interstellar medium (ISM), using the WILM cross section (Wilms et al. 2000) and the VERN cross-section (Verner et al. 1996). This is followed by two multiplicative absorption edges to manage features at ∼ 0 . 4 and ∼ 0 . 9 keV, which can interfere with the O VIII line. These features are potentially astrophysical in origin, owing to the ISM along the line of sight (Pinto et al. 2013). The inclusion or exclu-  \nof these edges may interfere with measurements of the O VIII feature. For the sake of consistency with Moutard et al. (2023), we continue this analysis with two edges, but we discuss the implications further in Section 4.4. The higher and lower energy edges roughly correspond to the Fe L and C K edge respectively, though the edge at lower energy has a centroid that sits near or below the band used, so it is less confidently measured. We are deferring to previous studies of this source and similar sources which use such edges (Madej et al. 2014; Ludlam et al. 2020; Moutard et al. 2023). After fully modeling the spectra (see Section 3.2 for more details), we find the edges have average energies at ∼ 0 . 40 and ∼ 0 . 87 keV, with respective optical depths of ∼ 0 . 85 and ∼ 0 . 18 . The lower energy edge is roughly consistent with a known nitrogen detector edge, while the higher energy edge sits in the region of a Fe-L or Ne-K edge. The remaining model components are BBODY to account for thermal emission from the NS, DISKBB to account for thermal emission from the disk, and CUTOFFPL to account for the nonthermal emission from the corona. We use CUTOFFPL rather than a standard POWERLAW because many observations in the highsoft state display a high energy cutoff within the bands observed (Degenaar et al. 2018). Using this while allowing the cutoff energy to increase beyond the band provided by NICER allows us to account for spectra both with and without this energy cutoff, allowing our model to be agnostic to the spectral state. The final model used to fit the continuum is TBABS*EDGE*EDGE*(BBODY+DISKBB+CUTOFFPL). All of the models are fitted in the 0 . 45 -9 keV range for consistency with Moutard et al. (2023).  \nWith the model in place, we proceed to fit all 51 NICER spectra. We use χ 2 statistics with Churazov weighting. Churazov weighting is similar to standard χ 2 statistics, but takes into account adjacent bins to help prevent overfitting (Churazov et al. 1996). However, this can also lead to higher values of χ 2 statistics. If we were using standard weighting, the reduced χ 2 values reported would often be above acceptable levels, but the values found in this analysis are consistent with those from Moutard et al. (2023). With the fit in place we then calculate the errors on continuum parameters by using the fitted model as a starting point for a Markov-chain Monte Carlo fit with 100 walkers, a burn-in of 100000, and a chain length of 10000. This is done using the CHAIN command in XSPEC. To demonstrate the quality of these fits, we show a histogram of reduced χ 2 values in Figure 3, and compare the continuum-only model to a model with a reflection component added (as discussed in Section 3.2). The inclusion of the reflection component leads to an improvement of in χ 2 values of ∼ 1566 for 6 degrees of freedom. A table of all continuum values can be found in Appendix B.', '3.2. NICER Reflection modeling': "With the continuum models in place we then proceed to add a reflection component. This reflection component is modeled using RELCONV*XILLVERCO. As mentioned previously, XILLVERCO is an adaptation of XILLVER with high carbon and oxygen abundances to account for UCXB fea-  \nFigure 2. Long term MAXI light curve in blue points with the start date of each observation overlaid as dashed lines. The top row of plots highlights specific 20-day regions of the light curve, demonstrating that the source is captured in various states. These zoomed-in plots also display the quasi-periodic nature of the flux of the source.  \n<!-- image -->  \ntures. RELCONV is a relativistic convolution kernel that accounts for the broadening of the features in XILLVERCO. In RELCONV we tie the two emissivity indices q 1 = q 2 = q . This makes the parameter R BR redundant, so we fix it to 500 R g . Since we are dealing primarily with the inner disk region, we also fix the outer radius to 990 R g , as we do not expect this parameter to affect our results. We fix the dimensionless spin parameter a to 0. The limb darkening parameter, used to control whether the model accounts for limb darkening in the source, is also set to 0. This leaves us with 3 free parameters in RELCONV: q , as mentioned above, the inclination i , and the location of the inner disk radius R in , displayed in this paper in units of R ISCO . For XILLVERCO and RELCONV, we fix the inclination to 55 · ; this value is an average of measured values from previous works (Madej et al. 2014; Ludlam et al. 2019), and serves to reduce the number of free parameters. The same was done in Moutard et al. (2023), which finds that the parameters do not vary greatly when the inclination is allowed to be free. We tie the powerlaw index and energy cutoff in XILLVERCO to those from CUTOFFPL.  \nThe remaining free parameters are the normalization 1 , the carbon-oxygen abundance A CO , the frac parameter, and the temperature kT xill . The frac parameter describes the ratio of the illuminating flux from the power law to that of the blackbody which emerges from the disk at the point of reprocessing. kT xill is the temperature of that emergent blackbody. The ionization log ξ is not a parameter in XILLVERCO, but can be calculated using  \nξ = 4 π n F x (1)  \nwhere F x = frac × σT 4 with frac being the same as the frac parameter, σ being the Stefan-Boltzmann constant, and T being the temperature from kT xill (Garc'ıa et al. 2013). n is the disk number density, which is fixed to 10 17 cm -3 in XILLVERCO.  \nFigure 3. A comparison of the reduced χ 2 values for the observations when fit only with a continuum, versus when fit with a reflection component. We see that a reflection component significantly improves the reduced χ 2 statistics for all observations. This should not be understood as a robust statistical comparison of the models, but rather a quick reference to determine whether a reflection component significantly improves the fit quality. The values shown here correspond to an average χ 2 improvement of ∼ 1566 for the additional 6 degrees of freedom. Note: the two models are binned differently for readability.  \n<!-- image -->  \nWe add the reflection component to the continuum model, and then fit using XSPEC. In order to calculate errors, we then repeat the process described for the continuum in Section 3.1. We also calculate the fluxes (both overall and for each component) using the multiplicative model component CFLUX once the models are fit. The results of fitting all of these models can be found in Appendix B. We measure the equivalent width (EW) of the O VIII feature by applying a DISKLINE model component to the continuum model. DISKLINE is a model designed to fit emission features in accretion disks (Fabian et al. 1989). We fit the spectrum and use the EQWIDTH command in XSPEC to extract the EW and 90% confidence limit errors. We fix the inclination and outer disk radius to the same values used in XILLVERCO, and allow the emissivity index, inner disk radius, line energy, and normalization to be free. We use the same chain method described above to calculate errors, and the results can be found in the Table B.4. These equivalent widths are consistently around ∼ 50 eV, which is similar to the EW of other measured O VIII lines, and is comparable to measured EWs for Fe lines in LMXBs (Madej & Jonker 2011; Cackett et al. 2010). The centroids of these lines are also consistently near 0.7 keV, as anticipated.", '3.3. NuSTAR Spectra': 'We use the 6 NuSTAR spectra primarily to test whether the high energy components of the continuum remain constant over time. These spectra unfortunately do not probe the  \nsoftest regions of the X-ray band. This means they can not be used to measure the O VIII feature, and so they are used primarily for constraining high energy continuum features. A similar power law index was seen in all observations of Moutard et al. (2023), so we see if this remains true for observations separated further in time. If that were the case, we could use this fact to fix certain parameters, which would be useful in dealing with the large number of spectra involved in this study. However, we see that the slope of the hard X-ray spectrum varies over time (see Figure 5). This means that we can not use this to constrain and fix the slope of the power law component, but does indicate that the spectral state is varying over time.  \nTo study the constraints we can place using the NuSTAR spectra only, we use a similar model to the one described above and model the spectra from 3 to 30 keV, with the spectra binned optimally again with a minimum of 30 counts per bin. This band is chosen to be consistent with Moutard et al. (2023). The main differences between the model used for the NuSTAR spectra and the NICER spectra is that we omit the EDGE components, as those model features in the low energy, and we add a constant component to account for differences between NuSTAR FPMA and FPMB. The constant is fixed to unity for FPMA and allowed to be free for FPMB. Since TBABS primarily affects the lowest energy bands, we fix the value of N H to 0 . 4 × 10 22 cm -2 , an approximate average of the reflection model. We find that the results of certain parameters, such as the power law and the blackbody, are comparable to NICER values and to those in literature (Moutard et al. 2023). However, other parameters vary more significantly, with predictions of extremely truncated inner disks, unrealistically low (and sometimes negative) emissivity indices, low A CO , and occasionally poorly constrained values from DISKBB. This highlights the utility of instruments with lower energy band passes for UCXB studies, and demonstrates the ability of NuSTAR to help constrain continuum parameters. The results of these fits can be found in table C.1 in Appendix C.', '4.1. Flux and Inner Disk Radius': "Moutard et al. (2023) discusses 4U 0614 + 091 in terms of the overall flux behavior. In that analysis, the flux of all model components scale with the overall flux except the flux of the powerlaw component, which seems to be anticorrelated. The decrease in flux of the XILLVERCO component as the illuminating powerlaw increases is explained by a truncated disk during the lowest flux observation. This can also be explained by the fact that the fraction of emitted coronal photons to those reflected decreases as the height of the corona increases. If the coronal height decreases as the source gets softer and the powerlaw becomes weaker, this means more photons illuminate the disk and are reflected, resulting in a greater contribution from XILLVERCO (Dauser et al. 2016). We find in this analysis that the flux patterns persist over longer time scales as well. We show in Figure 6 that, in general, the flux of each component is increasing  \nFigure 5. Shown here are the 6 archival NuSTAR spectra. The 5 most recent spectra have comparable slopes in the higher energy bands, but the earliest observation, 30363002002, is noticeably different. This means that we can not effectively use NuSTAR observations which are not simultaneous with NICER observations to constrain continuum parameters for the broader study.  \n<!-- image -->  \n<!-- image -->  \nFigure 4. An example of 3 unfolded NICER spectra with the relevant model components. We show here (a) the lowest flux observation, obsid = 4701010301, (b) an intermediate flux observation, obsid = 1050020236, and the highest flux observation, obsid = 1050020252. These display the range of fluxes we see across these observations, as well as the most relevant components for various fluxes. In (a), we see the power law is dominant, and remains relatively flat out to the highest energies, while in (c) the cutoff and preference for a stronger thermal component becomes more apparent. In each, the O VIII feature is visible around 0.7 keV.  \n<!-- image -->  \n<!-- image -->  \nTotal Flux (10  \nwith overall flux except the powerlaw component. This is consistent with the well-documented low-hard and high-soft states, where higher flux observations are more thermally dominated and lower flux states are more dominated by nonthermal coronal emission. These correlations are supported by a spearman rank ≥ 0 . 75 for all components except the anticorrelated power law, which is at a spearman rank of ∼ -0 . 5 . We also note that the anticorrelation of the powerlaw flux and XILLVERCO flux discussed above is supported by a spearman rank of -0 . 6 .  \nThe truncated inner disk at the low-flux state is pointed to as an analogy for BH-LMXB systems, which tend to show a  \nFigure 6. Shown are the fluxes of each component for each observation plotted against the overall flux of the system. In general, the component fluxes are well-correlated with the overall flux (spearman rank ≥ 0 . 75 ), with the exception of the powerlaw flux, which is anti-correlated (spearman rank ∼ -0 . 5 ).  \n<!-- image -->  \nlarger value of R in for low luminosities (which is a correlated with the mass accretion rate, Done et al. 2007; Tomsick et al. 2009; Garc'ıa et al. 2015). However, over longer timescales, 4U 0614 + 091 does not appear to strongly show any correlation between inner disk radius and flux, aside for a slight preference for truncation in the lowest flux states. This relationship is shown in Figure 7, and even though the lowest flux states display some disk truncation, we see a random assortment of higher flux observations with equally or more truncated disks. We test a handful of lower flux observations to see if the disk truncation can be accounted for by the magnetic field. In some NS LMXBs, the magnetic field is able to truncate the disk at the Alfv'en radius ( r A ), where the energy density of the magnetic field equals the energy density of the accreting material. To do this we calculate the strength of the  \nequatorial magnetic field during our most truncated observation using the following equation from Ludlam et al. (2020):  \nB = 3 . 5 × 10 5 k -7 / 4 A x 7 / 4 ( M 1 . 4 M ⊙ ) 2 ( 10 km R NS ) 3 × ( f ang η F bol 10 -9 ergs -1 cm -2 ) 1 / 2 D 3 . 5 kpc G (2)  \nwhere the efficiency η is assumed to be 0.2 (Sibgatullin & Sunyaev 2000), the conversion factor k A and the angular anisotropy f ang are set to unity (Ludlam et al. 2019), and the distance D to 4U 0614 + 091 is 3.3 kpc (determined using GAIA data in Arnason et al. 2021). We use canonical NS values M = 1 . 4 M ⊙ and R NS = 10 km, and the 0 . 5 -50 keV flux and the inner disk radius (in R g ) are used for F bol and x respectively. If we use the lowest flux observation to determine our magnetic field strength, this provides an upper limit on the magnetic field of B 1 ≤ 0 . 4 × 10 8 G. We convert this value to a magnetic dipole moment (again assuming the values above) and find µ 1 = 3 . 8 × 10 25 G cm -3 . We use this value and mass accretion rate ˙ m = L/ηc 2 where L is the luminosity during that observation to calculate  \nr A = ( µ 4 GM ˙ m 2 ) 1 / 7 (3)  \nwhere G is the gravitational constant. This radius is plotted as a red dashed line on Figure 7. We see that this radius is generally too high to explain the inner disk radius in this case. These radii are generally too large to assume that the disk is truncating due the magnetic field. If instead we use the next lowest flux observation with some indication of truncation to determine our magnetic field strength, we find B 2 ≤ 0 . 3 × 10 8 G and µ 2 = 2 . 6 × 10 25 G cm -3 , for which we plot r A as a blue dashed line in Figure 7.  \nIf instead we assume k A = 0 . 5 (Long et al. 2005), µ increases by roughly a factor of three, which can drastically increase r A , especially in the case of µ 1 . This would extend r A beyond the truncated inner disk, meaning truncation can not be caused by the magnetic field. If instead we assume k A = 1 . 05 (below the theoretical upper limit of 1.1, discussed by Ibragimov & Poutanen 2009), we find that in some cases, the magnetosphere is consistently at or below the disk truncation. In such a case, it is possible that disk truncation is caused by the magnetic field during the low flux states. We can see similar changes to r A by assuming different values of f ang , (though it is expected to be near or above unity Ibragimov & Poutanen 2009) or η (which can vary depending on the rotational frequency of the NS). Because of this, we can not confidently say that the low flux disk truncation is caused exclusively by the magnetic field.  \nAnother possibility for the cause of the disk truncation is disk depletion at low mass accretion rates. However, this was explored for 4U 0614 + 091 in (Moutard et al. 2023), which determined that the timescale required does not correspond with the cadence of observations. It should be noted that  \nFigure 7. The overall flux of the system is poorly correlated with the inner disk radius. We see the lowest flux states display some degree of truncation as stated in Moutard et al. (2023), but there are also higher flux states which are truncated just as heavily, if not moreso. To test whether the observed disk truncation is caused by magnetic fields surrounding the NS, we plot the magnetosphere radius for various limits. Assuming our magnetic field upper limit (red dashed line) is well characterized by the most lowest flux point does not yield a magnetosphere radius consistent with most measured R in values. If instead we use a different low flux, seemingly truncated observation to characterize the magnetic field strength (blue dashed line), we find a slightly better relation, though it still can not be concluded that the magnetic field is the sole source of inner disk truncation in this system. k A , f ang , and η are all degenerate, so it is possible to make these magnetic field radii very small or very large, even using reasonable values for this parameter.  \n<!-- image -->  \nthe degree of truncation is quite small, only extending up to ∼ 1 . 5 R ISCO (compare this to truncation up to ∼ 10 R ISCO as seen in Figure 5 of Ludlam 2024). This small degree of truncation, paired with the limitations of XILLVERCO at high flux (discussed in more detail in Section 4.3) makes drawing conclusions about the cause of the truncation difficult.", '4.2. Emissivity Index': 'The emissivity index q is a value that represents the power law index defining the illumination of the disk as a function of radius, r -q . In the classical limit, q = 3 , but this can be affected by multiple factors, such as the geometry of the illuminating corona and relativistic light bending (Wilkins 2018). One can use the emissivity index as a stand in for how heavily light is curved toward the disk; in BH systems, this index is often much higher as the gravitational effects curve more of the emission to the regions of the disk closest to the BH. This leads to a steeper fall-off in illumination near the inner disk ( q ≳ 6 in the case of rapidly spinning BHs), before returning to the classical limit at the outer disk (Wilkins & Fabian 2012).  \nFor NS systems, these values are typically closer to the classical limit. In this study we use only one emissivity index to minimize the number of free parameters, shown in Figure 8. By using only one value for q , we are essentially  \nFigure 8. Shown is the emissivity index q plotted against the total unabsorbed flux of the system. Flatter emissivity profiles are generally consistent with the system in a low-hard state, while steeper values of q correspond with the high/soft state. This is potentially explained by a corona which is extending spatially during low flux states.  \n<!-- image -->  \nadopting an average value of the inner and outer emissivity indices. Values of q < 3 , or a flattening of the emissivity index, may be indicative of a spatial extension of the illuminating corona. This may mean a vertical extension above the axis of the disk or a radial extension as the corona covers the inner region of the disk. We see in Figure 8 that for low flux sources, q is generally consistent with a slightly flatter profile, while in higher flux states, the index tends to increase. This is also correlated with the spectral state of the system, as shown by the color bar representing the index of the power law. A larger power law index is indicative of a softer spectrum, which is expected for higher flux states.  \nThis relationship is potentially analogous to what is seen by Kara et al. (2019) for the BH transient system MAXI J1820+070. In that study, they propose a geometry wherein the corona extends vertically above the disk as the spectrum becomes harder. Similarly, Buisson et al. (2019) finds a similar result for the same source, indicating that the height of the illuminating corona decreases as the spectrum softens. While NSand BH XRBs are not exactly the same, the accretion processes are similar, so an analogy between these systems is not unreasonable. 4U 0614 + 091 is known to contain a relativistic jet during its hard states, though the jet alone is likely not a viable mechanism to produce hard X-rays, and would also require some other source of Comptonized photons (Migliari et al. 2010). It is unlikely the jet is turning on and off on the time scale of weeks, following the quasi-periodic flux in the MAXI light curve. This provides evidence that the jet is not the driving source of the coronal photons, though it may contribute in some cases.  \nAnother possible geometry of the spatially extended corona is one that expands radially rather than vertically above the disk. If at some point the disk becomes highly magnetized, this could provide the fields needed to produce the Comptonized corona which would then illuminate the disk from directly above it. This too would lead to a flattening of the emissivity index over the region where the corona sits (Wilkins 2018). However, at least in AGN systems, this coronal geometry is predicted primarily for high-flux states (Wilkins et al. 2014). In 4U 0614 + 091, the flattening of q appears to occur in the lower luminosity states, opposite to what is seen for systems where this geometry is well studied. Because of this, as well as the known existence of a jet for this source, we deem that geometry to be less likely. Therefore, the most likely scenario is a vertically extended corona. However, as mentioned in Section 4.1, the limitations of XILLVERCO make drawing broad conclusions about the transition into high-soft states difficult.', '4.3. Limitations of XILLVERCO': "It is worth noting that XILLVERCO assumes disk illumination by a non thermal Comptonized component. This works well for UCXBs in the low-hard state, but may begin to falter when faced with disks that are illuminated by thermal emission from near the NS itself. This regime is represented by the the high-soft state, where we see q increase. This may serve to explain some of the discrepancies in the values of q for the few high flux but apparently soft points in Figure 8. These observations don't seem to agree with the overall trend of q increasing in the high flux states. These outlier observations are also among those with the worst reduced χ 2 value. The self consistent reflection model RELXILL also assumes illumination by a Comptonized corona, but the extension of the model, RELXILLNS, assumes thermal illumination. Similarly, a version of XILLVERCO with thermal illumination would be useful to describe the source in all states.  \nFigure 9 demonstrates one such limit of XILLVERCO. The bottom right region of the plot shows the coolest blackbody temperatures being associated with shallow photon indices and high energy cutoffs, as expected for harder spectral states. As the blackbody temperature increases, we see either an increase in the power law index or a decrease in the cutoff energy (the two parameters are somewhat degenerate, especially when the energy cutoff sits beyond the NICER band used in the study). However, in the top left region of the plot, we see that the observations with the highest blackbody temperatures demonstrate a very soft powerlaw spectrum. These represent a regime where the illuminating power law shifts to the lowest energies and cuts off quite steeply. A steep, absorbed power law with a very low energy cutoff mimics the shape of a blackbody. Figure 1 in Garc'ıa et al. (2022) demonstrates the difference between XILLVER and XILLVERNS. We see therein that even when the incident spectrum mimics the shape of a blackbody illuminating the disk, the resulting reflection spectrum depends heavily on the input model used. For this study, that means that using a model with the incor-  \nFigure 9. Shown here is the powerlaw photon index plotted against each spectrum's powerlaw cutoff energy, color coded by the blackbody temperature. While many observations demonstrate the expected behavior of a soft spectrum with a high blackbody temperature, some of the highest temperature blackbodies are coincident with steep illuminating powerlaws with very low energy cutoffs. These observations are at the edge of the limits of XILLVERCO, which assumes illumination by a powerlaw, even if the disk is being physically illuminated by the thermal NS boundary layer. In these cases, the model attempts to force the illuminating powerlaw into a shape resembling that of a blackbody.  \n<!-- image -->  \nct input spectrum can cause some of the modeled parameters to tend towards unphysical values.", '4.4. Absorption Edges': 'In the model we include two absorption edges in the lowest energies to account for apparent features at roughly 0.4 and 0.87 keV (corresponding roughly to a N detector edge and a Fe-L or Ne-K edge, respectively.) The inclusion of absorption edges surrounding the O VIII feature at ∼ 0 . 7 keV can cause some interference with the measurements of the feature itself. Between the two edges, the higher energy edge is much more well constrained, as shown in Figure 10. The lower edge is less constrained, but that is to be expected as it sits at or below the lower end of the NICER band in use. As a result we initially test fits using only one edge centered near 0.87 keV, while keeping the rest of the model the same. The reduced χ 2 of the fits worsen marginally when the lower energy edge is removed. We find that the distributions of most parameters do not vary greatly, with some exceptions. A CO tends to increase dramatically when one edge is removed, offering values up to 20 times higher than those from models with two edges. These carbon-oxygen abundances, reaching up to several hundred times solar is inconsistent with Madej et al. (2014) when accounting for the fact that the abundances in that paper are artificially increased in that version of the model. Moutard et al. (2023) finds the  \nFigure 10. The distribution of edge energies for the model containing both edges. We see that the higher energy edge near 0.87 keV (consistent with a Ne-K or Fe-L edge) is very well constrained, while the lower energy edge around 0.4 keV (consistent with a N detector edge) has a significant amount of spread. This is to be expected for a feature that sits at or below the given NICER band, but should be addressed more carefully in future studies with higher resolutions and improved collecting area at the softest energies.  \n<!-- image -->  \nvalue of A CO to be near 25, which is broadly consistent with Madej et al. (2014) (though it should be noted that 2 edges are also used in Moutard et al. 2023).  \nWe also find that the model with only one edge tends to fit the soft state using a harder powerlaw photon index than is physically reasonable for these states, although with a low energy cutoff. This model trends towards Γ ≤ 1 . 2 during the softest states, which is only observed in the extreme hard states for NS systems (Ludlam et al. 2016; Parikh et al. 2017). While the trends between the models do not match, the conclusions drawn from Figure 9 remain: in the soft state, the powerlaw component cutoff energy shifts to imitate a thermal component. With only one edge, the trend seen in Figures 8 and 7 is also less consistent, with more spread in these parameters. We further test the validity of the edges by fixing their energies at 0.4016 keV and 0.8669 keV, corresponding to the laboratory values of the N edge and the Ne-K edge. We choose the Ne-K edge value over the Fe-L edge as it is much closer to the average free value of 0.8677. Similarly to the test with one edge, while some outliers do appear, the overall distributions of parameters remain similar to the free energy fits. Again we see some increase in A CO , though not as dramatically as in our single edge test. We also see some discrepancies in the disk blackbody temperature, which is unsurprising as the central energy often coincides with the lower energy edge. However, other key parameters studied over the course of this paper, such as R in and the emissivity index, maintain very similar values after freezing the edge energies. For completeness, we also test a handful of observations with no edges at all. When not accounting for the  \nwell constrained higher energy edge, we found that the reduced χ 2 increases drastically for every observation checked, sometimes by a factor of greater than 5. Because of this, we do not pursue such a model further. Without a careful study of these absorption edges at low energy, it is difficult to say which model best explains the physical nature of the system. As a result we primarily focus on the model that is a more direct comparison with Moutard et al. (2023), but this does encourage a low energy, high resolution study of this source to better constrain the edge components.', '5. CONCLUSION': "In this study, we examine the long term behavior of 4U 0614 + 091 via 51 archival NICER spectra. These spectra cover a range of time spanning over 5 years, and capture the source in various different spectral and flux states. After using XILLVERCO, a reflection model designed for UCXBs, we are able to say the following about this source:  \n- · The flux behavior is consistent with Moutard et al. (2023), which follows along the conventional highsoft and low-hard states. That is, as the overall unabsorbed flux of the system increases, so do the fluxes of most of the individual model components. The exception is the flux of the powerlaw, which represents the comptonized corona. This component decreases as the source gets brighter.\n- · The relationship between the flux and R in is seemingly uncorrelated, but displays some agreement with Moutard et al. (2023), which finds that the lowest flux observations appear to show a truncated disk. This is true in some cases, but drawing any conclusion about the behavior of the inner disk is difficult with the given data.\n- · The emissivity index tends to increase during the highsoft state when fit with XILLVERCO. This is indicative of the coronal extent increasing during harder states. A jet may exist during these states, and could provide some contribution to the illumination of the disk from above the NS itself, but it is not likely the driving source. In softer states, the illuminating component may recede to be closer to the NS, causing the light to be curved more steeply toward the inner region of the disk. The picture of a corona extending further from the compact object as the spectrum hardens is also analagous to certain BH XRB systems.\n- · The reflection model used here, XILLVERCO, may begin to falter in high flux states, where disk illumination from a nonthermal corona gives way to thermal illumination from a NS boundary layer. In these cases, the model appears to shift the incident powerlaw to the absolute lowest energies, mimicking the shape of a thermal illuminating component. During such states, the reflection model depends heavily on the correct illuminating component, and therefore makes drawing conclusions about the source in the high flux state difficult. This combined with a poor characterization of the lowest energy absorption edges leads to a model that is not complete when discussing UCXBs.  \nThis long term study of 4U 0614 + 091 paints a picture of a jetted UCXB which behaves similarly to other NS or BH LMXBs. The hard component appears to decrease in flux as the source gets brighter, and the inner disk truncation is poorly correlated with luminosity. The emissivity index provides a clue that a jet may contribute to the corona in hard states, while the illumination geometry shifts to be closer to the NS in softer states. Future high resolution studies of reflection features with missions such as XRISM could help to better constrain the emissivity index and therefore the geometry of the illuminating component. This, paired with simultaneous radio measurements of the jet could also be enlightening to determine the extent of the jet's contribution to the corona as the spectrum changes between states. Polarimetry instruments such as IXPE can be used to probe coronal geometry in some systems, and could potentially help to constrain our predictions of the extent of the corona in the soft states (Krawczynski et al. 2022). Also useful for this analysis would be a UCXB model which accounts for disk illumination by the NS itself rather than a comptonized corona.  \nAcknowledgements: This work is supported by NASA under grant No. 80NSSC22K0054. This research has made use of MAXI data provided by RIKEN, JAXA and the MAXI team (Matsuoka et al. 2009). This research has made use of data and/or software provided by the High Energy Astrophysics Science Archive Research Center (HEASARC), which is a service of the Astrophysics Science Division at NASA/GSFC. This research has made use of the NuSTAR Data Analysis Software (NuSTARDAS) jointly developed by the ASI Science Data Center (ASDC, Italy) and the California Institute of Technology (USA).", 'REFERENCES': 'Arnason, R. M., Papei, H., Barmby, P., Bahramian, A., & Gorski, M. D. 2021, MNRAS, 502, 5455, doi: 10.1093/mnras/stab345  \nArnaud, K. A. 1996, in Astronomical Society of the Pacific Conference Series, Vol. 101, Astronomical Data Analysis Software and Systems V, ed. G. H. Jacoby & J. Barnes, 17  \nAvakyan, A., Neumann, M., Zainab, A., et al. 2023, A&A, 675, A199, doi: 10.1051/0004-6361/202346522  \nBahramian, A., & Degenaar, N. 2023, in Handbook of X-ray and Gamma-ray Astrophysics. Edited by Cosimo Bambi and Andrea Santangelo, 120, doi: 10.1007/978-981-16-4544-0 94-1  \nWilkins, D. R., Kara, E., Fabian, A. C., & Gallo, L. C. 2014, MNRAS, 443, 2746, doi: 10.1093/mnras/stu1273  \nWilms, J., Allen, A., & McCray, R. 2000, ApJ, 542, 914, doi: 10.1086/317016', 'APPENDIX': 'A. OBSERVATION INFORMATION  \nThis section of the appendix summarizes all observations from NICER and NuSTAR used in this analysis.  \nTable A.1. 4U 0614 + 091 Observation Information', 'B. NICER FIT PARAMETERS': 'This section of the appendix summarizes all of the final fit parameters used in both the initial continuum model as well as the reflection model for NICER.', 'C. NuSTAR FIT PARAMETERS': 'This section of the appendix summarizes the results of reflection modeling for the NuSTAR spectra. Note: † indicates the parameter reached an upper bound and had no constraints on uncertainty.  \n14', 'T able B.1. 4U 0614 + 091 Results of Continuum Modeling': '<!-- image -->  \nT able B.2. 4U 0614 + 091 Results of Reflection Modeling: Continuum P arameters  \n<!-- image -->  \n<!-- image -->  \n-  \n10  \n-  \n2  \n-  \n1  \nTable B.4. 4U 0614 + 091 Results of Reflection Modeling: Calculated Fluxes and Equivalent Widths  \n-  \n0  \n.  \n004', 'T able C.1. 4U 0614 + 091 Results of Reflection Modeling: NuST AR': '<!-- image -->'}