Examples with SAMRecord htsjdk.samtools.SAMRecord used on opensource projects

Example 1 with SAMRecord

use of htsjdk.samtools.SAMRecord 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;
}

Example 2 with SAMRecord

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

the class AlignedAssemblyOrExcuse method toSAMStreamForOneContig.

public static Stream<SAMRecord> toSAMStreamForOneContig(final SAMFileHeader header, final List<String> refNames, final int assemblyId, final int contigIdx, final byte[] contigSequence, final List<BwaMemAlignment> alignments) {
    if (alignments.isEmpty())
        return Stream.empty();
    final String readName = formatContigName(assemblyId, contigIdx);
    final Map<BwaMemAlignment, String> saTagMap = BwaMemAlignmentUtils.createSATags(alignments, refNames);
    return alignments.stream().map(alignment -> {
        final SAMRecord samRecord = BwaMemAlignmentUtils.applyAlignment(readName, contigSequence, null, null, alignment, refNames, header, false, false);
        final String saTag = saTagMap.get(alignment);
        if (saTag != null)
            samRecord.setAttribute("SA", saTag);
        return samRecord;
    });
}

Example 3 with SAMRecord

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

the class BestEndMapqPrimaryAlignmentStrategy method randomlySelectPrimaryFromBest.

/**
     * Randomly picks one of the best alignments and puts it into the 0th slot of the list.
     * @param recs List of alignments sorted in descending order of alignment quality.
     */
private void randomlySelectPrimaryFromBest(List<SAMRecord> recs) {
    if (recs.isEmpty())
        return;
    final int bestMapq = recs.get(0).getMappingQuality();
    int i;
    for (i = 1; i < recs.size() && recs.get(i).getMappingQuality() == bestMapq; ++i) {
    }
    final int bestIndex = random.nextInt(i);
    if (bestIndex == 0)
        return;
    final SAMRecord tmp = recs.get(0);
    recs.set(0, recs.get(bestIndex));
    recs.set(bestIndex, tmp);
}

Example 4 with SAMRecord

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

the class HitsForInsert method coordinateByHitIndex.

/**
     * Some alignment strategies expect to receive alignments for ends that are coordinated by
     * hit index (HI) tag.  This method lines up alignments for each end by HI tag value, and if there is
     * no corresponding alignment for an alignment, there is a null in the array at that slot.
     *
     * This method then renumbers the HI values so that they start at zero and have no gaps, because some
     * reads may have been filtered out.
     */
public void coordinateByHitIndex() {
    // Sort by HI value, with reads with no HI going at the end.
    Collections.sort(firstOfPairOrFragment, comparator);
    Collections.sort(secondOfPair, comparator);
    // and uncorrelated alignments have null in the other list at the corresponding index.
    for (int i = 0; i < Math.min(firstOfPairOrFragment.size(), secondOfPair.size()); ++i) {
        final Integer leftHi = firstOfPairOrFragment.get(i).getIntegerAttribute(SAMTag.HI.name());
        final Integer rightHi = secondOfPair.get(i).getIntegerAttribute(SAMTag.HI.name());
        if (leftHi != null) {
            if (rightHi != null) {
                if (leftHi < rightHi)
                    secondOfPair.add(i, null);
                else if (rightHi < leftHi)
                    firstOfPairOrFragment.add(i, null);
            // else the are correlated
            }
        } else if (rightHi != null) {
            firstOfPairOrFragment.add(i, null);
        } else {
            // Both alignments do not have HI, so push down the ones on the right.
            // Right is arbitrarily picked to push down.
            secondOfPair.add(i, null);
        }
    }
    // Now renumber any correlated alignments, and remove hit index if no correlated read.
    int hi = 0;
    for (int i = 0; i < numHits(); ++i) {
        final SAMRecord first = getFirstOfPair(i);
        final SAMRecord second = getSecondOfPair(i);
        if (first != null && second != null) {
            first.setAttribute(SAMTag.HI.name(), i);
            second.setAttribute(SAMTag.HI.name(), i);
            ++hi;
        } else if (first != null) {
            first.setAttribute(SAMTag.HI.name(), null);
        } else {
            second.setAttribute(SAMTag.HI.name(), null);
        }
    }
}

Example 5 with SAMRecord

use of htsjdk.samtools.SAMRecord 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);
}