Come analyze HEASARC, IRSA, and MAST data in the cloud! The Fornax Initiative is now welcoming all interested beta users.
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,
). 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
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
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
-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
/
/
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-
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
) 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.