Aped class

class Aped{
 public:

  vector<ApedTemperatureRecord> TemperatureRecord;

  // These store the information in the initial PARAMETERS extension

  RealArray Temperatures;
  IntegerVector NelementsLine;
  IntegerVector NelementsCoco;
  IntegerVector Nline;
  IntegerVector Ncont;

  string coconame;
  string linename;

  bool noLines;
  bool noResonanceLines;
  bool thermalBroadening;
  Real velocityBroadening;
  Real minimumLineFluxForBroadening;
  bool broadenPseudoContinuum;
  bool logTempInterpolation;

  Aped();     //default constructor
  ~Aped();    //destructor

  Aped& operator=(const Aped&);    // deep copy

The class provides methods to load the data from the AtomDB continuum and line files. The Read method just gets the information from the PARAMETERS extension, stores the filenames used and sizes the subsidiary objects correctly. To actually load the data use ReadTemperature.

  // Reads the continuum and line files and stores the names and
  // temperatures. Does not read the actual continuum and line data.

  int Read(string cocofilename, string linefilename);

  int ReadTemperature(const int TemperatureIndex);
  int ReadTemperature(const vector<int>& TemperatureIndex);

The noLines, noResonanceLines, thermalBroadening, velocityBroadening, minimumLinefluxForBroadening, broadenPseudoContinuum, logTempInterpolation elements can be set using the following methods. Each of these methods also checks for the corresponding xset variables: APECNOLINES, APECNORES, APECTHERMAL, APECVELOCITY, APECMINFLUX, APECBROADPSEUDO, APECLOGINTERP, respectively. For SetVelocityBroadening a warning will be issued if the input value is non-zero and differs from APECVELOCITY. If noLines is set then the spectrum will only be continuum. At the moment noResonanceLines has no effect but is there for a possible approach to resonance scattering. If thermalBroadening is set then any lines will be thermally broadened while if velocityBroadening is set then any lines will be broadened by the appropriate velocity. The line shape is assumed Gaussian and thermal and velocity broadening are added in quadrature. If minimumLinefluxForBroadening is set then all lines below that flux will not be broadened. If broadenPseudoContinuum is set then the pseudo-continuum will also be broadened. Note that line broadened spectra take longer to calculate.

  void SetNoLines(const bool qno);
  void SetNoResonanceLines(const bool qnores);
  void SetThermalBroadening(const bool qtherm);
  void SetVelocityBroadening(const Real velocity);
  void SetMinimumLinefluxForBroadening(const Real flux);
  void SetBroadenPseudoContinuum(const bool qpbroad);
  void SetLogTempInterpolation(const bool qltemp);

The following provide methods to extract some useful numbers.

  int NumberTemperatures();     // return number of tabulated temperatures
  int NumberElements();         // return number of tabulated elements
  int NumberIons(int Z);        // return number of ions across all temperatures
                                // for element Z.

Whether the data for a given temperature index has been loaded can be checked using IsTemperatureLoaded.

 bool IsTemperatureLoaded(const int TemperatureIndex);

The method to generate spectra for CEI plasmas is SumEqSpectra. This comes in several overloaded options depending on whether a single temperature of distribution of temperatures are required. It is also possible to use a different temperature for thermal broadening of lines than that used for the ionization balance. The energy bins on which the spectrum is to be calculated are specified by standard XSPEC energyArray. The elements required and their abundances are given by Zinput and abundance. Tinput and Dem give the temperature(s) and emission measure(s) required.

  void SumEqSpectra(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Dem, RealArray& fluxArray, 
                    RealArray& fluxErrArray);

  void SumEqSpectra(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Dem, RealArray& fluxArray, 
                    RealArray& fluxErrArray);

  // case where the temperature used for the thermal broadening differs from that
  // for the ionization

  void SumEqSpectra(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput, 
                    const Real& Tbinput, const Real& Dem,
                    RealArray& fluxArray, RealArray& fluxErrArray);

  void SumEqSpectra(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput, 
                    const RealArray& Tbinput, const RealArray& Dem, 
                    RealArray& fluxArray, RealArray& fluxErrArray);

For NEI plasmas, the method to calculate the spectrum is SumNeqSpectra. The inputs differ from SumEqSpectra only in including IonFrac, which specifies the ionization fractions required for each element, for each temperature.

  void SumNeqSpectra(const RealArray& energyArray, 
                     const IntegerVector& Zinput, const RealArray& abundance,
                     const Real Redshift, const Real& Tinput,
                     const vector<RealArray>& IonFrac, 
                     RealArray& fluxArray, RealArray& fluxErrArray);

  void SumNeqSpectra(const RealArray& energyArray, 
                     const IntegerVector& Zinput, const RealArray& abundance,
                     const Real Redshift, const RealArray& Tinput,
                     const vector<vector<RealArray> >& IonFrac, 
                     RealArray& fluxArray, RealArray& fluxErrArray);

  // case where the temperature used for the thermal broadening differs from that
  // for the ionization

  void SumNeqSpectra(const RealArray& energyArray, 
                     const IntegerVector& Zinput, const RealArray& abundance,
                     const Real Redshift, const Real& Tinput,
                     const Real& Tbinput, 
                     const vector<RealArray>& IonFrac, 
                     RealArray& fluxArray, RealArray& fluxErrArray);

  void SumNeqSpectra(const RealArray& energyArray, 
                     const IntegerVector& Zinput, const RealArray& abundance,
                     const Real Redshift, const RealArray& Tinput,
                     const RealArray& Tbinput, 
                     const vector<vector<RealArray> >& IonFrac, 
                     RealArray& fluxArray, RealArray& fluxErrArray);

All versions of SumEqSpectra and SumNeqSpectra operate through

  void SumSpectra(const RealArray& energyArray, 
                  const IntegerVector& Zinput, const RealArray& abundance,
                  const Real Redshift, const RealArray& Tinput,
                  const RealArray& Tbinput, 
                  const vector<vector<RealArray> >& IonFrac, const bool isCEI,
                  RealArray& fluxArray, RealArray& fluxErrArray);

Outside the class there are wrap-up routines which create an internal, static Aped object then call SumSpectra to calculate the output spectrum. There are overloaded versions corresponding to the options for SumEqSpectra and SumNeqSpectra. The int returned will be non-zero in the event of an error on reading the AtomDB files. For the CEI case these routines check the APECROOT variable to find the files to read. For the NEI case they check the NEIVERS and NEIAPECROOT variables. The calcRSSpectrum routines are for the old Raymond-Smith model and read the RS files in the AtomDB format.

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Dem, const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Dem, const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Tbinput, const Real& Dem, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Tbinput, const RealArray& Dem, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Dem, const bool qtherm, const Real velocity,
                    const bool noLines,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Dem, const bool qtherm, const Real velocity,
                    const bool noLines,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Tbinput, const Real& Dem, 
                    const bool qtherm, const Real velocity,
                    const bool noLines,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcCEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Tbinput, const RealArray& Dem, 
                    const bool qtherm, const Real velocity,
                    const bool noLines,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcRSSpectrum(const RealArray& energyArray, 
                   const IntegerVector& Zinput, const RealArray& abundance,
                   const Real Redshift, const Real& Tinput,
                   const Real& Dem, const bool qtherm, const Real velocity,
                   RealArray& fluxArray, RealArray& fluxErrArray);

int calcRSSpectrum(const RealArray& energyArray, 
                   const IntegerVector& Zinput, const RealArray& abundance,
                   const Real Redshift, const RealArray& Tinput,
                   const RealArray& Dem, const bool qtherm, const Real velocity,
                   RealArray& fluxArray, RealArray& fluxErrArray);

int calcRSSpectrum(const RealArray& energyArray, 
                   const IntegerVector& Zinput, const RealArray& abundance,
                   const Real Redshift, const Real& Tinput,
                   const Real& Tbinput, const Real& Dem, 
                   const bool qtherm, const Real velocity,
                   RealArray& fluxArray, RealArray& fluxErrArray);

int calcRSSpectrum(const RealArray& energyArray, 
                   const IntegerVector& Zinput, const RealArray& abundance,
                   const Real Redshift, const RealArray& Tinput,
                   const RealArray& Tbinput, const RealArray& Dem, 
                   const bool qtherm, const Real velocity,
                   RealArray& fluxArray, RealArray& fluxErrArray);
int calcNEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const vector<RealArray>& IonFrac, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcNEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const vector<vector<RealArray> >& IonFrac, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcNEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const Real& Tinput,
                    const Real& Tbinput,
                    const vector<RealArray>& IonFrac, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

int calcNEISpectrum(const RealArray& energyArray, 
                    const IntegerVector& Zinput, const RealArray& abundance,
                    const Real Redshift, const RealArray& Tinput,
                    const RealArray& Tbinput,
                    const vector<vector<RealArray> >& IonFrac, 
                    const bool qtherm, const Real velocity,
                    RealArray& fluxArray, RealArray& fluxErrArray);

A few other useful routines also live in Aped.h. getAtomicMass returns the atomic mass for the given atomic number. getApedFileNames returns the AtomDB version number string, and the continuum and line filenames. getNEIApedFileNames returns the AtomDB version number string, the NEI eigenvectors version number string, and the continuum and line filenames. apedInterpFlux is used to interpolate the continuum and pseudo-continuum arrays.

Real getAtomicMass(const int& AtomicNumber);
bool getApedFileNames(string& version, string& continuumFile, string& lineFile);
bool getNEIApedFileNames(string& version, string& eigenVersion,
                         string& continuumFile, string& lineFile);
void apedInterpFlux(const RealArray& inputEnergy, const RealArray& inputFlux, 
                    const Real& z15, const Real& coeff,
                    const RealArray& energyArray, RealArray& fluxArray);


HEASARC Home | Observatories | Archive | Calibration | Software | Tools | Students/Teachers/Public

Last modified: Wednesday, 23-Mar-2022 17:16:49 EDT

HEASARC Home | Observatories | Archive | Calibration | Software | Tools | Students/Teachers/Public

Last modified: 23-Mar-2022