Quick Tutorial
This assumes the user already has a basic familiarity with both XSPEC and Python. Everything in PyXspec is accessible by importing the package xspec into your Python script.
PyXspec can be utilized in a Python script or from the command line of the plain interactive Python interpreter. PyXspec does not implement its own command handler, so it is NOT intended to be run as the Python equivalent of a traditional interactive XSPEC session (which is really an enhanced interactive Tcl interpreter). In other words you launch an interactive PyXspec session with:
LINUX> python >>> import xspec >>>rather than:
LINUX> xspec XSPEC12>Note that in all the tutorial examples the xspec package name qualifier is left off. You must either include the xspec qualifier:
s = xspec.Spectrum("file1.pha")or use a variation of the Python
importorfrom...importcommands:from xspec import * s = Spectrum("file1.pha")
Jumping In: The Really Quick Tutorial
A simple Xspec load-fit-plot Python script may look something like this:
#!/usr/bin/python from xspec import * Spectrum("file1.pha") Model("phabs*pow") Fit.perform() Plot.device = "/xs" Plot("data")Keeping this template in mind, we'll proceed to fill in the details...
Terminology
This description uses the standard Python object-oriented terminology, distinguishing between classes and objects. Class is used when referring to the type or definition of an object. An object refers to a specific instance of a class and is normally assigned to a variable. For example a user may load 3 data files by creating 3 spectral data objects s1, s2, and s3, which are all instances of the class Spectrum.
The functions and stored data members that make up the definition of a class are referred to as methods and attributes respectively.
The term Standard XSPEC refers to the traditional ways of using XSPEC, either with a Tcl script or an interactive XSPEC session.
Getting Help
Documentation for PyXspec classes can be found in the Classes section of this manual. You can also access the same information from the Python shell command line using Python's built-in
help(<class-name>)function.
The 6 Global Objects
An XSPEC session fundamentally consists of loading data, fitting that data to a model, and plotting the results. To manage these operations, PyXspec offers the user 6 global objects: AllChains, AllData, AllModels, Fit, Plot, and Xset. Note that these are NOT the names of classes. They are instantiated objects of the class types shown below.
Note
Operations involving these 6 global ojbects should ALWAYS be performed through the object names and NEVER their class names. Their class names should never appear in your code.
Object Name
Class
Role
AllChains
Monte Carlo Markov Chain container
AllData
Container for all loaded data sets (objects of class Spectrum)
AllModels
Container for all Model objects
Fit
Manager class for setting properties and running a fit
Plot
Manager class for performing XSPEC plots
Xset
Storage class for XSPEC settings
PyXspec instantiates these objects immediately upon the importing of the xspec package. You cannot create any other objects of these class types, as they each allow only 1 instance of their type. (They are singletons in the language of design patterns.)
Loading And Removing Data
Spectral data files can be loaded in several ways. You can create an object of the Spectrum class by passing it the data file name:
s1 = Spectrum("file1.pha")which also adds the new object s1 to the AllData container. Or you can simply add the new file directly to the container without retrieving a Spectrum object:
AllData += "file1.pha"Later you can always obtain a Spectrum object reference to any of the loaded spectra by passing AllData an integer:
s2 = AllData(2) # s2 is a reference to the 2nd loaded spectrumFor more complicated data loading, you have access to the same functionality in Standard XSPEC's data command. Simply pass a string to the AllData object's __call__ method:
AllData("file1 file2 2:3 file3")Note that only the last example allows you to assign multiple data groups, the 3rd spectrum being assigned to data group 2. Also note that in the last example any previously loaded data sets are removed, thus reproducing the behavior of Standard XSPEC's data command.
Other ways of removing Spectrum objects (ie. data sets) from the container:
AllData -= 3 # Removes the 3rd Spectrum object (the spectrum with # index number 3) from the container. AllData -= s1 # Removes the Spectrum object s1. AllData -= "*" # Removes all Spectrum objects. AllData.clear() # Removes all Spectrum objects.You can check the current state of the AllData container at any time by doing:
AllData.show()Similarly, to view information about a single Spectrum object:
s2.show()
Defining Models
The basic way of defining an XSPEC model is to create an object of the PyXspec class Model. Simply pass in a string containing a combination of 1 or more XSPEC model components. Since this uses the same syntax as Standard XSPEC's model command, component abbreviations are allowed:
m1 = Model("phabs*po + ga")and to see a complete listing of available XSPEC model components, do:
Model.showList()When you define a model like this, PyXspec also automatically adds the new object to the global AllModels container. If the model is applied to multiple data groups, object copies are added to the container for each data group.
Similar to the case of spectral data, you can also load models directly into the global container:
# Another way to define a new model and create an object for each data group. AllModels += "phabs*po + ga" # Retrieve the model object assigned to data group 2. m2 = AllModels(2) # Various ways to remove all model objects from the container. AllModels.clear() AllModels -= "*"To display models and their parameters:
# This displays all parameters in all model objects: AllModels.show() # While this displays just parameters 1,2,3 and 5: AllModels.show("1-3, 5") # This displays a single model object: m2.show()For defining multiple (or named) models and assigning multiple sources, please see the Extended Tutorial section.
Component and Parameter Objects
Model objects contain Component objects and Component objects contain Parameter objects. There are several ways to access and set components and parameters individually (and if you want to change many parameter values at once, it may be faster to use the Model or AllModels setPars methods described in the next section). Examples of individual Component and Parameter object access:
# Component objects are accessible-by-name as Model object attributes: comp1 = m1.phabs comp2 = m1.powerlaw # Parameter objects are accessible-by-name as Component object attributes: par4 = m1.gaussian.LineE # ...and we can modify their values: par4.values = 3.895 m1.phabs.nH = 5.0 comp2.PhoIndex = 1.5 # Can also get a Parameter object directly from a Model, without going # through a Component. Just pass the Model the Parameter index number: par5 = m1(5) # Examples of numerical operations allowed with Parameter objects: par4 += 0.75 par4 *= 2.0 y1 = m1.phabs.nH*100.0 y2 = par4 + par5Note
For models with duplicate copies of components, see the Extended Tutorial for accessing Component objects by name.
Note that in the above examples, only the parameter's value is being accessed or modified. To change all or part of its FULL list of settings including auxiliary values: value, fit delta, min, bot, top, max, you can set its values attribute to a tuple or list of size 1-6:
par4.values = 4.3, .01, 1e-3 par4.values = [4.3, .01, 1e-3, 1e-2, 100, 200]Or for greater flexibility you can set it to a string using Standard XSPEC's newpar command syntax:
# This allows you to set new values non-consecutively. par4.values = "1.0, -.01,,,,150"A quick way to freeze or thaw a parameter is to toggle its frozen attribute:
par4.frozen = False par5.frozen = TrueTo link a parameter to one or more others, set its link attribute to a link expression string as you would have with the newpar command. Or if you simply want to link to a single parameter with no numerical operations, you can link it directly to that Parameter object. To remove the link, set link to an empty string or call the parameter's untie method:
par5.link = "2.3 * 4" # Link par 5 to par 4 with a multiplicative constant. par5.link = m.phabs.nH # Link par 5 to the specified par object. par5.link = "" # Removes the link. par5.untie() # Also removes the link.Also ALL linked parameters in a model object can be untied with a single call to the Model class untie method.
To display a parameter's full set of values (including auxiliary values), just print its values attribute:
>>> print(par4.values) [6.5, 0.05, 0.0, 0.0, 1000000.0, 1000000.0]
Setting Multiple Parameters At A Time
You can set multiple parameter values with a single call using the Model or AllModels (see ModelManager) setPars methods. This may be considerably faster than setting parameters one at a time through the individual Parameter objects as shown in the previous section. With setPars, the model will be recalculated just ONCE after all the changes have been made. But when setting through individual Parameter objects, the model will be recalculated after EACH parameter change:
# For Model object m1, supply 1 or more new parameter values in consecutive order: m1.setPars(2.5, 1.4, 1.0e3) # This changes pars 1, 2, and 3. # Can also change paramater auxiliary values by passing a string using the same # syntax as with Standard XSPEC's `newpar` command: m1.setPars(.95, "1.8,,-5,-4,10,10") # Now set parameters NON-CONSECUTIVELY by passing a Python dictionary object. # This example changes pars 1, 2, 4, and 6: m1.setPars(.95, 1.2, {4:9.8, 6:2.0})Parameters can also be initialized by passing values to the Model object constructor. You do this by setting the Model constructor's setPars keyword argument to a tuple, list, or dictionary (or just a single value or string if only setting the first parameter):
# Supply values for parameters 1 and 3, use defaults for the rest. m = Model("phabs*ga", setPars={1:1.5, 3:.2}) # Supply values for 1 and 2, use defaults for the rest. m = Model("phabs*ga", setPars=(1.5, 0.7)) # Supply value only for 1. m = Model("phabs*ga", setPars=1.5)Finally, if you wish to set multiple parameters that belong to different model objects, you must use the AllModels container's setPars method (see ModelManager). This follows the same syntax rules as the single Model setPars, except that you also supply the Model objects as arguments:
# Change pars 1 and 3 in m1, and pars 1 and 2 in m2: AllModels.setPars(m1, {1:6.4, 3:1.78}, m2, 3.5, 0.99)
Fitting
Once data and models are loaded, fitting is performed by calling the perform method of the Fit global object:
Fit.perform()Some of the more frequently modified fit settings are the type of statistic to minimize and the maximum number of fit iterations to perform. These settings are attributes of Fit:
Fit.nIterations = 100 Fit.statMethod = "cstat" Fit.statMethod = "chi"Please see FitManager and the Extended Tutorial for Fit's complete functionality.
To display the fit results at any time:
Fit.show()
Plotting
In Standard XSPEC, plot settings are adjusted using the setplot command while the plot is displayed through the plot command. In PyXspec, all plot settings and functionality is handled through the global Plot object (see PlotManager for full description). A device must be set before any plots can be displayed, and this done through the device attribute:
Plot.device = "/xs"The device can also be set to print to an output file in several formats. The list of possible devices is given by the cpd command in the Standard XSPEC manual.
A typical setting to adjust is the X-axis units. You can choose to plot channel numbers, or select from various energy and wavelength units. The strings can also be abbreviated. Examples:
Plot.xAxis = "channel" Plot.xAxis = "MeV" Plot.xAxis = "Hz" Plot.xAxis = "angstrom"The displays of individual additive components or background spectra is toggled by setting their attributes to a bool:
Plot.add = True Plot.background = FalseSimilarly log/linear settings for data plots (when using energy or wavelength units):
Plot.xLog = True Plot.yLog = FalseThe current plot settings are displayed with:
Plot.show()To actually display a plot, send 1 or more string arguments to the Plot __call__ method:
# Single panel plots Plot("data") Plot("model") Plot("ufspec") # Multi panel plots Plot("data chisq") Plot("data","model","resid") # Call Plot with no arguments to repeat the previously entered Plot command Plot()After displaying a plot, you can get an array of the plotted values by calling one of Plot's retrieval methods. All of these functions take an optional plot group number argument for the case of multiple plot groups, and all return the plot values in a Python list:
Plot("data") xVals = Plot.x() yVals = Plot.y() yVals2 = Plot.y(2) # Gets values for data in the second plot group modVals = Plot.model() # To get a background array, Plot.background must be set prior to plot Plot.background = True Plot("data") bkg = Plot.backgroundVals() # Retrieve error arrays xErrs = Plot.xErr() yErrs = Plot.yErr()