use of htsjdk.samtools.CigarElement in project ASCIIGenome by dariober.
the class TextReadTest method getNumAlignedBasesFromCigar.
@Test
public void getNumAlignedBasesFromCigar() {
// A reminder of how cigar works.
SAMRecord rec = new SAMRecord(null);
rec.setCigarString("2M2D4M");
int noRefBases = rec.getCigar().getReferenceLength();
assertEquals(8, noRefBases);
// N same as D
rec.setCigarString("2M2N4M");
noRefBases = rec.getCigar().getReferenceLength();
assertEquals(8, noRefBases);
rec.setCigarString("2M2I4M");
noRefBases = rec.getCigar().getReferenceLength();
assertEquals(6, noRefBases);
rec.setCigarString("2S10M4H");
noRefBases = rec.getCigar().getReferenceLength();
assertEquals(10, noRefBases);
// Cigar oprators consuming:
// Hard clipping doesn't consume bases on read or ref
rec.setCigarString("2H");
CigarElement el = rec.getCigar().getCigarElement(0);
assertFalse(el.getOperator().consumesReadBases());
assertFalse(el.getOperator().consumesReferenceBases());
// Clips the read not the ref.
rec.setCigarString("2S");
el = rec.getCigar().getCigarElement(0);
assertTrue(el.getOperator().consumesReadBases());
assertFalse(el.getOperator().consumesReferenceBases());
// This makes a gap in the read when aligned to ref
rec.setCigarString("2D");
el = rec.getCigar().getCigarElement(0);
assertFalse(el.getOperator().consumesReadBases());
assertTrue(el.getOperator().consumesReferenceBases());
// Opposite of D
rec.setCigarString("2I");
el = rec.getCigar().getCigarElement(0);
assertTrue(el.getOperator().consumesReadBases());
assertFalse(el.getOperator().consumesReferenceBases());
}
use of htsjdk.samtools.CigarElement in project ASCIIGenome by dariober.
the class TextRead method getSoftUnclippedAlignmentStart.
public int getSoftUnclippedAlignmentStart(SAMRecord rec) {
List<CigarElement> cigar = rec.getCigar().getCigarElements();
if (cigar.size() == 0) {
return rec.getUnclippedStart();
}
if (cigar.size() == 1 && cigar.get(0).getOperator().equals(CigarOperator.SOFT_CLIP)) {
return rec.getUnclippedStart();
}
int offset = 0;
for (int i = 0; i < cigar.size(); i++) {
CigarElement op = cigar.get(i);
if (op.getOperator().equals(CigarOperator.HARD_CLIP)) {
continue;
} else if (op.getOperator().equals(CigarOperator.SOFT_CLIP)) {
offset += op.getLength();
} else {
break;
}
}
int start = rec.getAlignmentStart() - offset;
// }
return start;
}
use of htsjdk.samtools.CigarElement in project ASCIIGenome by dariober.
the class TextRead method setTextPositionsOfSkippedBases.
/**
*List of positions on screen where the skipped bases (cigar op: N) start
* @throws IOException
* @throws InvalidGenomicCoordsException
*/
private void setTextPositionsOfSkippedBases() throws InvalidGenomicCoordsException, IOException {
int genomicPosition = this.getSamRecord().getAlignmentStart();
List<CigarElement> cigar = this.getSamRecord().getCigar().getCigarElements();
for (CigarElement el : cigar) {
if (el.getOperator().equals(CigarOperator.SKIPPED_REGION)) {
int[] textPositions = new int[2];
// +1 because textPosition is 1-based
textPositions[0] = Utils.getIndexOfclosestValue(genomicPosition, this.gc.getMapping()) + 1;
textPositions[1] = Utils.getIndexOfclosestValue(genomicPosition + el.getLength(), this.gc.getMapping()) + 1;
this.textPositionsOfSkippedBases.add(textPositions);
}
;
if (el.getOperator().consumesReferenceBases()) {
genomicPosition += el.getLength();
}
}
}
use of htsjdk.samtools.CigarElement in project ASCIIGenome by dariober.
the class Track method isSNVRead.
/**
*Return true if samrecord contains a mismatch or insertion/deletion in the target region.
*/
private boolean isSNVRead(SAMRecord rec, boolean variantOnly) {
boolean passed = false;
int varFrom = this.getFeatureFilter().getVariantFrom();
int varTo = this.getFeatureFilter().getVariantTo();
if (this.getFeatureFilter().getVariantChrom().equals(rec.getReferenceName()) && varFrom <= rec.getAlignmentEnd() && rec.getAlignmentStart() <= varTo) {
// Variant read filter is set and this read overlaps it.
if (!variantOnly) {
// No need to check whether read is variant.
return true;
}
int readPos = 0;
int refPos = rec.getAlignmentStart();
for (CigarElement cigar : rec.getCigar().getCigarElements()) {
if (cigar.getOperator().equals(CigarOperator.SOFT_CLIP)) {
readPos += cigar.getLength();
} else if (cigar.getOperator().equals(CigarOperator.M) || cigar.getOperator().equals(CigarOperator.EQ) || cigar.getOperator().equals(CigarOperator.X)) {
for (int i = 0; i < cigar.getLength(); i++) {
if (refPos >= varFrom && refPos <= varTo && rec.getReadLength() > 0) {
byte readBase = rec.getReadBases()[readPos];
byte refBase = this.getFeatureFilter().getFaSeq()[refPos - varFrom];
if (readBase != refBase) {
passed = true;
break;
}
}
readPos++;
refPos++;
}
} else if (cigar.getOperator().equals(CigarOperator.DELETION)) {
// ^^^^
for (int i = 0; i < cigar.getLength(); i++) {
if (refPos >= varFrom && refPos <= varTo) {
passed = true;
break;
}
refPos++;
}
} else if (cigar.getOperator().equals(CigarOperator.INSERTION)) {
// ^
for (int i = 0; i < cigar.getLength(); i++) {
if (refPos >= varFrom && refPos <= varTo) {
passed = true;
break;
}
readPos++;
}
} else if (cigar.getOperator().equals(CigarOperator.HARD_CLIP)) {
//
} else if (cigar.getOperator().equals(CigarOperator.SKIPPED_REGION)) {
// Same deletion but it's not a mismatch
refPos += cigar.getLength();
} else if (cigar.getOperator().equals(CigarOperator.PADDING)) {
// Not sure what to do with this...
}
if (passed) {
break;
}
}
} else {
// Variant read filter is set and this read does not overlap it.
passed = false;
}
return passed;
}
use of htsjdk.samtools.CigarElement in project gridss by PapenfussLab.
the class CigarUtil method convertOpenEndedOperations.
/**
* Converts indel operations that do not contain anchoring alignments on both sides to
* a soft cliped version of the operator
* @param list
*/
private static void convertOpenEndedOperations(List<CigarElement> list) {
for (int i = list.size() - 1; i >= 0; i--) {
CigarElement ce = list.get(i);
if (ce.getOperator() == CigarOperator.SOFT_CLIP || ce.getOperator() == CigarOperator.HARD_CLIP) {
} else if (ce.getOperator() == CigarOperator.INSERTION) {
list.set(i, new CigarElement(ce.getLength(), CigarOperator.SOFT_CLIP));
} else if (ce.getOperator() == CigarOperator.DELETION) {
list.remove(i);
} else {
break;
}
}
for (int i = 0; i < list.size(); i++) {
CigarElement ce = list.get(i);
if (ce.getOperator() == CigarOperator.SOFT_CLIP || ce.getOperator() == CigarOperator.HARD_CLIP) {
} else if (ce.getOperator() == CigarOperator.INSERTION) {
list.set(i, new CigarElement(ce.getLength(), CigarOperator.SOFT_CLIP));
} else if (ce.getOperator() == CigarOperator.DELETION) {
list.remove(i);
i--;
} else {
break;
}
}
}
Aggregations