4 * Copyright (c) 2007 Des Higgins, Julie Thompson and Toby Gibson.
9 #include "LowScoreSegProfile.h"
14 LowScoreSegProfile::LowScoreSegProfile(int prfLen, int firstS, int lastS)
19 profile.resize(prfLength + 2, vector<int>(LENCOL + 2));
22 void LowScoreSegProfile::calcLowScoreSegProfile(const SeqArray* seqArray,
23 int matrix[NUMRES][NUMRES], vector<int>* seqWeight)
25 vector<vector<int> > weighting;
29 int _gapPos1 = userParameters->getGapPos1();
30 int _gapPos2 = userParameters->getGapPos2();
31 int _maxAA = userParameters->getMaxAA();
33 weighting.resize(NUMRES + 2, vector<int>(prfLength + 2));
35 for (r = 0; r < prfLength; r++)
37 for (d = 0; d <= _maxAA; d++)
40 for (i = firstSeq; i < lastSeq; i++)
42 if (r + 1 < (int)(*seqArray)[i + 1].size() - 1)
44 if (d == (*seqArray)[i + 1][r + 1])
46 weighting[d][r] += (*seqWeight)[i];
52 weighting[_gapPos1][r] = 0;
54 for (i = firstSeq; i < lastSeq; i++)
56 if (r + 1 < (int)(*seqArray)[i + 1].size() - 1)
58 if (_gapPos1 == (*seqArray)[i + 1][r + 1])
60 weighting[_gapPos1][r] += (*seqWeight)[i];
65 weighting[_gapPos2][r] = 0;
67 for (i = firstSeq; i < lastSeq; i++)
69 if (r + 1 < (int)(*seqArray)[i + 1].size() - 1)
71 if (_gapPos2 == (*seqArray)[i + 1][r + 1])
73 weighting[_gapPos2][r] += (*seqWeight)[i];
79 for (pos = 0; pos < prfLength; pos++)
81 for (res = 0; res <= _maxAA; res++)
85 for (d = 0; d <= _maxAA; d++)
87 f += (weighting[d][pos] * matrix[d][res]);
90 f += (weighting[_gapPos1][pos] * matrix[_gapPos1][res]);
91 f += (weighting[_gapPos2][pos] * matrix[_gapPos2][res]);
92 profile[pos + 1][res] = f;
96 for (d = 0; d <= _maxAA; d++)
98 f += (weighting[d][pos] * matrix[d][_gapPos1]);
101 f += (weighting[_gapPos1][pos] * matrix[_gapPos1][_gapPos1]);
102 f += (weighting[_gapPos2][pos] * matrix[_gapPos2][_gapPos1]);
103 profile[pos + 1][_gapPos1] = f;
106 for (d = 0; d <= _maxAA; d++)
108 f += (weighting[d][pos] * matrix[d][_gapPos2]);
110 f += (weighting[_gapPos1][pos] * matrix[_gapPos1][_gapPos2]);
111 f += (weighting[_gapPos2][pos] * matrix[_gapPos2][_gapPos2]);
112 profile[pos + 1][_gapPos2] = f;