libatomprobe
Library for Atom Probe Tomography (APT) computation
|
Class to load abundance information for natural isotopes. More...
#include <abundance.h>
Public Member Functions | |
size_t | open (const char *file, bool strict=false) |
Attempt to open the abundance data file, return 0 on success, nonzero on failure. More... | |
size_t | numIsotopes () const |
Return the number of isotopes in the loaded table. More... | |
size_t | numElements () const |
Return the number of elements in the loaded table. More... | |
size_t | symbolIndex (const char *symbol, bool caseSensitive=true) const |
Return the element's position in table, starting from 0. More... | |
void | getSymbolIndices (const std::vector< std::string > &symbols, std::vector< size_t > &indices) const |
Return a vector of symbol indices. More... | |
void | getSymbolIndices (const std::vector< std::pair< std::string, size_t > > &symbols, std::vector< std::pair< size_t, unsigned int > > &indices) const |
Obtain the symbol indicies for a string,payload type. Resultant pair has first item in index, second is payload. More... | |
std::string | elementName (size_t elemIdx) const |
Obtain the name of an element from its index. More... | |
size_t | symbolIdxFromAtomicNumber (size_t atomicNumber) const |
Return the symbol index (offset in vector) for the given atomic number, -1 if it cannot be found. More... | |
std::string | elementNames (size_t start, size_t end, char separator=',') const |
Obtain a delimited list of element names, using the (optional) supplied separator, over the specified [start,end) range. More... | |
void | isotopeIndex (size_t elem, float mass, size_t &isotopeIdx) const |
Obtain the atom ID (first) and the isotope value(second). More... | |
void | isotopeIndex (size_t elem, float mass, float tolerance, size_t &isotopeIdx) const |
Return the first matching element in the table which has a given tolerance to the specified match. More... | |
unsigned int | getAtomicNumber (size_t elemIdx) const |
Obtain the atomic number for the given element, by element index. More... | |
const std::vector< ISOTOPE_ENTRY > & | isotopes (size_t offset) const |
Return the isotope at a given position in the table. Offset must exist. More... | |
void | generateIsotopeDist (const std::vector< size_t > &elementIdx, const std::vector< size_t > &frequency, std::vector< std::pair< float, float > > &massDist) const |
Compute the mass-probability distribution for the grouped set of ions. More... | |
void | generateGroupedIsotopeDist (const std::vector< size_t > &elementIdx, const std::vector< size_t > &frequency, std::vector< std::pair< float, float > > &massDist, float massTolerance) const |
As per generateIsotopeDist, however, this convenience groups the distribution to limit the effect of minute mass changes. More... | |
void | generateSingleAtomDist (size_t atomIdx, unsigned int repeatCount, std::vector< std::pair< float, float > > &massDist) const |
Obtain the mass distribution from a single isotope that repeats. Output is not grouped. More... | |
float | abundanceBetweenLimits (const std::vector< size_t > &elemIdx, const std::vector< size_t > &frequency, float massStart, float massEnd) const |
Obtain the fractional abundance between two limits for the given species [start,end) More... | |
const ISOTOPE_ENTRY & | isotope (size_t elementIdx, size_t isotopeIdx) const |
Obtain a reference to a particular isotope, using the element's index in the table, and the isotope index. More... | |
size_t | getMajorIsotopeFromElemIdx (size_t elementIdx) const |
Obtain the most prominent isotope's index from the element index. More... | |
float | getNominalMass (size_t elementIdx) const |
Obtain the isotope weighted mass for the given element. More... | |
size_t | getNearestCharge (const std::vector< std::pair< size_t, size_t > > &molecule, float targetMass, size_t maxCharge) const |
Obtain the closest charge state for the given mass and species group. More... | |
Static Public Member Functions | |
static const char * | getErrorText (size_t errorCode) |
Obtain a human-readable error message from the open( ) call. More... | |
Class to load abundance information for natural isotopes.
Definition at line 54 of file abundance.h.
float AtomProbe::AbundanceData::abundanceBetweenLimits | ( | const std::vector< size_t > & | elemIdx, |
const std::vector< size_t > & | frequency, | ||
float | massStart, | ||
float | massEnd | ||
) | const |
Obtain the fractional abundance between two limits for the given species [start,end)
Definition at line 463 of file abundance.cpp.
References ASSERT, generateIsotopeDist(), and isotope().
Referenced by AtomProbe::leastSquaresOverlapSolve().
|
inline |
Obtain the name of an element from its index.
Definition at line 107 of file abundance.h.
References AtomProbe::ISOTOPE_ENTRY::atomicNumber, and AtomProbe::ISOTOPE_ENTRY::mass.
Referenced by elementNames().
std::string AtomProbe::AbundanceData::elementNames | ( | size_t | start, |
size_t | end, | ||
char | separator = ',' |
||
) | const |
Obtain a delimited list of element names, using the (optional) supplied separator, over the specified [start,end) range.
Definition at line 593 of file abundance.cpp.
References ASSERT, elementName(), and symbolIdxFromAtomicNumber().
void AtomProbe::AbundanceData::generateGroupedIsotopeDist | ( | const std::vector< size_t > & | elementIdx, |
const std::vector< size_t > & | frequency, | ||
std::vector< std::pair< float, float > > & | massDist, | ||
float | massTolerance | ||
) | const |
As per generateIsotopeDist, however, this convenience groups the distribution to limit the effect of minute mass changes.
This groups the masses of nearby isotope combinations together. This is desirable due to species, such as Mo2, which have isotopes within <0.001 amu of one another. The mass centres are computed using a weighted average, and the intensities are summed
Definition at line 395 of file abundance.cpp.
References generateIsotopeDist(), and AtomProbe::weightedMean().
Referenced by AtomProbe::findOverlaps(), getMassDistributions(), getNearestCharge(), and AtomProbe::maxExplainedFraction().
void AtomProbe::AbundanceData::generateIsotopeDist | ( | const std::vector< size_t > & | elementIdx, |
const std::vector< size_t > & | frequency, | ||
std::vector< std::pair< float, float > > & | massDist | ||
) | const |
Compute the mass-probability distribution for the grouped set of ions.
E.g., for C2, this would be elementIdx = elementIdx("C"), and frequency 2. This would return probabilities for each unique mass {C-12 C-12},{C-12 C-13} and {C-13 C-13}.... note an internal "grouping" tolerance is used to group-together very close masses. This may be smaller than what is needed for your application, so if you wish to obtain a more specific distribution use "generateGroupedIsotopeDist"
Definition at line 277 of file abundance.cpp.
References ASSERT, and AtomProbe::vectorMultiErase().
Referenced by abundanceBetweenLimits(), AtomProbe::checkMassRangingCorrectness(), AtomProbe::computeIonDistAdjacency(), AtomProbe::filterPeakNeedBiggerObs(), generateGroupedIsotopeDist(), getNearestCharge(), and AtomProbe::RangeFile::guessChargeState().
void AtomProbe::AbundanceData::generateSingleAtomDist | ( | size_t | atomIdx, |
unsigned int | repeatCount, | ||
std::vector< std::pair< float, float > > & | massDist | ||
) | const |
Obtain the mass distribution from a single isotope that repeats. Output is not grouped.
Definition at line 517 of file abundance.cpp.
Referenced by AtomProbe::MultiRange::copyDataFromRange().
unsigned int AtomProbe::AbundanceData::getAtomicNumber | ( | size_t | elemIdx | ) | const |
Obtain the atomic number for the given element, by element index.
Definition at line 503 of file abundance.cpp.
References ASSERT.
Referenced by AtomProbe::computeRangeAdjacency(), AtomProbe::findOverlaps(), and AtomProbe::MultiRange::MultiRange().
|
static |
Obtain a human-readable error message from the open( ) call.
Definition at line 59 of file abundance.cpp.
References ASSERT.
size_t AtomProbe::AbundanceData::getMajorIsotopeFromElemIdx | ( | size_t | elementIdx | ) | const |
Obtain the most prominent isotope's index from the element index.
Definition at line 482 of file abundance.cpp.
References ASSERT.
Referenced by AtomProbe::maxExplainedFraction().
size_t AtomProbe::AbundanceData::getNearestCharge | ( | const std::vector< std::pair< size_t, size_t > > & | molecule, |
float | targetMass, | ||
size_t | maxCharge | ||
) | const |
Obtain the closest charge state for the given mass and species group.
In the "molecule" parameter, the pair is <element index, count>, the target mass is given in Da, and the maximum charge is the maximum the function will search for
Definition at line 642 of file abundance.cpp.
References ASSERT, generateGroupedIsotopeDist(), generateIsotopeDist(), open(), symbolIndex(), and TEST.
float AtomProbe::AbundanceData::getNominalMass | ( | size_t | elementIdx | ) | const |
Obtain the isotope weighted mass for the given element.
Definition at line 626 of file abundance.cpp.
Referenced by AtomProbe::convertMassToMol(), and AtomProbe::convertMolToMass().
void AtomProbe::AbundanceData::getSymbolIndices | ( | const std::vector< std::string > & | symbols, |
std::vector< size_t > & | indices | ||
) | const |
Return a vector of symbol indices.
Definition at line 559 of file abundance.cpp.
References symbolIndex().
Referenced by AtomProbe::MultiRange::copyDataFromRange(), AtomProbe::findOverlaps(), getMassDistributions(), and AtomProbe::MultiRange::MultiRange().
void AtomProbe::AbundanceData::getSymbolIndices | ( | const std::vector< std::pair< std::string, size_t > > & | symbols, |
std::vector< std::pair< size_t, unsigned int > > & | indices | ||
) | const |
Obtain the symbol indicies for a string,payload type. Resultant pair has first item in index, second is payload.
const ISOTOPE_ENTRY & AtomProbe::AbundanceData::isotope | ( | size_t | elementIdx, |
size_t | isotopeIdx | ||
) | const |
Obtain a reference to a particular isotope, using the element's index in the table, and the isotope index.
Definition at line 271 of file abundance.cpp.
References ASSERT.
Referenced by abundanceBetweenLimits(), AtomProbe::MultiRange::copyDataFromRange(), AtomProbe::filterBySolutionPPM(), and AtomProbe::maxExplainedFraction().
void AtomProbe::AbundanceData::isotopeIndex | ( | size_t | elem, |
float | mass, | ||
size_t & | isotopeIdx | ||
) | const |
Obtain the atom ID (first) and the isotope value(second).
The mass and element ID must match exactly NOTE: The element and isotope indices are NOT its isotopic number it is the offset to that isotope int the abundance table
Definition at line 238 of file abundance.cpp.
References ASSERT.
void AtomProbe::AbundanceData::isotopeIndex | ( | size_t | elem, |
float | mass, | ||
float | tolerance, | ||
size_t & | isotopeIdx | ||
) | const |
Return the first matching element in the table which has a given tolerance to the specified match.
Definition at line 254 of file abundance.cpp.
References ASSERT.
const std::vector< ISOTOPE_ENTRY > & AtomProbe::AbundanceData::isotopes | ( | size_t | offset | ) | const |
Return the isotope at a given position in the table. Offset must exist.
Definition at line 621 of file abundance.cpp.
References ASSERT.
Referenced by AtomProbe::leastSquaresOverlapSolve().
size_t AtomProbe::AbundanceData::numElements | ( | ) | const |
Return the number of elements in the loaded table.
Definition at line 74 of file abundance.cpp.
size_t AtomProbe::AbundanceData::numIsotopes | ( | ) | const |
Return the number of isotopes in the loaded table.
Definition at line 65 of file abundance.cpp.
size_t AtomProbe::AbundanceData::open | ( | const char * | file, |
bool | strict = false |
||
) |
Attempt to open the abundance data file, return 0 on success, nonzero on failure.
If nonzero, a human readable message can be obtained from getErrorText(...)
Definition at line 79 of file abundance.cpp.
References AtomProbe::ISOTOPE_ENTRY::abundance, AtomProbe::ISOTOPE_ENTRY::abundanceError, AtomProbe::ISOTOPE_ENTRY::atomicNumber, AtomProbe::ISOTOPE_ENTRY::mass, AtomProbe::ISOTOPE_ENTRY::massError, AtomProbe::ISOTOPE_ENTRY::massNumber, AtomProbe::ISOTOPE_ENTRY::symbol, AtomProbe::XMLHelpFwdToElem(), and AtomProbe::XMLHelpGetProp().
Referenced by AtomProbe::checkMassRangingCorrectness(), AtomProbe::computeRangeAdjacency(), AtomProbe::findOverlaps(), getNearestCharge(), AtomProbe::RangeFile::guessChargeState(), AtomProbe::leastSquaresOverlapSolve(), main(), AtomProbe::maxExplainedFraction(), AtomProbe::parseCompositionData(), and runTests().
size_t AtomProbe::AbundanceData::symbolIdxFromAtomicNumber | ( | size_t | atomicNumber | ) | const |
Return the symbol index (offset in vector) for the given atomic number, -1 if it cannot be found.
Definition at line 609 of file abundance.cpp.
Referenced by AtomProbe::computeIonDistAdjacency(), elementNames(), and AtomProbe::leastSquaresOverlapSolve().
size_t AtomProbe::AbundanceData::symbolIndex | ( | const char * | symbol, |
bool | caseSensitive = true |
||
) | const |
Return the element's position in table, starting from 0.
Notes
Definition at line 214 of file abundance.cpp.
References ASSERT, and AtomProbe::lowercase().
Referenced by AtomProbe::computeRangeAdjacency(), AtomProbe::MultiRange::copyDataFromRange(), getNearestCharge(), AtomProbe::getRangeMolecule(), getSymbolIndices(), AtomProbe::RangeFile::guessChargeState(), AtomProbe::leastSquaresOverlapSolve(), AtomProbe::maxExplainedFraction(), and AtomProbe::parseCompositionData().