Tuesday, April 14, 2009

XSPEC bug fixes update

Here are the bug fixes from the last month. They can be found in the usual place.

12.5.0w A new tclout option has been added to make it easier to retrieve the fit parameters' sigma values. The syntax for this option is:

tclout sigma [<>:]n

where n is the parameter number. If it is not a variable parameter or if the fit was unable to calculate its sigma, a value of -1.0 is returned. Report added on Mar 06, 2009.

12.5.0x The model.dat entry for the vequil model is missing the parameter for Ar abundance, causing the vequil parameters which follow it to be misinterpreted in the code. (This bug does not affect the equil model.) Report added on Mar 06, 2009.

12.5.0y The kerrconv convolution model code is still using the xspec11-only LMODDIR environment variable for locating the kerrtable.dat model data file. This should be modified for usage in xspec12, where LMODDIR doesn't exist. Report added on Mar 11, 2009.

12.5.0z This removes an ambiguous reference build error specific to gcc-4.1.1 on the Solaris-2.9 platform. Our thanks to Dacheng Lin for pointing this out. Report added on Mar 26, 2009.

12.5.0aa The diskir model goes to infinity when its rirr parameter = 1.0, which is also the default value of rirr's lower limit. Our thanks to Brian Refsdal for pointing this out. Report added on Mar 26, 2009.

12.5.0ab The nsmax model function is only able to find its auxiliary files when it is run directly from XSPEC's model data directory (heasoft-6.x/spectral/modelData). Our thanks to Stephen Doe for pointing this out. Report added on Mar 31, 2009.

12.5.0ac Fixes to a couple of obscure cases of gain parameter usage: If a gain parameter belongs to a response temporarily replaced by a dummy response, or a gain parameter is indirectly removed through the data command erasing its associated spectrum and response, the gain shift may still be applied to a response occupying its former spectrum and source number. This patch also improves the show response output. Report added on Apr 06, 2009.

12.5.0ad The nthcomp model has an uninitialized variable for the case of input energies less than kT_bb/10^4, which may cause a crash to occur. Our thanks to Gulab Dewangan for pointing this out. Report added on Apr 14, 2009.

12.5.0ae For those linking the models library into their own programs, the swind1 model crashes when outputing a warning message as it indirectly accesses a variable intended for use only in XSPEC. Our thanks to Stephen Doe for pointing this out. Report added on Apr 14, 2009.

Friday, April 03, 2009

Notes on using SWIG to import C++ in Python

As a test case using my heasp library which is being converted from C to C++. Created heasp.i file to define interface. Basically a concatenation of the heasp.h, PHA.h and PHAII.h files with the following at the top. Adding additional classes will just require cat'ing their *.h file to the end and adding the appropriate #include at the top. Note the %include of "std_string.i" which is required to make the string arguments work. Similar %include specifications will be required if other STL classes are arguments for methods.

%module heasp

%include "std_string.i"



#include "heasp.h"
#include "PHA.h"
#include "PHAII.h"

The interface file is then converted to C++ using

swig -python -c++ -classic heasp.i

which creates heasp_wrap.cxx and heasp.py. These must then be combined with the *.o files to build a shareable library called _heasp.so (note that the _ prefix and .so suffix are required). Compiling heasp_wrap.cxx will require include files from python. On my mac the relevant flag is
and on the lab Linux network

On the mac the command to make the shareable library is
ld -bundle -flat_namespace -undefined suppress -o _heasp.so ${OFILES} heasp_wrap.o -L${HEADAS}/lib -lCCfits_2.1 -lcfitsio_3.12
and on Linux
$(CXX) -shared -o _heasp.so ${OFILES} heasp_wrap.o -L${HEADAS}/lib -lCCfits_2.1 -lcfi

The module is loaded into python by
>>> import heasp
On the mac this generates an undefined symbol ___eprintf. This may require an update of Xcode to fix. On linux (running python2.5) there are no problems. Example use is

>>> import heasp
>>> spectrum = heasp.PHA()
>>> spectrum.read('file1.pha',1,1)
>>> spectrum.disp()
>>> spectrum.write('test.pha')

Thursday, April 02, 2009


Updated and much simplified chkarf so it actually lists the correct keywords as mandatory. This tool really should just be a perl script - the current program is overkill.

Wednesday, April 01, 2009

using openMP in gcc

I'm testing out parallelization using openMP on an 8-core Linux box. To build xspec to support openMP options requires the hmakerc to be edited to add -fopenmp to CFLAGS, FFLAGS, and CXXFLAGS and -lgomp to be added to F77LIBS4C (I added it immediately before -lgfortran).

First test is in sumape.f to parallelize over the individual elements when interpolating the continuum and pseudo-continuum...

C$OMP PARALLEL PRIVATE(ien,limdown,limup,ihigh,energy)


DO ielt = 1, nmelt





The variables defined as private are to avoid threads overwriting variables being set and used by other threads.

Timing tests show that parallelization doesn't win enough in this case. Including the parallel directives slows a newpar down from 0.02 seconds to 0.4 seconds due to the overhead.