From d976aae2b141567f8a91f53fe45b4925fff797c3 Mon Sep 17 00:00:00 2001
From: Eric Banks <ebanks@broadinstitute.org>
Date: Fri, 21 Jun 2013 16:59:22 -0400
Subject: [PATCH] Another fix for the Indel Realigner that arises because of
 secondary alignments.

This time we don't accidentally drop reads (phew), but this bug does cause us not to
update the alignment start of the mate.  Fixed and added unit test to cover it.
---
 .../indels/ConstrainedMateFixingManager.java  |  8 ++++--
 .../ConstrainedMateFixingManagerUnitTest.java | 28 ++++++++++++++++++-
 2 files changed, 32 insertions(+), 4 deletions(-)

diff --git a/protected/java/src/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManager.java b/protected/java/src/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManager.java
index c98fe4d3c..4d50ef951 100644
--- a/protected/java/src/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManager.java
+++ b/protected/java/src/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManager.java
@@ -130,7 +130,7 @@ public class ConstrainedMateFixingManager {
     private static final boolean DEBUG = false;
 
     /** How often do we check whether we want to emit reads? */
-    private final static int EMIT_FREQUENCY = 1000;
+    protected final static int EMIT_FREQUENCY = 1000;
 
     /**
      * How much could a single read move in position from its original position?
@@ -324,7 +324,8 @@ public class ConstrainedMateFixingManager {
                            || noReadCanMoveBefore(read.getMateAlignmentStart(), newRead ) ) ) { // we're already past where the mate started
 
                     // remove reads from the map that we have emitted -- useful for case where the mate never showed up
-                    forMateMatching.remove(read.getReadName());
+                    if ( !read.getNotPrimaryAlignmentFlag() )
+                        forMateMatching.remove(read.getReadName());
 
                     if ( DEBUG )
                         logger.warn(String.format("EMIT!  At %d: read %s at %d with isize %d, mate start %d, op = %s",
@@ -346,7 +347,8 @@ public class ConstrainedMateFixingManager {
 
     private void writeRead(SAMRecord read) {
         try {
-            writer.addAlignment(read);
+            if ( writer != null )
+                writer.addAlignment(read);
         } catch (IllegalArgumentException e) {
             throw new UserException("If the maximum allowable reads in memory is too small, it may cause reads to be written out of order when trying to write the BAM; please see the --maxReadsInMemory argument for details.  " + e.getMessage(), e);
         }
diff --git a/protected/java/test/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManagerUnitTest.java b/protected/java/test/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManagerUnitTest.java
index 9bcd7a3a3..0f910507e 100644
--- a/protected/java/test/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManagerUnitTest.java
+++ b/protected/java/test/org/broadinstitute/sting/gatk/walkers/indels/ConstrainedMateFixingManagerUnitTest.java
@@ -66,9 +66,10 @@ public class ConstrainedMateFixingManagerUnitTest extends BaseTest {
 
     @BeforeClass
     public void beforeClass() {
-        header = ArtificialSAMUtils.createArtificialSamHeader(3, 1, 100);
+        header = ArtificialSAMUtils.createArtificialSamHeader(3, 1, 10000);
         genomeLocParser = new GenomeLocParser(header.getSequenceDictionary());
     }
+
     @Test
     public void testSecondaryAlignmentsDoNotInterfere() {
         final List<GATKSAMRecord> properReads = ArtificialSAMUtils.createPair(header, "foo", 1, 10, 30, true, false);
@@ -105,4 +106,29 @@ public class ConstrainedMateFixingManagerUnitTest extends BaseTest {
         }
     }
 
+    @Test
+    public void testSecondaryAlignmentsDoNotCauseAccidentalRemovalOfMate() {
+        final List<GATKSAMRecord> properReads = ArtificialSAMUtils.createPair(header, "foo", 1, 530, 1594, true, false);
+        final GATKSAMRecord read1 = properReads.get(0);
+        read1.setFlags(99);   // first in proper pair, mate negative strand
+
+        final GATKSAMRecord read2Primary = properReads.get(1);
+        read2Primary.setFlags(147);   // second in pair, mate unmapped, not primary alignment
+        read2Primary.setAlignmentStart(1596); // move the read
+
+        final GATKSAMRecord read2NonPrimary = new GATKSAMRecord(read2Primary);
+        read2NonPrimary.setReadName("foo");
+        read2NonPrimary.setFlags(393);   // second in proper pair, on reverse strand
+        read2NonPrimary.setAlignmentStart(451);
+        read2NonPrimary.setMateAlignmentStart(451);
+
+        final ConstrainedMateFixingManager manager = new ConstrainedMateFixingManager(null, genomeLocParser, 10000, 200, 10000);
+        manager.addRead(read2NonPrimary, false, false);
+        manager.addRead(read1, false, false);
+
+        for ( int i = 0; i < ConstrainedMateFixingManager.EMIT_FREQUENCY; i++ )
+            manager.addRead(ArtificialSAMUtils.createArtificialRead(header, "foo" + i, 0, 1500, 10), false, false);
+
+        Assert.assertTrue(manager.forMateMatching.containsKey("foo"));
+    }
 }
\ No newline at end of file