Picard modifications to support flow based sequencing

.gitignore

@@ -15,3 +15,4 @@ report
 jacoco.data
 .gradle
 build
+*.swp

src/main/java/picard/analysis/AlignmentSummaryMetrics.java

@@ -139,7 +139,12 @@ public enum Category {UNPAIRED, FIRST_OF_PAIR, SECOND_OF_PAIR, PAIR}
     /** The standard deviation of the read lengths. Computed using all read lengths including clipped bases. */
     public double SD_READ_LENGTH;
 
-    /** The median read length. Computed using all read lengths including clipped bases. */
+    /**
+     * The median read length of the set of reads examined.  When looking at the data for a single lane with
+     * equal length reads this number is just the read length.  When looking at data for merged lanes with
+     * differing read lengths this is the median read length of all reads. Computed using all bases in reads,
+     * including clipped bases.
+     */
     public double MEDIAN_READ_LENGTH;
 
     /**
@@ -155,6 +160,13 @@ public enum Category {UNPAIRED, FIRST_OF_PAIR, SECOND_OF_PAIR, PAIR}
     /** The maximum read length. Computed using all read lengths including clipped bases. */
     public double MAX_READ_LENGTH;
 
+    /**
+     * The mean aligned read length of the set of reads examined.  When looking at the data for a single lane with
+     * equal length reads this number is just the read length.  When looking at data for merged lanes with
+     * differing read lengths this is the mean read length of all reads. Clipped bases are not counted.
+     */
+    public double MEAN_ALIGNED_READ_LENGTH;
+
     /**
      * The number of aligned reads whose mate pair was also aligned to the reference.
      */

src/main/java/picard/analysis/AlignmentSummaryMetricsCollector.java

@@ -285,6 +285,7 @@ public void finish() {
                     metrics.PCT_PF_READS_ALIGNED = MathUtil.divide(metrics.PF_READS_ALIGNED, (double) metrics.PF_READS);
                     metrics.PCT_READS_ALIGNED_IN_PAIRS = MathUtil.divide(metrics.READS_ALIGNED_IN_PAIRS, (double) metrics.PF_READS_ALIGNED);
                     metrics.PCT_PF_READS_IMPROPER_PAIRS = MathUtil.divide(metrics.PF_READS_IMPROPER_PAIRS, (double) metrics.PF_READS_ALIGNED);
+                    metrics.MEAN_ALIGNED_READ_LENGTH = alignedReadLengthHistogram.getMean();
                     metrics.STRAND_BALANCE = MathUtil.divide(numPositiveStrand, (double) metrics.PF_READS_ALIGNED);
                     metrics.PCT_CHIMERAS = MathUtil.divide(chimeras, (double) chimerasDenominator);
                     metrics.PF_INDEL_RATE = MathUtil.divide(indels, (double) metrics.PF_ALIGNED_BASES);

src/main/java/picard/analysis/CollectQualityYieldMetrics.java

@@ -26,13 +26,13 @@
 
 import htsjdk.samtools.SAMFileHeader;
 import htsjdk.samtools.SAMRecord;
-import htsjdk.samtools.metrics.MetricBase;
 import htsjdk.samtools.metrics.MetricsFile;
 import htsjdk.samtools.reference.ReferenceSequence;
 import htsjdk.samtools.util.IOUtil;
 import org.broadinstitute.barclay.argparser.Argument;
 import org.broadinstitute.barclay.argparser.CommandLineProgramProperties;
 import org.broadinstitute.barclay.help.DocumentedFeature;
+import picard.PicardException;
 import picard.cmdline.StandardOptionDefinitions;
 import picard.cmdline.programgroups.DiagnosticsAndQCProgramGroup;
 import picard.util.help.HelpConstants;
@@ -92,6 +92,9 @@ public class CollectQualityYieldMetrics extends SinglePassSamProgram {
             "of bases if there are supplemental alignments in the input file.")
     public boolean INCLUDE_SUPPLEMENTAL_ALIGNMENTS = false;
 
+    @Argument(doc = "If true, calculates flow-specific READ_LENGTH_AVG_Q metrics.")
+    public boolean FLOW_MODE = false;
+
     /**
      * Ensure that we get all reads regardless of alignment status.
      */
@@ -103,7 +106,7 @@ protected boolean usesNoRefReads() {
     @Override
     protected void setup(final SAMFileHeader header, final File samFile) {
         IOUtil.assertFileIsWritable(OUTPUT);
-        this.collector = new QualityYieldMetricsCollector(USE_ORIGINAL_QUALITIES, INCLUDE_SECONDARY_ALIGNMENTS, INCLUDE_SUPPLEMENTAL_ALIGNMENTS);
+        this.collector = new QualityYieldMetricsCollector(USE_ORIGINAL_QUALITIES, INCLUDE_SECONDARY_ALIGNMENTS, INCLUDE_SUPPLEMENTAL_ALIGNMENTS, FLOW_MODE);
     }
 
     @Override
@@ -132,15 +135,30 @@ public static class QualityYieldMetricsCollector {
         // of bases if there are supplemental alignments in the input file.
         public final boolean includeSupplementalAlignments;
 
+        // If true collects RLQ25/RLQ30
+        private final boolean flowMode;
         // The metrics to be accumulated
-        private final QualityYieldMetrics metrics = new QualityYieldMetrics();
+        private final QualityYieldMetrics metrics;
+
+        public QualityYieldMetricsCollector(final boolean useOriginalQualities,
+                                            final boolean includeSecondaryAlignments,
+                                            final boolean includeSupplementalAlignments){
+            this(useOriginalQualities, includeSecondaryAlignments, includeSupplementalAlignments, false);
+        }
 
         public QualityYieldMetricsCollector(final boolean useOriginalQualities,
                                             final boolean includeSecondaryAlignments,
-                                            final boolean includeSupplementalAlignments) {
+                                            final boolean includeSupplementalAlignments,
+                                            final boolean flowMode) {
             this.useOriginalQualities = useOriginalQualities;
             this.includeSecondaryAlignments = includeSecondaryAlignments;
             this.includeSupplementalAlignments = includeSupplementalAlignments;
+            this.flowMode = flowMode;
+            if (flowMode){
+                this.metrics = new QualityYieldMetricsFlow(useOriginalQualities);
+            } else {
+                this.metrics = new QualityYieldMetrics(useOriginalQualities);
+            }
         }
 
         public void acceptRecord(final SAMRecord rec, final ReferenceSequence ref) {
@@ -187,12 +205,15 @@ public void acceptRecord(final SAMRecord rec, final ReferenceSequence ref) {
                     }
                 }
             }
+
+            if (flowMode) {
+                ((QualityYieldMetricsFlow)metrics).addRecordToHistogramGenerator(rec);
+            }
         }
 
         public void finish() {
             metrics.Q20_EQUIVALENT_YIELD = metrics.Q20_EQUIVALENT_YIELD / 20;
             metrics.PF_Q20_EQUIVALENT_YIELD = metrics.PF_Q20_EQUIVALENT_YIELD / 20;
-
             metrics.calculateDerivedFields();
         }
 
@@ -201,12 +222,69 @@ public void addMetricsToFile(final MetricsFile<QualityYieldMetrics, Integer> met
         }
     }
 
+    public static class QualityYieldMetricsFlow extends QualityYieldMetrics{
+        /** The length of the longest interval on the reads where the average quality per-base is above (Q30) */
+        @NoMergingIsDerived
+        public long READ_LENGTH_AVG_Q_ABOVE_30 = 0;
+
+        /** The length of the longest interval on the reads where the average quality per-base is above (Q25) */
+        @NoMergingIsDerived
+        public long READ_LENGTH_AVG_Q_ABOVE_25 = 0;
+
+        @MergingIsManual
+        protected final HistogramGenerator histogramGenerator;
+
+        public QualityYieldMetricsFlow(){
+            this(false);
+        }
+
+        public QualityYieldMetricsFlow(final boolean useOriginalBaseQualities){
+
+            super(useOriginalBaseQualities);
+            histogramGenerator=new HistogramGenerator(useOriginalQualities);
+        }
+
+        public QualityYieldMetricsFlow(final boolean useOriginalBaseQualities, final HistogramGenerator hg) {
+            histogramGenerator=hg;
+        }
+
+        @Override
+        public void calculateDerivedFields() {
+            super.calculateDerivedFields();
+            this.READ_LENGTH_AVG_Q_ABOVE_25 = histogramGenerator.calculateLQ(25, 1, 5);
+            this.READ_LENGTH_AVG_Q_ABOVE_30 = histogramGenerator.calculateLQ(30, 1, 5);
+        }
+
+        @Override
+        public MergeableMetricBase merge(final MergeableMetricBase other) {
+            if (!(other instanceof QualityYieldMetricsFlow)){
+                throw new PicardException("Only objects of the same type can be merged");
+            }
+            this.histogramGenerator.addOtherHistogramGenerator(((QualityYieldMetricsFlow)other).histogramGenerator);
+            super.merge(other);
+            return this;
+        }
+
+        protected void addRecordToHistogramGenerator(final SAMRecord rec) {
+            histogramGenerator.addRecord(rec);
+        }
+
+    }
     /**
      * A set of metrics used to describe the general quality of a BAM file
      */
     @DocumentedFeature(groupName = HelpConstants.DOC_CAT_METRICS, summary = HelpConstants.DOC_CAT_METRICS_SUMMARY)
     public static class QualityYieldMetrics extends MergeableMetricBase {
 
+        public QualityYieldMetrics() {
+            this(false);
+        }
+
+        public QualityYieldMetrics(final boolean useOriginalQualities) {
+            super();
+            this.useOriginalQualities = useOriginalQualities;
+        }
+
         /**
          * The total number of reads in the input file
          */
@@ -220,7 +298,7 @@ public static class QualityYieldMetrics extends MergeableMetricBase {
         public long PF_READS = 0;
 
         /**
-         * The average read length of all the reads (will be fixed for a lane)
+         * The average read length of all the reads
          */
         @NoMergingIsDerived
         public int READ_LENGTH = 0;
@@ -273,12 +351,26 @@ public static class QualityYieldMetrics extends MergeableMetricBase {
         @MergeByAdding
         public long PF_Q20_EQUIVALENT_YIELD = 0;
 
+        @MergeByAssertEquals
+        protected final boolean useOriginalQualities;
+
         @Override
         public void calculateDerivedFields() {
             super.calculateDerivedFields();
-
             this.READ_LENGTH = this.TOTAL_READS == 0 ? 0 : (int) (this.TOTAL_BASES / this.TOTAL_READS);
         }
-    }
 
+        @Override
+        public MergeableMetricBase merge(final MergeableMetricBase other) {
+            if (!(other instanceof QualityYieldMetrics)){
+                throw new PicardException("Only objects of the same type can be merged");
+            }
+
+            final QualityYieldMetrics otherMetric = (QualityYieldMetrics) other;
+
+            super.merge(otherMetric);
+            calculateDerivedFields();
+            return this;
+        }
+    }
 }