Skip to main content

Come analyze HEASARC, IRSA, and MAST data in the cloud! The Fornax Initiative is now welcoming all interested beta users.

Xspec Home Page

New in version 13.0.0

Version 13.0 is a major release whose principal theme is Bayesian inference and global optimization. A new chapter on Bayesian methods (Chapter 3) and a new Algorithms appendix (Appendix H) describe the analytic-gradient and response-compression machinery that underpins many of the features below.

Bayesian inference and optimization

  • Closed-form analytic derivatives (Jacobians and their reverse-mode vector-Jacobian products) are now available for the great majority of the model library. These provide a fast Levenberg-Marquardt path, in which the gradient is computed analytically in a single pass rather than by finite differences, and supply the posterior gradient used by the new Hamiltonian Monte Carlo sampler. The pipeline is described in Appendix H.1.

  • New hmc command: a No-U-Turn-Sampler variant of Hamiltonian Monte Carlo over the thawed parameters, using the analytic posterior gradient. It supports warmup adaptation of the leap-frog step size and diagonal mass matrix, multiple independent chains (parallelized with parallel hmc), and checkpoint/resume. Output is written in a format compatible with chain load.

  • New nest command: a multi-ellipsoidal nested sampler that produces weighted posterior samples together with an estimate of the marginal likelihood (model evidence, $\log Z$). It supports checkpoint/resume and parallel sampling (parallel nest).

  • New compare command: automated model selection. Alternative models are evaluated against the current data and ranked by the Akaike and Bayesian information criteria, by Bayesian evidence (compare nest), or by fitting each alternative with the global optimizer (compare global).

  • New fit global option: a derivative-free global optimization by Differential Evolution followed by a polish with the current fit method. It is never-regress (the original fit is restored if no improvement is found), reproducible under xset seed, parallelizable with parallel global, and reports any competitive minima within $\Delta\mathrm{stat} < 9$ of the best fit. The improve command now uses the same Differential Evolution machinery as a local warm-restart and works with all fit methods, not just the Minuit methods.

  • The Minuit fit methods (migrad, simplex) were overhauled: Minuit now works directly in the parameter's stored value (so a fit started between a parameter's soft and hard limits reaches the same minimum as one started inside them), hard limits are passed to Minuit only when they span a range of 100 or less, and a fit that stops short of the minimum now says so instead of reporting convergence — with different warnings for running out of function evaluations (fit again to continue) and for Minuit finding no direction of improvement. See the method command.

  • New bayes smooth subcommand: a multi-parameter smoothness prior that couples a chosen parameter across datagroups, penalizing shell-to-shell roughness. It is the in-engine regularizer for projct de-projection, with automatic selection of the regularization strength by the L-curve or discrepancy principle. Several new single-parameter prior families were also added (cauchy, beta, studentst, halfnormal, powerlaw, isotropic, lognormal, gamma, and shiftedlognormal).

  • Added an MCMC diagnostics suite, including the Rubin-Gelman $\hat{R}$ convergence statistic, and gave chain, hmc and nest a common sample format.

  • PyXspec gained a unified Bayesian-inference module (runHMC, runChain, shared result objects, and arviz interoperability), together with Fit.compare and Fit.globalFit.

New fit statistics

  • chistokes: a chi-squared variant that uses the per-bin Stokes Q/U cross-covariance when fitting polarization spectra, rather than treating Q and U as independent. It has been generalized to an arbitrary $n$-member cross-spectrum covariance group (with the chicov alias), so any set of bin-by-bin correlated spectra — overlapping extraction regions, grating orders, or the full Stokes I/Q/U triple — can be fit jointly with a multivariate-Gaussian likelihood, per-bin Cholesky solve and graceful shrinkage of non-positive-definite covariances. The covariances are carried in the spectrum files (COVFRMT='XCOV-SET-1'; legacy Q–U pair files are auto-promoted), produced by the ftextractor tool (stokescov=iqu) and summed under grouping by ftgrouppha. See Appendix B.

New models

  • Four X-ray spectropolarimetric reflection models based on precomputed STOKES tables, designed for joint fitting of the Stokes $I$/$Q$/$U$ spectra: stdisc (reflection of a power-law source of arbitrary incident polarization from the distant, nearly neutral equatorial regions of an accretion disc), and sttorus, stcone and stbowl (reflection from an optically thick, axially symmetric torus, double cone, or double bowl – representing e.g. an opaque AGN torus, a broad-line region, or a super-Eddington accretion funnel – with partially ionized, neutral, or Chandrasekhar electron-scattering surfaces).

Modifications to models

  • The broadening of the APEC/SPEX pseudo-continuum (set by the BROADPSEUDO xset option) now uses an FFT and accepts fast and slow options. The fast option (equivalent to on) uses an abundance-weighted atomic mass so that a single FFT broadens all elements at once, about five times faster than the per-element slow calculation.

  • The gsmooth and lsmooth convolution models use an FFT for the index 0 and index 1 cases, giving order-of-magnitude speed-ups, and cache the kernel FFT when it has not changed.

  • Added a check that warns (yellow alert) when the fiducial energy at which a line width is specified falls outside the input energy range.

  • Table models (atable, mtable, etable) may now carry a neural network (NN) emulation of the tabulated spectra, indicated by an NNPTFILE keyword in the table file pointing to a PyTorch (.pt) weights file. When such a table is loaded XSPEC automatically evaluates the model with the neural network instead of interpolating on the parameter grid, with no change to how the model is set up. This requires XSPEC to have been built with libtorch (the –with-libtorch configure option).

  • The crossarf mixing model can now be stacked, one copy per source, so that a single sky region may host several emission components with different spatial distributions – for example an extended source and a point source, which spill into neighbouring regions by different fractions and so require separate arfs. Each source's cross arfs are supplied by source-tagged XFLT keywords of the form “CrossArfSrcSFromNToM” (source S, contribution to region M from region N); the original “CrossArfFromNToM” form is unchanged and applies to source 1.

  • The recorn mixing model is once again applied separately within each data group, so the correction norms of several spectra can be fitted at the same time. Its first parameter names a single spectrum, so unlike the other mixing models – whose one transformation spans every data group – recorn needs one copy per corrected spectrum. Each data group's copy now has its own spectrum number and cornorm, all of which appear in show par and can be set and fitted; in the intervening releases only the lowest data group's copy did anything.

Other

  • Spectral model data files are no longer shipped with HEAsoft. The new ftgetmodeldata tool downloads, on demand, the data files required by the models in use, over the network (via cURL), with version lookup (including a latest option), a manifest file (modelDataFiles.csv), and clobber and dry-run options. See
    https://heasarc.gsfc.nasa.gov/docs/software/xspec/modeldata.html.

  • New group command for runtime (re)grouping of spectra without re-reading the files. Its most important use is group back, which lays a coarser “super-bin” grouping on top of the current source binning — each super-bin spanning a whole number of source bins — to reduce a bias in the C statistic when a background is present (the W statistic, see Appendix B). That statistic profiles out one background-rate nuisance parameter per source bin; when the background has few counts per bin these estimates are noisy and bias the fitted source parameters — an incidental-parameters problem that can reach tens of percent for faint sources on a low-count background. Sharing a single background rate across the source bins of each super-bin shrinks the number of nuisance parameters and can bring the bias down to a few percent. The pooled-background path is selected automatically whenever a background super-bin grouping has been set up at runtime with group back; a GROUPING column carried by the background file (e.g. from ftgrouppha) is not used for this, since on input XSPEC applies the source grouping to the background. The command also regroups the source spectrum itself (group data with optbin, minsn, mincounts, const, or file); all grouping is persisted in save files, recoverable with tclout group, and available in PyXspec.

  • New sim command for bulk posterior-predictive simulation: it draws many parameter sets (from the loaded chains, or from the fit covariance scaled by sim fsigma), realizes each as a simulated dataset in the loaded spectra, optionally refits it (a parametric bootstrap), and accumulates the results — the in-program replacement for the scripted tclout simpars + newpar + fakeit loop (the multifake script). Because the realization happens in place no files are written and the original state is restored on completion. An optional command file (sim @cmdfile) is run on every realization and may return a row of numbers or append to a per-realization output file; the products are recovered with tclout sim. The whole run is reproducible under xset seed and parallelizable with parallel sim.

  • The plot command allows finer control over plot stacks by using the pipe (|) character to delimit them, e.g. plot da res | mo | ld rat draws three stacks. At present up to four stacks may be drawn, each holding up to four panes.

  • New plot options dspec, despec and deespec: the model-independent counterparts of ufspec, eufspec and eeufspec. The flux density is obtained by a regularized truncated-SVD inversion of the detector response rather than by scaling the data by the model-dependent unfolded/folded ratio, so the plotted data points do not move when the model is changed.

  • Operational switches of the form xset CONTROL_STRING yes control the numerical paths used in fitting and model evaluation and enable cross-check validation modes. For example, DISABLE_ANALYTIC_GRAD reverts fits to finite-difference derivatives. They are listed in a new table in the xset entry, separate from the model control strings.

  • chain info and margin accept a parameter number of 0 to report and plot the distribution of the statistic values stored in a chain.

  • The identify command gained a spex option to search the SPEX line list.

  • Python local models are now supported in the standalone (Tcl) version of XSPEC, not only in PyXspec. A Python model function may also opt into a faster vectorized interface with the @xspec.vectorized decorator, receiving zero-copy numpy arrays instead of Python tuples and lists — avoiding the per-element Python-object overhead (Appendix C).

  • Rank-$K$ SVD compression of response matrices is now supported: a compressed side-car file (produced by the ftsvdcmprmf tool) can be loaded with the response command and tuned with the svdrank and svdeps options for fast folding (Appendix H.2).

  • The XSPEC libraries were reorganized so they can be linked into third-party C++ programs without Tcl or PGPLOT: core start-up was split into a new XSCore library, and the model and data classes were split into separate XSModel and XSData libraries. The PLT plotting package was modernized to call PLT and PGPLOT directly from C++.

Fixes

Many user-visible bugs present in earlier versions have been fixed. The more significant fixes are grouped by area below; a complete list of the user-facing fixes in this release is given in Appendix I.

Models

  • The btapec model overwrote the helium abundance with the metal abundance, producing a wrong continuum whenever the abundance parameter was not 1.

  • The zvgabs absorption-line width was applied in km/s instead of the equivalent keV value.

  • The velocity-broadened Lorentzian and Voigt models (the vlorentz and vvoigt families) shipped with a width-parameter hard maximum of 20 km/s, below their own default of 100 km/s, so a fit could never raise the width above 20 km/s. The maxima now match the vgauss family.

  • The xmmpsf and suzpsf PSF mixing models set the model flux to zero above 6 keV.

  • Out-of-bounds memory writes in the diskline model (when a line's blue edge fell in the first energy bin) and the pileup model (when max_ph$\,=1$) could corrupt memory or crash XSPEC.

  • The vvnpshock and vvsedov models (and their bvvnpshock/bvvsedov variants) reused stale ionization fractions when only the electron temperature was changed, giving a wrong spectrum.

Data, responses and setup

  • Defining a model with more than one data group loaded read out of bounds, aborting XSPEC on builds compiled with array-bounds assertions.

  • An ignore or notice command whose channel range lay entirely above the spectrum's channel count could crash.

  • Deleting several spectra in a single data command could remove the wrong spectra.

  • A table-model file path containing a slash or a hyphen was corrupted when the component was added with addcomp.

  • Channels re-noticed after a dummyrsp command were incorrectly folded to zero.

  • Variable-length QUALITY and GROUPING columns were incorrectly assumed to be scalar when their repeat count was 1.

Fitting and statistics

  • A segmentation fault in the cstat statistic when all channels are ignored.

  • The default Levenberg-Marquardt (leven) method could leave parameters stuck pegged at a hard limit.

  • A parameter resting on one of its limits was given a numerical derivative half its true value, so the Levenberg-Marquardt step for it was twice as long as it should be. Additive norms met this routinely. On the paths that differentiate the statistic directly (xset usechainrule no, and the recorn cornorm) the same parameter could hang fit outright.

  • fit could abort the whole session, or read out of bounds with a Minuit method, when the covariance matrix was ill-conditioned or not positive-definite; it now issues a recoverable error.

  • The empirical distribution function test statistics could crash or return NaN: with every channel ignored, with a model integrating to zero over the noticed range, or with a spectrum that was folded but not fitted.

  • The recorn cornorm parameter could not be fitted: because it rescales the data rather than the model, the chain-rule derivative was identically zero, so it was pegged as insensitive and fit stopped with “No variable parameters for fit”.

  • untie could free a mixing-model parameter outside the lowest data group, producing a parameter that could not be displayed, set or fitted but still consumed a degree of freedom.

Plotting

  • plot sensitivity and insensitivity drew the first detector's curve for every response of a multiple-response spectrum.

  • plot delchi, ufspec, fitstat and icounts could crash for a plot group whose spectrum has no response for an active model source.

  • The polarization plots polfrac and polangle returned NaN or zero for multiple-source models and could contain spurious or infinite bins.

Commands and scripting

  • autosave 0 crashed XSPEC with a divide-by-zero on the next command.

  • Running a script by a path containing a directory component (for example @../model) failed to append the default .xcm suffix.

  • energies with the log option and a low energy of 0 silently built a degenerate energy grid for every model.