Merge branch 'fnothaft-440-fix-indelrealigner' into ADAM-439

* fnothaft-440-fix-indelrealigner: Fixed flipped index in indel realigner. [ADAM-439] Fix ADAM to account for BDG-FORMATS-35
bigdatagenomics · Nov 17, 2014 · e12c48d · e12c48d
2 parents 65d18a9 + 16655b1
commit e12c48d
Show file tree

Hide file tree

Showing 3 changed files with 18 additions and 15 deletions.
diff --git a/adam-core/src/main/scala/org/bdgenomics/adam/converters/VCFLineConverter.scala b/adam-core/src/main/scala/org/bdgenomics/adam/converters/VCFLineConverter.scala
@@ -128,7 +128,7 @@ object VCFLineConverter {
         None
       } else {
         val gts = sampleFieldMap("GT").split("\\||/").map(_.toInt)
-	val genotypes: Seq[String] = gts.map(idx => line.alleleArray(idx))
+        val genotypes: Seq[String] = gts.map(idx => line.alleleArray(idx))
         val sampleId = line.samples(i)
 
         val flatGenotype = FlatGenotype.newBuilder()

diff --git a/adam-core/src/main/scala/org/bdgenomics/adam/models/SequenceDictionary.scala b/adam-core/src/main/scala/org/bdgenomics/adam/models/SequenceDictionary.scala
@@ -186,12 +186,12 @@ object SequenceRecord {
 
   def apply(name: String,
             length: Long,
-	    md5: String = null,
-	    url: String = null,
-	    refseq: String = null,
-	    genbank: String = null,
-	    assembly: String = null,
-	    species: String = null): SequenceRecord = {
+            md5: String = null,
+            url: String = null,
+            refseq: String = null,
+            genbank: String = null,
+            assembly: String = null,
+            species: String = null): SequenceRecord = {
     new SequenceRecord(
       name,
       length,

diff --git a/adam-core/src/main/scala/org/bdgenomics/adam/rdd/read/realignment/RealignIndels.scala b/adam-core/src/main/scala/org/bdgenomics/adam/rdd/read/realignment/RealignIndels.scala
@@ -326,6 +326,7 @@ private[rdd] class RealignIndels(val consensusModel: ConsensusGenerator = new Co
 
             // if read alignment is improved by aligning against new consensus, realign
             if (remapping != -1) {
+
               realignedReadCount += 1
 
               // bump up mapping quality by 10
@@ -335,24 +336,26 @@ private[rdd] class RealignIndels(val consensusModel: ConsensusGenerator = new Co
               builder.setStart(refStart + remapping)
 
               // recompute cigar
-              val newCigar: Cigar = if (refStart + remapping >= bestConsensus.index.start && refStart + remapping <= bestConsensus.index.end - 1) {
+              val newCigar: Cigar = {
                 // if element overlaps with consensus indel, modify cigar with indel
-                val (idElement, endLength) = if (bestConsensus.index.start == bestConsensus.index.end - 1) {
+                val (idElement, endLength, endPenalty) = if (bestConsensus.index.start == bestConsensus.index.end - 1) {
                   (new CigarElement(bestConsensus.consensus.length, CigarOperator.I),
-                    r.getSequence.length - bestConsensus.consensus.length - (bestConsensus.index.start - (refStart + remapping)))
+                    r.getSequence.length - bestConsensus.consensus.length - (bestConsensus.index.start - (refStart + remapping)),
+                    -bestConsensus.consensus.length)
                 } else {
                   (new CigarElement((bestConsensus.index.end - 1 - bestConsensus.index.start).toInt, CigarOperator.D),
-                    r.getSequence.length - (bestConsensus.index.start - (refStart + remapping)))
+                    r.getSequence.length - (bestConsensus.index.start - (refStart + remapping)),
+                    bestConsensus.consensus.length)
                 }
 
-                val cigarElements = List[CigarElement](new CigarElement((refStart + remapping - bestConsensus.index.start).toInt, CigarOperator.M),
+                // compensate the end
+                builder.setEnd(refStart + remapping + r.getSequence.length + endPenalty)
+
+                val cigarElements = List[CigarElement](new CigarElement((bestConsensus.index.start - (refStart + remapping)).toInt, CigarOperator.M),
                   idElement,
                   new CigarElement(endLength.toInt, CigarOperator.M))
 
                 new Cigar(cigarElements)
-              } else {
-                // else, new cigar is all matches
-                new Cigar(List[CigarElement](new CigarElement(r.getSequence.length, CigarOperator.M)))
               }
 
               // update mdtag and cigar