在 C++ 中读取 QUASAR 评分矩阵

/**
 * QUASARReader.cpp
 * QUASAR 替换矩阵格式的高效读取器。
 * 参见 http://www.bio.ifi.lmu.de/QUASAR/
 *
 * 需要 Boost String 算法（仅头文件）
 *
 * 编译方式：clang++ -o quasar quasar.cpp
 *
 * 版本 1.0.3
 *
 * 版权所有 (c) 2013 Uli Koehler
 * 此文件作为公共领域发布。
 * 如果在你的软件中使用，请提及原作者。
 */
#include <fstream>
#include <iostream>
#include <cassert>
#include <string>
#include <cstring>
#include <map>
#include <sstream>
#include <algorithm>
#include <vector>
#include <boost/algorithm/string.hpp>

/**
 * 方阵替换矩阵。
 * 此类只支持映射相同字符集的
 * 对称替换矩阵。
 */
class SubstitutionMatrix {
public:
    /**
     * 通过字符索引和默认分数构造新的替换矩阵。
     */
    SubstitutionMatrix(const std::string& index, double defaultScore = 0) : index(index) {
        size = index.size();
        matrix = new double[size * size];
        for(size_t i = 0; i < size * size; i++) {
            matrix[i] = defaultScore;
        }
    }
    ~SubstitutionMatrix() {
        if(matrix != nullptr) {
            delete[] matrix;
        }
    }
    /**
     * 对于给定字符，确定矩阵中的索引。
     * @return 索引，或 string::npos 如果未找到
     */
    int getIndexForChar(char c) {
        return index.find(c);
    }
    /**
     * 设置用字符 b 替换字符 a 的分数。
     */
    void setScore(char a, char b, double score) {
        setScore(getIndexForChar(a), getIndexForChar(b), score);
    }
    /**
     * 按索引设置分数。
     * @param a 字符 a 的索引
     * @param b 字符 b 的索引
     * @param score 用 b 替换 a（或用 a 替换 b）的分数
     */
    void setScore(int a, int b, double score) {
        assert(a >= 0);
        assert(b >= 0);
        assert(a < size);
        assert(b < size);
        //对称矩阵 - a->b == b->a
        matrix[a + size*b] = score;
        matrix[b + size*a] = score;
    }
    /**
     * 获取用字符 b 替换字符 a 的分数。
     * @return 分数
     */
    double getScore(char a, char b) {
        //替换 a 和 b 无关紧要（对称矩阵）
        return getScore(getIndexForChar(a), getIndexForChar(b));
    }
    /**
     * 获取用位置 b 的字符替换位置 a 的字符的分数。
     * @return 分数
     */
    double getScore(int a, int b) {
        assert(a >= 0);
        assert(b >= 0);
        assert(a < size);
        assert(b < size);
        //替换 a 和 b 无关紧要（对称矩阵）
        return matrix[b + size * a];
    }
    /**
     * 获取矩阵大小。
     * 矩阵有 getMatrixSize()**2 个条目。
     */
    size_t getMatrixSize() {
        return size * size;
    }
    /**
     * 获取矩阵字母表中所有字符的字符串
     */
    std::string getAlphabet() {
        return index;
    }
    /**
     * 获取原始矩阵缓冲区。
     * 有效长度为 sizeof(double) * size * size
     */
    double* getMatrix() {
        return matrix;
    }
private:
    /**
     * 包含分数的二维矩阵。
     * 行优先或列优先索引在此无关紧要
     * 因为矩阵是对称的。
     */
    double* matrix;
    size_t size;
    /**
     * 字母字符列表，
     * 分配给矩阵的数组索引。
     * 第一个字符分配给索引 0
     */
    std::string index;
};

SubstitutionMatrix* readSubstitutionMatrix(const char* filename) {
    static const char rowIndexPrefix[] = "ROWINDEX";
    static const char colIndexPrefix[] = "COLINDEX";
    static const char matrixPrefix[] = "MATRIX";
    std::ifstream infile(filename);
    std::string line = "", rowIndex = "", colIndex = "";
    int matrixRowIndex = 0; //已处理的 MATRIX 行数
    SubstitutionMatrix* matrix = nullptr;
    while(std::getline(infile, line)) {
        if(line.compare(0, sizeof(rowIndexPrefix) - 1, rowIndexPrefix) == 0) {
            rowIndex = line.substr(sizeof(rowIndexPrefix));
            boost::algorithm::trim(rowIndex);
        } else if(line.compare(0, sizeof(colIndexPrefix) - 1, colIndexPrefix) == 0) {
            colIndex = line.substr(sizeof(colIndexPrefix));
            boost::algorithm::trim(colIndex);
        } else if(line.compare(0, sizeof(matrixPrefix) - 1, matrixPrefix) == 0) {
            //如果尚未完成则初始化矩阵
            if(matrix == nullptr) {
                //目前不支持具有不同行和列索引的矩阵
                assert(!rowIndex.empty());
                assert(!colIndex.empty());
                assert(rowIndex == colIndex);
                matrix = new SubstitutionMatrix(rowIndex);
            }
            //分割元素
            int matrixColIndex = 0; //当前行中已处理的列数
            std::vector<std::string> elements;
            std::string lineContent = boost::trim_copy(line.substr(sizeof(matrixPrefix)));
            boost::split(elements, lineContent, boost::is_any_of("\t "), boost::token_compress_on);
            for(std::string elem : elements) {
                boost::trim(elem);
                double score = atof(elem.c_str());
                //行和列索引可互换
                matrix->setScore(matrixRowIndex, matrixColIndex, score);
                matrixColIndex++;
            }
            matrixRowIndex++;
        }
    }
    infile.close();
    return matrix;
}

using std::cout;
using std::cerr;
using std::endl;

/**
 * 主函数仅测试矩阵读取器对于给定矩阵是否不会失败。
 */
int main(int argc, char** argv) {
    if(argc < 2) {
        cerr << "Usage: " << argv[0] << " <QUASAR input file>" << endl;
        return 1;
    }
    SubstitutionMatrix* matrix = readSubstitutionMatrix(argv[1]);
    delete matrix;
    return 0;
}