diff --git a/c/bwa/bwa_gateway.cpp b/c/bwa/bwa_gateway.cpp index cdaa93c7d..538b8b16a 100644 --- a/c/bwa/bwa_gateway.cpp +++ b/c/bwa/bwa_gateway.cpp @@ -105,6 +105,8 @@ void BWA::copy_bases_into_sequence(bwa_seq_t* sequence, const char* bases, const } void BWA::create_alignments(bwa_seq_t* sequence, Alignment*& alignments, unsigned& num_alignments) { + bool debug = false; + num_alignments = 0; for(unsigned i = 0; i < (unsigned)sequence->n_aln; i++) num_alignments += (sequence->aln + i)->l - (sequence->aln + i)->k + 1; @@ -127,8 +129,22 @@ void BWA::create_alignments(bwa_seq_t* sequence, Alignment*& alignments, unsigne sequence->aln = &working_alignment_block; sequence->n_aln = 1; + sequence->sa = sa_idx; + sequence->strand = alignment_block->a; + sequence->score = alignment_block->score; + + // Each time through bwa_refine_gapped, seq gets reversed. Revert the reverse. + // TODO: Fix the interface to bwa_refine_gapped so its easier to work with. + if(alignment_idx > 0) + seq_reverse(sequence->len, sequence->seq, 0); + // Calculate the local alignment. bwa_cal_pac_pos_core(bwts[0],bwts[1],sequence,options.max_diff,options.fnr); + if(debug) { + printf("alignment_idx = %d, k = %d, l = %d, sa_idx = %d\n", alignment_idx, alignment_block->k, alignment_block->l, sa_idx); + printf("sequence->pos = %d\n",sequence->pos); + } + bwa_refine_gapped(bns, 1, sequence, 0, NULL); // Copy the local alignment data into the alignment object. diff --git a/java/src/org/broadinstitute/sting/alignment/bwa/AlignmentValidationWalker.java b/java/src/org/broadinstitute/sting/alignment/bwa/AlignmentValidationWalker.java index a2f8731fc..72c566f6a 100644 --- a/java/src/org/broadinstitute/sting/alignment/bwa/AlignmentValidationWalker.java +++ b/java/src/org/broadinstitute/sting/alignment/bwa/AlignmentValidationWalker.java @@ -25,6 +25,11 @@ public class AlignmentValidationWalker extends ReadWalker { */ private BWACAligner aligner = null; + /** + * Number of reads processed. + */ + private int count = 0; + /** * Create an aligner object. The aligner object will load and hold the BWT until close() is called. */ @@ -39,6 +44,15 @@ public class AlignmentValidationWalker extends ReadWalker { prefix + ".rsa"); } + /** Must return true for reads that need to be processed. Reads, for which this method return false will + * be skipped by the engine and never passed to the walker. + */ + @Override + public boolean filter(char[] ref, SAMRecord read) { + return true; + //return read.getReadName().contains("SL-XBC:1:76:604:397#0"); + } + /** * Aligns a read to the given reference. * @param ref Reference over the read. Read will most likely be unmapped, so ref will be null. @@ -47,6 +61,8 @@ public class AlignmentValidationWalker extends ReadWalker { */ @Override public Integer map(char[] ref, SAMRecord read) { + //logger.info(String.format("examining read %s", read.getReadName())); + byte[] bases = read.getReadBases(); if(read.getReadNegativeStrandFlag()) bases = BaseUtils.simpleReverseComplement(bases); @@ -85,6 +101,10 @@ public class AlignmentValidationWalker extends ReadWalker { throw new StingException(String.format("Read %s mismatches!", read.getReadName())); } + if(++count % 10000 == 0) logger.info(String.format("Processed %d reads", count)); + + //logger.info(String.format("validated read %s", read.getReadName())); + return 1; }