Skip to content

Commit feb1c2f

Browse files
SilinPavelZaal Lyanov
authored andcommitted
Overlapping reads handling implemented
1 parent a69332f commit feb1c2f

File tree

1 file changed

+23
-1
lines changed

1 file changed

+23
-1
lines changed

src/main/java/com/astrazeneca/vardict/VarDict.java

Lines changed: 23 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1596,6 +1596,10 @@ static Tuple4<Map<Integer, Map<String, Variation>>, Map<Integer, Map<String, Var
15961596
rlen = rlen2;
15971597
}
15981598

1599+
int mateAlignmentStart = record.getMateAlignmentStart();
1600+
boolean readsOverlap = isReadsOverlap(record);
1601+
1602+
processReadCigar:
15991603
//Loop over CIGAR records
16001604
for (int ci = 0; ci < cigar.numCigarElements(); ci++) {
16011605

@@ -2208,13 +2212,14 @@ && isHasAndEquals(querySequence.charAt(n - 1), ref, start - 1)) {
22082212
//adjust read position by m (CIGAR segment length) + offset + multoffp
22092213
n += offset + multoffp;
22102214
p += offset + multoffp;
2215+
2216+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
22112217
continue;
22122218
}
22132219
default:
22142220
break;
22152221
}
22162222
//End branch on CIGAR segment type
2217-
22182223
// Now dealing with matching part
22192224
int nmoff = 0;
22202225
int moffset = 0;
@@ -2243,6 +2248,7 @@ && isHasAndEquals(querySequence.charAt(n - 1), ref, start - 1)) {
22432248
start++;
22442249
n++;
22452250
p++;
2251+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
22462252
continue;
22472253
}
22482254

@@ -2295,6 +2301,7 @@ && isNotEquals('N', ref.get(start))) {
22952301
//shift reference position by 1
22962302
start++;
22972303
nmoff++;
2304+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
22982305
} else { //if bases match, exit loop
22992306
break;
23002307
}
@@ -2334,6 +2341,7 @@ && isNotEquals('N', ref.get(start))) {
23342341
n++;
23352342
p++;
23362343
start++;
2344+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
23372345
}
23382346

23392347
//remove '&' delimiter from s
@@ -2444,11 +2452,13 @@ && isNotEquals('N', ref.get(start))) {
24442452
//If variation starts with a deletion ('-' character)
24452453
if (startWithDelition) {
24462454
start += ddlen;
2455+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
24472456
}
24482457

24492458
//Shift reference position by 1 if CIGAR segment is not insertion
24502459
if (operator != CigarOperator.I) {
24512460
start++;
2461+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
24522462
}
24532463
//Shift read position by 1 if CIGAR segment is not deletion
24542464
if (operator != CigarOperator.D) {
@@ -2461,6 +2471,7 @@ && isNotEquals('N', ref.get(start))) {
24612471
n += moffset;
24622472
start += moffset;
24632473
p += moffset;
2474+
if (isPositionOverlapMate(start, mateAlignmentStart, readsOverlap)) break processReadCigar;
24642475
}
24652476
if (start > region.end) { //end if reference position is outside region of interest
24662477
break;
@@ -2500,6 +2511,17 @@ && isNotEquals('N', ref.get(start))) {
25002511
return tuple(hash, iHash, cov, rlen);
25012512
}
25022513

2514+
private static boolean isPositionOverlapMate(int start, int mateAlignmentStart, boolean readsOverlap) {
2515+
return readsOverlap && start >= mateAlignmentStart;
2516+
}
2517+
2518+
private static boolean isReadsOverlap(SAMRecord record) {
2519+
return !record.getReadNegativeStrandFlag()
2520+
&& record.getMateNegativeStrandFlag()
2521+
&& record.getMateReferenceName().equals(record.getReferenceName())
2522+
&& record.getMateAlignmentStart() <= record.getAlignmentEnd();
2523+
}
2524+
25032525
private static CigarOperator getCigarOperator(Cigar cigar, int ci) {
25042526
CigarOperator operator = cigar.getCigarElement(ci).getOperator();
25052527
// Treat insertions at the edge as soft-clipping

0 commit comments

Comments
 (0)