Matrix

Author: kromberg

General/General.vcxproj

@@ -116,11 +116,15 @@
  </ItemDefinitionGroup>
  <ItemGroup>
  <ClInclude Include="Knapsack.h" />
+ <ClInclude Include="Matrix.h" />
  <ClInclude Include="MaxWIS.h" />
  <ClInclude Include="SequenceAlignment.h" />
  </ItemGroup>
  <ItemGroup>
  <ClCompile Include="Knapsack.cpp" />
+ <ClInclude Include="Matrix.hpp">
+ <FileType>CppCode</FileType>
+ </ClInclude>
  <ClCompile Include="MaxWIS.cpp" />
  <ClCompile Include="SequenceAlignment.cpp" />
  </ItemGroup>

General/General.vcxproj.filters

@@ -24,6 +24,9 @@
  <ClInclude Include="SequenceAlignment.h">
  <Filter>Header Files</Filter>
  </ClInclude>
+ <ClInclude Include="Matrix.h">
+ <Filter>Header Files</Filter>
+ </ClInclude>
  </ItemGroup>
  <ItemGroup>
  <ClCompile Include="MaxWIS.cpp">
@@ -35,5 +38,8 @@
  <ClCompile Include="SequenceAlignment.cpp">
  <Filter>Source Files</Filter>
  </ClCompile>
+ <ClCompile Include="Matrix.hpp">
+ <Filter>Header Files</Filter>
+ </ClCompile>
  </ItemGroup>
 </Project>

General/Matrix.h

@@ -0,0 +1,50 @@
+#ifndef MATRIX_H
+#define MATRIX_H
+
+#include <memory>
+#include <cstdint>
+#include <exception>
+#include <algorithm>
+
+template<typename T>
+class Matrix
+{
+public:
+ typedef T ValueType;
+private:
+ uint64_t m_rowsCount;
+ uint64_t m_colsCount;
+
+ std::unique_ptr<ValueType> m_matrix;
+
+public:
+ Matrix();
+ Matrix(uint64_t rowsCount, uint64_t colsCount) throw (std::bad_alloc);
+ Matrix(const Matrix<ValueType>& matrix) throw (std::bad_alloc);
+ Matrix(Matrix<ValueType>&& matrix);
+ ~Matrix();
+ 
+ ValueType* operator[] (uint64_t i);
+
+ void reset();
+ void resize(uint64_t rowsCount, uint64_t colsCount) throw (std::bad_alloc);
+
+ const uint64_t getRowsCount() const;
+ const uint64_t getColsCount() const;
+};
+
+template<typename T>
+inline const uint64_t Matrix<T>::getRowsCount() const
+{
+ return m_rowsCount;
+}
+
+template<typename T>
+inline const uint64_t Matrix<T>::getColsCount() const
+{
+ return m_colsCount;
+}
+
+#include "Matrix.hpp"
+
+#endif

General/Matrix.hpp

@@ -0,0 +1,72 @@
+template<typename T>
+Matrix<T>::Matrix():
+ m_rowsCount(0),
+ m_colsCount(0),
+ m_matrix()
+{
+}
+
+template<typename T>
+Matrix<T>::Matrix(uint64_t rowsCount, uint64_t colsCount) throw (std::bad_alloc) :
+ m_rowsCount(rowsCount),
+ m_colsCount(colsCount),
+ m_matrix()
+{
+ m_matrix.reset(new ValueType[m_rowsCount * m_colsCount]);
+}
+
+template<typename T>
+Matrix<T>::Matrix(const Matrix<T>& matrix) throw (std::bad_alloc) :
+ m_rowsCount(matrix.m_rowsCount),
+ m_colsCount(matrix.m_colsCount),
+ m_matrix()
+{
+ m_matrix.reset(new ValueType[m_rowsCount * m_colsCount]);
+ memcpy(m_matrix.get(), matrix.m_matrix.get(), sizeof(ValueType) * m_rowsCount * m_colsCount);
+}
+
+template<typename T>
+Matrix<T>::Matrix(Matrix&& matrix):
+ m_rowsCount(matrix.m_rowsCount),
+ m_colsCount(matrix.m_colsCount)
+{
+ m_matrix = std::move(matrix.m_matrix);
+ matrix.m_rowsCount = 0;
+ matrix.m_colsCount = 0;
+}
+
+template<typename T>
+Matrix<T>::~Matrix()
+{
+ m_matrix.reset();
+}
+
+template<typename T>
+T* Matrix<T>::operator[] (uint64_t i)
+{
+ return &m_matrix.get()[i * m_colsCount];
+}
+
+template<typename T>
+void Matrix<T>::reset()
+{
+ m_matrix.reset();
+ m_rowsCount = 0;
+ m_colsCount = 0;
+}
+
+template<typename T>
+void Matrix<T>::resize(uint64_t rowsCount, uint64_t colsCount) throw (std::bad_alloc)
+{
+ std::unique_ptr<ValueType> matrix(new ValueType[rowsCount * colsCount]);
+ uint64_t minColsCount = std::min(m_colsCount, colsCount);
+ uint64_t minRowsCount = std::min(m_rowsCount, rowsCount);
+ for (uint64_t i = 0; i < minRowsCount; ++i)
+ {
+ memcpy(&matrix.get()[i * colsCount], &m_matrix.get()[i * m_colsCount], sizeof(ValueType) * minColsCount);
+ }
+ m_matrix.reset();
+ m_rowsCount = rowsCount;
+ m_colsCount = colsCount;
+ m_matrix = std::move(matrix);
+}

General/SequenceAlignment.cpp

@@ -17,7 +17,7 @@ void SequenceAlignment::clear()
  m_charToIndexMap.clear();
  m_numScores = 0;
  m_gapScore = 0;
- m_scores.clear();
+ m_scores.reset();
 }
 
 bool SequenceAlignment::read(std::istream& in)
@@ -28,20 +28,20 @@ bool SequenceAlignment::read(std::istream& in)
  std::string alphabet;
  in >> m_alphabet;
 
- uint32_t idx = 0;
+ uint64_t idx = 0;
  for (auto c : m_alphabet)
  {
  m_charToIndexMap.insert(std::make_pair(c, idx++));
  }
 
  m_numScores = m_alphabet.size();
- m_scores.resize(m_numScores * m_numScores, 0);
+ m_scores.resize(m_numScores, m_numScores);
 
- for (uint32_t i = 0; i < m_numScores; ++i)
+ for (uint64_t i = 0; i < m_numScores; ++i)
  {
- for (uint32_t j = 0; j < m_numScores; ++j)
+ for (uint64_t j = 0; j < m_numScores; ++j)
  {
- in >> m_scores[i * m_numScores + j];
+ in >> m_scores[i][j];
  }
  }
  in >> m_gapScore;
@@ -55,53 +55,53 @@ bool SequenceAlignment::run(
  uint32_t &scoreNeedlemanWunsch,
  std::string(&seq)[2])
 {
- uint32_t size[2] = { m_sequnces[0].size() + 1, m_sequnces[1].size() + 1 };
- std::vector<char> nwScore(size[0] * size[1], 0);
- for (uint32_t i = 0; i < size[1]; ++i)
+ uint64_t size[2] = { m_sequnces[0].size() + 1, m_sequnces[1].size() + 1 };
+ Matrix<uint32_t> nwScore(size[0], size[1]);
+ for (uint64_t i = 0; i < size[0]; ++i)
  {
- nwScore[i] = (i + 1) * m_gapScore;
+ nwScore[i][0] = (i + 1) * m_gapScore;
  }
- for (uint32_t i = 0; i < size[0]; ++i)
+ for (uint64_t i = 0; i < size[1]; ++i)
  {
- nwScore[i * size[1]] = (i + 1) * m_gapScore;
+ nwScore[0][i] = (i + 1) * m_gapScore;
  }
 
  uint32_t tmpScore[3];
- for (uint32_t i = 1; i < size[0]; ++i)
+ for (uint64_t i = 1; i < size[0]; ++i)
  {
- for (uint32_t j = 1; j < size[1]; ++j)
+ for (uint64_t j = 1; j < size[1]; ++j)
  {
- uint32_t idx1 = m_charToIndexMap[m_sequnces[0][i - 1]];
- uint32_t idx2 = m_charToIndexMap[m_sequnces[1][j - 1]];
- tmpScore[0] = nwScore[(i - 1) * size[1] + (j - 1)] + m_scores[idx1 * m_numScores + idx2];
- tmpScore[1] = nwScore[i * size[1] + (j - 1)] + m_gapScore;
- tmpScore[2] = nwScore[(i - 1) * size[1] + j] + m_gapScore;
- nwScore[i * size[1] + j] = std::min(tmpScore[0], std::min(tmpScore[1], tmpScore[2]));
+ uint64_t idx1 = m_charToIndexMap[m_sequnces[0][i - 1]];
+ uint64_t idx2 = m_charToIndexMap[m_sequnces[1][j - 1]];
+ tmpScore[0] = nwScore[i - 1][j - 1] + m_scores[idx1][idx2];
+ tmpScore[1] = nwScore[i][j - 1] + m_gapScore;
+ tmpScore[2] = nwScore[i - 1][j] + m_gapScore;
+ nwScore[i][j] = std::min(tmpScore[0], std::min(tmpScore[1], tmpScore[2]));
  }
  }
 
- scoreNeedlemanWunsch = nwScore[size[0] * size[1] - 1];
+ scoreNeedlemanWunsch = nwScore[size[0] - 1][size[1] - 1];
 
- for (uint32_t i = size[0] - 1, j = size[1] - 1; (i > 0) && (j > 0);)
+ for (uint64_t i = size[0] - 1, j = size[1] - 1; (i > 0) && (j > 0);)
  {
- uint32_t idx1 = m_charToIndexMap[m_sequnces[0][i - 1]];
- uint32_t idx2 = m_charToIndexMap[m_sequnces[1][j - 1]];
- tmpScore[0] = nwScore[(i - 1) * size[1] + (j - 1)] + m_scores[idx1 * m_numScores + idx2];
- tmpScore[1] = nwScore[i * size[1] + (j - 1)] + m_gapScore;
- tmpScore[2] = nwScore[(i - 1) * size[1] + j] + m_gapScore;
- if (nwScore[i * size[1] + j] == tmpScore[0])
+ uint64_t idx1 = m_charToIndexMap[m_sequnces[0][i - 1]];
+ uint64_t idx2 = m_charToIndexMap[m_sequnces[1][j - 1]];
+ tmpScore[0] = nwScore[i - 1][j - 1] + m_scores[idx1][idx2];
+ tmpScore[1] = nwScore[i][j - 1] + m_gapScore;
+ tmpScore[2] = nwScore[i - 1][j] + m_gapScore;
+ if (nwScore[i][j] == tmpScore[0])
  {
  seq[0].push_back(m_sequnces[0][i - 1]);
  seq[1].push_back(m_sequnces[1][j - 1]);
  --i, --j;
  }
- else if (nwScore[i * size[1] + j] == tmpScore[1])
+ else if (nwScore[i][j] == tmpScore[1])
  {
  seq[0].push_back(' ');
  seq[1].push_back(m_sequnces[1][j - 1]);
  --j;
  }
- else if (nwScore[i * size[1] + j] == tmpScore[2])
+ else if (nwScore[i][j] == tmpScore[2])
  {
  seq[0].push_back(m_sequnces[0][i - 1]);
  seq[1].push_back(' ');

General/SequenceAlignment.h

@@ -6,14 +6,18 @@
 #include <string>
 #include <cstdint>
 
+#include "Matrix.h"
+
 class SequenceAlignment
 {
 private:
- typedef std::unordered_map<char, uint32_t> CharToIndexMap;
+ typedef std::unordered_map<char, uint64_t> CharToIndexMap;
  CharToIndexMap m_charToIndexMap;
  std::string m_alphabet;
- uint32_t m_numScores;
- std::vector<uint32_t> m_scores;
+ uint64_t m_numScores;
+ //std::vector<uint32_t> m_scores;
+ Matrix<uint32_t> m_scores;
+
  uint32_t m_gapScore;
  std::string m_sequnces[2];