diff --git a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManager.java b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManager.java index 4ac4f5739..cd4f359b4 100644 --- a/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManager.java +++ b/protected/gatk-tools-protected/src/main/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManager.java @@ -279,7 +279,7 @@ public class ConstrainedMateFixingManager { // fix mates, as needed // Since setMateInfo can move reads, we potentially need to remove the mate, and requeue // it to ensure proper sorting - if ( newRead.getReadPairedFlag() && !newRead.getNotPrimaryAlignmentFlag() ) { + if ( isMateFixableRead(newRead) ) { SAMRecordHashObject mate = forMateMatching.get(newRead.getReadName()); if ( mate != null ) { // 1. Frustratingly, Picard's setMateInfo() method unaligns (by setting the reference contig @@ -307,13 +307,9 @@ public class ConstrainedMateFixingManager { reQueueMate = false; } - // we've already seen our mate -- set the mate info and remove it from the map - // Via Nils Homer: - // There will be two SamPairUtil.setMateInfo functions. The default will not update the mate - // cigar tag; in fact, it will remove it if it is present. An alternative SamPairUtil.setMateInfo - // function takes a boolean as an argument ("addMateCigar") and will add/update the mate cigar if - // set to true. This is the one you want to use. - SamPairUtil.setMateInfo(mate.record, newRead, null, true); + // we've already seen our mate -- set the mate info and remove it from the map; + // add/update the mate cigar if appropriate + SamPairUtil.setMateInfo(mate.record, newRead, true); if ( reQueueMate ) waitingReads.add(mate.record); } @@ -364,6 +360,16 @@ public class ConstrainedMateFixingManager { } } + /** + * Is the given read one for which we can fix its mate? + * + * @param read the read + * @return true if we could fix its mate, false otherwise + */ + protected boolean isMateFixableRead(final SAMRecord read) { + return read.getReadPairedFlag() && !read.isSecondaryOrSupplementary(); + } + /** * @param read the read * @return true if the read shouldn't be moved given the constraints of this SAMFileWriter diff --git a/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManagerUnitTest.java b/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManagerUnitTest.java index 4f5453fdb..d35536df3 100644 --- a/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManagerUnitTest.java +++ b/protected/gatk-tools-protected/src/test/java/org/broadinstitute/gatk/tools/walkers/indels/ConstrainedMateFixingManagerUnitTest.java @@ -52,7 +52,9 @@ package org.broadinstitute.gatk.tools.walkers.indels; import htsjdk.samtools.SAMFileHeader; +import htsjdk.samtools.SAMFileWriter; import htsjdk.samtools.SAMRecord; +import htsjdk.samtools.util.ProgressLoggerInterface; import org.broadinstitute.gatk.utils.BaseTest; import org.broadinstitute.gatk.utils.GenomeLocParser; import org.broadinstitute.gatk.utils.sam.ArtificialSAMUtils; @@ -61,6 +63,7 @@ import org.testng.Assert; import org.testng.annotations.BeforeClass; import org.testng.annotations.Test; +import java.util.ArrayList; import java.util.List; @@ -136,4 +139,60 @@ public class ConstrainedMateFixingManagerUnitTest extends BaseTest { Assert.assertTrue(manager.forMateMatching.containsKey("foo")); } + + @Test + public void testSupplementaryAlignmentsDoNotCauseBadMateFixing() { + final List properReads = ArtificialSAMUtils.createPair(header, "foo", 1, 1000, 2000, true, false); + final GATKSAMRecord read1 = properReads.get(0); + read1.setFlags(99); // first in pair, negative strand + + final GATKSAMRecord read2 = properReads.get(1); + read2.setFlags(161); // second in pair, mate negative strand + + final GATKSAMRecord read2Supp = new GATKSAMRecord(read2); + read2Supp.setReadName("foo"); + read2Supp.setFlags(2209); // second in pair, mate negative strand, supplementary + read2Supp.setAlignmentStart(100); + read2Supp.setMateAlignmentStart(1000); + + final DummyWriter writer = new DummyWriter(); + final ConstrainedMateFixingManager manager = new ConstrainedMateFixingManager(writer, genomeLocParser, 10000, 200, 10000); + manager.addRead(read2Supp, false, false); + manager.addRead(read1, false, false); + manager.addRead(read2, false, false); + manager.close(); // "write" the reads to our dummy writer + + // check to make sure that none of the mate locations were changed, which is the problem brought to us by a user + for ( final SAMRecord read : writer.reads ) { + final int start = read.getAlignmentStart(); + switch (start) { + case 100: + Assert.assertEquals(read.getMateAlignmentStart(), 1000); + break; + case 1000: + Assert.assertEquals(read.getMateAlignmentStart(), 2000); + break; + case 2000: + Assert.assertEquals(read.getMateAlignmentStart(), 1000); + break; + default: + Assert.assertTrue(false, "We saw a read located at the wrong position"); + } + } + } + + private class DummyWriter implements SAMFileWriter { + + public List reads; + + public DummyWriter() { reads = new ArrayList<>(10); } + + public void addAlignment(final SAMRecord alignment) { reads.add(alignment);} + + public SAMFileHeader getFileHeader() { return null; } + + public void setProgressLogger(final ProgressLoggerInterface progress) {} + + public void close() {} + } } \ No newline at end of file