Search in sources :

Example 6 with SamReader

use of htsjdk.samtools.SamReader in project gatk by broadinstitute.

the class UnmarkDuplicatesIntegrationTest method getDuplicateCountForBam.

private static int getDuplicateCountForBam(final File bam, final File referenceFile) throws IOException {
    int duplicateCount = 0;
    final SamReaderFactory factory = SamReaderFactory.makeDefault();
    if (referenceFile != null) {
        factory.referenceSequence(referenceFile);
    }
    try (final SamReader reader = factory.open(bam)) {
        for (SAMRecord read : reader) {
            if (read.getDuplicateReadFlag()) {
                ++duplicateCount;
            }
        }
    }
    return duplicateCount;
}
Also used : SamReader(htsjdk.samtools.SamReader) SamReaderFactory(htsjdk.samtools.SamReaderFactory) SAMRecord(htsjdk.samtools.SAMRecord)

Example 7 with SamReader

use of htsjdk.samtools.SamReader in project gatk-protected by broadinstitute.

the class BayesianHetPulldownCalculatorUnitTest method initHeaders.

@BeforeClass
public void initHeaders() throws IOException {
    try (final SamReader normalBamReader = SamReaderFactory.makeDefault().open(NORMAL_BAM_FILE);
        final SamReader tumorBamReader = SamReaderFactory.makeDefault().open(TUMOR_BAM_FILE)) {
        normalHeader = normalBamReader.getFileHeader();
        tumorHeader = tumorBamReader.getFileHeader();
    }
}
Also used : SamReader(htsjdk.samtools.SamReader) BeforeClass(org.testng.annotations.BeforeClass)

Example 8 with SamReader

use of htsjdk.samtools.SamReader in project gatk by broadinstitute.

the class CompareBaseQualities method doWork.

@Override
protected Object doWork() {
    if (roundDown && (staticQuantizationQuals == null || staticQuantizationQuals.isEmpty())) {
        throw new CommandLineException.BadArgumentValue("round_down_quantized", "true", "This option can only be used if static_quantized_quals is also used");
    }
    staticQuantizedMapping = constructStaticQuantizedMapping(staticQuantizationQuals, roundDown);
    IOUtil.assertFileIsReadable(samFiles.get(0));
    IOUtil.assertFileIsReadable(samFiles.get(1));
    final SamReaderFactory factory = SamReaderFactory.makeDefault().validationStringency(VALIDATION_STRINGENCY);
    final SamReader reader1 = factory.referenceSequence(REFERENCE_SEQUENCE).open(samFiles.get(0));
    final SamReader reader2 = factory.referenceSequence(REFERENCE_SEQUENCE).open(samFiles.get(1));
    final SecondaryOrSupplementarySkippingIterator it1 = new SecondaryOrSupplementarySkippingIterator(reader1.iterator());
    final SecondaryOrSupplementarySkippingIterator it2 = new SecondaryOrSupplementarySkippingIterator(reader2.iterator());
    final CompareMatrix finalMatrix = new CompareMatrix(staticQuantizedMapping);
    final ProgressMeter progressMeter = new ProgressMeter(1.0);
    progressMeter.start();
    while (it1.hasCurrent() && it2.hasCurrent()) {
        final SAMRecord read1 = it1.getCurrent();
        final SAMRecord read2 = it2.getCurrent();
        progressMeter.update(read1);
        if (!read1.getReadName().equals(read2.getReadName())) {
            throw new UserException.BadInput("files do not have the same exact order of reads:" + read1.getReadName() + " vs " + read2.getReadName());
        }
        finalMatrix.add(read1.getBaseQualities(), read2.getBaseQualities());
        it1.advance();
        it2.advance();
    }
    progressMeter.stop();
    if (it1.hasCurrent() || it2.hasCurrent()) {
        throw new UserException.BadInput("files do not have the same exact number of reads");
    }
    CloserUtil.close(reader1);
    CloserUtil.close(reader2);
    finalMatrix.printOutResults(outputFilename);
    if (throwOnDiff && finalMatrix.hasNonDiagonalElements()) {
        throw new UserException("Quality scores from the two BAMs do not match");
    }
    return finalMatrix.hasNonDiagonalElements() ? 1 : 0;
}
Also used : SamReader(htsjdk.samtools.SamReader) SamReaderFactory(htsjdk.samtools.SamReaderFactory) ProgressMeter(org.broadinstitute.hellbender.engine.ProgressMeter) SAMRecord(htsjdk.samtools.SAMRecord) SecondaryOrSupplementarySkippingIterator(htsjdk.samtools.SecondaryOrSupplementarySkippingIterator) UserException(org.broadinstitute.hellbender.exceptions.UserException)

Example 9 with SamReader

use of htsjdk.samtools.SamReader in project gatk by broadinstitute.

the class PrintReadsIntegrationTest method testReadFilters.

@Test(dataProvider = "readFilterTestData")
public void testReadFilters(final String input, final String reference, final String extOut, final List<String> inputArgs, final int expectedCount) throws IOException {
    final File outFile = createTempFile("testReadFilter", extOut);
    final ArgumentsBuilder args = new ArgumentsBuilder();
    args.add("-I");
    args.add(new File(TEST_DATA_DIR, input).getAbsolutePath());
    args.add("-O");
    args.add(outFile.getAbsolutePath());
    if (reference != null) {
        args.add("-R");
        args.add(new File(TEST_DATA_DIR, reference).getAbsolutePath());
    }
    for (final String filter : inputArgs) {
        args.add(filter);
    }
    runCommandLine(args);
    SamReaderFactory factory = SamReaderFactory.makeDefault();
    if (reference != null) {
        factory = factory.referenceSequence(new File(TEST_DATA_DIR, reference));
    }
    int count = 0;
    try (final SamReader reader = factory.open(outFile)) {
        Iterator<SAMRecord> it = reader.iterator();
        while (it.hasNext()) {
            SAMRecord rec = it.next();
            count++;
        }
    }
    Assert.assertEquals(count, expectedCount);
}
Also used : SamReader(htsjdk.samtools.SamReader) SamReaderFactory(htsjdk.samtools.SamReaderFactory) SAMRecord(htsjdk.samtools.SAMRecord) ArgumentsBuilder(org.broadinstitute.hellbender.utils.test.ArgumentsBuilder) File(java.io.File) BaseTest(org.broadinstitute.hellbender.utils.test.BaseTest) Test(org.testng.annotations.Test) CommandLineProgramTest(org.broadinstitute.hellbender.CommandLineProgramTest)

Example 10 with SamReader

use of htsjdk.samtools.SamReader in project gatk by broadinstitute.

the class SinglePassSamProgram method makeItSo.

public static void makeItSo(final File input, final File referenceSequence, final boolean assumeSorted, final long stopAfter, final Collection<SinglePassSamProgram> programs) {
    // Setup the standard inputs
    IOUtil.assertFileIsReadable(input);
    final SamReader in = SamReaderFactory.makeDefault().referenceSequence(referenceSequence).open(input);
    // Optionally load up the reference sequence and double check sequence dictionaries
    final ReferenceSequenceFileWalker walker;
    if (referenceSequence == null) {
        walker = null;
    } else {
        IOUtil.assertFileIsReadable(referenceSequence);
        walker = new ReferenceSequenceFileWalker(referenceSequence);
        if (!in.getFileHeader().getSequenceDictionary().isEmpty()) {
            SequenceUtil.assertSequenceDictionariesEqual(in.getFileHeader().getSequenceDictionary(), walker.getSequenceDictionary());
        }
    }
    // Check on the sort order of the BAM file
    {
        final SortOrder sort = in.getFileHeader().getSortOrder();
        if (sort != SortOrder.coordinate) {
            if (assumeSorted) {
                logger.warn("File reports sort order '" + sort + "', assuming it's coordinate sorted anyway.");
            } else {
                throw new UserException("File " + input.getAbsolutePath() + " should be coordinate sorted but " + "the header says the sort order is " + sort + ". If you believe the file " + "to be coordinate sorted you may pass ASSUME_SORTED=true");
            }
        }
    }
    // Call the abstract setup method!
    boolean anyUseNoRefReads = false;
    for (final SinglePassSamProgram program : programs) {
        program.setup(in.getFileHeader(), input);
        anyUseNoRefReads = anyUseNoRefReads || program.usesNoRefReads();
    }
    final ProgressLogger progress = new ProgressLogger(logger);
    for (final SAMRecord rec : in) {
        final ReferenceSequence ref;
        if (walker == null || rec.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
            ref = null;
        } else {
            ref = walker.get(rec.getReferenceIndex());
        }
        for (final SinglePassSamProgram program : programs) {
            program.acceptRead(rec, ref);
        }
        progress.record(rec);
        // See if we need to terminate early?
        if (stopAfter > 0 && progress.getCount() >= stopAfter) {
            break;
        }
        // And see if we're into the unmapped reads at the end
        if (!anyUseNoRefReads && rec.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
            break;
        }
    }
    CloserUtil.close(in);
    for (final SinglePassSamProgram program : programs) {
        program.finish();
    }
}
Also used : SamReader(htsjdk.samtools.SamReader) SAMRecord(htsjdk.samtools.SAMRecord) SortOrder(htsjdk.samtools.SAMFileHeader.SortOrder) ProgressLogger(org.broadinstitute.hellbender.utils.runtime.ProgressLogger) UserException(org.broadinstitute.hellbender.exceptions.UserException) ReferenceSequence(htsjdk.samtools.reference.ReferenceSequence) ReferenceSequenceFileWalker(htsjdk.samtools.reference.ReferenceSequenceFileWalker)

Aggregations

SamReader (htsjdk.samtools.SamReader)211 SAMRecord (htsjdk.samtools.SAMRecord)137 File (java.io.File)111 SAMRecordIterator (htsjdk.samtools.SAMRecordIterator)89 SAMFileHeader (htsjdk.samtools.SAMFileHeader)83 IOException (java.io.IOException)71 SamReaderFactory (htsjdk.samtools.SamReaderFactory)65 ArrayList (java.util.ArrayList)63 SAMFileWriter (htsjdk.samtools.SAMFileWriter)58 SAMSequenceDictionary (htsjdk.samtools.SAMSequenceDictionary)44 SAMSequenceDictionaryProgress (com.github.lindenb.jvarkit.util.picard.SAMSequenceDictionaryProgress)42 List (java.util.List)39 CigarElement (htsjdk.samtools.CigarElement)32 IndexedFastaSequenceFile (htsjdk.samtools.reference.IndexedFastaSequenceFile)32 HashMap (java.util.HashMap)31 Cigar (htsjdk.samtools.Cigar)30 SAMReadGroupRecord (htsjdk.samtools.SAMReadGroupRecord)30 PrintWriter (java.io.PrintWriter)27 Interval (htsjdk.samtools.util.Interval)26 HashSet (java.util.HashSet)26