Recreating the contour plots (e.g. from reg_proj)

Contour visualizations can be created by commands like contour_data or by the "region" error routines like reg_proj. If the visuzalization data is retrieved - for the contour or region error cases - then they can be used to create your own visualization of the data.

Using the example from the region error case:

sherpa> rdata = get_reg_proj()
sherpa> print(rdata)
x0      = [1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ,
 1.5989,1.6397,1.6805,1.7213,1.7621,1.8029,1.8437,1.8844,1.9252,1.966 ]
x1      = [0.    ,0.    ,0.    ,0.    ,0.    ,0.    ,0.    ,0.    ,0.    ,0.    ,
 0.0044,0.0044,0.0044,0.0044,0.0044,0.0044,0.0044,0.0044,0.0044,0.0044,
 0.0088,0.0088,0.0088,0.0088,0.0088,0.0088,0.0088,0.0088,0.0088,0.0088,
 0.0131,0.0131,0.0131,0.0131,0.0131,0.0131,0.0131,0.0131,0.0131,0.0131,
 0.0175,0.0175,0.0175,0.0175,0.0175,0.0175,0.0175,0.0175,0.0175,0.0175,
 0.0219,0.0219,0.0219,0.0219,0.0219,0.0219,0.0219,0.0219,0.0219,0.0219,
 0.0263,0.0263,0.0263,0.0263,0.0263,0.0263,0.0263,0.0263,0.0263,0.0263,
 0.0307,0.0307,0.0307,0.0307,0.0307,0.0307,0.0307,0.0307,0.0307,0.0307,
 0.0351,0.0351,0.0351,0.0351,0.0351,0.0351,0.0351,0.0351,0.0351,0.0351,
 0.0394,0.0394,0.0394,0.0394,0.0394,0.0394,0.0394,0.0394,0.0394,0.0394]
y       = [120.4008,111.1954,105.9672,104.6679,107.2309,113.5727,123.5955,137.1891,
 154.2326,174.5966,124.2947,113.1681,105.8671,102.3463,102.5419,106.3749,
 113.7521,124.5691,138.7115,156.057 ,129.7901,117.002 ,107.9083,102.4659,
 100.6141,102.2773,107.3664,115.7815,127.4135,142.1456,136.5014,122.2577,
 111.5935,104.468 ,100.8234,100.5874,103.6746,109.989 ,119.4257,131.8729,
 144.1443,128.6101,116.5533,107.9356,102.7023,100.784 ,102.099 ,106.5548,
 114.0505,124.4784,152.5089,135.8178,122.5124,112.5571,105.8999,102.4747,
 102.2031,104.9961,110.7565,119.3804,161.4389,143.701 ,129.2652,118.0987,
 110.1525,105.3636,103.6567,104.9465,109.1391,116.1344,170.8176,152.1246,
 136.6568,124.3843,115.2607,109.2264,106.2092,106.1271,108.8898,114.4004,
 180.5566,160.9863,144.5699,131.2797,121.0728,113.8921,109.669 ,108.3246,
 109.772 ,113.9176,190.5884,170.2079,152.9141,138.6822,127.4719,119.2292,
 113.8887,111.3747,111.6034,114.4846]
min     = [1.5989155 0.       ]
max     = [1.96602156 0.03943551]
nloop   = (10, 10)
fac     = 4
delv    = None
log     = (False, False)
sigma   = (1, 2, 3)
parval0 = 1.7824685288808277
parval1 = 0.012577007195034283
levels  = [102.56510032 106.44942569 112.09850947]

As with the plot objects, the print display does not include all possible attributes. In this case we have (dropping the "internal" attributes that start with an underscore):

sherpa> print([n for n in dir(rdata) if not n.startswith("_")])
['calc', 'conf_type', 'contour', 'contour_prefs', 'delv', 'fac', 'fast', 'levels', 'log', 'max', 'min', 'nloop', 'numcores', 'overcontour', 'parval0', 'parval1', 'point', 'point_prefs', 'prepare', 'sigma', 'stat', 'title', 'x0', 'x1', 'xlabel', 'y', 'ylabel']

This data can be used with the matplotlib contour command to create a visualization. Note that Sherpa returns the axes as 1D arrays, so they need to be reshaped to match the 2D data that contour expects:

from matplotlib import pyplot as pt

nx0, nx1 = rdata.nloop
x0 = rdata.x0.reshape((nx1, nx0))
x1 = rdata.x1.reshape((nx1, nx0))
y = rdata.y.reshape((nx1, nx0))

plt.contour(x0, x1, y, levels=rdata.levels,
            linestyles=["dotted", "dashed", "solid"], colors="brown")

plt.axvline(rdata.parval0, linestyle="dotted", alpha=0.5, color="black")
plt.axhline(rdata.parval1, linestyle="dotted", alpha=0.5, color="black")

plt.xlabel(rdata.xlabel)
plt.ylabel(rdata.ylabel)

Figure 1: manually displaying the error region

[The contours are in brown, with the inner contour as a dotted line, the middle contour as a dashed line, and the outer contour is solid. The "best-fit" location is indicated by two dotted-black lines that cover the whole plot and intersect at the location. The axes labels are "pl.gamma" for the horizontal axis and "gal.nH" for the vertical axis.]

[Print media version: The contours are in brown, with the inner contour as a dotted line, the middle contour as a dashed line, and the outer contour is solid. The "best-fit" location is indicated by two dotted-black lines that cover the whole plot and intersect at the location. The axes labels are "pl.gamma" for the horizontal axis and "gal.nH" for the vertical axis.]

Figure 1: manually displaying the error region

The data can also be displayed as an image, with the imshow command, but it requires more work to set the axes up correctly. In the following we also restrict the display range to exclude the points furthest from the best fit, invert the colormap, and add lines to the color bar to show the contour levels:

from matplotlib import pyplot as plt
from matplotlib.lines import Line2D

x0min = rdata.x0.min()
x0max = rdata.x0.max()
x1min = rdata.x1.min()
x1max = rdata.x1.max()

# Calculate the half-width of each pixel
hx0 = 0.5 * (x0max - x0min) / (nx0 - 1)
hx1 = 0.5 * (x1max - x1min) / (nx1 - 1)

# The edge of the image
extent = (x0min - hx0, x0max + hx0, x1min - hx1, x1max + hx1)

# The maximum value to display
vmax = rdata.levels[-1] + 10

plt.imshow(y, origin="lower", cmap="magma_r", vmax=vmax, extent=extent, aspect='auto')
cbar = plt.colorbar()

# Overlay the contours
plt.contour(x0, x1, y, levels=rdata.levels, linestyles=["dotted", "dashed", "solid"], colors="white")

# Add the contour lines to the colorbar
def add_line(lvl, ls):
    l = Line2D([0, 1], [lvl, lvl], linestyle=ls, color="white")
    cbar.ax.add_line(l)

add_line(rdata.levels[0], "dotted")
add_line(rdata.levels[1], "dashed")
add_line(rdata.levels[2], "solid")

plt.xlabel(rdata.xlabel)
plt.ylabel(rdata.ylabel)

Figure 2: Overlaying the contours on the statistic image

[The image pixels match the contours (as they should), but outside the solid line the image pixels are set to black, as they exceed the maximum-displayed value. The center of the contour (within the dotted line) is the lightest color. The colorbar has vertical lines indicating the values of the three contours (solid, dashed, and dotted).]

[Print media version: The image pixels match the contours (as they should), but outside the solid line the image pixels are set to black, as they exceed the maximum-displayed value. The center of the contour (within the dotted line) is the lightest color. The colorbar has vertical lines indicating the values of the three contours (solid, dashed, and dotted).]