Known Issues and Limitations for Sherpa 4.17

IPython and Sherpa

The sherpa application is a specialized version of IPython, in that several of the settings have been changed from their default values (such as the banner and prompt). To avoid colliding with existing IPython setup scripts the Sherpa-related settings are written to the $HOME/.ipython-ciao/ directory. If IPython users want any personal customizations to be available when running CIAO, they will have to copy them from $HOME/.ipython/ to $HOME/.ipython-ciao/.

IPython version

CIAO 4.17 only supports Python 3.11, but the version of IPython provided with CIAO depends on whether it was installed with the conda package manager or ciao-install.

RMF and ARF energy range must be the same.

If an RMF and ARF with different energy ranges are assigned to a data set in Sherpa, Sherpa will print an error message when trying to fit the data: ValueError: RMF data is invalid or inconsistent. See the CIAO Imaging Spectroscopy threads to learn how to create an ARF and RMF with the same energy range.

RMF with the energy grid starting with 0

If the RMF energy grid starts with 0, this first energ_lo bin will be set to 1e-10 to avoid issues evaluation of some models.

sherpa In [2]: rmf = unpack_rmf("mos1_p012_src_gti035.rmf")
/ciao4.14/lib/python3.8/site-packages/sherpa/astro/data.py:952: 
UserWarning: The minimum ENERG_LO in the RMF 'mos1_p012_src_gti035.rmf' was 0 and has been replaced by 1e-10
warnings.warn(wmsg)

sherpa In [3]: print(rmf.energ_lo[0])
1e-10

set_psf should not be used when a PSF model has been applied manually.

The Sherpa set_full_model functionality provides users with the flexibility to define model expressions that apply instrument responses or PSFs to some model components, while not applying the response or the PSF to other model components in the model expression. When a PSF has been applied to a model expression in this way - e.g., "set_full_model(psf1(gauss2d.g1) + beta2d.b1+const2d.c1)"- it is not necessary to apply the PSF with set_psf in the usual way when defining model expressions automatically with set_model; doing so will erroneously apply the PSF twice in the final model expression used for fitting.

# Manual application of a PSF

sherpa> load_image("image.fits")
sherpa> load_psf("psf1","psf.fits")
sherpa> set_full_model(psf1(gauss2d.g1) + beta2d.b1+const2d.c1)
sherpa> show_model()
Model: 1
((psfmodel.psf1(gauss2d.g1) + beta2d.b1) + const2d.c1)
...

# Automatic application of a PSF with set_psf()

sherpa> load_image("image.fits")
sherpa> load_psf("psf1","psf.fits") 
sherpa> set_model(gauss2d.g1)
sherpa> set_psf(psf1)
sherpa> show_model()
...
Model: 1
psfmodel.psf1(gauss2d.g1)
...

User-specified PSF 'origin' position is off by 1 pixel in x and y directions.

If the user chooses to set the optional, hidden 'origin' parameter of a PSF model to the (x,y) coordinates of a chosen pixel in that kernel, the position should be decremented by one in both x and y directions. For example, if DS9 was used to find that the brightest pixel is at (122, 131), 'psf.origin' should be set equal to '(121, 130)'. This is necessary because Sherpa internally stores data, PSF, and kernel information as NumPy arrays - which start counting from element 0 - whereas DS9 uses pixel coordinates, with the first pixel having a value of (1,1).

sherpa> print get_psf()
psfmodel.psf
   Param        Type          Value          Min          Max      Units
   -----        ----          -----          ---          ---      -----
   psf.kernel  frozen psf_f1_norm_0.25pix.fits
   psf.size    frozen   (256, 256)   (256, 256)   (256, 256)
   psf.center  frozen   (128, 128)   (128, 128)   (128, 128)
   psf.radial  frozen            0            0            1           
   psf.norm    frozen            1            0            1    

sherpa> psf.origin = [121,130]

sherpa> print get_psf()
psfmodel.psf
   Param        Type          Value          Min          Max      Units
   -----        ----          -----          ---          ---      -----
   psf.kernel  frozen psf_f1_norm_0.25pix.fits
   psf.size    frozen   (256, 256)   (256, 256)   (256, 256)
   psf.center  frozen   (128, 128)   (128, 128)   (128, 128)
   psf.origin  frozen   (121, 130)   (121, 130)   (121, 130)
   psf.radial  frozen            0            0            1           
   psf.norm    frozen            1            0            1 

If psf.origin is not set by the user - leaving Sherpa to automatically locate the position of the brightest pixel in the kernel and set this as the origin - Sherpa will find the correct location, provided that a) the brightest pixel should coincide with the origin in the analysis, and b) that there is only one pixel that has that brightest value.

PSF Images must have even-numbered axis lengths.

Caveat: PSF images should have axis lengths with even number of pixel. An odd number of pixels may result in the data being fitted to appear to be shifted by 1 pixel in either direction.

Users can crop the PSF to an even number of pixels using DM filtering syntax

sherpa> load_psf("mypsf", "psf.fits[#1=1:100,#2=1:100]")

Coordinate systems

Image fitting in sherpa should only be done in either physical or image (logical) coordinates. Sherpa does not take size of pixels into account when set_coord is used to change to world (wcs).

Narrow 2D Gaussians

Sherpa evaluates 2D models at the center of the pixel. If the model changes rapidly within a pixel, it can lead to unexpected results. This can happen if a Gaussian, or other 2D model, is given a small FWHM (much less than a pixel). The example below illustrates the problem

sherpa> dataspace2d(100,100)
sherpa> set_model(gauss2d.g1)
sherpa> g1.xpos = 50
sherpa> g1.ypos = 50
sherpa> g1.fwhm = 0.01
sherpa> g1.ampl = 1.0
sherpa> calc_model_sum2d()
           1.0
sherpa> g1.xpos = 50.1
sherpa> calc_model_sum2d()
           3.8725919148189143e-121

Here a 2D Gaussian with a FWHM=0.01 pixel is created. When the Gaussian is center on the pixel center, the sum of the model is 1.0. However, shifting the narrow Gaussian by a small amount, 0.1 pixels in the X direction, results in the model now essentially summing to 0.

This leads to very steep slopes in the error surface separated by large areas where the fit statistic barely changes with positon. The error surface looks like a stair case rather than say the classic parabolic shape.

sherpa> load_data("psf.fits")
sherpa> load_psf("p", "psf.fits")
sherpa> set_psf("p")
sherpa> set_model(gauss2d.g1)
sherpa> g1.fwhm=0.1
sherpa> freeze(g1.fwhm)
sherpa> set_method("neldermead")
sherpa> g1.xpos = 196/2.0
sherpa> g1.ypos = 138/2.0
sherpa> fit()
...
   g1.xpos        97.9973     
   g1.ypos        68.9982     
   g1.ampl        1.00023     

sherpa> int_unc( g1.xpos, min=96, max=100, nloop=400)

Here the same data are loaded as both the dataset and the PSF. A 2D Gaussian is used to approximate a delta function and the data are fit. Trying to run conf will yeild tiny errors on the xpos and ypos values. It is best to look at the error surface using int_unc.

The stair-step error surface is a result of the Gaussian only being evaluated at the pixel center. The large region between the pixels is when the model evaluates to approximately zero.

This is not restricted just to Gaussians. All 2D models behave the same way. It can also be a problem for 1D datasets which are evaluated over each bin - such as with PHA data.