Java SAMRecord.getMappingQuality方法代码示例

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
@Override
public void run() {

    theslog.log(true, "Pass one: " + bamfile.getAbsolutePath());

    // set up input/output objects using SAM library
    SAMFileReader inputSam = new SAMFileReader(bamfile);
    inputSam.setValidationStringency(SAMFileReader.ValidationStringency.SILENT);

    for (final SAMRecord record : inputSam) {

        if (record.getReadUnmappedFlag()) {
            // if read is unaligned, output coordinates into a bed file
            processOneUnalignedRecord(record);
        } else {
            // if aligned, then make thesaurus entries
            int mapqual = record.getMappingQuality();
            // only consider reads with a minimum mapping quality
            // (this here is not a "real" mapping quality, but a proxy for number of mismatches)
            if (mapqual >= minmapqual) {
                // check that read is aligned                
                processOneRecord(record);
            }
        }
    }

    // and close the input file stream
    inputSam.close();
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
@Override
public void processSequence(SAMRecord read) {
	int quality = read.getMappingQuality();

	//log.debug("quality = " + quality);

	distribution[quality]++;
	readNumber++;
	
	//log.debug("quality count = " + distribution[quality]);

	if (distribution[quality] > maxCount) {
		maxCount = distribution[quality];
	}
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
public void processSAMRecord(SAMRecord r){
	totalReads++;
	if(!r.getReadUnmappedFlag()){
		totalHits++;
		int count = 1;  //Have to figure out something for BWA when reporting multiple alignments
		if(r.getIntegerAttribute("NH")!=null)
			count = r.getIntegerAttribute("NH");
		if(count==1 && r.getMappingQuality()!=0) //Second clause for BWA
			uniquelyMapped++;
		
		weight += 1/(float)count;

		if(r.getReadPairedFlag()){
			if(r.getMateUnmappedFlag()){
				singleEnd++;
			}else{
				pairMapped++;
				if(r.getMateReferenceName().equals(r.getReferenceName())){
					pairedEndSameChr++;
				}else{
					pairedEndDiffChr++;
				}
			}
		}else{
			singleEnd++;
		}
		
		List<AlignmentBlock> blocks = r.getAlignmentBlocks();
   		if(blocks.size()>=2){
   			junctions+=blocks.size()-1;
   		}
		
		if(!r.getNotPrimaryAlignmentFlag()){
			if(!r.getReadPairedFlag() || r.getFirstOfPairFlag()){
				LHits++;
				if(r.getReadPairedFlag() && r.getProperPairFlag()){
					properPair++;
					if(!r.getReadNegativeStrandFlag() && r.getMateNegativeStrandFlag()){
						double dist = (r.getMateAlignmentStart()+r.getReadLength())-r.getAlignmentStart();
						histo.addValue(dist);
					}
				}
			}else if(r.getSecondOfPairFlag()){
				RHits++;
			}
		}else{
			notPrimary++;
		}
	}
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
private boolean acceptMappingQuality(SAMRecord record) {
	return record.getMappingQuality() >= MINIMUM_MAPPING_QUALITY;
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
/**
 * Get sam line for record
 * @param rec The record
 * @return Sam formatted line ending in newline
 */
public static String getSamString(SAMRecord rec) {
	
	// Make string representation of record in sam format
	String rtrn = rec.getReadName() + "\t";
	
	// Get sam flags
	int flags = 0;
	// 0x1 template having multiple segments in sequencing
	if(rec.getReadPairedFlag()) flags += 1;
	// 0x4 segment unmapped
	if(rec.getReadUnmappedFlag()) flags += 4;
	// 0x10 SEQ being reverse complemented
	if(rec.getReadNegativeStrandFlag()) flags += 16;
	// 0x100 secondary alignment
	if(rec.getNotPrimaryAlignmentFlag()) flags += 256;
	// 0x200 not passing quality controls
	if(rec.getReadFailsVendorQualityCheckFlag()) flags += 512;
	// 0x400 PCR or optical duplicate
	if(rec.getDuplicateReadFlag()) flags += 1024;
	
	if(rec.getReadPairedFlag()) {
		// 0x2 each segment properly aligned according to the aligner
		if(rec.getProperPairFlag()) flags += 2;
		// 0x8 next segment in the template unmapped
		if(rec.getMateUnmappedFlag()) flags += 8;
		// 0x20 SEQ of the next segment in the template being reversed
		if(rec.getMateNegativeStrandFlag()) flags += 32;
		// 0x40 the first segment in the template
		if(rec.getFirstOfPairFlag()) flags += 64;
		// 0x80 the last segment in the template
		if(rec.getSecondOfPairFlag()) flags += 128;
	}
	
	rtrn += Integer.valueOf(flags).toString() + "\t";
	
	// The rest of the sam fields
	rtrn += rec.getReferenceName() + "\t";
	rtrn += rec.getAlignmentStart() + "\t";
	rtrn += rec.getMappingQuality() + "\t";
	rtrn += rec.getCigarString() + "\t";
	rtrn += rec.getMateReferenceName() + "\t";
	rtrn += rec.getMateAlignmentStart() + "\t";
	rtrn += rec.getInferredInsertSize() + "\t";
	rtrn += rec.getReadString() + "\t";
	rtrn += rec.getBaseQualityString() + "\n";
	
	return rtrn;
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
public WrapSamRecord(SAMRecord record){
	
	this.mappingQuality=record.getMappingQuality();
	this.isDuplicate=record.getDuplicateReadFlag();
	this.readSequence=record.getReadString();
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
private boolean passFilter(SAMRecord alignment) {
	return !alignment.getReadUnmappedFlag() && !alignment.getDuplicateReadFlag() && alignment.getMappingQuality()>mappingQuality;
}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
/**
 * The map function runs once per single-base locus, and accepts a 'context', a
 * data structure consisting of the reads which overlap the locus, the sites over
 * which they fall, and the base from the reference that overlaps.
 *
 * @param tracker The accessor for reference metadata.
 * @param ref     The reference base that lines up with this locus.
 * @param context Information about reads aligning to this locus.
 * @return In this case, returns a count of how many loci were seen at thisb site (1).
 */
@Override
public Integer map(RefMetaDataTracker tracker, ReferenceContext ref, AlignmentContext context) {

    char refBase = (char) ref.getBase();
    double ref_count = 0;
    double nonref_count = 0;
    double ar;


    for (PileupElement p : context.getBasePileup()) {
        SAMRecord read = p.getRead();

        //filter
        if (read.getMappingQuality() < MIN_MAPPING_QUALTY_SCORE || p.getQual() < MIN_BASE_QUALTY_SCORE)
            continue;

        char base = (char) p.getBase();

        if (base != refBase)
            nonref_count++;
        else
            ref_count++;

    }

    out.printf("%s\t%d\t%.0f\t%.0f\t%.3f\n", context.getLocation().getContig(), context.getLocation().getStart(), nonref_count, ref_count, nonref_count / (nonref_count + ref_count));

    return 0;


}
 

import net.sf.samtools.SAMRecord; //导入方法依赖的package包/类
private void setOptionsFromFile (File file) {

		// This just reads the first few thousand lines from the first file and
		// tries to set the preferences options to the correct defaults.
		
		SAMFileReader.setDefaultValidationStringency(SAMFileReader.ValidationStringency.SILENT);

		SAMFileReader inputSam = new SAMFileReader(file); 

		int lineCount = 0;
		// Now process the file

		boolean pairedEnd = false;
		boolean spliced = false;
		
		int maxQ = 0;
		
		
		for (SAMRecord samRecord : inputSam) {
			++lineCount;
			
			if (lineCount == 100000) break;
			
			if (samRecord.getMappingQuality() > maxQ) {
				maxQ = samRecord.getMappingQuality();
			}
			
			if (samRecord.getCigarString().contains("N")) {
				spliced = true;
			}
			if (samRecord.getReadPairedFlag()) {
				pairedEnd = true;
			}
			
		}
		
		inputSam.close();
		
		prefs.setPairedEnd(pairedEnd);
		prefs.setSpliced(spliced);
		prefs.setMinMappingQuality(Math.min(maxQ,20));
		

		
	}