The 6-s X-ray pulsar 1E1048.1-5937 was observed by *EXOSAT* in June
1985 for 20 ks. In this example, we'll conduct a general investigation
of the spectrum from the Medium Energy (ME) instrument, i.e. follow the
same sort of steps as the original investigators (Seward, Charles &
Smale,1986). The ME spectrum and corresponding response matrix were
obtained from the HEASARC On-line service.

Once installed, XSPEC is invoked by typing `xspec`, as in this
example.

%xspec Xspec 11.1.0 21:50:14 08-Jul-2001 For documentation, notes, and fixes see http://xspec.gsfc.nasa.gov/ Plot device not set, use "cpd" to set it Type "help" or "?" for further information XSPEC>data s54405.pha Net count rate (cts/s) for file 1 3.783 +/- 0.1367 using response (RMF) file... s54405.rsp 1 data set is in use

The `data` command tells the program to read the data as well as the
response file that is named in the header of the data file. In general,
XSPEC commands can be truncated, provided they remain unambiguous.
Since the default extension of a data file is `.pha`,
and the abbreviation of `data` to the first two letters is unambiguous, the above
command can be abbreviated to `da s54405`, if desired. To obtain
help on the `data` command, or on any other command, type `help command`
followed by the name of the command.

One of the first things most users will want to do at this stage--even
before fitting models--is to look at their data. The plotting device
should be set first, with the command `cpd` (*change plotting
device*). Here, we use the pgplot X-Window server, /xw.

XSPEC> cpd /xw

To see a list of alternative plot devices, type `cpd ?` There are
more than 20 different things that can be plotted, all related in
some way to the data, the model, the fit and the instrument. To see
them, type:

XSPEC> plot ? Choose from the following `plot' sub-commands: data counts ldata residuals chisq delchi ratio summary model emodel eemodel contour efficien ufspec eufspec eeufspec dem insensitv sensitvty genetic foldmodel icounts Insert selection: (data)

The most fundamental is the data plotted against instrument channel
(`data`); next most fundamental, and more informative, is the data
plotted against channel energy. To do this plot, use the XSPEC command
`setplot energy`. Figure A shows the result of the commands:

XSPEC> setplot energy XSPEC> plot data

People familiar with *EXOSAT* ME data will recognize the spectrum to
be soft, absorbed and without an obvious bright iron emission line.

We are now ready to fit the data with a model. Models in XSPEC are
specified using the `model` command, followed by an algebraic
expression of a combination of model components. There
are two basic kinds of model components: *additive* which represent
X-Ray sources of different kinds. After being
convolved with the instrument response, the components prescribe the number of counts
per energy bin (e.g., a bremsstrahlung continuum); and *multiplicative* models components, which represent phenomena that
modify the observed X-Radiation (e.g. reddening or an absorption edge).
They apply an energy-dependent multiplicative
factor to the source radiation before the result is convolved with the instrumental
response.
More generally, XSPEC allows three types of modifying components: convolutions
and mixing models in addition to the multiplicative type.
Since there must be a source, there must be
least one additive component in a model, but there
is no restriction on the number of modifying components. To see what
components are available, type `model ?`:

XSPEC>model ? Possible additive models are : apec bbody bbodyrad bexrav bexriv bknpower bmc bremss c6mekl c6pmekl c6pvmkl c6vmekl cemekl cevmkl cflow compbb compLS compST compTT cutoffpl disk diskbb diskline diskm disko diskpn equil gaussian gnei grad grbm laor lorentz meka mekal mkcflow nei npshock nteea pegpwrlw pexrav pexriv plcabs powerlaw posm pshock raymond redge refsch sedov srcut sresc step vapec vbremss vequil vgnei vmeka vmekal vmcflow vnei vnpshock vpshock vraymond vsedov zbbody zbremss zgauss zpowerlw atable Possible multiplicative models are : absori constant cabs cyclabs dust edge expabs expfac highecut hrefl notch pcfabs phabs plabs redden smedge spline SSS_ice TBabs TBgrain TBvarabs uvred varabs vphabs wabs wndabs xion zedge zhighect zpcfabs zphabs zTBabs zvarabs zvfeabs zvphabs zwabs zwndabs mtable etable Possible mixing models are : ascac projct Possible convolution models are : gsmooth lsmooth reflect Possible pile-up models are : pileup

For information about a specific component, type `help Models` followed by the name of the
component):

XSPEC>help Models raymond XSPEC_11.1_commands Models raymond An emission spectrum from hot diffuse gas based on the model calculations of Raymond and Smith including line emission from several elements. This model interpolates on a grid of spectra for different temperatures. The grid is logarithmically spaced with 80 temperatures ranging from 0.008 to 80 keV. par1 = plasma temperature in keV par2 = Metal abundances (He fixed at cosmic) The elements included are C, N, O, Ne, Mg, Si, S, Ar, Ca, Fe, Ni and their relative abundances are set by the abund command. par3 = redshift, z K = 10**-14 / (4 pi (D_L/(1+z))**2) Int n_e n_H dV, where D_L is the luminosity distance to the source (cm), n_e is the electron density (cm**-3), and n_H is the hydrogen density (cm**-3) help>

Given the quality of our data, as shown by the plot, we'll choose an absorbed power law, specified as follows :

XSPEC> model phabs(powerlaw)

Or, abbreviating unambiguously:

XSPEC> mo pha(po)

The user is then prompted for the initial values of
the parameters. Entering <`return`> or `/` in response to a prompt
uses the default values.
We could also have set all parameters to their default values by entering `/*` at
the first prompt.
As well as the parameter values themselves,
users also may
specify step sizes and ranges (`value`, `delta`, `min`, `bot`,
`top`, and `max values`), but here, we'll enter the defaults:

XSPEC>mo pha(po) Model: phabs[1]( powerlaw[2] ) Input parameter value, delta, min, bot, top, and max values for ... Current: 1 0.001 0 0 1E+05 1E+06 phabs:nH>/* --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 1.000 +/- 0. 2 2 2 powerlaw PhoIndex 1.000 +/- 0. 3 3 2 powerlaw norm 1.000 +/- 0. --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 4.8645994E+08 using 125 PHA bins. Reduced chi-squared = 3987376. for 122 degrees of freedom Null hypothesis probability = 0.

Note that most of the other numerical values in this section
are have changed from those produced by earlier versions.
This is because the default photoelectric absorption
cross-sections have changed since XSPEC V.10. To recover identical
results to earlier versions, use
`xsect obcm`. `bcmc` is the default; see `help xsect` for details).

The `Current statistic` is
and is huge for the initial,
default values--mostly because the power law normalization is two
orders of magnitude too large. This is easily fixed using the
`renorm` command :

XSPEC> renorm Chi-Squared = 852.1660 using 125 PHA bins. Reduced chi-squared = 6.984967 for 122 degrees of freedom Null hypothesis probability = 0.

We are not quite ready to fit the data (and obtain a better
), because not all of the 125 PHA
bins should be included in the fitting: some are below the lower
discriminator of the instrument and therefore do not contain valid data;
some have imperfect background
subtraction at the margins of the pass band; and some may not contain
enough counts for
to be strictly meaningful. To find out
which channels to discard (*ignore* in `XSPEC` terminology),
users must consult mission-specific documentation, which will inform them
of discriminator settings, background subtraction problems and other
issues. For the mature missions in the HEASARC archives, this
information already has been encoded in the headers of the spectral
files as a list of ``bad'' channels. Simply issue the command:

XSPEC> ignore bad Chi-Squared = 799.7109 using 85 PHA bins. Reduced chi-squared = 9.752572 for 82 degrees of freedom Null hypothesis probability = 0. XSPEC> setplot chan XSPEC> plot data

We can see that 40 channels are bad--but do we need to ignore any more?
These channels are bad because of certain instrument properties: other
channels still may need to be ignored because of the shape and
brightness of the spectrum itself. Figure B shows the result of plotting
the data in channels (using the commands `setplot channel` and
`plot data`). We see that above about channel 33 the S/N becomes
small. We also see, comparing Figure B with Figure A, which bad channels were
ignored. Although
visual inspection is not the most rigorous method for deciding which
channels to ignore (more on this subject later), it's good enough for
now, and will at least prevent us from getting grossly misleading
results from the fitting. To ignore channels above 33:

XSPEC> ignore 34-** Chi-Squared = 677.6218 using 31 PHA bins. Reduced chi-squared = 24.20078 for 28 degrees of freedom Null hypothesis probability = 0.

The same result can be achieved with the command `ignore 34-125`.
The inverse of `ignore` is `notice`, which has the same syntax.

We are now ready to fit the data. Fitting is initiated by the command
`fit`. As the fit proceeds, the screen displays the status of the
fit for each iteration until either the fit converges to the minimum
,
or the user is asked whether the fit is to go through
another set of iterations to find the minimum. The default number of
iterations is ten.

XSPEC>fit Chi-Squared Lvl Fit param # 1 2 3 204.136 -3 7.9869E-02 1.564 4.4539E-03 84.5658 -4 0.3331 2.234 1.0977E-02 30.2511 -5 0.4422 2.174 1.1965E-02 30.1202 -6 0.4648 2.196 1.2264E-02 30.1189 -7 0.4624 2.195 1.2244E-02 --------------------------------------------------------------------------- Variances and Principal axes : 1 2 3 4.14E-08 | 0.00 -0.01 1.00 8.70E-02 | -0.91 -0.41 -0.01 2.32E-03 | -0.41 0.91 0.01 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 0.4624 +/- 0.2698 2 2 2 powerlaw PhoIndex 2.195 +/- 0.1288 3 3 2 powerlaw norm 1.2244E-02 +/- 0.2415E-02 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 30.11890 using 31 PHA bins. Reduced chi-squared = 1.075675 for 28 degrees of freedom Null hypothesis probability = 0.358

The fit is good: reduced is 1.075 for 31 - 3 = 28 degrees of freedom. The null hypothesis probability is the probability of getting a value of as large or larger than observed if the model is correct. If this probability is small then the model is not a good fit. The matrix of principal axes printed out at the end of a fit provides an indication of whether parameters are correlated (at least local to the best fit). In this example the powerlaw norm is not correlated with any other parameter while the column and powerlaw index are slightly correlated.

To see the fit and the residuals, we use the command

XSPEC>plot data resid

The result is shown in Figure C.

The screen output shows the best-fitting parameter values, as well as
approximations to their errors. These errors should be regarded as
indications of the uncertainties in the parameters and should *not*
be quoted in publications. The true errors, i.e. the confidence ranges,
are obtained using the `error` command:

XSPEC>error 1 2 3 Parameter Confidence Range ( 2.706) 1 3.254765E-02 0.932778 2 1.99397 2.40950 *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization *WARNING*:RENORM: No variable model to allow renormalization 3 8.942987E-03 1.711637E-02

Here, the numbers 1, 2, 3 refer to the parameter numbers in the `Model
par` column of the screen output. For the first parameter, the column of
absorbing hydrogen atoms, the 90%
confidence range is
cm^{-2}. This corresponds to an
excursion in
of 2.706. The reason these ``better'' errors are
not given automatically as part of the `fit` output is that they
entail further fitting. When the model is simple, this does not require
much CPU, but for complicated models the extra time can be
considerable. The warning message is generated because there are
no free normalizations in the model while the error is being
calculated on the normalization itself. In this case, the warning may
safely be ignored.

What else can we do with the fit? One thing is to derive the flux of the
model. The data by themselves only give the instrument-dependent count
rate. The model, on the other hand, is an estimate of the true spectrum
emitted. In XSPEC, the model is defined in physical units
independent of the instrument. The command `flux` integrates the
current model over the range specified by the user:

XSPEC> flux 2 10 Model flux 3.5496E-03 photons ( 2.2492E-11 ergs)cm**-2 s**-1 ( 2.000- 10.000)

Here, we have chosen the standard X-ray range of 2-10 keV and find that
the energy flux is
erg cm^{-2} s^{-1}. Note
that `flux` will integrate only within the energy range of the
current response matrix. If the model flux outside this range is
desired--in effect, an extrapolation beyond the data--then the command
`dummyrsp` should be used. This command sets up a dummy response
that covers the range required. For example, if we want to know the flux
of our model in the *ROSAT* PSPC band of 0.2-2 keV, we enter:

XSPEC>dummy 0.2 2. Chi-Squared = 3583.779 using 31 PHA bins. Reduced chi-squared = 127.9921 for 28 degrees of freedom Null hypothesis probability = 0. XSPEC>flux 0.2 2. Model flux 4.5306E-03 photons ( 9.1030E-12 ergs)cm**-2 s**-1 ( 0.200- 2.000)

The energy flux, at
erg cm^{-2} s^{-1}, is
lower in this band but the photon flux is higher. To get our original
response matrix back we enter :

XSPEC> response Chi-Squared = 30.11890 using 31 PHA bins. Reduced chi-squared = 1.075675 for 28 degrees of freedom Null hypothesis probability = 0.358

The fit, as we've remarked, is good, and the parameters are constrained. But unless the purpose of our investigation is merely to measure a photon index, it's a good idea to check whether alternative models can fit the data just as well. We also should derive upper limits to components such as iron emission lines and additional continua, which, although not evident in the data nor required for a good fit, are nevertheless important to constrain. First, let's try an absorbed black body:

XSPEC>mo pha(bb) Model: phabs[1]( bbody[2] ) Input parameter value, delta, min, bot, top, and max values for ... Current: 1 0.001 0 0 1E+05 1E+06 phabs:nH>/* --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( bbody[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 1.000 +/- 0. 2 2 2 bbody kT keV 3.000 +/- 0. 3 3 2 bbody norm 1.000 +/- 0. --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 3.3142067E+09 using 31 PHA bins. Reduced chi-squared = 1.1836453E+08 for 28 degrees of freedom Null hypothesis probability = 0. XSPEC>fit Chi-Squared Lvl Fit param # 1 2 3 1420.96 0 0.2116 2.987 7.7666E-04 1387.72 0 4.4419E-02 2.975 7.7003E-04 1376.39 0 4.3009E-03 2.963 7.6354E-04 1371.67 0 1.8192E-03 2.951 7.5730E-04 1367.04 0 5.8100E-04 2.939 7.5130E-04 1362.47 0 1.9059E-04 2.926 7.4548E-04 1357.84 0 6.1455E-05 2.913 7.3982E-04 1353.14 0 2.5356E-05 2.900 7.3429E-04 1348.36 0 6.7582E-06 2.887 7.2885E-04 1343.50 0 2.0479E-06 2.874 7.2350E-04 Number of trials exceeded - last iteration delta = 4.861 Continue fitting? (Y)y ... 113.954 0 0. 0.8907 2.7865E-04 113.954 -1 0. 0.8905 2.7859E-04 113.954 4 0. 0.8905 2.7859E-04 --------------------------------------------------------------------------- Variances and Principal axes : 2 3 2.88E-04 | -1.00 0.00 8.45E-11 | 0.00 -1.00 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( bbody[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 0. +/- -1.000 2 2 2 bbody kT keV 0.8905 +/- 0.1696E-01 3 3 2 bbody norm 2.7859E-04 +/- 0.9268E-05 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 113.9542 using 31 PHA bins. Reduced chi-squared = 4.069792 for 28 degrees of freedom Null hypothesis probability = 2.481E-12

Note that when more than 10 iterations are required for convergence
the user is asked whether or not to continue at the end of each set
of 10. Saying *no* at these prompts is a good idea if the fit is not
converging quickly. Conversely, to avoid having to keep answering the
question, i.e., to increase the number of iterations before the prompting
question appears, begin the fit with, say `fit 100`. This command will put
the fit through 100 iterations before pausing.

Plotting the data and residuals again with

XSPEC> plot data resid

we obtain Figure D.

The black body fit is obviously not a good one. Not only is
large, but the best-fitting *N*_{H} is rather low. Inspection of the
residuals confirms this: the pronounced wave-like shape is indicative of
a bad choice of overall continuum (see Figure D). Let's try thermal
bremsstrahlung next:

XSPEC>mo pha(br) Model: phabs[1]( bremss[2] ) Input parameter value, delta, min, bot, top, and max values for ... Current: 1 0.001 0 0 1E+05 1E+06 phabs:nH>/* --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( bremss[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 1.000 +/- 0. 2 2 2 bremss kT keV 7.000 +/- 0. 3 3 2 bremss norm 1.000 +/- 0. --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 4.5311800E+07 using 31 PHA bins. Reduced chi-squared = 1618279. for 28 degrees of freedom Null hypothesis probability = 0. XSPEC>fit Chi-Squared Lvl Fit param # 1 2 3 113.305 -3 0.2441 6.557 6.8962E-03 40.4519 -4 0.1173 5.816 7.7944E-03 36.0549 -5 4.4750E-02 5.880 7.7201E-03 33.4168 -6 1.8882E-02 5.868 7.7476E-03 32.6766 -7 7.8376E-03 5.864 7.7495E-03 32.3192 -8 2.7059E-03 5.862 7.7515E-03 32.1512 -9 2.3222E-04 5.861 7.7525E-03 32.1471 -10 1.0881E-04 5.861 7.7523E-03 --------------------------------------------------------------------------- Variances and Principal axes : 1 2 3 2.29E-08 | 0.00 0.00 1.00 3.18E-02 | 0.95 0.31 0.00 8.25E-01 | 0.31 -0.95 0.00 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( bremss[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 1.0881E-04 +/- 0.3290 2 2 2 bremss kT keV 5.861 +/- 0.8651 3 3 2 bremss norm 7.7523E-03 +/- 0.8122E-03 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 32.14705 using 31 PHA bins. Reduced chi-squared = 1.148109 for 28 degrees of freedom Null hypothesis probability = 0.269

Bremsstrahlung is a better fit than the black body--and is as good as
the power law--although it shares the low *N*_{H}. With two good fits,
the power law and the bremsstrahlung, it's time to scrutinize their
parameters in more detail.

First, we reset our fit to the powerlaw (output omitted):

XSPEC>mo pha(po)

From the *EXOSAT* database on HEASARC, we know that the target in
question, 1E1048.1-5937, has a Galactic latitude of 24 arcmin, i.e.,
almost on the plane of the Galaxy. In fact, the data base also provides
the value of the Galactic *N*_{H} based on 21-cm radio observations. At
cm^{-2}, it is higher than the 90 percent-confidence upper
limit from the power-law fit. Perhaps, then, the power-law fit is not so
good after all. What we can do is fix (*freeze* in XSPEC
terminology) the value of *N*_{H} at the Galactic value and refit the
power law. Although we won't get a good fit, the shape of the residuals
might give us a clue to what is missing. To freeze a parameter in
XSPEC, use the command `freeze` followed by the parameter number,
like this:

XSPEC> freeze 1 Number of variable fit parameters = 2

The inverse of `freeze` is `thaw`:

XSPEC> thaw 1 Number of variable fit parameters = 3

Alternatively, parameters can be frozen using the `newpar` command,
which allows all the quantities associated with a parameter to be
changed. The second quantity, `delta`, is the step size used to
calculate the derivative in the fitting, and, if set to a negative
number, will cause the parameter to
be frozen. In our case, we want *N*_{H} frozen at
cm^{-2}, so we go back to the power law best fit
and do the following :

XSPEC>newpar 1 Current: 0.463 0.001 0 0 1E+05 1E+06 phabs:nH>4,0 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 4.000 frozen 2 2 2 powerlaw PhoIndex 2.195 +/- 0.1287 3 3 2 powerlaw norm 1.2247E-02 +/- 0.2412E-02 --------------------------------------------------------------------------- --------------------------------------------------------------------------- 2 variable fit parameters Chi-Squared = 829.3545 using 31 PHA bins. Reduced chi-squared = 28.59843 for 29 degrees of freedom Null hypothesis probability = 0.

Note the useful trick of giving a value of zero for `delta` in the
`newpar` command. This has the effect of changing `delta` to the
negative of its current value. If the parameter is free, it will be
frozen, and if frozen, thawed. The same result can be obtained by
putting everything onto the command line, i.e., `newpar 1 4, 0`, or
by issuing the two commands, `newpar 1 4` followed by `freeze
1`. Now, if we fit and plot again, we get the following model (Fig. E).

XSPEC>fit ... --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 4.000 frozen 2 2 2 powerlaw PhoIndex 3.594 +/- 0.6867E-01 3 3 2 powerlaw norm 0.1161 +/- 0.9412E-02 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 125.5134 using 31 PHA bins. Reduced chi-squared = 4.328048 for 29 degrees of freedom Null hypothesis probability = 5.662E-14 XSPEC>plot data resid

The fit is not good. In Figure E we can see why: there appears to be a
surplus of softer photons, perhaps indicating a second continuum
component. To investigate this possibility we can add a component to our
model. The `editmod` command lets us do this without having to
respecify the model from scratch. Here, we'll add a black body
component.

XSPEC>editmod pha(po+bb) Model: phabs[1]( powerlaw[2] + bbody[3] ) Input parameter value, delta, min, bot, top, and max values for ... Current: 3 0.01 0.0001 0.01 100 200 bbody:kT>2,0 Current: 1 0.01 0 0 1E+24 1E+24 bbody:norm>1.e-5 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] + bbody[3] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 4.000 frozen 2 2 2 powerlaw PhoIndex 3.594 +/- 0.6867E-01 3 3 2 powerlaw norm 0.1161 +/- 0.9412E-02 4 4 3 bbody kT keV 2.000 frozen 5 5 3 bbody norm 1.0000E-05 +/- 0. --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 122.1538 using 31 PHA bins. Reduced chi-squared = 4.362635 for 28 degrees of freedom Null hypothesis probability = 9.963E-14

Notice that in specifying the initial values of the black body, we have frozen kT at 2 keV (the canonical temperature for nuclear burning on the surface of a neutron star in a low-mass X-ray binary) and started the normalization at zero. Without these measures, the fit might ``lose its way''. Now, if we fit, we get (not showing all the iterations this time):

--------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] + bbody[3] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 4.000 frozen 2 2 2 powerlaw PhoIndex 4.932 +/- 0.1618 3 3 2 powerlaw norm 0.3761 +/- 0.5449E-01 4 4 3 bbody kT keV 2.000 frozen 5 5 3 bbody norm 2.3212E-04 +/- 0.3966E-04 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 55.63374 using 31 PHA bins. Reduced chi-squared = 1.986919 for 28 degrees of freedom Null hypothesis probability = 1.425E-03

The fit is better than the one with just a power law *and* the fixed
Galactic column, but it is still not good. Thawing the black body
temperature and fitting gives us:

XSPEC>thaw 4 Number of variable fit parameters = 4 XSPEC>fit ... --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] + bbody[3] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 4.000 frozen 2 2 2 powerlaw PhoIndex 6.401 +/- 0.3873 3 3 2 powerlaw norm 1.086 +/- 0.3032 4 4 3 bbody kT keV 1.199 +/- 0.8082E-01 5 5 3 bbody norm 2.6530E-04 +/- 0.3371E-04 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 37.21207 using 31 PHA bins. Reduced chi-squared = 1.378225 for 27 degrees of freedom Null hypothesis probability = 9.118E-02

This, of course, is a better fit, but the photon index of the power law has ended up extremely and implausibly steep. Looking at this odd model with the command

XSPEC> plot model

we see, in Figure F, that the black body and the power law have changed places, in that the power law provides the soft photons required by the high absorption, while the black body provides the harder photons.

We could continue to search for a plausible, well-fitting model, but the
data, with their limited signal-to-noise and energy resolution, probably
don't warrant it (the original investigators published only the power
law fit). There is, however, one final, useful thing to do with the
data: derive an upper limit to the presence of a fluorescent iron
emission line. First we delete the black body component using `delcomp`:

XSPEC>delcomp 3 Model: phabs[1]( powerlaw[2] ) Chi-Squared = 1285.487 using 31 PHA bins. Reduced chi-squared = 44.32712 for 29 degrees of freedom Null hypothesis probability = 0.

Then we thaw *N*_{H} and refit to recover our original, best fit:

XSPEC>thaw 1 Number of variable fit parameters = 3 XSPEC>fit Chi-Squared Lvl Fit param # 1 2 3 924.178 -2 5.087 5.076 0.4056 305.507 -2 4.525 3.791 0.1249 140.460 -2 2.930 3.367 6.5553E-02 64.4275 -3 0.6068 2.244 1.4635E-02 30.3738 -4 0.4837 2.201 1.2279E-02 30.1189 -5 0.4641 2.195 1.2258E-02 30.1189 -6 0.4637 2.195 1.2255E-02 --------------------------------------------------------------------------- Variances and Principal axes : 1 2 3 4.13E-08 | 0.00 -0.01 1.00 8.69E-02 | -0.91 -0.41 -0.01 2.31E-03 | -0.41 0.91 0.01 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 0.4637 +/- 0.2696 2 2 2 powerlaw PhoIndex 2.195 +/- 0.1287 3 3 2 powerlaw norm 1.2255E-02 +/- 0.2412E-02 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 30.11890 using 31 PHA bins. Reduced chi-squared = 1.075675 for 28 degrees of freedom Null hypothesis probability = 0.358

Now, we add a gaussian emission line of fixed energy and width:

XSPEC>editmod pha(po+ga) Model: phabs[1]( powerlaw[2] + gaussian[3] ) Input parameter value, delta, min, bot, top, and max values for ... Current: 6.5 0.05 0 0 1E+06 1E+06 gaussian:LineE>6.4 0 Current: 0.1 0.05 0 0 10 20 gaussian:Sigma>0.1 0 Current: 1 0.01 0 0 1E+24 1E+24 gaussian:norm>1.e-4 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] + gaussian[3] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 0.4637 +/- 0.2696 2 2 2 powerlaw PhoIndex 2.195 +/- 0.1287 3 3 2 powerlaw norm 1.2255E-02 +/- 0.2412E-02 4 4 3 gaussian LineE keV 6.400 frozen 5 5 3 gaussian Sigma keV 0.1000 frozen 6 6 3 gaussian norm 1.0000E-04 +/- 0. --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 32.75002 using 31 PHA bins. Reduced chi-squared = 1.212964 for 27 degrees of freedom Null hypothesis probability = 0.205 XSPEC>fit ... --------------------------------------------------------------------------- --------------------------------------------------------------------------- Model: phabs[1]( powerlaw[2] + gaussian[3] ) Model Fit Model Component Parameter Unit Value par par comp 1 1 1 phabs nH 10^22 0.6562 +/- 0.3193 2 2 2 powerlaw PhoIndex 2.324 +/- 0.1700 3 3 2 powerlaw norm 1.4636E-02 +/- 0.3642E-02 4 4 3 gaussian LineE keV 6.400 frozen 5 5 3 gaussian Sigma keV 0.1000 frozen 6 6 3 gaussian norm 9.6462E-05 +/- 0.9542E-04 --------------------------------------------------------------------------- --------------------------------------------------------------------------- Chi-Squared = 29.18509 using 31 PHA bins. Reduced chi-squared = 1.080929 for 27 degrees of freedom Null hypothesis probability = 0.352

The energy and width have to be frozen because,in the absence of an
obvious line in the data, the fit would be completely unable to converge
on meaningful values. Besides, our aim is to see how bright a line at
6.4 keV can be and still not ruin the fit. To do this, we fit first and
then use the `error` command to derive the maximum allowable iron
line normalization. We then set the normalization at this maximum value
with `newpar` and, finally, derive the equivalent width using the
`eqwidth` command. That is:

XSPEC>err 6 Parameter Confidence Range ( 2.706) Parameter pegged at hard limit 0. with delta ftstat= 0.9338 6 0. 1.530722E-04 XSPEC>new 6 1.530722E-04 4 variable fit parameters Chi-Squared = 34.91923 using 31 PHA bins. Reduced chi-squared = 1.293305 for 27 degrees of freedom Null hypothesis probability = 0.141 XSPEC>eqwidth 3 Additive group equiv width for model 3 (gaussian): 839. eV

Things to note:

- The true minimum value of the gaussian normalization is less than
zero, but the
`error`command stopped searching for a of 2.706 when the minimum value hit zero, the ``hard" lower limit of the parameter. Hard limits can be adjusted with the`newpar`command, and they correspond to the quantities`min`and`max`associated with the parameter values. In fact, according to the screen output, the value of corresponding to zero normalization is 0.934. - The command
`eqwidth`takes the component number as its argument. - The upper limit on the equivalent width of a 6.4 keV emission line is high (839 eV)!