% sherpa
-----------------------------------------------------
Welcome to Sherpa: CXC's Modeling and Fitting Package
-----------------------------------------------------
CIAO 4.2 Monday, November 30, 2009
sherpa-1> def flatmod(params, x, y):
print("Lower left = {0},{1}".format(x[0],y[0]))
print("Upper right = {0},{1}".format(x[-1],y[-1]))
return x*0 + params[0]
sherpa-2> load_user_model(flatmod, "reallyflat")
sherpa-3> add_user_pars("reallyflat", ["ampl"])
sherpa-4> print(reallyflat)
usermodel.reallyflat
Param Type Value Min Max Units
----- ---- ----- --- --- -----
reallyflat.ampl thawed 0 -3.40282e+38 3.40282e+38
sherpa-5> reallyflat.ampl = 2
So here I have created the model reallyflat which sets every pixel to the value of 2 (i.e. it ignores the coordinates of the input). As a side effect it prints out the first and last coordinate sent to it (which will be the bottom-left and top-right pixels if there are no image filters). I created the model on the fly from within Sherpa for demonstration purposes; in normal use I would write the model definition in an external file and load it in using execfile() or import.
Now we load in an image and use the model:
sherpa-6> load_data("img.fits")
sherpa-7> set_source(reallyflat)
sherpa-8> image_data()
sherpa-9> image_model()
Lower left = 1.0,1.0
Upper right = 421.0,367.0
Displaying the model causes it to be evaluated, as shown above. If we use image_source() we see the same result (we have no PSF set up so there is no actual difference between the source and model views in this case).
sherpa-10> image_source()
Lower left = 1.0,1.0
Upper right = 421.0,367.0
The default coordinate system is "logical", which matches the FITS convention, where (1,1) is the center of the bottom-left pixel and each pixel is square with a width of 1. Below we change to the "physical" system, which uses the ACIS SKY coordinates (as this is an ACIS image). Note that the coordinate values sent to the model have changed.
sherpa-11> set_coord("physical")
sherpa-12> image_source()
Lower left = 4064.0,3511.0
Upper right = 4484.0,3877.0
sherpa-13> image_model()
Lower left = 4064.0,3511.0
Upper right = 4484.0,3877.0
sherpa-14> image_fit()
Lower left = 4064.0,3511.0
Upper right = 4484.0,3877.0
Lower left = 4064.0,3511.0
Upper right = 4484.0,3877.0
The image_fit() output above shows that the model is evaluated twice; once for the model display and once to create the residual image.
An example of a slightly-more complex 2D model would be:
def mygauss2d(pars,x,y):
"My version of a gaussian"
(sigma, ampl, xpos, ypos) = pars
r2 = (x - xpos)**2 + (y - ypos)**2
return ampl * np.exp(-r2/(sigma*sigma))
load_user_model(mygauss2d, "g2d")
add_user_pars("g2d",
["sigma", "ampl", "xpos", "ypos"],
[5,1,1,1], [0.1,0,0,0])
where I use the optional third and fourth arguments of add_user_pars() to set up both the starting and minimum values of the parameters.