zzh
/
gatk-3.8
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Fix some special-case cleaner issues. 

We now do the same as brute force in all examples to date.


git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@626 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
ebanks 2009-05-07 18:16:35 +00:00
parent 93211c1cd8
commit 283a4d1b54
1 changed files with 15 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

26
java/src/org/broadinstitute/sting/playground/gatk/walkers/IntervalCleanerWalker.java
View File

@ -135,24 +135,30 @@ public class IntervalCleanerWalker extends LocusWindowWalker<Integer, Integer>
// do a pairwise alignment against the reference // do a pairwise alignment against the reference
AlignedRead aRead = altReads.get(index); AlignedRead aRead = altReads.get(index);
SWPairwiseAlignment swConsensus = new SWPairwiseAlignment(reference, aRead.getReadString()); SWPairwiseAlignment swConsensus = new SWPairwiseAlignment(reference, aRead.getReadString());
int idx = swConsensus.getAlignmentStart2wrt1(); int refIdx = swConsensus.getAlignmentStart2wrt1();
// create the new consensus // create the new consensus
StringBuffer sb = new StringBuffer(); StringBuffer sb = new StringBuffer();
sb.append(reference.substring(0, idx)); sb.append(reference.substring(0, refIdx));
sb.append(aRead.getReadString());
Cigar c = swConsensus.getCigar(); Cigar c = swConsensus.getCigar();
int indelCount = 0; int indelCount = 0;
int altIdx = 0;
for ( int i = 0 ; i < c.numCigarElements() ; i++ ) { for ( int i = 0 ; i < c.numCigarElements() ; i++ ) {
CigarElement ce = c.getCigarElement(i); CigarElement ce = c.getCigarElement(i);
switch( ce.getOperator() ) { switch( ce.getOperator() ) {
case D: case D:
indelCount++; indelCount++;
refIdx += ce.getLength();
break;
case M: case M:
idx += ce.getLength(); sb.append(reference.substring(refIdx, refIdx+ce.getLength()));
refIdx += ce.getLength();
altIdx += ce.getLength();
break; break;
case I: case I:
sb.append(aRead.getReadString().substring(altIdx, altIdx+ce.getLength()));
altIdx += ce.getLength();
indelCount++; indelCount++;
break; break;
} }
@ -161,16 +167,12 @@ public class IntervalCleanerWalker extends LocusWindowWalker<Integer, Integer>
if ( indelCount > 1 ) if ( indelCount > 1 )
continue; continue;
sb.append(reference.substring(idx)); sb.append(reference.substring(refIdx));
String altConsensus = sb.toString(); String altConsensus = sb.toString();
// for each imperfect match to the reference, score it against this alternative // for each imperfect match to the reference, score it against this alternative
Consensus consensus = new Consensus(altConsensus, c, swConsensus.getAlignmentStart2wrt1()); Consensus consensus = new Consensus(altConsensus, c, swConsensus.getAlignmentStart2wrt1());
for ( int j = 0; j < altReads.size(); j++ ) { for ( int j = 0; j < altReads.size(); j++ ) {
if (j == index) {
consensus.readIndexes.add(new Pair<Integer, Integer>(j, swConsensus.getAlignmentStart2wrt1()));
continue;
}
AlignedRead toTest = altReads.get(j); AlignedRead toTest = altReads.get(j);
Pair<Integer, Integer> altAlignment = findBestOffset(altConsensus, toTest); Pair<Integer, Integer> altAlignment = findBestOffset(altConsensus, toTest);
@ -282,7 +284,8 @@ public class IntervalCleanerWalker extends LocusWindowWalker<Integer, Integer>
indelOffsetOnRead = altCE2.getLength(); indelOffsetOnRead = altCE2.getLength();
} }
} else if ( altCE2.getOperator() == CigarOperator.D ) { } else if ( altCE2.getOperator() == CigarOperator.D ) {
readCigar.add(altCE2); if ( sawAlignmentStart )
readCigar.add(altCE2);
indelOffsetOnRef = altCE2.getLength(); indelOffsetOnRef = altCE2.getLength();
} else { } else {
throw new RuntimeException("Operator of middle block is not I or D: " + altCE2.getOperator()); throw new RuntimeException("Operator of middle block is not I or D: " + altCE2.getOperator());
@ -301,7 +304,8 @@ public class IntervalCleanerWalker extends LocusWindowWalker<Integer, Integer>
if ( ce.getOperator() != CigarOperator.D ) if ( ce.getOperator() != CigarOperator.D )
readRemaining -= ce.getLength(); readRemaining -= ce.getLength();
} }
readCigar.add(new CigarElement(readRemaining, CigarOperator.M)); if ( readRemaining > 0 )
readCigar.add(new CigarElement(readRemaining, CigarOperator.M));
aRead.getRead().setCigar(readCigar); aRead.getRead().setCigar(readCigar);
} }