Java SAMTag类代码示例

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
public int callElPrep(String input, String output, String rg, int threads, 
        SAMRecordIterator SAMit,
        SAMFileHeader header, String dictFile, boolean updateRG, boolean keepDups, String RGID) throws InterruptedException, QualityException {
    
    SAMRecord sam;
    SAMFileWriterFactory factory = new SAMFileWriterFactory();
    SAMFileWriter Swriter = factory.makeSAMWriter(header, true, new File(input));
    
    int reads = 0;
    while(SAMit.hasNext()) {
        sam = SAMit.next();
        if(updateRG)
            sam.setAttribute(SAMTag.RG.name(), RGID);
        Swriter.addAlignment(sam);
        reads++;
    }
    Swriter.close();
    
    String customArgs = HalvadeConf.getCustomArgs(context.getConfiguration(), "elprep", "");  
    String[] command = CommandGenerator.elPrep(bin, input, output, threads, true, rg, null, !keepDups, customArgs);
    long estimatedTime = runProcessAndWait("elPrep", command);
    if(context != null)
        context.getCounter(HalvadeCounters.TIME_ELPREP).increment(estimatedTime);
    
    return reads;
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
/**
 * Creates a GATKRead from a FastqRecord and a header.
 *
 * @param header the header for the record.
 * @param record the record to use as GATKRead.
 */
public FastqGATKRead(final SAMFileHeader header, final FastqRecord record) {
    super(new SAMRecord(header));
    Utils.nonNull(record, "null record");
    // update the record with the read name information
    FastqReadNameEncoding.updateReadFromReadName(this, record.getReadName());
    // set the bases and the qualities
    this.setBases(record.getReadBases());
    this.setBaseQualities(record.getBaseQualities());
    // add the comments in the quality header to the comment if present
    final String baseQualHeader = record.getBaseQualityHeader();
    if (baseQualHeader != null) {
        // the default tag in the specs is CO
        this.setAttribute(SAMTag.CO.toString(), baseQualHeader);
    }
    this.setIsUnmapped();
    if (this.isPaired()) {
        this.setMateIsUnmapped();
    }
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@Override
public boolean test( final GATKRead read ) {
    final SAMReadGroupRecord readGroup = ReadUtils.getSAMReadGroupRecord(read, samHeader);
    if ( readGroup == null ) {
        return true;
    }

    for (final String attributeType : blacklistEntries.keySet()) {

        final String attribute;
        if (SAMReadGroupRecord.READ_GROUP_ID_TAG.equals(attributeType) || SAMTag.RG.name().equals(attributeType)) {
            attribute = readGroup.getId();
        } else {
            attribute = readGroup.getAttribute(attributeType);
        }
        if (attribute != null && blacklistEntries.get(attributeType).contains(attribute)) {
            return false;
        }
    }

    return true;
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
/**
 * Write out a representation of this haplotype as a read
 *
 * @param haplotype a haplotype to write out, must not be null
 * @param paddedRefLoc the reference location, must not be null
 * @param isAmongBestHaplotypes true if among the best haplotypes, false if it was just one possible haplotype
 */
private void writeHaplotype(final Haplotype haplotype,
                            final Locatable paddedRefLoc,
                            final boolean isAmongBestHaplotypes) {
    Utils.nonNull(haplotype, "haplotype cannot be null");
    Utils.nonNull(paddedRefLoc, "paddedRefLoc cannot be null");

    final SAMRecord record = new SAMRecord(output.getBAMOutputHeader());
    record.setReadBases(haplotype.getBases());
    record.setAlignmentStart(paddedRefLoc.getStart() + haplotype.getAlignmentStartHapwrtRef());
    // Use a base quality value "!" for it's display value (quality value is not meaningful)
    record.setBaseQualities(Utils.dupBytes((byte) '!', haplotype.getBases().length));
    record.setCigar(AlignmentUtils.consolidateCigar(haplotype.getCigar()));
    record.setMappingQuality(isAmongBestHaplotypes ? bestHaplotypeMQ : otherMQ);
    record.setReadName(output.getHaplotypeSampleTag() + uniqueNameCounter++);
    record.setAttribute(output.getHaplotypeSampleTag(), haplotype.hashCode());
    record.setReadUnmappedFlag(false);
    record.setReferenceIndex(output.getBAMOutputHeader().getSequenceIndex(paddedRefLoc.getContig()));
    record.setAttribute(SAMTag.RG.toString(), output.getHaplotypeReadGroupID());
    record.setFlags(SAMFlag.READ_REVERSE_STRAND.intValue());

    output.add(new SAMRecordToGATKReadAdapter(record));
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
/**
 * Returns input read with alignment-related info cleared
 */
private static GATKRead clearReadAlignment(final GATKRead read, final SAMFileHeader header) {
    final GATKRead newRead = new SAMRecordToGATKReadAdapter(new SAMRecord(header));
    newRead.setName(read.getName());
    newRead.setBases(read.getBases());
    newRead.setBaseQualities(read.getBaseQualities());
    if (read.isReverseStrand()) {
        SequenceUtil.reverseComplement(newRead.getBases());
        SequenceUtil.reverseQualities(newRead.getBaseQualities());
    }
    newRead.setIsUnmapped();
    newRead.setIsPaired(read.isPaired());
    if (read.isPaired()) {
        newRead.setMateIsUnmapped();
        if (read.isFirstOfPair()) {
            newRead.setIsFirstOfPair();
        } else if (read.isSecondOfPair()) {
            newRead.setIsSecondOfPair();
        }
    }
    final String readGroup = read.getReadGroup();
    if (readGroup != null) {
        newRead.setAttribute(SAMTag.RG.name(), readGroup);
    }
    return newRead;
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@Test
public void testUsingOriginalQualities() throws Exception {
    final MeanQualityByCycleSpark.HistogramGenerator hg = new MeanQualityByCycleSpark.HistogramGenerator(true);
    Assert.assertEquals(hg.useOriginalQualities, true);

    GATKRead read1 = ArtificialReadUtils.createArtificialRead("aa".getBytes(), new byte[]{50, 50}, "2M");
    hg.addRead(read1);
    assertEqualsLongArray(hg.firstReadCountsByCycle, new long[0]);
    assertEqualsDoubleArray(hg.firstReadTotalsByCycle, new double[0], 1e-05);
    assertEqualsLongArray(hg.secondReadCountsByCycle, new long[0]);
    assertEqualsDoubleArray(hg.secondReadTotalsByCycle, new double[0], 1e-05);

    GATKRead read2 = ArtificialReadUtils.createArtificialRead("aa".getBytes(), new byte[]{50, 50}, "2M");
    read2.setAttribute(SAMTag.OQ.name(), SAMUtils.phredToFastq(new byte[]{30, 40}));
    hg.addRead(read2);

    assertEqualsLongArray(hg.firstReadCountsByCycle, new long[]{0, 1, 1});
    assertEqualsDoubleArray(hg.firstReadTotalsByCycle, new double[]{0, 30, 40}, 1e-05);
    assertEqualsLongArray(hg.secondReadCountsByCycle, new long[]{0, 0, 0});
    assertEqualsDoubleArray(hg.secondReadTotalsByCycle, new double[]{0, 0, 0}, 1e-05);
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@Override
protected void reduce(ChromosomeRegion key, Iterable<SAMRecordWritable> values, Context context) throws IOException, InterruptedException {
    Iterator<SAMRecordWritable> it = values.iterator();
    SAMRecord sam = null;
    while(it.hasNext()) {
        sam = it.next().get();
        if(!inputIsBam  || updateRG) {
            sam.setAttribute(SAMTag.RG.name(), RGID);
        }
        samWritable.set(sam);
        recordWriter.write(outKey, samWritable);
                
    }
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
public int streamElPrep(Reducer.Context context, String output, String rg, 
            int threads, SAMRecordIterator SAMit, 
            SAMFileHeader header, String dictFile, boolean updateRG, boolean keepDups, String RGID) throws InterruptedException, IOException, QualityException {
        long startTime = System.currentTimeMillis();
        String customArgs = HalvadeConf.getCustomArgs(context.getConfiguration(), "elprep", "");  
        String[] command = CommandGenerator.elPrep(bin, "/dev/stdin", output, threads, true, rg, null, !keepDups, customArgs);
//        runProcessAndWait(command);
        ProcessBuilderWrapper builder = new ProcessBuilderWrapper(command, null);
        builder.startProcess(true);        
        BufferedWriter localWriter = builder.getSTDINWriter();
        
        // write header
        final StringWriter headerTextBuffer = new StringWriter();
        new SAMTextHeaderCodec().encode(headerTextBuffer, header);
        final String headerText = headerTextBuffer.toString();
        localWriter.write(headerText, 0, headerText.length());
        
        
        SAMRecord sam;
        int reads = 0;
        while(SAMit.hasNext()) {
            sam = SAMit.next();
            if(updateRG)
                sam.setAttribute(SAMTag.RG.name(), RGID);
            String samString = sam.getSAMString();
            localWriter.write(samString, 0, samString.length());
            reads++;
        }
        localWriter.flush();
        localWriter.close();
                
        int error = builder.waitForCompletion();
        if(error != 0)
            throw new ProcessException("elPrep", error);
        long estimatedTime = System.currentTimeMillis() - startTime;
        Logger.DEBUG("estimated time: " + estimatedTime / 1000);
        if(context != null)
            context.getCounter(HalvadeCounters.TIME_ELPREP).increment(estimatedTime);
        return reads;
    }
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@Override
public void addRead(final GATKRead read) {
    // adding the raw barcode information if found
    String readName = RTReadUtils.getReadNameWithIlluminaBarcode(read);
    // adding the pair information
    if (read.isPaired()) {
        readName += (read.isFirstOfPair())
                ? FastqConstants.FIRST_OF_PAIR : FastqConstants.SECOND_OF_PAIR;
    }
    writer.write(new FastqRecord(readName,
            read.getBasesString(),
            read.getAttributeAsString(SAMTag.CO.name()),
            ReadUtils.getBaseQualityString(read)));
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@Override
public boolean test(final GATKRead read) {
    metric.TOTAL++;
    // trimming function modify in place the read
    final boolean pass = delegate.test(read);
    // update the metrics
    if (pass) {
        metric.PASSED++;
    } else {
        // if it does not pass, add a FT tag with the name of the filter
        read.setAttribute(SAMTag.FT.name(), metric.FILTER);
    }
    return pass;
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
@DataProvider(name = "fastqRecordData")
public Iterator<Object[]> fastqRecordDataProvider() {
    final String baseQualities = "FFGCHI5";
    final String bases = "ACTGTTAG";
    final GATKRead baseRecord = ArtificialReadUtils
            .createArtificialUnmappedRead(null,
                    new byte[] {'A', 'C', 'T', 'G', 'T', 'T', 'A', 'G'},
                    new byte[] {37, 37, 38, 34, 39, 40, 20});
    baseRecord.setName("baseRecord");
    final List<Object[]> data = new ArrayList<>();
    // simple case test
    data.add(new Object[] {new FastqRecord(baseRecord.getName(), bases, null, baseQualities),
            baseRecord.deepCopy()});
    // case with comment information
    baseRecord.setAttribute(SAMTag.CO.name(), "quality comment");
    data.add(new Object[] {
            new FastqRecord(baseRecord.getName(), bases, "quality comment", baseQualities),
            baseRecord.deepCopy()});
    // case with a read name with pair-end information
    baseRecord.setIsSecondOfPair();
    baseRecord.setIsUnmapped();
    data.add(new Object[] {
            new FastqRecord(baseRecord.getName() + "/2", bases, "quality comment",
                    baseQualities),
            baseRecord.deepCopy()});
    // case with read name as CASAVA format
    baseRecord.setName("baseRecord");
    baseRecord.setAttribute("BC", "ATCG");
    data.add(new Object[] {
            new FastqRecord("baseRecord 2:N:3:ATCG", bases, "quality comment", baseQualities),
            baseRecord.deepCopy()});
    // case with PF flag
    baseRecord.setFailsVendorQualityCheck(true);
    data.add(new Object[] {
            new FastqRecord("baseRecord 2:Y:3:ATCG", bases, "quality comment", baseQualities),
            baseRecord.deepCopy()});
    return data.iterator();
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
public int compare(final SAMRecord rec1, final SAMRecord rec2) {
    final Integer hi1 = rec1.getIntegerAttribute(SAMTag.HI.name());
    final Integer hi2 = rec2.getIntegerAttribute(SAMTag.HI.name());
    if (hi1 == null) {
        if (hi2 == null) return 0;
        else return 1;
    } else if (hi2 == null) {
        return -1;
    } else {
        return hi1.compareTo(hi2);
    }
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
/**
 * Helper function that writes reads from iterator it into writer out, updating each SAMRecord along the way
 * according to the newOrder mapping from dictionary index -> index.  Name is used for printing only.
 */
private void writeReads(final SAMFileWriter out,
                        final SAMRecordIterator it,
                        final Map<Integer, Integer> newOrder,
                        final String name) {
    long counter = 0;
    log.info("  Processing " + name);

    while (it.hasNext()) {
        counter++;
        final SAMRecord read = it.next();
        final int oldRefIndex = read.getReferenceIndex();
        final int oldMateIndex = read.getMateReferenceIndex();
        final int newRefIndex = newOrderIndex(read, oldRefIndex, newOrder);

        read.setHeader(out.getFileHeader());
        read.setReferenceIndex(newRefIndex);

        final int newMateIndex = newOrderIndex(read, oldMateIndex, newOrder);
        if (oldMateIndex != -1 && newMateIndex == -1) { // becoming unmapped
            read.setMateAlignmentStart(0);
            read.setMateUnmappedFlag(true);
            read.setAttribute(SAMTag.MC.name(), null);      // Set the Mate Cigar String to null
        }
        read.setMateReferenceIndex(newMateIndex);

        out.addAlignment(read);
    }

    it.close();
    log.info("Wrote " + counter + " reads");
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
/**
 * Calculates the mode for outward-facing pairs, using the first SAMPLE_FOR_MODE
 * outward-facing pairs found in INPUT
 */
private double getOutieMode() {

    int samplePerFile = SAMPLE_FOR_MODE / INPUT.size();

    Histogram<Integer> histo = new Histogram<Integer>();

    for (File f : INPUT) {
        SamReader reader = SamReaderFactory.makeDefault().open(f);
        int sampled = 0;
        for (Iterator<SAMRecord> it = reader.iterator(); it.hasNext() && sampled < samplePerFile; ) {
            SAMRecord sam = it.next();
            if (!sam.getFirstOfPairFlag()) {
                continue;
            }
            // If we get here we've hit the end of the aligned reads
            if (sam.getReadUnmappedFlag() && sam.getReferenceIndex() == SAMRecord.NO_ALIGNMENT_REFERENCE_INDEX) {
                break;
            } else if (sam.getReadUnmappedFlag() || sam.getMateUnmappedFlag()) {
                continue;
            } else if ((sam.getAttribute(SAMTag.MQ.name()) == null ||
                    sam.getIntegerAttribute(SAMTag.MQ.name()) >= MINIMUM_MAPPING_QUALITY) &&
                    sam.getMappingQuality() >= MINIMUM_MAPPING_QUALITY &&
                    sam.getMateNegativeStrandFlag() != sam.getReadNegativeStrandFlag() &&
                    sam.getMateReferenceIndex().equals(sam.getReferenceIndex()) &&
                    SamPairUtil.getPairOrientation(sam) == PairOrientation.RF) {
                histo.increment(Math.abs(sam.getInferredInsertSize()));
                sampled++;
            }
        }
        CloserUtil.close(reader);
    }

    return histo.size() > 0 ? histo.getMode() : 0;
}
 

import htsjdk.samtools.SAMTag; //导入依赖的package包/类
private String getPlatformUnit( final GATKRead read ) {
    final String pu_attr = read.getAttributeAsString(SAMTag.PU.name());
    if ( pu_attr != null ) {
        return pu_attr;
    }

    return ReadUtils.getPlatformUnit(read, samHeader);
}