Java SAMRecord.getCigar方法代码示例

import htsjdk.samtools.SAMRecord; //导入方法依赖的package包/类
private int getNumClippedBasesAtStart(final SAMRecord read) {
    // compute total number of clipped bases (soft or hard clipped)
    // check for hard clips (never consider these bases):
    final Cigar c = read.getCigar();
    final CigarElement first = c.getCigarElement(0);

    int numStartClippedBases = 0;
    if (first.getOperator() == CigarOperator.H) {
        numStartClippedBases = first.getLength();
    }
    final byte[] unclippedReadBases = read.getReadBases();
    final byte[] unclippedReadQuals = read.getBaseQualities();

    // Do a stricter base clipping than provided by CIGAR string, since this one may be too conservative,
    // and may leave a string of Q2 bases still hanging off the reads.
    for (int i = numStartClippedBases; i < unclippedReadBases.length; i++) {
        if (unclippedReadQuals[i] < PairHMMIndelErrorModel.BASE_QUAL_THRESHOLD)
            numStartClippedBases++;
        else
            break;

    }

    return numStartClippedBases;
}
 

import htsjdk.samtools.SAMRecord; //导入方法依赖的package包/类
public boolean startWIns(SAMRecord sr){
Cigar cigar = sr.getCigar();
if(cigar == null){
    return true;
}else{
    CigarOperator op = cigar.getCigarElements().get(0).getOperator();
    if(op == CigarOperator.I){
	if(HLA.DEBUG)
	    HLA.log.appendln("SKIPPING(Start with Insertion):\t" + sr.getReadName());
	return true;
	
    }
}
return false;
   }
 

import htsjdk.samtools.SAMRecord; //导入方法依赖的package包/类
/**
 * Returns number of alignment blocks (continuous stretches of aligned bases) in the specified alignment.
 * This method follows closely the SAMRecord::getAlignmentBlocks() implemented in samtools library, but
 * it only counts blocks without actually allocating and filling the list of blocks themselves. Hence, this method is
 * a much more efficient alternative to r.getAlignmentBlocks.size() in the situations when this number is all that is needed.
 * Formally, this method simply returns the number of M elements in the cigar.
 *
 * @param r alignment
 * @return number of continuous alignment blocks (i.e. 'M' elements of the cigar; all indel and clipping elements are ignored).
 */
public static int getNumAlignmentBlocks(final SAMRecord r) {
    if ( r == null ) throw new IllegalArgumentException("read cannot be null");
    final Cigar cigar = r.getCigar();
    if (cigar == null) return 0;

    int n = 0;
    for (final CigarElement e : cigar.getCigarElements()) {
        if (ALIGNED_TO_GENOME_OPERATORS.contains(e.getOperator()))
            n++;
    }

    return n;
}
 

import htsjdk.samtools.SAMRecord; //导入方法依赖的package包/类
/**
 * Count the number of bases hard clipped from read
 *
 * If read's cigar is null, return 0
 *
 * @param r a non-null read
 * @return a positive integer
 */
public static int getNumHardClippedBases(final SAMRecord r) {
    if ( r == null ) throw new IllegalArgumentException("Read cannot be null");

    int n = 0;
    final Cigar cigar = r.getCigar();
    if (cigar == null) return 0;

    for (final CigarElement e : cigar.getCigarElements())
        if (e.getOperator() == CigarOperator.H)
            n += e.getLength();

    return n;
}
 

import htsjdk.samtools.SAMRecord; //导入方法依赖的package包/类
public boolean qcCheck(SAMRecord sr){
Cigar cigar = sr.getCigar();
int rLen = sr.getReadLength();
int effectiveLen = 0;
if(cigar==null) 
    return false;
else{
    for(final CigarElement ce : cigar.getCigarElements()){
	CigarOperator op = ce.getOperator();
	int cigarLen = ce.getLength();
	switch(op)
	    {
	    case M:
		{
		    effectiveLen += cigarLen;
		    break;
		}
	    case I:
		{
		    effectiveLen += cigarLen;
		    break;
		}
	    default: 
		break;
	    }
    }
}
boolean readdebug = false;
if(readdebug){
    HLA.log.appendln(sr.getSAMString());
    HLA.log.appendln("EffectiveLen:\t" + effectiveLen);
    HLA.log.appendln("ReadLen:\t" + rLen);
}
Integer i = sr.getIntegerAttribute("NM");
int nm = 0;
if(i!=null)
    nm = i.intValue();
if(readdebug)
    HLA.log.appendln("NM=\t" + nm);
if(nm < 16){
    if(readdebug)
	HLA.log.appendln("PASSWED QC");
    return true;
}
if(readdebug){
    HLA.log.appendln("FAILED QC");
    HLA.log.appendln(sr.getSAMString());
}
return false;
   }