|
libStatGen Software
1
|
|
libStatGen Software
1
|
Class for helping to filter a SAM/BAM record. More...
#include <CigarHelper.h>
Static Public Member Functions | |
| static int32_t | softClipBeginByRefPos (SamRecord &record, int32_t refPosition0Based, CigarRoller &newCigar, int32_t &new0BasedPosition) |
| Soft clip the cigar from the beginning of the read at the specified reference position. | |
| static int32_t | softClipEndByRefPos (SamRecord &record, int32_t refPosition0Based, CigarRoller &newCigar) |
| Soft clip the cigar from the back of the read at the specified reference position. | |
Static Public Attributes | |
| static const int32_t | NO_CLIP = -1 |
Class for helping to filter a SAM/BAM record.
Definition at line 24 of file CigarHelper.h.
| int32_t CigarHelper::softClipBeginByRefPos | ( | SamRecord & | record, |
| int32_t | refPosition0Based, | ||
| CigarRoller & | newCigar, | ||
| int32_t & | new0BasedPosition | ||
| ) | [static] |
Soft clip the cigar from the beginning of the read at the specified reference position.
If the clip position is deleted/skipped or is immediately followed by a deletion/skip/pad/insert, that entire CIGAR operation is also removed. Nothing is clipped if the reference position is before the read starts, everything is clipped if the reference position is after the read ends.
| record | record to calculate the clip for. |
| refPosition0Based | 0-based reference position to end the clip at |
| newCigar | cigar object to set with the updated cigar. |
| new0BasedPosition | new 0-based reference position of the read. |
| read | position where the clip ends (last clipped position) or |
Definition at line 23 of file CigarHelper.cpp.
References CigarRoller::Add(), CigarRoller::clear(), Cigar::foundInQuery(), Cigar::foundInReference(), SamRecord::get0BasedPosition(), SamRecord::getCigar(), SamRecord::getCigarInfo(), Cigar::getOperator(), ErrorHandler::handleError(), Cigar::hardClip, CigarRoller::Set(), Cigar::size(), and Cigar::softClip.
{
newCigar.clear();
Cigar* cigar = record.getCigarInfo();
if(cigar == NULL)
{
// Failed to get the cigar.
ErrorHandler::handleError("Soft clipping, but failed to read the cigar");
return(NO_CLIP);
}
// No cigar or position in the record, so return no clip.
if((cigar->size() == 0) || (record.get0BasedPosition() == -1))
{
return(NO_CLIP);
}
// Check to see if the reference position occurs before the record starts,
// if it does, do no clipping.
if(refPosition0Based < record.get0BasedPosition())
{
// Not within this read, so nothing to clip.
newCigar.Set(record.getCigar());
return(NO_CLIP);
}
// The position falls after the read starts, so loop through until the
// position or the end of the read is found.
int32_t readClipPosition = 0;
bool clipWritten = false;
new0BasedPosition = record.get0BasedPosition();
for(int i = 0; i < cigar->size(); i++)
{
const Cigar::CigarOperator* op = &(cigar->getOperator(i));
if(clipWritten)
{
// Clip point has been found, so just add everything.
newCigar += *op;
// Go to the next operation.
continue;
}
// The clip point has not yet been found, so check to see if we found
// it now.
// Not a clip, check to see if the operation is found in the
// reference.
if(Cigar::foundInReference(*op))
{
// match, mismatch, deletion, skip
// increment the current reference position to just past this
// operation.
new0BasedPosition += op->count;
// Check to see if this is also in the query, because otherwise
// the operation is still being consumed.
if(Cigar::foundInQuery(*op))
{
// Also in the query, determine if the entire thing should
// be clipped or just part of it.
uint32_t numKeep = 0;
// Check to see if we have hit our clip position.
if(refPosition0Based < new0BasedPosition)
{
// The specified clip position is in this cigar operation.
numKeep = new0BasedPosition - refPosition0Based - 1;
if(numKeep > op->count)
{
// Keep the entire read. This happens because
// we keep reading until the first match/mismatch
// after the clip.
numKeep = op->count;
}
}
// Add the part of this operation that is being clipped
// to the clip count.
readClipPosition += (op->count - numKeep);
// Only write the clip if we found a match/mismatch
// to write. Otherwise we will keep accumulating clips
// for the case of insertions.
if(numKeep > 0)
{
new0BasedPosition -= numKeep;
newCigar.Add(Cigar::softClip, readClipPosition);
// Add the clipped part of this cigar to the clip
// position.
newCigar.Add(op->operation, numKeep);
// Found a match after the clip point, so stop
// consuming cigar operations.
clipWritten = true;
continue;
}
}
}
else
{
// Only add hard clips. The softclips will be added in
// when the total number is found.
if(op->operation == Cigar::hardClip)
{
// Check if this is the first operation, if so, just write it.
if(i == 0)
{
newCigar += *op;
}
// Check if it is the last operation (otherwise skip it).
else if(i == (cigar->size() - 1))
{
// Check whether or not the clip was ever written, and if
// not, write it.
if(clipWritten == false)
{
newCigar.Add(Cigar::softClip, readClipPosition);
// Since no match/mismatch was ever found, set
// the new ref position to the original one.
new0BasedPosition = record.get0BasedPosition();
clipWritten = true;
}
// Add the hard clip.
newCigar += *op;
}
}
// Not yet to the clip position, so do not add this operation.
if(Cigar::foundInQuery(*op))
{
// Found in the query, so update the read clip position.
readClipPosition += op->count;
}
}
} // End loop through cigar.
// Check whether or not the clip was ever written, and if
// not, write it.
if(clipWritten == false)
{
newCigar.Add(Cigar::softClip, readClipPosition);
// Since no match/mismatch was ever found, set
// the new ref position to the original one.
new0BasedPosition = record.get0BasedPosition();
}
// Subtract 1 since readClipPosition atually contains the first 0based
// position that is not clipped.
return(readClipPosition - 1);
}
| int32_t CigarHelper::softClipEndByRefPos | ( | SamRecord & | record, |
| int32_t | refPosition0Based, | ||
| CigarRoller & | newCigar | ||
| ) | [static] |
Soft clip the cigar from the back of the read at the specified reference position.
If the clip position is deleted/skipped or is immediately preceded by a deletion/skip/pad, that entire CIGAR operation is also removed. If the clip position is immediately preceded by an insertion, the insertion is left in the CIGAR. Nothing is clipped if the reference position is after the read ends, everything is clipped if the reference position is before the read starts (including insertions).
| record | record to calculate the clip for. |
| refPosition0Based | 0-based reference position to start clip at |
| newCigar | cigar object to set with the updated cigar. |
| read | position where the clip starts or |
Definition at line 184 of file CigarHelper.cpp.
References CigarRoller::Add(), CigarRoller::clear(), Cigar::foundInQuery(), Cigar::foundInReference(), SamRecord::get0BasedAlignmentEnd(), SamRecord::get0BasedPosition(), SamRecord::getCigar(), SamRecord::getCigarInfo(), Cigar::getOperator(), SamRecord::getReadLength(), ErrorHandler::handleError(), Cigar::hardClip, Cigar::isClip(), CigarRoller::Remove(), CigarRoller::Set(), Cigar::size(), and Cigar::softClip.
{
newCigar.clear();
Cigar* cigar = record.getCigarInfo();
if(cigar == NULL)
{
// Failed to get the cigar.
ErrorHandler::handleError("Soft clipping, but failed to read the cigar");
return(NO_CLIP);
}
// No cigar or position in the record, so return no clip.
if((cigar->size() == 0) || (record.get0BasedPosition() == -1))
{
return(NO_CLIP);
}
// Check to see if the reference position occurs after the record ends,
// if so, do no clipping.
if(refPosition0Based > record.get0BasedAlignmentEnd())
{
// Not within this read, so nothing to clip.
newCigar.Set(record.getCigar());
return(NO_CLIP);
}
// The position falls before the read ends, so loop through until the
// position is found.
int32_t currentRefPosition = record.get0BasedPosition();
int32_t readClipPosition = 0;
for(int i = 0; i < cigar->size(); i++)
{
const Cigar::CigarOperator* op = &(cigar->getOperator(i));
// If the operation is found in the reference, increase the
// reference position.
if(Cigar::foundInReference(*op))
{
// match, mismatch, deletion, skip
// increment the current reference position to just past
// this operation.
currentRefPosition += op->count;
}
// Check to see if we have hit our clip position.
if(refPosition0Based < currentRefPosition)
{
// If this read is also in the query (match/mismatch),
// write the partial op to the new cigar.
int32_t numKeep = 0;
if(Cigar::foundInQuery(*op))
{
numKeep = op->count - (currentRefPosition - refPosition0Based);
if(numKeep > 0)
{
newCigar.Add(op->operation, numKeep);
readClipPosition += numKeep;
}
}
else if(Cigar::isClip(*op))
{
// This is a hard clip, so write it.
newCigar.Add(op->operation, op->count);
}
else
{
// Not found in the query (skip/deletion),
// so don't write any of the operation.
}
// Found the clip point, so break.
break;
}
else if(refPosition0Based == currentRefPosition)
{
newCigar += *op;
if(Cigar::foundInQuery(*op))
{
readClipPosition += op->count;
}
}
else
{
// Not yet to the clip position, so add this operation/size to
// the new cigar.
newCigar += *op;
if(Cigar::foundInQuery(*op))
{
// Found in the query, so update the read clip position.
readClipPosition += op->count;
}
}
} // End loop through cigar.
// Before adding the softclip, read from the end of the cigar checking to
// see if the operations are in the query, removing operations that are
// not (pad/delete/skip) until a hardclip or an operation in the query is
// found. We do not want a pad/delete/skip right before a softclip.
for(int j = newCigar.size() - 1; j >= 0; j--)
{
const Cigar::CigarOperator* op = &(newCigar.getOperator(j));
if(!Cigar::foundInQuery(*op) && !Cigar::isClip(*op))
{
// pad/delete/skip
newCigar.Remove(j);
}
else if(Cigar::foundInQuery(*op) & Cigar::isClip(*op))
{
// Soft clip, so increment the clip position for the return value.
// Remove the softclip since the readClipPosition is used to
// calculate teh size of the soft clip added.
readClipPosition -= op->count;
newCigar.Remove(j);
}
else
{
// Found a cigar operation that should not be deleted, so stop deleting.
break;
}
}
// Determine the number of soft clips.
int32_t numSoftClips = record.getReadLength() - readClipPosition;
// NOTE that if the previous operation is a softclip, the CigarRoller logic
// will merge this with that one.
newCigar.Add(Cigar::softClip, numSoftClips);
// Check if an ending hard clip needs to be added.
if(cigar->size() != 0)
{
const Cigar::CigarOperator* lastOp =
&(cigar->getOperator(cigar->size() - 1));
if(lastOp->operation == Cigar::hardClip)
{
newCigar += *lastOp;
}
}
return(readClipPosition);
}
1.7.6.1