1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
|
#include "IncludeDefine.h"
#include "Parameters.h"
#include "Transcript.h"
#include "ReadAlign.h"
#include "ErrorWarning.h"
#include <random>
void ReadAlign::multMapSelect() {//select multiple mappers from all transcripts of all windows
nTr=0;
if (nW==0) {//no good windows
return;
};
maxScore=-10*Lread;
for (uint iW=0; iW<nW; iW++) {//scan windows
if (maxScore < trAll[iW][0]->maxScore) maxScore = trAll[iW][0]->maxScore;
};
if (maxScore!=trBest->maxScore) {
ostringstream errOut;
errOut << "BUG: maxScore!=trBest->maxScore in multMapSelect";
exitWithError(errOut.str(), std::cerr, P.inOut->logMain, EXIT_CODE_BUG, P);
};
for (uint iW=0; iW<nW; iW++) {//scan windows
for (uint iTr=0; iTr<nWinTr[iW]; iTr++) {//scan transcripts
if ( (trAll[iW][iTr]->maxScore + P.outFilterMultimapScoreRange) >= maxScore ) {//record this alignment
// if paired-end, record alignments from ALL windows
if (nTr==MAX_N_MULTMAP) {//too many alignments for this read, do not record it
ostringstream errOut;
errOut << "EXITING: Fatal ERROR: number of alignments exceeds MAX_N_MULTMAP, increase it and re-compile STAR";
exitWithError(errOut.str(), std::cerr, P.inOut->logMain, EXIT_CODE_PARAMETER, P);
};
trMult[nTr]=trAll[iW][iTr];
trMult[nTr]->Chr = trAll[iW][0]->Chr;
trMult[nTr]->Str = trAll[iW][0]->Str;
trMult[nTr]->roStr = trAll[iW][0]->roStr;
nTr++;
};
};
};
if (nTr > P.outFilterMultimapNmax || nTr==0){
//too multi OR no alignments, no need for further processing, since it will be considered unmapped
return;
};
for (uint iTr=0; iTr<nTr; iTr++) {
trMult[iTr]->roStart = trMult[iTr]->roStr==0 ? trMult[iTr]->rStart : Lread - trMult[iTr]->rStart - trMult[iTr]->rLength;
trMult[iTr]->cStart=trMult[iTr]->gStart - mapGen.chrStart[trMult[iTr]->Chr];
};
if (nTr==1){//unique mappers
trMult[0]->primaryFlag=true;
} else {//multimappers
int nbest=0;
if (P.outMultimapperOrder.random || P.outSAMmultNmax != (uint) -1 ) {//bring the best alignment to the top of the list. TODO sort alignments by the score?
for (uint itr=0; itr<nTr; itr++) {//move the best aligns to the top of the list
if ( trMult[itr]->maxScore == maxScore ) {
swap(trMult[itr],trMult[nbest]);
++nbest;
};
};
};
if (P.outMultimapperOrder.random) {//shuffle separately the best aligns, and the rest
for (int itr=nbest-1; itr>=1; itr--) {//Fisher-Yates-Durstenfeld-Knuth shuffle
int rand1=int (rngUniformReal0to1(rngMultOrder)*itr+0.5);
swap(trMult[itr],trMult[rand1]);
};
for (int itr=nTr-nbest-1; itr>=1; itr--) {//Fisher-Yates-Durstenfeld-Knuth shuffle
int rand1=int (rngUniformReal0to1(rngMultOrder)*itr+0.5);
swap(trMult[nbest+itr],trMult[nbest+rand1]);
};
};
if ( P.outSAMprimaryFlag=="AllBestScore" ) {
for (uint itr=0; itr<nTr; itr++)
{//mark all best score aligns as primary
if ( trMult[itr]->maxScore == maxScore ) trMult[itr]->primaryFlag=true;
};
} else if (P.outMultimapperOrder.random || P.outSAMmultNmax != (uint) -1) {
trMult[0]->primaryFlag=true;//mark as primary the first one in the random ordered list: best scoring aligns are already in front of the list
} else {//old way
trBest->primaryFlag=true;
};
};
//TODO re-point trBest to the primary one. This should not make a difference in the output
};
|
