From ae4c339094396ee85bb2c70b03dafd3a60d1e965 Mon Sep 17 00:00:00 2001 From: Doug Burke Date: Thu, 8 Jun 2023 16:16:28 -0400 Subject: [PATCH 1/6] Tests: xspec edge-handling cases Check some error conditions - eg. using an unknown table name - to ensure that we catch the error. These are often coming from XSPEC C++ code so it would be nice if we knew we were creating any extra string output, but it's not clear how best to do that. Three unrelated changes - use the tmp_path fixture for providing access to a temporary directly fo a test rather than using routines from the tempfile module - a test of set_xsabund() with a file containing < 30 elements - remove the use of the clean_astro_ui fixture as not needed --- sherpa/astro/xspec/tests/test_xspec_unit.py | 80 ++++++++++++++++++++- 1 file changed, 78 insertions(+), 2 deletions(-) diff --git a/sherpa/astro/xspec/tests/test_xspec_unit.py b/sherpa/astro/xspec/tests/test_xspec_unit.py index 8e174b7376..e3f0340c4c 100644 --- a/sherpa/astro/xspec/tests/test_xspec_unit.py +++ b/sherpa/astro/xspec/tests/test_xspec_unit.py @@ -269,6 +269,49 @@ def test_abund_element(): assert fe == pytest.approx(2.69e-05) +@requires_xspec +def test_abund_get_invalid_element(caplog): + """Check what happens if sent the wrong element name""" + + from sherpa.astro import xspec + + # TODO: TypeError is not the best error type here. + with pytest.raises(TypeError, match="^could not find element 'O3'$"): + xspec.get_xsabund("O3") + + assert len(caplog.records) == 0 + + +@requires_xspec +def test_abund_set_invalid_name(caplog): + """Check what happens if sent an unknown table + + It is unlikely that the name "foo-foo" will become valid. + """ + + from sherpa.astro import xspec + + with pytest.raises(ValueError, match="^Cannot read file 'foo-foo'. It may not exist or contains invalid data$"): + xspec.set_xsabund("foo-foo") + + assert len(caplog.records) == 0 + + +@requires_xspec +def test_xsect_set_invalid_name(caplog): + """Check what happens if sent an unknown table + + It is unlikely that the name "foo-foo" will become valid. + """ + + from sherpa.astro import xspec + + with pytest.raises(ValueError, match="^could not set XSPEC photoelectric cross-section to 'foo-foo'$"): + xspec.set_xsxsect("foo-foo") + + assert len(caplog.records) == 0 + + def validate_xspec_setting(getfunc, setfunc, newval, altval): """Check we can change an XSPEC setting. @@ -397,6 +440,39 @@ def test_abund_change_file(tmp_path): assert out[n] == pytest.approx(elems[n]) +@requires_xspec +def test_abund_change_file_subset(tmp_path): + """What happens if send in too-few elements?""" + + from sherpa.astro import xspec + + elems = {n: i * 0.1 for i, n in enumerate(ELEMENT_NAMES) + if i < 10} + + tmpname = tmp_path / "abunds.xspec" + with open(tmpname, "w") as tfh: + for v in elems.values(): + tfh.write(f"{v}\n") + + oval = xspec.get_xsabund() + try: + xspec.set_xsabund(str(tmpname)) + + abund = xspec.get_xsabund() + out = {n: xspec.get_xsabund(n) + for n in ELEMENT_NAMES} + + finally: + xspec.set_xsabund(oval) + + assert abund == 'file' + for i, n in enumerate(ELEMENT_NAMES): + if i < 10: + assert out[n] == pytest.approx(elems[n]) + else: + assert out[n] == pytest.approx(0) + + @requires_xspec def test_xset_change(): """Can we change the xset setting. @@ -781,7 +857,7 @@ def test_xspec_model_requires_bins_very_low_level(clsname, xsmodel): @requires_fits @requires_data @requires_xspec -def test_xspec_tablemodel_requires_bin_edges(make_data_path, clean_astro_ui): +def test_xspec_tablemodel_requires_bin_edges(make_data_path): """Check we can not call a table model with a single grid. This used to be supported in Sherpa 4.13 and before. @@ -800,7 +876,7 @@ def test_xspec_tablemodel_requires_bin_edges(make_data_path, clean_astro_ui): @requires_fits @requires_data @requires_xspec -def test_xspec_tablemodel_requires_bin_edges_low_level(make_data_path, clean_astro_ui): +def test_xspec_tablemodel_requires_bin_edges_low_level(make_data_path): """Check we can not call a table model with a single grid (calc). This used to be supported in Sherpa 4.13 and before. From 45ab2ed4bbdd3ce29b88f1def97ddbd4e81877f0 Mon Sep 17 00:00:00 2001 From: Douglas Burke Date: Wed, 28 Jul 2021 12:58:22 -0400 Subject: [PATCH 2/6] xspec: move from xsFortran to FunctionUtility interface Rather than using the xsFortran routines - like FGSOLR - we now directly use the FunctionUtility C++ interface. Note that the code is based on the xsFortran.cxx file as of HEASOFT 6.28 (although as it is not the era of 6.31.1 it's a bit hard to check that everything is supported in XSPEC 12.12.0). Additional changes: - the removal of a lot of try/except handling since these routines should not be throwing errors. - use IosHolder::setStreams() to temporarily replace the cout or cerr channel used by XSPEC rather than manually redirecting the std::cout/cerr stream, which is a lot cleaner This includes adding this logic around the call to the FunctionUtility::checkXsect call, which - without this change - will cause output to the stderr channel if an invalid name is used (which is what we are trying to validate in this check so it's annoying to have XSPEC tell us it's invalid as well as return false so we can tell the user). - code calling PyErr_Format has been changed to return it, rather than have a separate "return NULL;" following statement, as the Python C-API guarantees that PyErr_Format returns NULL. --- sherpa/astro/xspec/src/_xspec.cc | 631 +++++------------- .../include/sherpa/astro/xspec_extension.hh | 2 +- 2 files changed, 167 insertions(+), 466 deletions(-) diff --git a/sherpa/astro/xspec/src/_xspec.cc b/sherpa/astro/xspec/src/_xspec.cc index 9959027166..e6e84c5390 100644 --- a/sherpa/astro/xspec/src/_xspec.cc +++ b/sherpa/astro/xspec/src/_xspec.cc @@ -23,91 +23,20 @@ #include #include -// The symbols listed in XSPEC version 12.9.1 -// at https://heasarc.gsfc.nasa.gov/xanadu/xspec/manual/XSappendixExternal.html -// are given below. Note that this is the C/FORTRAN interface, not the -// more-featureful FunctionUtility module. -// -// Functions which are used below: -// FNINIT Initializes data directory locations needed by the models. See below for a fuller description. -// FGABND Get an element abundance. -// FGCHAT Get current chatter level setting for model functions' output verbosity. -// FPCHAT Set the chatter level. Default is 10, higher chatter levels produce more output. -// FGMSTR Get a model string value (see XSPEC xset command). -// FPMSTR Set a model string value. -// FGDATD Get the model .dat files path. -// FPDATD Set the model .dat files path. -// FGMODF Get the model ion data path. -// FPSLFL Load values of a file solar abundance table (see abund command). -// FGSOLR Get the solar abundance table setting. -// FPSOLR Set the solar abundance table. -// FGXSCT Get the cross section table setting. -// FPXSCT Set the cross section table. -// csmgh0 Get the cosmology H$_0$ setting (see the cosmo command). -// csmph0 Set H$_0$. -// csmgl0 Get $\Lambda_0$. -// csmpl0 Set $\Lambda_0$. -// csmgq0 Get q$_0$. -// csmpq0 Put q$_0$. -// xs_getVersion (or xgvers) Retrieve XSPEC's version string. -// -// Functions which are not wrapped as their functionality is available: -// RFLABD Read abundance data from a file, then load and set this to be the current abundance table. (Essentially this combines a file read with the FPSLFL and FPSOLR functions.) -// -// Functions not wrapped as not felt to be that useful: -// fzsq Computes the luminosity distance, (c/H$_0$)*fzsq. The function is valid for small values of q$_0$*z for the case of no cosmological constant and uses the approximation of Pen (1999 ApJS 120, 49) for the case of a cosmological constant and a flat Universe. The function is not valid for non-zero cosmological constant if the Universe is not flat. -// -// Functions not wrapped since they are not useful as is (they need -// functionality from 12.9.1 to set the XFLT keywords): -// DGFILT Get a particular XFLT keyword value from a data file. -// DGNFLT Get the number of XFLT keywords in a data file. -// -// Other symbols in xsFortran.h are: -// DGQFLT Does a XFLT keyword exist? -// PDBVAL Set a database value -// -// Symbols in 12.9.1/HEASOFT 6.22 but not in 12.9.0/HEASOFT 6.19 -// FGABNZ -// FGTABN -// FGTABZ -// FGELTI -// FGNELT -// FGABFL -// FPABFL -// FGAPTH -// FPAPTH -// csmpall -// DPFILT -// DCLFLT -// GDBVAL -// CDBASE -// FGATDV -// FPATDV -// -// These seem unlikely to be useful for Sherpa -// xs_getChat -// xs_write -// xs_read -// -// These are numeric functions which we should have available elsewhere -// xs_erf -// xs_erfc -// gammap -// gammq -// -// TODO:: switch to C++ FuntionUtility interface rather than use xsFortran.h +// xsfortran is only needed to support FNINIT as there's (as of XSPEC 12.13.0) +// no version of this functionality in FunctionUtility. // +#include +#include +#include +#include + // funcWrappers: C_ are declared here; the other models are defined in // functionMap.h // -#include #include #include -// TODO: is this defined in an XSPEC header file? -#define ABUND_SIZE (30) // number of elements in Solar Abundance table - - // The XSPEC initialization used to be done lazily - that is, only // when the first routine from XSPEC was about to be called - but // the module is now set up so that we need to know the version @@ -133,262 +62,118 @@ static int _sherpa_init_xspec_library() return EXIT_FAILURE; } - // Stream buffer for redirected stdout - std::streambuf* cout_sbuf = NULL; - std::ofstream fout; + // Redirect the stdout channel for the duration of the FNINIT call. + // The replacement is set to ostringstream so that we can find out + // what was sent to it in case of emergency (but would need code + // changes to do this). This used to be done by manually replacing + // the std::cout buffer. + // + // Note that only FNINIT needs this, as the FunctionUtility calls do + // not create string output (at least for 12.12.x/12.13.0). + // + std::ostream* outStream = IosHolder::outHolder(); + std::ostringstream tmpStream; + IosHolder::setStreams(IosHolder::inHolder(), + &tmpStream, + IosHolder::errHolder()); try { - - // Redirect std::cout during call to FNINIT; this causes - // us to swallow annoying messages about having set the - // XSPEC abundance and cross section by default. - - // We have to do it this way because XSPEC code sends message - // to XSPEC stream buffers, that in turn talk to std::cin, - // std::cout, and std::cerr. When we are not given a toggle - // to turn off output to output streams, we can instead redirect - // the stream itself. - - // In this case, the only stream we want to redirect is stdout. - - // If creating any stream buffer fails for any reason, - // make sure we can still use original std::cout stream buffer. - - // Note: since we are *not* redirecting std::cerr, any error - // messages are still seen; if, for example, HEADAS is not set, - // user will still see error message printed to stderr (i.e., to - // screen by default). We still want such error messages to be - // seen! So do *not* redirect std::cerr. - // - // Unfortunately it appears that XSPEC does not use stderr for its - // error messages, but stdout, so the following code will hide - // any error messages from FNINIT (XSPEC versions 12.8.2 and 12.9.0). - // Perhaps the screen output could be saved, rather than redirected - // to /dev/null, and checked for the string '***Error: ', which - // appears to be used in the error messages from XSPEC. - - cout_sbuf = std::cout.rdbuf(); - fout.open("/dev/null"); - if (cout_sbuf != NULL && fout.is_open()) - std::cout.rdbuf(fout.rdbuf()); // temporary redirect stdout to /dev/null - // Initialize XSPEC model library FNINIT(); - // Get back original std::cout - if (cout_sbuf != NULL) { - std::cout.clear(); - std::cout.rdbuf(cout_sbuf); - } - fout.clear(); - fout.close(); - - // We used to set the chatter to 0 but this misses useful info - // (like XSPEC can't find the data files) that would be reported - // by XSPEC, so use the default XSPEC setting. It does not appear - // that this value is read in from ~/.xspec/Xspec.init so set - // it here so we have repeatable behavior. - // - FPCHAT( 10 ); - - // Set cosmology initial values to XSPEC initial values - csmph0( 70.0 ); - csmpq0( 0.0 ); - csmpl0( 0.73 ); - } catch(...) { - - // Get back original std::cout - if (cout_sbuf != NULL) { - std::cout.clear(); - std::cout.rdbuf(cout_sbuf); - } - fout.clear(); - fout.close(); + IosHolder::setStreams(IosHolder::inHolder(), + outStream, + IosHolder::errHolder()); // Raise appropriate error message that XSPEC initialization failed. PyErr_SetString( PyExc_ImportError, (char*)"XSPEC initialization failed; " "check HEADAS environment variable" ); return EXIT_FAILURE; - } - init = true; - return EXIT_SUCCESS; - -} - -static PyObject* get_version( PyObject *self ) -{ - - char version[256]; - int retval; - - try { + IosHolder::setStreams(IosHolder::inHolder(), + outStream, + IosHolder::errHolder()); - retval = xs_getVersion(version, 256); + // We used to set the chatter to 0 but this misses useful info + // (like XSPEC can't find the data files) that would be reported + // by XSPEC, so use the default XSPEC setting. It does not appear + // that this value is read in from ~/.xspec/Xspec.init so set + // it here so we have repeatable behavior. + // + FunctionUtility::xwriteChatter( 10 ); - } catch(...) { + // Set cosmology initial values to XSPEC initial values + FunctionUtility::setH0( 70.0 ); + FunctionUtility::setq0( 0.0 ); + FunctionUtility::setlambda0( 0.73 ); - PyErr_SetString( PyExc_LookupError, - (char*)"could not get XSPEC version string" ); - return NULL; - - } - - if( retval < 0 ) { - PyErr_SetString( PyExc_LookupError, - (char*)"XSPEC version string was truncated" ); - return NULL; - } - - return Py_BuildValue( (char*)"s", version ); + init = true; + return EXIT_SUCCESS; } static PyObject* get_chatter( PyObject *self ) { - - int chatter = 0; - - try { - - chatter = FGCHAT(); - - } catch(...) { - - PyErr_SetString( PyExc_LookupError, - (char*)"could not get XSPEC chatter level" ); - return NULL; - - } - - return Py_BuildValue( (char*)"i", chatter ); - + return Py_BuildValue( (char*)"i", FunctionUtility::xwriteChatter() ); } +// TODO: we could send in an integer for the Z number (ie either name or number) static PyObject* get_abund( PyObject *self, PyObject *args ) { - char* abund = NULL; char* element = NULL; - PyObject *retval = NULL; - if ( !PyArg_ParseTuple( args, (char*)"|s", &element ) ) return NULL; - try { - - abund = FGSOLR(); - - } catch(...) { - - PyErr_SetString( PyExc_LookupError, - (char*)"could not get XSPEC solar abundance" ); - return NULL; - + // Not asked for an element so return the table name + if ( !element ) { + return (PyObject*) Py_BuildValue( (char*)"s", FunctionUtility::ABUND().c_str() ); } - if( !element ) { - - retval = (PyObject*) Py_BuildValue( (char*)"s", abund ); - - } else { - - float abundVal = 0.0; - std::streambuf *cerr_sbuf = NULL; - std::ostringstream fout; - - try { - - cerr_sbuf = std::cerr.rdbuf(); - - if (cerr_sbuf != NULL) - std::cerr.rdbuf(fout.rdbuf()); - - abundVal = FGABND(element); - - // Get back original std::cerr - if (cerr_sbuf != NULL) - std::cerr.rdbuf(cerr_sbuf); - - - } catch(...) { - - // Get back original std::cerr - if (cerr_sbuf != NULL) - std::cerr.rdbuf(cerr_sbuf); - - PyErr_Format( PyExc_ValueError, - (char*)"could not get XSPEC abundance for '%s'", - element); - return NULL; - - } - - if( fout.str().size() > 0 ) { - PyErr_Format( PyExc_TypeError, - (char*)"could not find element '%s'", element); - return NULL; - } - - retval = (PyObject*) Py_BuildValue( (char*)"f", abundVal ); + // Get the specific abundance. Unfortunately getAbundance reports an + // error to stderr when an invalid element is used, so we need to + // hide this. + // + std::ostream* errStream = IosHolder::errHolder(); + std::ostringstream tmpStream; + IosHolder::setStreams(IosHolder::inHolder(), + IosHolder::outHolder(), + &tmpStream); + + float abundVal = FunctionUtility::getAbundance(string(element)); + + IosHolder::setStreams(IosHolder::inHolder(), + IosHolder::outHolder(), + errStream); + + // Was there an error? + // + if( tmpStream.str().size() > 0 ) { + return PyErr_Format( PyExc_TypeError, // TODO: change from TypeError to ValueError? + (char*)"could not find element '%s'", element); } - return retval; + return (PyObject*) Py_BuildValue( (char*)"f", abundVal ); + } static PyObject* get_cosmo( PyObject *self ) { - - float h0; - float l0; - float q0; - - try { - - h0 = csmgh0(); - l0 = csmgl0(); - q0 = csmgq0(); - - } catch(...) { - - PyErr_SetString( PyExc_LookupError, - (char*)"could not get XSPEC cosmology settings" ); - return NULL; - - } + // Assume these can not throw errors + float h0 = FunctionUtility::getH0(); + float l0 = FunctionUtility::getlambda0(); + float q0 = FunctionUtility::getq0(); return Py_BuildValue( (char*)"fff", h0, q0, l0 ); } -static PyObject* get_cross( PyObject *self ) -{ - - char* cross = NULL; - - try { - - cross = FGXSCT(); - - } catch(...) { - - PyErr_SetString( PyExc_LookupError, - (char*)"could not get XSPEC cross-section" ); - return NULL; - - } - - return Py_BuildValue( (char*)"s", cross ); - -} - - static PyObject* set_chatter( PyObject *self, PyObject *args ) { @@ -397,92 +182,70 @@ static PyObject* set_chatter( PyObject *self, PyObject *args ) if ( !PyArg_ParseTuple( args, (char*)"i", &chatter ) ) return NULL; - try { - - FPCHAT( chatter ); - - } catch(...) { - - PyErr_Format( PyExc_ValueError, - (char*)"could not set XSPEC chatter level to %d", - chatter); - return NULL; - - } - + FunctionUtility::xwriteChatter(chatter); Py_RETURN_NONE; } +// Based on xsFortran::FPSOLR +// +// TODO: add a version where we can send in an array of numbers +// static PyObject* set_abund( PyObject *self, PyObject *args ) { char* table = NULL; - int status = 0; - if ( !PyArg_ParseTuple( args, (char*)"s", &table ) ) return NULL; - try { - - FPSOLR( table, &status ); - - } catch(...) { - - status = 1; + string tableName = string(table); + tableName = XSutility::lowerCase(tableName); + if (tableName == "file") { + FunctionUtility::ABUND(tableName); + Py_RETURN_NONE; } - // if abundance table name fails, try it as a filename - if( status ) { - - std::ifstream fileStream(table); - std::vector vals(ABUND_SIZE, 0); - size_t count(0); - - try { - - float element; - fileStream.exceptions(std::ios_base::failbit); + if (FunctionUtility::checkAbund(tableName)) { + FunctionUtility::ABUND(tableName); + Py_RETURN_NONE; + } - while (count < ABUND_SIZE && fileStream >> element) { - vals[count] = element; - ++count; - } + // If we've got here then try to read the data from a file + // + const size_t nelems = FunctionUtility::NELEMS(); + std::vector vals(nelems, 0); + size_t count(0); - status = 0; - } - catch ( std::exception& ) { + std::ifstream fileStream(table); - if( !fileStream.eof() ) { - PyErr_Format( PyExc_ValueError, - (char*)"Cannot read file '%s'. It may not exist or contains invalid data", - table); - return NULL; - } + try { + float element; + fileStream.exceptions(std::ios_base::failbit); - status = 1; + while (count < nelems && fileStream >> element) { + vals[count] = element; + ++count; } + } + catch ( std::exception& ) { - try { - - FPSOLR((char*)"file", &status); - FPSLFL( &vals[0], ABUND_SIZE, &status ); - - } catch(...) { - - status = 1; - + // We do not error out if it was an eofbit failure, as we assume + // this just means that the file contains < NELEMS elements, and + // the missing elements will be set to 0 abundance. We could try + // and explicitly handle this case from the exception, but let's + // follow XSPEC and check the eof status. + // + if (!fileStream.eof()) { + return PyErr_Format( PyExc_ValueError, + (char*)"Cannot read file '%s'. It may not exist or contains invalid data", + table); } } - if ( 0 != status ) { - PyErr_Format( PyExc_ValueError, - (char*)"could not set XSPEC abundance to %s", - table ); - return NULL; - } + FunctionUtility::ABUND("file"); + FunctionUtility::abundanceVectors("file", vals); Py_RETURN_NONE; @@ -499,22 +262,9 @@ static PyObject* set_cosmo( PyObject *self, PyObject *args ) if ( !PyArg_ParseTuple( args, (char*)"fff", &h0, &q0, &l0 ) ) return NULL; - try { - - csmph0( h0 ); - csmpl0( l0 ); - csmpq0( q0 ); - - } catch(...) { - - PyErr_Format( PyExc_ValueError, - (char*)"could not set XSPEC cosmology settings to " - "H0=%g q0=%g Lambda0=%g", - h0, l0, q0); - return NULL; - - } - + FunctionUtility::setH0(h0); + FunctionUtility::setq0(q0); + FunctionUtility::setlambda0(l0); Py_RETURN_NONE; } @@ -524,61 +274,56 @@ static PyObject* set_cross( PyObject *self, PyObject *args ) { char* csection = NULL; - int status = 0; if ( !PyArg_ParseTuple( args, (char*)"s", &csection ) ) return NULL; - try { - - FPXSCT( csection, &status ); - - } catch(...) { - - status = 1; - - } - - if ( 0 != status ) { - PyErr_Format( PyExc_ValueError, - (char*)"could not set XSPEC photoelectric " - "cross-section to '%s'", - csection); - return NULL; + string tableName = string(csection); + tableName = XSutility::lowerCase(tableName); + + // On failure, checkXsect catches a YellowAlert which creates output + // to the error channel, so we over-ride it for this call. + // + std::ostream* errStream = IosHolder::errHolder(); + std::ostringstream tmpStream; + IosHolder::setStreams(IosHolder::inHolder(), + IosHolder::outHolder(), + &tmpStream); + const bool known = FunctionUtility::checkXsect(tableName); + IosHolder::setStreams(IosHolder::inHolder(), + IosHolder::outHolder(), + errStream); + + if (!known) { + return PyErr_Format( PyExc_ValueError, + (char*)"could not set XSPEC photoelectric " + "cross-section to '%s'", + csection); } + FunctionUtility::XSECT(tableName); Py_RETURN_NONE; } +// TODO: We could have a seperate "reset" command +// static PyObject* set_xset( PyObject *self, PyObject *args ) { char* str_name = NULL; char* str_value = NULL; - int status = 0; if ( !PyArg_ParseTuple( args, (char*)"ss", &str_name, &str_value ) ) return NULL; - try { - - FPMSTR( str_name, str_value ); - - } catch(...) { - - status = 1; - - } - - if ( 0 != status ) { - PyErr_Format( PyExc_ValueError, - (char*)"could not set XSPEC model strings '%s: %s'", - str_name, str_value); - return NULL; + string name = XSutility::upperCase(string(str_name)); + if (name == "INITIALIZE") { + FunctionUtility::eraseModelStringDataBase(); + } else { + FunctionUtility::setModelString(name, string(str_value)); } - Py_RETURN_NONE; } @@ -587,102 +332,58 @@ static PyObject* get_xset( PyObject *self, PyObject *args ) { char* str_name = NULL; - char* str_value = NULL; if ( !PyArg_ParseTuple( args, (char*)"s", &str_name ) ) return NULL; - try { - - str_value = FGMSTR( str_name ); - - } catch(...) { - - PyErr_Format( PyExc_KeyError, - (char*)"could not get XSPEC model string '%s'", - str_name); - return NULL; - + static string value; + value = FunctionUtility::getModelString(string(str_name)); + if (value == FunctionUtility::NOT_A_KEY()) { + value.erase(); } - return Py_BuildValue( (char*)"s", str_value ); + return Py_BuildValue( (char*)"s", value.c_str() ); } -// Perhaps this should be expanded to some of the other routines -// that return a string, rather than just the paths? -// -static PyObject* get_xspec_path( const char *label, char *getfunc() ) -{ - - char* str_value = NULL; - try { - str_value = getfunc(); - } catch(...) { - - std::ostringstream emsg; - emsg << "could not get XSPEC " << label << " path"; - PyErr_SetString( PyExc_LookupError, - emsg.str().c_str() ); - return NULL; - - } - - return Py_BuildValue( (char*)"s", str_value ); - -} - -static PyObject* get_manager_data_path( PyObject *self ) -{ - return get_xspec_path("manager", FGDATD); -} - -static PyObject* get_model_data_path( PyObject *self ) -{ - return get_xspec_path("model", FGMODF); +template +static PyObject* get_xspec_string( PyObject *self ) { + return Py_BuildValue( (char*)"s", get().c_str() ); } static PyObject* set_manager_data_path( PyObject *self, PyObject *args ) { char* path = NULL; - if ( !PyArg_ParseTuple( args, (char*)"s", &path ) ) return NULL; - try { - - FPDATD( path ); - - } catch(...) { - - std::ostringstream emsg; - emsg << "could not set XSPEC manager path to '" << path << "'"; - PyErr_SetString( PyExc_ValueError, - emsg.str().c_str() ); - return NULL; - } - + FunctionUtility::managerPath(string(path)); Py_RETURN_NONE; } +#define NOARGSPEC(name, func) \ + { (char *)#name, (PyCFunction)func, METH_NOARGS, NULL } + + static PyMethodDef XSpecMethods[] = { - { (char*)"get_xsversion", (PyCFunction)get_version, METH_NOARGS, NULL }, - { (char*)"get_xschatter", (PyCFunction)get_chatter, METH_NOARGS, NULL }, + // Utility routines + // + NOARGSPEC(get_xsversion, get_xspec_string), + NOARGSPEC(get_xschatter, get_chatter), FCTSPEC(set_xschatter, set_chatter), FCTSPEC(get_xsabund, get_abund), FCTSPEC(set_xsabund, set_abund), FCTSPEC(set_xscosmo, set_cosmo), - { (char*)"get_xscosmo", (PyCFunction)get_cosmo, METH_NOARGS, NULL }, - { (char*)"get_xsxsect", (PyCFunction)get_cross, METH_NOARGS, NULL }, + NOARGSPEC(get_xscosmo, get_cosmo), + NOARGSPEC(get_xsxsect, get_xspec_string), + FCTSPEC(set_xsxsect, set_cross), FCTSPEC(set_xsxset, set_xset), FCTSPEC(get_xsxset, get_xset), - { (char*)"get_xspath_manager", - (PyCFunction)get_manager_data_path, METH_NOARGS, NULL }, - { (char*)"get_xspath_model", - (PyCFunction)get_model_data_path, METH_NOARGS, NULL }, + NOARGSPEC(get_xspath_manager, get_xspec_string), + NOARGSPEC(get_xspath_model, get_xspec_string), FCTSPEC(set_xspath_manager, set_manager_data_path), // Start model definitions diff --git a/sherpa/include/sherpa/astro/xspec_extension.hh b/sherpa/include/sherpa/astro/xspec_extension.hh index 8ce82d5ace..d19457940e 100644 --- a/sherpa/include/sherpa/astro/xspec_extension.hh +++ b/sherpa/include/sherpa/astro/xspec_extension.hh @@ -989,7 +989,7 @@ PyObject* xspecmodelfct_con_f77( PyObject* self, PyObject* args, PyObject* kwarg // As of XSPEC 12.10.1, the table-model routines have been // consolidated into one routine, so there is no need for -// a template. A templace could be used to allow compile-time +// a template. A template could be used to allow compile-time // specialization over additive versus multiplicative, but // for now have a run-time check rather than multiple versions // of this routine. From 4ab8f9d894953bedaa179784c2c98eedc647afa5 Mon Sep 17 00:00:00 2001 From: Doug Burke Date: Fri, 9 Jun 2023 08:37:40 -0400 Subject: [PATCH 3/6] xspec: remove support for old compilers Remove an ifdef check that allowed code to run with old compilers (that didn't support the 1997 C++ standard). This lets us slightly simplify the macro as we do not need to send in the length. --- .../include/sherpa/astro/xspec_extension.hh | 19 +++++-------------- 1 file changed, 5 insertions(+), 14 deletions(-) diff --git a/sherpa/include/sherpa/astro/xspec_extension.hh b/sherpa/include/sherpa/astro/xspec_extension.hh index d19457940e..725b1b8f8c 100644 --- a/sherpa/include/sherpa/astro/xspec_extension.hh +++ b/sherpa/include/sherpa/astro/xspec_extension.hh @@ -79,19 +79,10 @@ namespace sherpa { namespace astro { namespace xspec { typedef sherpa::Array< float, NPY_FLOAT > FloatArray; typedef float FloatArrayType; -// Try and support the use of std::transform while still building -// against C++-98 compilers. -// -#if __cplusplus > 199711L -#define CONVERTARRAY(orig, out, npts) \ +// Assume std::transform is available on any system we support +#define CONVERTARRAY(orig, out) \ std::transform(std::begin(orig), std::end(orig), std::begin(out), \ [](const double val) -> FloatArrayType { return static_cast(val); }); -#else -#define CONVERTARRAY(orig, out, npts) \ - for (int i = 0; i < npts; i++) { \ - out[i] = static_cast(orig[i]); \ - } -#endif // XSpec models can be called from Sherpa using either @@ -634,7 +625,7 @@ static void create_output(int nbins, T &a, T &b) { void call_xspec( RealArray& result ) { // convert to 32-byte float std::vector fear(this->ngrid); - CONVERTARRAY(this->ear, fear, this->ngrid); + CONVERTARRAY(this->ear, fear); XSpecFunc( &fear[0], this->npts, &this->pars[0], this->ifl, &result[0], &this->error[0] ); return; @@ -751,7 +742,7 @@ static void create_output(int nbins, T &a, T &b) { void call_xspec( RealArray& result ) { // convert to 32-byte float std::vector fear(this->ngrid); - CONVERTARRAY(this->ear, fear, this->ngrid); + CONVERTARRAY(this->ear, fear); XSpecFunc( &fear[0], this->npts, &this->pars[0], this->ifl, &result[0], &this->error[0] ); return; @@ -787,7 +778,7 @@ static void create_output(int nbins, T &a, T &b) { // convert to 32-byte float std::vector fear(ngrid); - CONVERTARRAY(ear, fear, ngrid); + CONVERTARRAY(ear, fear); // Number of bins to send to XSPEC nout = ngrid; From 762f9203e424054841f9ac428ead1d9f7ceaf531 Mon Sep 17 00:00:00 2001 From: Douglas Burke Date: Wed, 28 Aug 2024 09:17:04 -0400 Subject: [PATCH 4/6] address review comments Improve the documentation on the initialization of the XSPEC module. --- sherpa/astro/xspec/src/_xspec.cc | 28 ++++++++++++---------------- 1 file changed, 12 insertions(+), 16 deletions(-) diff --git a/sherpa/astro/xspec/src/_xspec.cc b/sherpa/astro/xspec/src/_xspec.cc index e6e84c5390..deecd186d1 100644 --- a/sherpa/astro/xspec/src/_xspec.cc +++ b/sherpa/astro/xspec/src/_xspec.cc @@ -23,6 +23,9 @@ #include #include +// Documentation for the "XSPEC internal functions" is at +// https://heasarc.gsfc.nasa.gov/docs/xanadu/xspec/internal/XspecInternalFunctionsGuide.html +// // xsfortran is only needed to support FNINIT as there's (as of XSPEC 12.13.0) // no version of this functionality in FunctionUtility. // @@ -63,13 +66,6 @@ static int _sherpa_init_xspec_library() } // Redirect the stdout channel for the duration of the FNINIT call. - // The replacement is set to ostringstream so that we can find out - // what was sent to it in case of emergency (but would need code - // changes to do this). This used to be done by manually replacing - // the std::cout buffer. - // - // Note that only FNINIT needs this, as the FunctionUtility calls do - // not create string output (at least for 12.12.x/12.13.0). // std::ostream* outStream = IosHolder::outHolder(); std::ostringstream tmpStream; @@ -78,7 +74,7 @@ static int _sherpa_init_xspec_library() IosHolder::errHolder()); try { - // Initialize XSPEC model library + // Initialize XSPEC model library. FNINIT(); } catch(...) { @@ -86,7 +82,11 @@ static int _sherpa_init_xspec_library() outStream, IosHolder::errHolder()); - // Raise appropriate error message that XSPEC initialization failed. + // The contents of tmpStream could be inspected to see if it + // contains useful information for the user, but at this point of + // the initialization it is not obvious that it would provide any + // extra information. + // PyErr_SetString( PyExc_ImportError, (char*)"XSPEC initialization failed; " "check HEADAS environment variable" ); @@ -97,15 +97,11 @@ static int _sherpa_init_xspec_library() outStream, IosHolder::errHolder()); - // We used to set the chatter to 0 but this misses useful info - // (like XSPEC can't find the data files) that would be reported - // by XSPEC, so use the default XSPEC setting. It does not appear - // that this value is read in from ~/.xspec/Xspec.init so set - // it here so we have repeatable behavior. + // Set a number of values to their XSPEC defaults (as of XSPEC + // 12.14.1, but they have not changed for a long time). It appears + // these are not read in from the user's ~/.xspec/Xspec.init file. // FunctionUtility::xwriteChatter( 10 ); - - // Set cosmology initial values to XSPEC initial values FunctionUtility::setH0( 70.0 ); FunctionUtility::setq0( 0.0 ); FunctionUtility::setlambda0( 0.73 ); From cfb1ea9204647791423ec19074bc49ec53d54325 Mon Sep 17 00:00:00 2001 From: Douglas Burke Date: Thu, 29 Aug 2024 14:12:58 -0400 Subject: [PATCH 5/6] Revert "xspec: remove support for old compilers" This reverts commit 4ab8f9d894953bedaa179784c2c98eedc647afa5. --- .../include/sherpa/astro/xspec_extension.hh | 19 ++++++++++++++----- 1 file changed, 14 insertions(+), 5 deletions(-) diff --git a/sherpa/include/sherpa/astro/xspec_extension.hh b/sherpa/include/sherpa/astro/xspec_extension.hh index 725b1b8f8c..d19457940e 100644 --- a/sherpa/include/sherpa/astro/xspec_extension.hh +++ b/sherpa/include/sherpa/astro/xspec_extension.hh @@ -79,10 +79,19 @@ namespace sherpa { namespace astro { namespace xspec { typedef sherpa::Array< float, NPY_FLOAT > FloatArray; typedef float FloatArrayType; -// Assume std::transform is available on any system we support -#define CONVERTARRAY(orig, out) \ +// Try and support the use of std::transform while still building +// against C++-98 compilers. +// +#if __cplusplus > 199711L +#define CONVERTARRAY(orig, out, npts) \ std::transform(std::begin(orig), std::end(orig), std::begin(out), \ [](const double val) -> FloatArrayType { return static_cast(val); }); +#else +#define CONVERTARRAY(orig, out, npts) \ + for (int i = 0; i < npts; i++) { \ + out[i] = static_cast(orig[i]); \ + } +#endif // XSpec models can be called from Sherpa using either @@ -625,7 +634,7 @@ static void create_output(int nbins, T &a, T &b) { void call_xspec( RealArray& result ) { // convert to 32-byte float std::vector fear(this->ngrid); - CONVERTARRAY(this->ear, fear); + CONVERTARRAY(this->ear, fear, this->ngrid); XSpecFunc( &fear[0], this->npts, &this->pars[0], this->ifl, &result[0], &this->error[0] ); return; @@ -742,7 +751,7 @@ static void create_output(int nbins, T &a, T &b) { void call_xspec( RealArray& result ) { // convert to 32-byte float std::vector fear(this->ngrid); - CONVERTARRAY(this->ear, fear); + CONVERTARRAY(this->ear, fear, this->ngrid); XSpecFunc( &fear[0], this->npts, &this->pars[0], this->ifl, &result[0], &this->error[0] ); return; @@ -778,7 +787,7 @@ static void create_output(int nbins, T &a, T &b) { // convert to 32-byte float std::vector fear(ngrid); - CONVERTARRAY(ear, fear); + CONVERTARRAY(ear, fear, ngrid); // Number of bins to send to XSPEC nout = ngrid; From e2f11f32979502ec1ac9bc8426e6bd0fed98d276 Mon Sep 17 00:00:00 2001 From: Douglas Burke Date: Thu, 29 Aug 2024 14:13:28 -0400 Subject: [PATCH 6/6] XSPEC: update copyright year This code was written a long time ago. --- sherpa/include/sherpa/astro/xspec_extension.hh | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/sherpa/include/sherpa/astro/xspec_extension.hh b/sherpa/include/sherpa/astro/xspec_extension.hh index d19457940e..f17e5cd021 100644 --- a/sherpa/include/sherpa/astro/xspec_extension.hh +++ b/sherpa/include/sherpa/astro/xspec_extension.hh @@ -1,5 +1,5 @@ // -// Copyright (C) 2009, 2015, 2017, 2020, 2021, 2022, 2023 +// Copyright (C) 2009, 2015, 2017, 2020 - 2024 // Smithsonian Astrophysical Observatory // //