Java Code Examples for htsjdk.samtools.SAMRecordSetBuilder#addPair()

private List<SAMRecord> getPairedRead (final String name, final int contig, final int start1, final int start2) {
	List<SAMRecord> result = new ArrayList<> ();

	SAMRecordSetBuilder builder = new SAMRecordSetBuilder();
	builder.addPair(name, contig, start1, start2);

	Collection<SAMRecord> recs = builder.getRecords();

	for (SAMRecord r: recs) {
		if (r.getFirstOfPairFlag()) result.add(0, r);
		if (r.getSecondOfPairFlag()) result.add(1, r);
		r.setMappingQuality(10);
	}

	return (result);

}
 

@Test
public void testBasicPairedRead() {
    SAMRecordSetBuilder builder = new SAMRecordSetBuilder(false, SAMFileHeader.SortOrder.queryname);
    builder.setReadLength(READ_LENGTH);
    builder.addPair("mapped_paired", 1, 1, 31);
    PeekableIterator<SAMRecord> iterator = new PeekableIterator<SAMRecord>(builder.iterator());

    QuerySortedReadPairIteratorUtil.ReadPair pair = QuerySortedReadPairIteratorUtil.getNextReadPair(iterator);
    Assert.assertNotNull(pair);
    Assert.assertNotNull(pair.read1);
    Assert.assertNotNull(pair.read2);
    Assert.assertEquals("mapped_paired", pair.read1.getReadName());
    Assert.assertEquals("mapped_paired", pair.read2.getReadName());

    pair = QuerySortedReadPairIteratorUtil.getNextReadPair(iterator);
    Assert.assertNull(pair);
}
 

@Test
public void testBasicHalfmappedReadPair() {
    SAMRecordSetBuilder builder = new SAMRecordSetBuilder(false, SAMFileHeader.SortOrder.queryname);
    builder.setReadLength(READ_LENGTH);
    builder.addPair("halfmapped_paired", 1, 1, 31, false, true, "20M", "20M", true, false, 20);
    PeekableIterator<SAMRecord> iterator = new PeekableIterator<SAMRecord>(builder.iterator());

    QuerySortedReadPairIteratorUtil.ReadPair pair = QuerySortedReadPairIteratorUtil.getNextReadPair(iterator);
    Assert.assertNotNull(pair);
    Assert.assertNotNull(pair.read1);
    Assert.assertNotNull(pair.read2);
    Assert.assertEquals("halfmapped_paired", pair.read1.getReadName());
    Assert.assertEquals("halfmapped_paired", pair.read2.getReadName());

    pair = QuerySortedReadPairIteratorUtil.getNextReadPair(iterator);
    Assert.assertNull(pair);
}
 

@Test
public void testBaseErrorAggregation() {
    final SAMSequenceRecord samSequenceRecord = new SAMSequenceRecord("chr1", 200);
    final SAMFileHeader samFileHeader = new SAMFileHeader();
    samFileHeader.addSequence(samSequenceRecord);

    final SAMRecordSetBuilder builder = new SAMRecordSetBuilder();

    final List<SAMRecord> samRecords = builder.addPair("Read1234", 0, 1, 1,
            false, false, "16M", "16M", true, false, 20);
    final SAMRecord samRecord1 = samRecords.get(0);
    final SAMRecord samRecord2 = samRecords.get(1);

    samRecord1.setReadBases("CgTGtGGAcAAAgAAA".getBytes());
    samRecord2.setReadBases("CcTGGtGAcAAAgAAA".getBytes());
    final byte[] refBases = "CATGGGGAAAAAAAAA".getBytes();

    BaseErrorAggregation<?> baseErrorAggregation = new BaseErrorAggregation<>(OverlappingReadsErrorCalculator::new, ReadBaseStratification.readOrdinalityStratifier);

    final int length = getLengthAndAddBases(samSequenceRecord, samRecord1, samRecord2, refBases, baseErrorAggregation);

    final ErrorMetric[] metrics = baseErrorAggregation.getMetrics();
    final OverlappingErrorMetric metric1 = (OverlappingErrorMetric) metrics[0];
    metric1.calculateDerivedFields();
    Assert.assertEquals(metric1.COVARIATE, ReadBaseStratification.ReadOrdinality.FIRST.name());
    Assert.assertEquals(metric1.TOTAL_BASES, length);
    Assert.assertEquals(metric1.NUM_BASES_WITH_OVERLAPPING_READS, length);
    Assert.assertEquals(metric1.NUM_DISAGREES_WITH_REFERENCE_ONLY, 2L);
    Assert.assertEquals(metric1.NUM_DISAGREES_WITH_REF_AND_MATE, 1L);
    Assert.assertEquals(metric1.NUM_THREE_WAYS_DISAGREEMENT, 1L);

    final OverlappingErrorMetric metric2 = (OverlappingErrorMetric) metrics[1];
    metric2.calculateDerivedFields();
    Assert.assertEquals(metric2.COVARIATE, ReadBaseStratification.ReadOrdinality.SECOND.name());
    Assert.assertEquals(metric2.TOTAL_BASES, length);
    Assert.assertEquals(metric2.NUM_BASES_WITH_OVERLAPPING_READS, length);
    Assert.assertEquals(metric2.NUM_DISAGREES_WITH_REFERENCE_ONLY, 2L);
    Assert.assertEquals(metric2.NUM_DISAGREES_WITH_REF_AND_MATE, 1L);
    Assert.assertEquals(metric2.NUM_THREE_WAYS_DISAGREEMENT, 1L);
}
 

@Test
public void testBaseErrorAggregation2() {
    final SAMSequenceRecord samSequenceRecord = new SAMSequenceRecord("chr1", 200);
    final SAMFileHeader samFileHeader = new SAMFileHeader();
    samFileHeader.addSequence(samSequenceRecord);

    final SAMRecordSetBuilder builder = new SAMRecordSetBuilder();

    final List<SAMRecord> samRecords = builder.addPair("Read1234", 0, 1, 1,
            false, false, "16M", "16M", true, false, 20);
    final SAMRecord samRecord1 = samRecords.get(0);
    final SAMRecord samRecord2 = samRecords.get(1);

    samRecord1.setReadBases("CgTGtGGAcAAAgAAA".getBytes());
    samRecord2.setReadBases("CcTGGtGAcAAAgAAA".getBytes());
    final byte[] refBases = "CATGGGGAAAAAAAAA".getBytes();

    BaseErrorAggregation<?> baseErrorAggregation = new BaseErrorAggregation<>(OverlappingReadsErrorCalculator::new, ReadBaseStratification.readDirectionStratifier);

    final int length = getLengthAndAddBases(samSequenceRecord, samRecord1, samRecord2, refBases, baseErrorAggregation);

    final ErrorMetric[] metrics = baseErrorAggregation.getMetrics();
    final OverlappingErrorMetric metric1 = (OverlappingErrorMetric) metrics[0];
    metric1.calculateDerivedFields();
    Assert.assertEquals(metric1.COVARIATE, ReadBaseStratification.ReadDirection.POSITIVE.toString());
    Assert.assertEquals(metric1.TOTAL_BASES, length);
    Assert.assertEquals(metric1.NUM_BASES_WITH_OVERLAPPING_READS, length);
    Assert.assertEquals(metric1.NUM_DISAGREES_WITH_REFERENCE_ONLY, 2L);
    Assert.assertEquals(metric1.NUM_DISAGREES_WITH_REF_AND_MATE, 1L);
    Assert.assertEquals(metric1.NUM_THREE_WAYS_DISAGREEMENT, 1L);

    final OverlappingErrorMetric metric2 = (OverlappingErrorMetric) metrics[1];
    metric2.calculateDerivedFields();
    Assert.assertEquals(metric2.COVARIATE, ReadBaseStratification.ReadDirection.NEGATIVE.toString());
    Assert.assertEquals(metric2.TOTAL_BASES, length);
    Assert.assertEquals(metric2.NUM_BASES_WITH_OVERLAPPING_READS, length);
    Assert.assertEquals(metric2.NUM_DISAGREES_WITH_REFERENCE_ONLY, 2L);
    Assert.assertEquals(metric2.NUM_DISAGREES_WITH_REF_AND_MATE, 1L);
    Assert.assertEquals(metric2.NUM_THREE_WAYS_DISAGREEMENT, 1L);
}
 

public void setupTest1(final int ID, final String readGroupId, final SAMReadGroupRecord readGroupRecord, final String sample,
                  final String library, final SAMFileHeader header, final SAMRecordSetBuilder setBuilder)
        throws IOException {

    final String separator = ":";
    final int contig1 = 0;
    final int contig2 = 1;
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    header.addReadGroup(readGroupRecord);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);

    setBuilder.setHeader(header);

    final int max = 800;
    final int min = 1;
    final Random rg = new Random(5);

    //add records that align to chrM and O but not N
    for (int i = 0; i < NUM_READS; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = READ_NAME + separator + ID + separator + i;

        if (i != NUM_READS - 1) {
            setBuilder.addPair(newReadName, contig1, start + ID, start + ID + LENGTH);
        } else {
            setBuilder.addPair(newReadName, contig2, start + ID, start + ID + LENGTH);
        }
    }
}
 

public void setupTest2(final int ID, final String readGroupId, final SAMReadGroupRecord readGroupRecord, final String sample,
                       final String library, final SAMFileHeader header, final SAMRecordSetBuilder setBuilder)
        throws IOException {

    final String separator = ":";
    final int contig1 = 0;
    final int contig2 = 1;
    final int contig3 = 2;
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);

    setBuilder.setHeader(header);

    final int max = 800;
    final int min = 1;
    final Random rg = new Random(5);

    //add records that align to all 3 chr in reference file
    for (int i = 0; i < NUM_READS; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = READ_NAME + separator + ID + separator + i;

        if (i<=NUM_READS/3) {
            setBuilder.addPair(newReadName, contig1, start + ID, start + ID + LENGTH);
        } else if (i< (NUM_READS - (NUM_READS/3))) {
            setBuilder.addPair(newReadName, contig2, start + ID, start + ID + LENGTH);
        } else {
            setBuilder.addPair(newReadName, contig3, start + ID, start + ID + LENGTH);
        }
    }
}
 

void setup(final int numReads,
           final String readName,
           final int ID,
           final String readGroupId,
           final SAMReadGroupRecord readGroupRecord,
           final String sample,
           final String library,
           final SAMFileHeader header,
           final SAMRecordSetBuilder setBuilder) throws IOException {
    final String separator = ":";
    readGroupRecord.setSample(sample);
    readGroupRecord.setPlatform(platform);
    readGroupRecord.setLibrary(library);
    readGroupRecord.setPlatformUnit(readGroupId);
    header.addReadGroup(readGroupRecord);
    setBuilder.setReadGroup(readGroupRecord);
    setBuilder.setUseNmFlag(true);
    setBuilder.setHeader(header);

    final int max = 15000;
    final int min = 1;
    final Random rg = new Random(5);

    for (int i = 0; i < numReads; i++) {
        final int start = rg.nextInt(max) + min;
        final String newReadName = readName + separator + ID + separator + i;
        setBuilder.addPair(newReadName, 0, start + ID, start + ID + 99);
    }
}
 

public static File writeBamFileWithLargeHeader() throws IOException {
  SAMRecordSetBuilder samRecordSetBuilder =
      new SAMRecordSetBuilder(true, SAMFileHeader.SortOrder.queryname);
  for (int i = 0; i < 1000; i++) {
    int chr = 20;
    int start1 = (i + 1) * 1000;
    int start2 = start1 + 100;
    samRecordSetBuilder.addPair(String.format("test-read-%03d", i), chr, start1,
        start2);
  }

  final File bamFile = File.createTempFile("test", ".bam");
  bamFile.deleteOnExit();
  SAMFileHeader samHeader = samRecordSetBuilder.getHeader();
  StringBuffer sb = new StringBuffer();
  for (int i = 0; i < 1000000; i++) {
    sb.append("0123456789012345678901234567890123456789012345678901234567890123456789012345678901234567890123456789");
  }
  samHeader.addComment(sb.toString());
  final SAMFileWriter bamWriter = new SAMFileWriterFactory()
      .makeSAMOrBAMWriter(samHeader, true, bamFile);
  for (final SAMRecord rec : samRecordSetBuilder.getRecords()) {
    bamWriter.addAlignment(rec);
  }
  bamWriter.close();

  return bamFile;
}
 

@Test(dataProvider = "overlapReadData")
public void testOverlappedReadClipping(final int readLength, final int start1, final int start2, final String cigar1, final String cigar2,
                                       final boolean strand1, final boolean strand2,
                                       final int r1ExpectedAlignmentStart, final int r2ExpectedAlignmentStart,
                                       final String expectedR1Cigar, final String expectedR2Cigar, final boolean hardClipOverlappingReads,
                                       final String read1Bases, final String read2Bases, final String expectedR1Bases, final String expectedR2Bases,
                                       final String expectedR1ClippedBases, final String expectedR2ClippedBases, final String read1Qualities,
                                       final String read2Qualities, final String expectedR1Qualities, final String expectedR2Qualities,
                                       final String expectedR1ClippedQualities, final String expectedR2ClippedQualities) {

    final SAMRecordSetBuilder set = new SAMRecordSetBuilder();
    set.setReadLength(readLength);
    final List<SAMRecord> recs = set.addPair("q1", 0, start1, start2, false, false, cigar1, cigar2, strand1, strand2, 30);
    final SAMRecord r1 = recs.get(0);
    final SAMRecord r2 = recs.get(1);

    r1.setReadBases(StringUtil.stringToBytes(read1Bases));
    r2.setReadBases(StringUtil.stringToBytes(read2Bases));

    r1.setBaseQualities(SAMUtils.fastqToPhred(read1Qualities));
    r2.setBaseQualities(SAMUtils.fastqToPhred(read2Qualities));

    AbstractAlignmentMerger.clipForOverlappingReads(r1, r2, hardClipOverlappingReads);
    Assert.assertEquals(r1.getAlignmentStart(), r1ExpectedAlignmentStart);
    Assert.assertEquals(r1.getCigarString(), expectedR1Cigar);
    Assert.assertEquals(r2.getAlignmentStart(), r2ExpectedAlignmentStart);
    Assert.assertEquals(r2.getCigarString(), expectedR2Cigar);
    Assert.assertEquals(r1.getReadString(), expectedR1Bases);
    Assert.assertEquals(r2.getReadString(), expectedR2Bases);
    Assert.assertEquals(SAMUtils.phredToFastq(r1.getBaseQualities()), expectedR1Qualities);
    Assert.assertEquals(SAMUtils.phredToFastq(r2.getBaseQualities()), expectedR2Qualities);

    Assert.assertEquals(r1.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASES_TAG), expectedR1ClippedBases);
    Assert.assertEquals(r2.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASES_TAG), expectedR2ClippedBases);

    Assert.assertEquals(r1.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASE_QUALITIES_TAG), expectedR1ClippedQualities);
    Assert.assertEquals(r2.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASE_QUALITIES_TAG), expectedR2ClippedQualities);


    // Swap first and second read to ensure logic is correct for both F1R2 and F2R1
    final SAMRecordSetBuilder setSwapped = new SAMRecordSetBuilder();
    setSwapped.setReadLength(readLength);
    final List<SAMRecord> recsSwapped = set.addPair("q1", 0, start2, start1, false, false, cigar2, cigar1, strand2, strand1, 30);
    final SAMRecord r1Swapped = recsSwapped.get(0);
    final SAMRecord r2Swapped = recsSwapped.get(1);

    r1Swapped.setReadBases(StringUtil.stringToBytes(read2Bases));
    r2Swapped.setReadBases(StringUtil.stringToBytes(read1Bases));

    r1Swapped.setBaseQualities(SAMUtils.fastqToPhred(read2Qualities));
    r2Swapped.setBaseQualities(SAMUtils.fastqToPhred(read1Qualities));

    AbstractAlignmentMerger.clipForOverlappingReads(r1Swapped, r2Swapped, hardClipOverlappingReads);
    Assert.assertEquals(r1Swapped.getAlignmentStart(), r2ExpectedAlignmentStart);
    Assert.assertEquals(r1Swapped.getCigarString(), expectedR2Cigar);
    Assert.assertEquals(r2Swapped.getAlignmentStart(), r1ExpectedAlignmentStart);
    Assert.assertEquals(r2Swapped.getCigarString(), expectedR1Cigar);
    Assert.assertEquals(r1Swapped.getReadString(), expectedR2Bases);
    Assert.assertEquals(r2Swapped.getReadString(), expectedR1Bases);
    Assert.assertEquals(SAMUtils.phredToFastq(r1Swapped.getBaseQualities()), expectedR2Qualities);
    Assert.assertEquals(SAMUtils.phredToFastq(r2Swapped.getBaseQualities()), expectedR1Qualities);

    Assert.assertEquals(r1Swapped.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASES_TAG), expectedR2ClippedBases);
    Assert.assertEquals(r2Swapped.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASES_TAG), expectedR1ClippedBases);

    Assert.assertEquals(r1Swapped.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASE_QUALITIES_TAG), expectedR2ClippedQualities);
    Assert.assertEquals(r2Swapped.getAttribute(AbstractAlignmentMerger.HARD_CLIPPED_BASE_QUALITIES_TAG), expectedR1ClippedQualities);

}
 

@Test
    public void testUnmapBacterialContamination() throws IOException {
        final SAMRecordSetBuilder builder = new SAMRecordSetBuilder(true, SAMFileHeader.SortOrder.queryname);
        final SAMFileHeader header = builder.getHeader();
        final SAMFileHeader.SortOrder sortOrder = header.getSortOrder();
        final SAMFileHeader newHeader = SAMRecordSetBuilder.makeDefaultHeader(sortOrder, 100000,true);
        builder.setHeader(newHeader);

        final File reference = File.createTempFile("reference",".fasta");
        reference.deleteOnExit();

        builder.writeRandomReference(reference.toPath());

        builder.addPair("overlappingpair", 0,500,500, false,false,"20S20M60S","20S20M60M",true,false,45);
        builder.addPair("overlappingpairFirstAligned", 0,500,500, false,true,"20S20M60S",null,true,false,45);
        builder.addPair("overlappingpairSecondAligned", 0,500,500, true,false,null,"20S20M60S",true,false,45);
        builder.addPair("overlappingpairFirstAlignedB", 0,500,500, false,true,"20S20M60S",null,false,true,45);
        builder.addPair("overlappingpairSecondAlignedB", 0,500,500, true,false,null,"20S20M60S",false,true,45);

//        builder.addFrag("frag",1,500,false,false,"20S20M60S",null, 45);
//        builder.addFrag("frag2",1,500,true,false,"20S20M60S",null, 45);
//        builder.addFrag("frag3",1,500,false,false,"20S20M60S",null, 45);
//        builder.addFrag("frag4",1,500,true,false,"20S20M60S",null, 45);

        final File file = newTempSamFile("aligned");

        try (SAMFileWriter writer = new SAMFileWriterFactory().makeWriter(builder.getHeader(), true, file, null)) {
            builder.getRecords().forEach(writer::addAlignment);
        }

        final RevertSam revertSam = new RevertSam();

        revertSam.INPUT = file;
        final File fileUnaligned = newTempSamFile("unaligned");
        revertSam.OUTPUT = fileUnaligned;

        revertSam.SANITIZE = false;
        revertSam.REMOVE_ALIGNMENT_INFORMATION=true;
        revertSam.REMOVE_DUPLICATE_INFORMATION=true;

        revertSam.SORT_ORDER = SAMFileHeader.SortOrder.queryname;

        Assert.assertEquals(revertSam.doWork(),0);

        MergeBamAlignment mergeBamAlignment = new MergeBamAlignment();

        mergeBamAlignment.ALIGNED_BAM = Collections.singletonList(file);
        mergeBamAlignment.UNMAPPED_BAM = fileUnaligned;
        mergeBamAlignment.UNMAP_CONTAMINANT_READS = true;

        //yuck!
        final RequiredReferenceArgumentCollection requiredReferenceArgumentCollection = new RequiredReferenceArgumentCollection();
        requiredReferenceArgumentCollection.REFERENCE_SEQUENCE = reference;
        mergeBamAlignment.referenceSequence = requiredReferenceArgumentCollection;

        final File fileMerged = newTempSamFile("merged");

        mergeBamAlignment.OUTPUT = fileMerged;

        // merge file with itself.
        Assert.assertEquals(mergeBamAlignment.doWork(),0);

        // check that all reads have been unmapped due to bacterial contamination as needed.
        try (SamReader mergedReader = SamReaderFactory.makeDefault().open(fileMerged)) {
            for (SAMRecord mergedRecord : mergedReader) {
                Assert.assertTrue(mergedRecord.getReadUnmappedFlag(), mergedRecord.toString());
                Assert.assertTrue(!mergedRecord.getReadPairedFlag() || mergedRecord.getMateUnmappedFlag(), mergedRecord.toString());
            }
        }
    }
 

@Test
public void testBaseErrorAggregation3() {
    final SAMSequenceRecord samSequenceRecord = new SAMSequenceRecord("chr1", 200);
    final SAMFileHeader samFileHeader = new SAMFileHeader();
    samFileHeader.addSequence(samSequenceRecord);

    final SAMRecordSetBuilder builder = new SAMRecordSetBuilder();

    final List<SAMRecord> samRecords = builder.addPair("Read1234", 0, 1, 1,
            false, false, "16M", "16M", true, false, 20);
    final SAMRecord samRecord1 = samRecords.get(0);
    final SAMRecord samRecord2 = samRecords.get(1);

    samRecord1.setReadBases("CgTGtGGAcAAAgAAA".getBytes());
    samRecord2.setReadBases("CcTGGtGAcAAAgAAA".getBytes());
    final byte[] refBases = "CATGGGGAAAAAAAAA".getBytes();

    BaseErrorAggregation<?> baseErrorAggregation = new BaseErrorAggregation<>(SimpleErrorCalculator::new, ReadBaseStratification.referenceBaseStratifier);

    final int length = getLengthAndAddBases(samSequenceRecord, samRecord1, samRecord2, refBases, baseErrorAggregation);

    final ErrorMetric[] metrics = baseErrorAggregation.getMetrics();
    final BaseErrorMetric metricA = Arrays.stream(metrics).map(a -> (BaseErrorMetric) a)
            .filter(a -> a.COVARIATE.equals("A")).findFirst().get();
    metricA.calculateDerivedFields();
    Assert.assertEquals(metricA.COVARIATE, "A");
    Assert.assertEquals(metricA.TOTAL_BASES, 11);
    Assert.assertEquals(metricA.ERROR_BASES, 3);

    final BaseErrorMetric metricC = Arrays.stream(metrics).map(a -> (BaseErrorMetric) a)
            .filter(a -> a.COVARIATE.equals("C")).findFirst().get();
    metricA.calculateDerivedFields();
    Assert.assertEquals(metricC.COVARIATE, "C");
    Assert.assertEquals(metricC.TOTAL_BASES, 5);
    Assert.assertEquals(metricC.ERROR_BASES, 1);

    final BaseErrorMetric metricG = Arrays.stream(metrics).map(a -> (BaseErrorMetric) a)
            .filter(a -> a.COVARIATE.equals("G")).findFirst().get();
    metricA.calculateDerivedFields();
    Assert.assertEquals(metricG.COVARIATE, "G");
    Assert.assertEquals(metricG.TOTAL_BASES, 5);
    Assert.assertEquals(metricG.ERROR_BASES, 1);

    final BaseErrorMetric metricT = Arrays.stream(metrics).map(a -> (BaseErrorMetric) a)
            .filter(a -> a.COVARIATE.equals("T")).findFirst().get();
    metricA.calculateDerivedFields();
    Assert.assertEquals(metricT.COVARIATE, "T");
    Assert.assertEquals(metricT.TOTAL_BASES, 11);
    Assert.assertEquals(metricT.ERROR_BASES, 3);

}
 

@Test
public void testWidthOfMetrics() throws IOException {
    final File testSamFile = File.createTempFile("CollectInsertSizeMetrics", ".bam");
    testSamFile.deleteOnExit();
    final File testSamFileIndex = new File(testSamFile.getParentFile(),testSamFile.getName().replaceAll("bam$","bai"));
    testSamFileIndex.deleteOnExit();


    new File(testSamFile.getAbsolutePath().replace(".bam$",".bai")).deleteOnExit();

    final SAMRecordSetBuilder setBuilder = new SAMRecordSetBuilder(true, SAMFileHeader.SortOrder.coordinate, true);
    setBuilder.setReadLength(10);

    final int insertBy = 3; // the # of bases to increase the insert by in the records below.
    int queryIndex = 0;

    // Create records such that we have 10 records in the 10th through the 90th percentiles, and 5 for the 95th and 99th percentiles.
    // WIDTH_OF_10_PERCENT through WIDTH_OF_90_PERCENT (90 pairs total))
    for (int j = 0; j < 9; j++) {
        for (int i = 0; i < 5; i++, queryIndex++) {
            setBuilder.addPair("query:" + queryIndex, 0, 1, 50 + j*insertBy, false, false, "10M", "10M", false, true, 60);
        }
        for (int i = 0; i < 5; i++, queryIndex++) {
            setBuilder.addPair("query:" + queryIndex, 0, 1, 50 - j*insertBy, false, false, "10M", "10M", false, true, 60);
        }
    }
    // WIDTH_OF_95_PERCENT through WIDTH_OF_99_PERCENT (10 pairs total)
    for (int j = 9; j < 11; j++) {
        for (int i = 0; i < 3; i++, queryIndex++) {
            setBuilder.addPair("query:" + queryIndex, 0, 1, 50 + j*insertBy, false, false, "10M", "10M", false, true, 60);
        }
        for (int i = 0; i < 2; i++, queryIndex++) {
            setBuilder.addPair("query:" + queryIndex, 0, 1, 50 - j*insertBy, false, false, "10M", "10M", false, true, 60);
        }
    }

    // Add one to make the an odd # of pairs for the median
    setBuilder.addPair("query:" + queryIndex, 0, 1, 50, false, false, "10M", "10M", false, true, 60);

    final SAMFileWriter writer = new SAMFileWriterFactory().setCreateIndex(true).makeBAMWriter(setBuilder.getHeader(), false, testSamFile);
    setBuilder.forEach(writer::addAlignment);
    writer.close();


    final File outfile = File.createTempFile("test", ".insert_size_metrics");
    final File pdf     = File.createTempFile("test", ".pdf");
    outfile.deleteOnExit();
    pdf.deleteOnExit();
    final String[] args = new String[] {
            "INPUT="  + testSamFile.getAbsolutePath(),
            "OUTPUT=" + outfile.getAbsolutePath(),
            "Histogram_FILE=" + pdf.getAbsolutePath()
    };
    Assert.assertEquals(runPicardCommandLine(args), 0);

    final MetricsFile<InsertSizeMetrics, Comparable<?>> output = new MetricsFile<InsertSizeMetrics, Comparable<?>>();
    output.read(new FileReader(outfile));
    
    final List<InsertSizeMetrics> metrics = output.getMetrics();
    
    Assert.assertEquals(metrics.size(), 1);
    
    final InsertSizeMetrics metric = metrics.get(0);

    Assert.assertEquals(metric.PAIR_ORIENTATION.name(), "FR");
    Assert.assertEquals((int) metric.MEDIAN_INSERT_SIZE, 59);
    Assert.assertEquals((int) metric.MODE_INSERT_SIZE, 59);
    Assert.assertEquals(metric.MIN_INSERT_SIZE, 29);
    Assert.assertEquals(metric.MAX_INSERT_SIZE, 89);
    Assert.assertEquals(metric.READ_PAIRS, 101);
    Assert.assertEquals(metric.WIDTH_OF_10_PERCENT, 1);
    Assert.assertEquals(metric.WIDTH_OF_20_PERCENT, 1 + insertBy*2);
    Assert.assertEquals(metric.WIDTH_OF_30_PERCENT, 1 + insertBy*4);
    Assert.assertEquals(metric.WIDTH_OF_40_PERCENT, 1 + insertBy*6);
    Assert.assertEquals(metric.WIDTH_OF_50_PERCENT, 1 + insertBy*8);
    Assert.assertEquals(metric.WIDTH_OF_60_PERCENT, 1 + insertBy*10);
    Assert.assertEquals(metric.WIDTH_OF_70_PERCENT, 1 + insertBy*12);
    Assert.assertEquals(metric.WIDTH_OF_80_PERCENT, 1 + insertBy*14);
    Assert.assertEquals(metric.WIDTH_OF_90_PERCENT, 1 + insertBy*16);
    Assert.assertEquals(metric.WIDTH_OF_95_PERCENT, 1 + insertBy*18);
    Assert.assertEquals(metric.WIDTH_OF_99_PERCENT, 1 + insertBy*20);
}
 

public static File writeBamFile(int numPairs, SAMFileHeader.SortOrder sortOrder)
    throws IOException {
  // file will be both queryname and coordinate sorted, so use one or the other
  SAMRecordSetBuilder samRecordSetBuilder = new SAMRecordSetBuilder(true, sortOrder);
  for (int i = 0; i < numPairs; i++) {
    int chr = 20;
    int start1 = (i + 1) * 1000;
    int start2 = start1 + 100;
    if (i == 5) { // add two unmapped fragments instead of a mapped pair
      samRecordSetBuilder.addFrag(String.format("test-read-%03d-1", i), chr, start1,
          false, true, null,
          null,
          -1, false);
      samRecordSetBuilder.addFrag(String.format("test-read-%03d-2", i), chr, start2,
          false, true, null,
          null,
          -1, false);
    } else {
      samRecordSetBuilder.addPair(String.format("test-read-%03d", i), chr, start1,
          start2);
    }
  }
  if (numPairs > 0) { // add two unplaced unmapped fragments if non-empty
    samRecordSetBuilder.addUnmappedFragment(String.format
        ("test-read-%03d-unplaced-unmapped", numPairs++));
    samRecordSetBuilder.addUnmappedFragment(String.format
        ("test-read-%03d-unplaced-unmapped", numPairs++));
  }

  final File bamFile = File.createTempFile("test", ".bam");
  bamFile.deleteOnExit();
  SAMFileHeader samHeader = samRecordSetBuilder.getHeader();
  final SAMFileWriter bamWriter = new SAMFileWriterFactory()
      .makeSAMOrBAMWriter(samHeader, true, bamFile);
  for (final SAMRecord rec : samRecordSetBuilder.getRecords()) {
    bamWriter.addAlignment(rec);
  }
  bamWriter.close();

  // create BAM index
  if (sortOrder.equals(SAMFileHeader.SortOrder.coordinate)) {
    SamReader samReader = SamReaderFactory.makeDefault()
        .enable(SamReaderFactory.Option.INCLUDE_SOURCE_IN_RECORDS)
        .open(bamFile);
    BAMIndexer.createIndex(samReader, new File(bamFile.getAbsolutePath()
        .replaceFirst("\\.bam$", BAMIndex.BAMIndexSuffix)));
  }

  return bamFile;
}