Plotting in Sherpa Using Common Options

Introduction

Sherpa uses the Matplotlib plotting package for 1-D (and some 2-D) visualization. This thread presents two examples of creating default Sherpa plots with commonly used Sherpa plotting commands, and modifying the plots using Matplotlib commands from the Sherpa prompt.

New in CIAO 4.16 is support for using the bokeh library for plotting rather than matplotlib. This thread will only use matplotlib for its examples: please see the Starting the bokeh backend for advice on using it.

Getting Started

The sample data files used in this thread are available in sherpa.tar.gz, as explained in the Sherpa Getting Started thread, in the setplot_manual sub-directory directory:

case1_1.pha
core1.arf
core1.rmf

In the first example, we create X-ray spectral plots, altering plotting preferences to create custom labels, colors, and axes ranges. In the second example, we use the same X-ray spectral profile, but re-plot it as a line plot instead of a scatter plot.

Plotting X-ray Spectra

We begin by loading a PHA data set and its associated instrument response using load_pha, and filter the data to include only the range between 0.3-7.0 keV with notice_id. Then plotting the data is simply done with plot_data (Figure 1):

sherpa> load_pha("case1_1.pha")
read ARF file core1.arf
read RMF file core1.rmf
sherpa> notice_id(1, 0.3, 7.)
dataset 1: 0.00146:14.9504 -> 0.292:7.008 Energy (keV)

sherpa> plot_data()

Figure 1: Screenshot of PHA Dataset

Figure 1: Screenshot of PHA Dataset

Plot of count-rate spectrum with the default plotting parameters.

A "hardcopy" version of the plot can be generated via the "save figure" button in the window or with the plt.savefig command (Figure 2):

sherpa> plt.savefig('pha_data.png')

Figure 2: Plot of PHA Dataset

Figure 2: Plot of PHA Dataset

The PNG version of Figure 1 created with plt.savefig.

As we can see, Sherpa plots data in linear-scale by default, but sometimes, it is more useful to plot a logarithmic-scale. To switch to log-scale, we can use Matplotlib commands to change the existing plot; for instance

sherpa> plt.xscale('log')
sherpa> plt.yscale('log')
sherpa> plt.ylim(0.001, 0.1)
(0.001, 0.1)

You can also set the xlog or ylog plot options when you call one of the plot functions (Figure 3):

sherpa> plot_data(xlog=True, ylog=True)
sherpa> plt.ylim(0.001, 0.1)
(0.001, 0.1)

Figure 3: Plot of PHA Dataset on Log-Scale

Figure 3: Plot of PHA Dataset on Log-Scale

Plot of count-rate spectrum on logarithmic-scale with the y-axis range set to 0.001-0.1 cts/sec/keV.

There are a number of ways to change the scale of plots, including:

When the plot is created

Set the xlog and ylog attributes of the plot call, which will over-ride any plot preference. This is shown above.

After the plot has been created

This is shown above with the plt.xscale() and plt.yscale() commands.

Changing the plot preferences

There are several ways to change the plot preferences:

1. The easiest is with the commands: set_xlog, set_ylog, set_xlinear, and set_ylinear. These can be used to change all plots, or a single type. For instance, we can set all plots to use a linear y-axis apart from residual-style plots (this does not include "ratio" or "delchi" style plots):

   sherpa> set_ylog()
   sherpa> set_ylinear('resid')
   

2. The plot preferences can be changed directly (this is what the set_* commands do). The preferences can be retrieved with the get_plot_prefs function, and changing an element of the returned dictionay will change the preferences. For example, to change data plots to always be displayed with a logarithmic scale for both axes you would say:

   sherpa> p = get_plot_prefs("data") 
   sherpa> p["xlog"] = True
   sherpa> p["ylog"] = True
   

   Changes in CIAO 4.16

   The get_plot_prefs routine is new in CIAO 4.16 and can be used to access the preferences of all the plot types (e.g. "data", "model", "arf"). The existing routines - get_data_plot_prefs and get_model_plot_prefs - can still be used.

The default labelling for tickmarks with a logarithmic scale does not work well for Chandra-style data, but we can adjust this if required using the matplotlib ticker API. An example of this is shown below. First we create the plot (matching Figure 3):

sherpa> plot_data(xlog=True, ylog=True)
sherpa> plt.ylim(0.001, 0.1)
(0.001, 0.1)

The tick formatting is controlled with (in this case) a simple string formatter:

sherpa> import matplotlib.ticker as mtick
sherpa> tformatter = mtick.FormatStrFormatter("%g")

We can apply this new formatter using the gca function to obtain the current axes. This can then be used to change the formatting of both axes (Figure 4):

sherpa> ax = plt.gca()
sherpa> ax.xaxis.set_major_formatter(tformatter)
sherpa> ax.yaxis.set_major_formatter(tformatter)

Figure 4: Adjusting the tick labels for a logarithmic scale

Figure 4: Adjusting the tick labels for a logarithmic scale

Grouping Data

In order to use Gaussian statistics to fit a model to a data set, it is often necessary to "group" the data—i.e., combine channels until you have enough counts—before use. It is possible to set and change the grouping of a file after it has been read into Sherpa by using the commands: set_grouping, group, group_counts, group_snr, group_adapt, group_adapt_snr, group_bins, and group_width. See the Sherpa thread Changing the grouping scheme of a data set within Sherpa for details.

We use the group_counts function in this example to force each bin to have a minimum number of counts and apply the grouping after the filter ( the behavior here has changed in CIAO 4.16 so that you no-longer have to set the tabStops argument):

sherpa> notice_id(1, 0.3, 7.)
dataset 1: 0.292:7.008 Energy (keV) (unchanged)

sherpa> group_counts(15)
dataset 1: 0.292:7.008 Energy (keV) (unchanged)

The resulting plot is shown in Figure 5.

sherpa> plot_data(xlog=True, ylog=True, linestyle="solid", marker="none")

Figure 5: Plot of Grouped PHA Dataset on Log-Scale

Figure 5: Plot of Grouped PHA Dataset on Log-Scale

Plot of the grouped, filtered count-rate spectrum on logarithmic-scale. The linestyle and marker options have been over-ridden so that the individual bins are more obvious.

Notice that the noise/error-bars have been reduced by the grouping.

Grouping and Filters

Note that after using a Sherpa group function, any notice or ignore filter previously applied will be reset if the data has been pre-grouped, e.g. with dmgroup, and must be re-applied, as the group functions restore use of the full range of the data in the analysis; however, if the data file is ungrouped, any applied notice/ignore filters will be maintained after grouping.

---

Fitting a Model to Data

Next we will fit a simple, 1-D absorbed power-law model to the data and plot the fit using plot_fit (Figure 6). Note that by default, plot_fit creates a plot with data points marked by circles and associated error bars, and the model over-plotted as a solid orange line.

sherpa> set_source(xsphabs.abs1 * powlaw1d.p1)
sherpa> fit()
Dataset               = 1
Method                = levmar
Statistic             = chi2gehrels
Initial fit statistic = 7.30455e+11
Final fit statistic   = 59.8111 at function evaluation 81
Data points           = 88
Degrees of freedom    = 85
Probability [Q-value] = 0.982656
Reduced statistic     = 0.70366
Change in statistic   = 7.30455e+11
   abs1.nH        0.0343182    +/- 0.014151
   p1.gamma       1.75149      +/- 0.0824047
   p1.ampl        4.01875e-05  +/- 2.59962e-06

sherpa> plot_fit(xlog=True, ylog=True)

Figure 6: Plot of Grouped PHA Dataset and Fitted Model

Figure 6: Plot of Grouped PHA Dataset and Fitted Model

The grouped, filtered count-rate spectrum on logarithmic-scale, with an absorbed power-law model fitted to the data (orange line).

The plot_model and plot_source commands are available for plotting the total convolved and unconvolved source model, respectively; and plot_model_component and plot_source_component plot one or a combination of individual model components contributing to a fit, where a multi-component source model is used (see the thread "Fitting a PHA Data Set with Multiple Responses" for a demonstration of this functionality).

The plot_source command supports the factor setting of the set_analysis command. This means that calling plot_source while the set_analysis 'factor' setting is 1 will plot E F(E) versus E in keV, or λ f(λ) versus λ in Angstroms, depending on which units are set for the spectral analysis. A set_analysis 'factor=2' setting will plot E² F(E) versus E, or λ² f(λ) versus λ.

sherpa> plt.clf()
sherpa> plt.subplot(2, 1, 1)
sherpa> set_analysis("energy", factor=1)
dataset 1: 0.292:7.008 Energy (keV)
sherpa> plot_source(xlog=True, ylog=True, clearwindow=False)

sherpa> plt.subplot(2, 1, 2)
sherpa> set_analysis("wavelength", factor=2)
dataset 1: 1.76918:42.4603 Wavelength (Angstrom)
sherpa> plot_source(xlog=True, ylog=True, clearwindow=False)

sherpa> plt.subplots_adjust(left=0.2, hspace=0.5)

Figure 7: Changing the axis units

Figure 7: Changing the axis units

In the top plot the unconvolved source model for data set 1 is plotted in units of photons sec^-1 cm^-2 versus keV. In the bottom plot, the units are be Angstroms photons sec^-1 cm^-2 versus Angstroms.

Note the use of Matplotlib commands to create the two plots, and setting clearwindow=False in the plot_source calls to make sure the newly-created plot areas were not automatically erased.

Let's return to the default setting before going any further (and ensure we use logarithmic scales for the plot axes):

sherpa> set_analysis('energy', factor=0)
dataset 1: 0.292:7.008 Energy (keV)
sherpa> set_xlog()
sherpa> set_ylog()

---

Creating Multi-Plots

We will now create a multi-plot with the plot_fit_delchi function, which plots the fitted spectrum and the δχ residuals of the fit.

sherpa> plot_fit_delchi()

Figure 8: Model Fit and δχ Residuals Multi-Plot

Figure 8: Model Fit and δχ Residuals Multi-Plot

The upper plot shows the data (blue points) with fitted model (orange), and the δχ residuals plot is shown on the bottom.

The plot can be "inspected" and modified, using a variety of Matplotlib commands. The plt.gcf and plt.gca commands return objects that represent the current figure and axis-pair (respectively). The axes field of the figure lists all axis-pairs shown on the figure—in this case there are two of them (top plot and bottom plot):

sherpa> fig = plt.gcf()
sherpa> ax1, ax2 = fig.axes
sherpa> print(ax1.get_ylabel())
Counts/sec/keV
sherpa> print(ax2.get_ylabel())
Sigma

By default the second (lower) plot is the "active" axis. That is shown by:

sherpa> print(plt.gca().get_ylabel())
Sigma

The order of the plots can be shown by using the set_facecolor method on each axis, as shown in Figure 9:

sherpa> ax1.set_facecolor('lightgray')
sherpa> ax2.set_facecolor('lightgreen')

Figure 9: Changing the plots

Figure 9: Changing the plots

Many things can be changed—e.g. the axis fonts, sizes, point styles and colors—but here we will concentrate on adding a label to the upper plot (Figure 10):

sherpa> plt.sca(ax1)
sherpa> lbl = plt.text(0.3, 1e-3, 'Chandra', fontsize=14)
sherpa> lbl.set_position((1, 1e-3))
sherpa> lbl.set_horizontalalignment('center')

Figure 10: Adding a label to a plot

Figure 10: Adding a label to a plot

As we stored the object returned by plt.text we could use it to change the appearance of the label, in this case it's position. Tab-completion in an interactive environment like Sherpa is very helpful in situations like this (saving you from having to type out all of set_horizontalalignment!).

---

Creating Plots of Multiple Datasets

As an aside, the Sherpa FAQ includes guides to creating many different types of plots, including those that display multiple datasets. The two following linked examples may be of interest to users:

• How can I create multiple plots in the same window?
• How can I include multiple datasets in the same plot?

Changing a Scatter Plot to a Line Plot

We continue on with an example of how to modify the current scatter plot into a line plot.

First we over-ride the yerrorbars and xerrorbars preferences so that error bars are not shown (Figure 11):

sherpa> plot_fit_delchi(yerrorbars=False, xerrorbars=False)

Figure 11: Removing Error Bars

Figure 11: Removing Error Bars

There are multiple ways to modify the plot. For example, we could say:

sherpa> ax1, ax2 = plt.gcf().axes
sherpa> print(ax1.lines)
<Axes.ArtistList of 4 lines>

sherpa> ax1.lines[0].set_linestyle('-')
sherpa> ax1.lines[1].set_marker(None)

but in this particular case we can use the plot preferences to generate Figure 12:

sherpa> dprefs = get_plot_prefs("data")
sherpa> dprefs["linestyle"] = "solid"
sherpa> dprefs["marker"] = "none"
sherpa> plot_fit_delchi(yerrorbars=False, xerrorbars=False)
sherpa> dprefs["linestyle"] = "none"
sherpa> dprefs["marker"] = "."

Figure 12: Removing Error Bars for a Spectral Profile

Figure 12: Removing Error Bars for a Spectral Profile

If a plot is manipulate with one or more matplotlib commands then they may not cause the display to update. The plot can be redrawn with a call to the redraw_in_frame method of the axis. For example, if we were to change the color of the markers in the residual display:

sherpa> ax1 = plt.gcf().axes[1]
sherpa> ax1.set_color("k")

we may need to add

sherpa> ax1.redraw_in_frame()

Changes to Sherpa plots

The plots created by Sherpa functions like plot_data and plot_fit_resid are not guaranteed to remain the same between CIAO versions, as they are updated to improve the display and take advantage of capabilities of the plotting library. This is of note for people who want to adjust the plots created by these commands, as it may depend on the version of CIAO in use as to what changes to make.

Setting Default Plot Preferences

In order to change the default plot preferences for Sherpa plotting commands—e.g., to have plot_data create data plots with a logarithmic scaling by default, as opposed to linear—you should add the following to a Sherpa Python customization file, located in ~/.ipython-ciao/profile_sherpa/startup, as shown in the following example. The file name should start with a two-digit integer and generally you would want it to be loaded last, so a value like 90 should suffice; in this case we have created the file ~/.ipython-ciao/profile_sherpa/startup/90-plot_customize.py:

dplot = get_plot_prefs("data")
dplot['xlog'] = True
dplot['ylog'] = True
dplot['linestyle'] = 'solid'
dplot['marker'] = none'

mplot = get_plot_prefs("model")
mplot['xlog'] = True
mplot['ylog'] = True

These settings specify that the data pointed are connected with a solid-line when plotted, and both the x- and y-axes of data and model plots created in Sherpa should be drawn using a logarithmic scale, each time a function like plot_data, plot_model, or plot_fit is called.

You must create a Python script that begins with a number, which is the order they will be called on, beginning with 00.

Writing Plots to File

Before exiting Sherpa, you may want to save your plot to file. In this example, we will create encapsulated Postscript, PDF, and PNG versions of the plot.

sherpa> plt.savefig("line.eps", orientation="landscape")
sherpa> plt.savefig("line.pdf", orientation="landscape")
sherpa> plt.savefig("line.png")

Scripting It

The file fit.py is a Python script which performs the primary commands used above; it can be executed by typing %run -i fit.py on the Sherpa command line.

The Sherpa script command may be used to save everything typed on the command line in a Sherpa session:

sherpa> script(filename="sherpa.log", clobber=False)

(Note that restoring a Sherpa session from such a file could be problematic since it may include syntax errors, unwanted fitting trials, et cetera.)

The CXC is committed to helping Sherpa users transition to Matplotlib syntax from ChIPS. Please see the ChIPS to Matplotlib conversion guide for help, and contact the CXC Helpdesk if you still have problems.