-
Notifications
You must be signed in to change notification settings - Fork 67
/
HybridAligner.cpp
107 lines (82 loc) · 2.7 KB
/
HybridAligner.cpp
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
95
96
97
98
99
100
101
102
103
104
105
106
107
/*++
Module Name:
HybridAligner.h
Abstract:
Header for SNAP genome aligner
Authors:
Bill Bolosky, October, 2011
Environment:
User mode service.
This class is NOT thread safe. It's the caller's responsibility to ensure that
at most one thread uses an instance at any time.
Revision History:
Adapted from Matei Zaharia's Scala implementation.
--*/
#include "stdafx.h"
#include "HybridAligner.h"
#if 0
HybridAligner::HybridAligner(
GenomeIndex *genomeIndex,
unsigned confDiff,
unsigned maxHitsToConsider,
unsigned maxK,
unsigned maxReadSize,
unsigned maxSeedsToUse,
unsigned lvCutoff)
{
baseAligner = new BaseAligner(genomeIndex,confDiff,maxHitsToConsider,maxK,maxReadSize,maxSeedsToUse,lvCutoff,100000);
intersectingAligner = new IntersectingAligner(genomeIndex,confDiff,280000,5,maxReadSize,20,lvCutoff);
if (NULL == baseAligner || NULL == intersectingAligner) {
fprintf(stderr,"HybridAligner: Unable to allocate underlying aligner\n");
exit(1);
}
}
AlignmentResult
HybridAligner::AlignRead(
Read *read,
unsigned *genomeLocation,
bool *hitIsRC,
int *score)
{
unsigned genomeLocation1;
bool hitIsRC1;
int score1;
AlignmentResult result1 = baseAligner->AlignRead(read,&genomeLocation1,&hitIsRC1,&score1);
if (CertainHit == result1 /*|| SingleHit == result1 && score1 <= 4*/) {
*genomeLocation = genomeLocation1;
*hitIsRC = hitIsRC1;
if (NULL != score) {
*score = score1;
}
return result1;
}
unsigned genomeLocation2;
bool hitIsRC2;
int score2;
AlignmentResult result2 = intersectingAligner->AlignRead(read,&genomeLocation2,&hitIsRC2,&score2);
if (SingleHit == result1 && SingleHit == result2 && genomeLocation1 != genomeLocation2 &&
(score1 <= score2 && score1 + baseAligner->getConfDiff() > score2 ||
score2 <= score1 && score2 + baseAligner->getConfDiff() > score1)) {
//
// Both aligners had single hits at different locations with scores within
// confDiff of one another. Declare a MultiHit.
//
return MultipleHits;
}
if (CertainHit == result2 ||
SingleHit == result2 && (SingleHit == result1 && score2 < score1 || SingleHit != result1)) {
*genomeLocation = genomeLocation2;
*hitIsRC = hitIsRC2;
if (NULL != score) {
*score = score2;
}
return result2;
}
*genomeLocation = genomeLocation1;
*hitIsRC = hitIsRC1;
if (NULL != score) {
*score = score1;
}
return result1;
}
#endif // 0