Refactor
authorMartin C. Frith
Fri Mar 13 06:22:27 2020 +0900 (2 months ago)
changeset 10553887935792ab
parent 1054 8621c8e512ff
child 1056 6f4ec587fe25
Refactor
src/split/cbrc_split_aligner.cc
src/split/cbrc_split_aligner.hh
     1.1 --- a/src/split/cbrc_split_aligner.cc	Wed Mar 11 21:32:19 2020 +0900
     1.2 +++ b/src/split/cbrc_split_aligner.cc	Fri Mar 13 06:22:27 2020 +0900
     1.3 @@ -319,9 +319,6 @@
     1.4  }
     1.5  
     1.6  long SplitAligner::viterbiSplit() {
     1.7 -  resizeMatrix(Vmat);
     1.8 -  resizeVector(Vvec);
     1.9 -
    1.10    unsigned *inplayAlnBeg = &newInplayAlnIndices[0];
    1.11    unsigned *inplayAlnEnd = inplayAlnBeg;
    1.12    unsigned *sortedAlnPtr = &sortedAlnIndices[0];
    1.13 @@ -360,9 +357,6 @@
    1.14  }
    1.15  
    1.16  long SplitAligner::viterbiSplice() {
    1.17 -    resizeMatrix(Vmat);
    1.18 -    resizeVector(Vvec);
    1.19 -
    1.20      unsigned sortedAlnPos = 0;
    1.21      unsigned oldNumInplay = 0;
    1.22      unsigned newNumInplay = 0;
    1.23 @@ -548,9 +542,6 @@
    1.24  }
    1.25  
    1.26  void SplitAligner::forwardSplit() {
    1.27 -  resizeVector(rescales);
    1.28 -  resizeMatrix(Fmat);
    1.29 -
    1.30    unsigned *inplayAlnBeg = &newInplayAlnIndices[0];
    1.31    unsigned *inplayAlnEnd = inplayAlnBeg;
    1.32    unsigned *sortedAlnPtr = &sortedAlnIndices[0];
    1.33 @@ -592,9 +583,6 @@
    1.34  }
    1.35  
    1.36  void SplitAligner::forwardSplice() {
    1.37 -    resizeVector(rescales);
    1.38 -    resizeMatrix(Fmat);
    1.39 -
    1.40      unsigned sortedAlnPos = 0;
    1.41      unsigned oldNumInplay = 0;
    1.42      unsigned newNumInplay = 0;
    1.43 @@ -641,8 +629,6 @@
    1.44  }
    1.45  
    1.46  void SplitAligner::backwardSplit() {
    1.47 -  resizeMatrix(Bmat);
    1.48 -
    1.49    unsigned *inplayAlnBeg = &newInplayAlnIndices[0];
    1.50    unsigned *inplayAlnEnd = inplayAlnBeg;
    1.51    unsigned *sortedAlnPtr = &sortedAlnIndices[0];
    1.52 @@ -679,8 +665,6 @@
    1.53  }
    1.54  
    1.55  void SplitAligner::backwardSplice() {
    1.56 -    resizeMatrix(Bmat);
    1.57 -
    1.58      unsigned sortedAlnPos = 0;
    1.59      unsigned oldNumInplay = 0;
    1.60      unsigned newNumInplay = 0;
    1.61 @@ -723,22 +707,22 @@
    1.62  std::vector<double>
    1.63  SplitAligner::marginalProbs(unsigned queryBeg, unsigned alnNum,
    1.64  			    unsigned alnBeg, unsigned alnEnd) const {
    1.65 -    std::vector<double> output;
    1.66 -    unsigned i = alnNum;
    1.67 -    unsigned j = queryBeg;
    1.68 -    for (unsigned pos = alnBeg; pos < alnEnd; ++pos) {
    1.69 -	size_t ij = matrixRowOrigins[i] + j;
    1.70 -        if (alns[i].qalign[pos] == '-') {
    1.71 -            double value = Fmat[ij] * Bmat[ij] * Dexp[ij] * cell(rescales, j);
    1.72 -            output.push_back(value);
    1.73 -        } else {
    1.74 -            double value = Fmat[ij + 1] * Bmat[ij] / Aexp[ij];
    1.75 -            if (value != value) value = 0.0;
    1.76 -            output.push_back(value);
    1.77 -            j++;
    1.78 -        }
    1.79 +  std::vector<double> output;
    1.80 +  unsigned i = alnNum;
    1.81 +  unsigned j = queryBeg;
    1.82 +  for (unsigned pos = alnBeg; pos < alnEnd; ++pos) {
    1.83 +    size_t ij = matrixRowOrigins[i] + j;
    1.84 +    if (alns[i].qalign[pos] == '-') {
    1.85 +      double value = Fmat[ij] * Bmat[ij] * Dexp[ij] * cell(rescales, j);
    1.86 +      output.push_back(value);
    1.87 +    } else {
    1.88 +      double value = Fmat[ij + 1] * Bmat[ij] / Aexp[ij];
    1.89 +      if (value != value) value = 0.0;
    1.90 +      output.push_back(value);
    1.91 +      j++;
    1.92      }
    1.93 -    return output;
    1.94 +  }
    1.95 +  return output;
    1.96  }
    1.97  
    1.98  // The next routine represents affine gap scores in a cunning way.
     2.1 --- a/src/split/cbrc_split_aligner.hh	Wed Mar 11 21:32:19 2020 +0900
     2.2 +++ b/src/split/cbrc_split_aligner.hh	Fri Mar 13 06:22:27 2020 +0900
     2.3 @@ -74,10 +74,9 @@
     2.4      // Call this before viterbi/forward/backward, and after layout
     2.5      void initMatricesForOneQuery();
     2.6  
     2.7 -    long viterbiSplit();
     2.8 -    long viterbiSplice();
     2.9 -
    2.10      long viterbi() {  // returns the optimal split-alignment score
    2.11 +      resizeMatrix(Vmat);
    2.12 +      resizeVector(Vvec);
    2.13        return (restartProb > 0) ? viterbiSplit() : viterbiSplice();
    2.14      }
    2.15  
    2.16 @@ -96,12 +95,10 @@
    2.17      int segmentScore(unsigned alnNum,
    2.18  		     unsigned queryBeg, unsigned queryEnd) const;
    2.19  
    2.20 -    void forwardSplit();
    2.21 -    void backwardSplit();
    2.22 -    void forwardSplice();
    2.23 -    void backwardSplice();
    2.24 -
    2.25      void forwardBackward() {
    2.26 +      resizeVector(rescales);
    2.27 +      resizeMatrix(Fmat);
    2.28 +      resizeMatrix(Bmat);
    2.29        if (restartProb > 0) {
    2.30  	forwardSplit();
    2.31  	backwardSplit();
    2.32 @@ -243,6 +240,14 @@
    2.33  				 unsigned& oldNumInplay,
    2.34  				 unsigned& newNumInplay, unsigned j);
    2.35  
    2.36 +    long viterbiSplit();
    2.37 +    long viterbiSplice();
    2.38 +
    2.39 +    void forwardSplit();
    2.40 +    void backwardSplit();
    2.41 +    void forwardSplice();
    2.42 +    void backwardSplice();
    2.43 +
    2.44      unsigned findScore(unsigned j, long score) const;
    2.45      unsigned findSpliceScore(unsigned i, unsigned j, long score) const;
    2.46      long scoreFromSplice(unsigned i, unsigned j, unsigned oldNumInplay,