IonBalNei

class IonBalNei{
 public:

  vector<IonBalTemperatureRecord> TemperatureRecord;
  RealArray Temperature;
  string Version;

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

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

  void Clear();  // clear out the contents

The Version string is used to control which version of the eigenvector data files are to be used. If the version is changed then the object is reset using Clear and setVersion returns true. The setVersion method with no input checks NEIVERS and uses that if it has been set.

  string getVersion();
  bool setVersion();
  bool setVersion(string version);

This class provides a Read method to load data from the eigenvector files. This can either be done using a vector of atomic numbers and the directory name and version string, so that the relevant filenames can be constructed, or using a single atomic number and filename.

  //  Read loads eigenvector data for a set of elements

  int Read(vector<int> Z, string dirname, string versionname);

  //  ReadElement loads eigenvector data for that element

  int ReadElement(int Z, string filename);

  void LoadTemperatureRecord(IonBalTemperatureRecord input);

The following provide methods to extract some useful numbers.

  RealArray Temperatures();     // return tabulated temperatures
  int NumberTemperatures();     // return number of tabulated temperatures
  int NumberElements();         // return number of tabulated elements
  bool ContainsElement(const int& Z); // returns true if data for Z

The CEI method returns the ionization fractions for the specified element for a collisional ionization equilibrium plasma with the given electron temperature.

  RealArray CEI(const Real& Te, const int& Z);

The method to calculate ion fractions is Calc. It is overloaded to give a number of different options. If the initIonFrac array is used then this gives the initial ion fractions for the element Z, if not the element is assumed to start un-ionized.

  // calculate the NEI ion fractions for electron temperature Te, ionization
  // parameter tau and element Z.

  RealArray Calc(const Real& Te, const Real& tau, const int& Z);
  RealArray Calc(const Real& Te, const Real& tau, const int& Z, 
                 const RealArray& initIonFrac);

  // Calculates ionization fractions at electron temperatures
  // Te and a set of ionization parameters tau(i), i=1,..,n,
  // where each tau is given weight(i). Electron temperature is
  // assumed to be linear function of tau.
  // Based on the old noneq.f.

  RealArray Calc(const RealArray& Te, const RealArray& tau, 
                 const RealArray& weight, const int& Z);
  RealArray Calc(const RealArray& Te, const RealArray& tau, 
                 const RealArray& weight, const int& Z, 
                 const RealArray& initIonFrac);

  // Calculates ionization fractions at electron temperature
  // Te(n) and ionization parameter tau(n), for electron 
  // temperatures Te given in a tabular form as a function of 
  // ionization parameter tau.
  // Based on the old noneqr.f.
  // Does initIonFrac make sense in this case.

  RealArray Calc(const RealArray& Te, const RealArray& tau, 
                 const int& Z);
  RealArray Calc(const RealArray& Te, const RealArray& tau, 
                 const int& Z, const RealArray& initIonFrac);

};

Outside the class there are wrap-up functions to read (if necessary) the eigenvector files and calculate ion fractions. They use the NEISPEC xset variable to choose which version of the files to use. The various overloaded versions match to the Calc methods.

void calcNEIfractions(const Real& Te, const Real& tau, const int& Z, 
                      RealArray& IonFrac);
void calcNEIfractions(const Real& Te, const Real& tau, const IntegerVector& Z, 
                      vector<RealArray>& IonFrac);

void calcNEIfractions(const Real& Te, const Real& tau, const int& Z, 
                      const RealArray& initIonFrac, RealArray& IonFrac);
void calcNEIfractions(const Real& Te, const Real& tau, const IntegerVector& Z, 
                      const vector<RealArray>& initIonFrac, 
                      vector<RealArray>& IonFrac);

void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const RealArray& weight, const int& Z, RealArray& IonFrac);
void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const RealArray& weight, const IntegerVector& Z, 
                      vector<RealArray>& IonFrac);

void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const RealArray& weight, const int& Z, 
                      const RealArray& initIonFrac, RealArray& IonFrac);
void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const RealArray& weight, const IntegerVector& Z, 
                      const vector<RealArray>& initIonFrac, 
                      vector<RealArray>& IonFrac);

void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const int& Z, RealArray& IonFrac);
void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const IntegerVector& Z, vector<RealArray>& IonFrac);

void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const int& Z, const RealArray& initIonFrac, 
                      RealArray& IonFrac);
void calcNEIfractions(const RealArray& Te, const RealArray& tau, 
                      const IntegerVector& Z, const RealArray& initIonFrac, 
                      vector<RealArray>& IonFrac);

The wrapper routine to return the collisional ionization equilibrium fractions is

void calcCEIfractions(const Real Te, const IntegerVector& Z, 
                      vector<RealArray>& IonFrac);

There are also handy routines to return and index into the temperatures and the number of temperatures

int getNEItempIndex(const Real& tkeV);
int getNEInumbTemp();

to providing debugging information

string writeIonFrac(const IntegerVector& Zarray, 
                    const vector<RealArray>& IonFrac);
string writeIonFrac(const int& Z, const IntegerVector& Zarray, 
                    const vector<RealArray>& IonFrac);

to do a binary search on a RealArray and return the index of the element immediately less than the input target

int locateIndex(const RealArray& xx, const Real x);

and to check whether arrays are identical (note that the C++11 standard includes this as part of the valarray class but this is not yet implemented in all compilers.

bool identicalArrays(const vector<RealArray>& a, const vector<RealArray>& b);
bool identicalArrays(const RealArray& a, const RealArray& b);
bool identicalArrays(const IntegerVector& a, const IntegerVector& b);


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