abundance.h

Go to the documentation of this file.

/* 
 * Copyright (C) 2015  D Haley
 * 
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 * 
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 * 
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

#ifndef ATOMPROBE_ABUNDANCE_H
#define ATOMPROBE_ABUNDANCE_H

#include <vector>
#include <string>
#include <utility>

namespace AtomProbe
{ 

//Structure containing information relating to a specific isotope
// such as say, U-238.
// Example (2^H, Deuterium)
//======
//  |------|
//  | 2    |<- mass number
//  |    H |<- Symbol
//  | 1    |<- atomic number
//  |------|
//    2.014101 <- mass
//    0.00000000006 <- mass error
//    0.000115 <- abundance
//    0.000070 <- abundance error
//======
struct ISOTOPE_ENTRY
{
    std::string symbol; //Isotope symbol
    size_t massNumber; //Mass number (for C, eg 12 or 13)
    size_t atomicNumber; //Atomic number (for C, this can only be 12)
    float mass; // in units specified in as-loaded XML file.
    float massError; //positive if known, zero if unknown
    float abundance; //Natural abundance. Zero is not naturally occurring. Between 0 and 1
    float abundanceError; //positive if known, zero if unknown
};

class AbundanceData
{
    enum
    {
        ABUNDANCE_ERR_BAD_DOC=1,
        ABUNDANCE_ERR_NO_CONTEXT,
        ABUNDANCE_ERROR_BAD_VALUE,
        ABUNDANCE_ERROR_FAILED_VALIDATION,
        ABUNDANCE_ERROR_MISSING_NODE,
        ABUNDANCE_ERROR_MISSING_ROOT_NODE,
        ABUNDANCE_ERROR_WRONG_ROOT_NODE,
        ABUNDANCE_ERROR_ENUM_END
    };
    
    // First vector is indexed by position, not atomic number because 
    //elements may be missing from abundance file, eg Tc
    std::vector<std::vector<ISOTOPE_ENTRY> >  isotopeData;

    //Check the abundance table for inconsistencies
    void checkErrors() const; 
    
public:

        size_t open(const char *file, bool strict=false);   

        static const char *getErrorText(size_t errorCode);

        size_t numIsotopes() const;
        size_t numElements() const;


        size_t symbolIndex(const char *symbol, bool caseSensitive=true) const;

        void getSymbolIndices(const std::vector<std::string> &symbols,
                        std::vector<size_t> &indices) const;
        
        void getSymbolIndices(const std::vector<std::pair<std::string,size_t> > &symbols,
                            std::vector<std::pair<size_t,unsigned int> > &indices) const;

        std::string elementName(size_t elemIdx ) const { return isotopeData[elemIdx][0].symbol;}

        size_t symbolIdxFromAtomicNumber(size_t atomicNumber) const;
    
        std::string elementNames(size_t start, size_t end, char separator=',') const; 
    

        void isotopeIndex(size_t elem, float mass, size_t &isotopeIdx) const;

        void isotopeIndex(size_t elem, float mass, float tolerance, size_t &isotopeIdx) const;

        unsigned int getAtomicNumber(size_t elemIdx) const;

        const std::vector<ISOTOPE_ENTRY> &isotopes(size_t offset) const; 


        void generateIsotopeDist(const std::vector<size_t> &elementIdx,
                    const std::vector<size_t> &frequency,
                    std::vector<std::pair<float,float> > &massDist) const;


        void generateGroupedIsotopeDist(const std::vector<size_t> &elementIdx,
                        const std::vector<size_t> &frequency, std::vector<std::pair<float,float> > &massDist, 
                        float massTolerance) const;

        void generateSingleAtomDist(size_t atomIdx, unsigned int repeatCount, 
                    std::vector<std::pair<float,float> > &massDist) const;

        float abundanceBetweenLimits(const std::vector<size_t> &elemIdx,
                        const std::vector<size_t> &frequency,
                        float massStart, float massEnd) const;

        const ISOTOPE_ENTRY &isotope(size_t elementIdx, size_t isotopeIdx) const;

        size_t getMajorIsotopeFromElemIdx(size_t elementIdx) const;

        float getNominalMass(size_t elementIdx) const;


        size_t getNearestCharge(const std::vector<std::pair<size_t,size_t> > &molecule,float targetMass,size_t maxCharge) const;

#ifdef DEBUG
        //Run the unit testing code
        static bool runUnitTests();
#endif
};

};

#endif