public Cigar getCigar() { return (newCigar != null ? newCigar : read.getCigar()); }
@Test(enabled = !DEBUG) public void testRevertSoftClippedBasesWithThreshold() { for (Cigar cigar : cigarList) { final int leadingSoftClips = leadingCigarElementLength(cigar, CigarOperator.SOFT_CLIP); final int tailSoftClips = leadingCigarElementLength(CigarUtils.invertCigar(cigar), CigarOperator.SOFT_CLIP); final GATKSAMRecord read = ReadClipperTestUtils.makeReadFromCigar(cigar, 0); final GATKSAMRecord unclipped = ReadClipper.revertSoftClippedBases(read); assertUnclippedLimits(read, unclipped); // Make sure limits haven't changed Assert.assertNull(read.getCigar().isValid(null, -1)); Assert.assertNull(unclipped.getCigar().isValid(null, -1)); if (!(leadingSoftClips > 0 || tailSoftClips > 0)) Assert.assertEquals(read.getCigarString(), unclipped.getCigarString()); } }
@Override protected Double getElementForRead(final GATKSAMRecord read, final int refLoc) { final int offset = ReadUtils.getReadCoordinateForReferenceCoordinate( read.getSoftStart(), read.getCigar(), refLoc, ReadUtils.ClippingTail.RIGHT_TAIL, true ); if ( offset == ReadUtils.CLIPPING_GOAL_NOT_REACHED ) return null; // If the offset inside a deletion, it does not lie on a read. if ( AlignmentUtils.isInsideDeletion(read.getCigar(), offset) ) { return INVALID_ELEMENT_FROM_READ; } int readPos = AlignmentUtils.calcAlignmentByteArrayOffset( read.getCigar(), offset, false, 0, 0 ); final int numAlignedBases = AlignmentUtils.getNumAlignedBasesCountingSoftClips( read ); if (readPos > numAlignedBases / 2) readPos = numAlignedBases - (readPos + 1); return (double)readPos; }
@Test(dataProvider = "ComplexReadAlignedToRef", enabled = !DEBUG) public void testReadAlignedToRefComplexAlignment(final int testIndex, final GATKSAMRecord read, final String reference, final Haplotype haplotype, final int expectedMaxMismatches) throws Exception { final HaplotypeBAMWriter writer = new CalledHaplotypeBAMWriter(new MockDestination()); final GATKSAMRecord alignedRead = AlignmentUtils.createReadAlignedToRef(read, haplotype, new Haplotype(reference.getBytes(),true), 1, true); if ( alignedRead != null ) { final int mismatches = AlignmentUtils.getMismatchCount(alignedRead, reference.getBytes(), alignedRead.getAlignmentStart() - 1).numMismatches; Assert.assertTrue(mismatches <= expectedMaxMismatches, "Alignment of read to ref looks broken. Expected at most " + expectedMaxMismatches + " but saw " + mismatches + " for readBases " + new String(read.getReadBases()) + " with cigar " + read.getCigar() + " reference " + reference + " haplotype " + haplotype + " with cigar " + haplotype.getCigar() + " aligned read cigar " + alignedRead.getCigarString() + " @ " + alignedRead.getAlignmentStart()); } } }
/** * Is this read all insertion? * * @param read * @return whether or not the only element in the cigar string is an Insertion */ public static boolean readIsEntirelyInsertion(GATKSAMRecord read) { for (CigarElement cigarElement : read.getCigar().getCigarElements()) { if (cigarElement.getOperator() != CigarOperator.INSERTION) return false; } return true; }
@Override protected boolean isUsableBase(final PileupElement p) { return super.isUsableBase(p) && p.getRead().getCigar() != null; }
public static int getReadCoordinateForReferenceCoordinateUpToEndOfRead(GATKSAMRecord read, int refCoord, ClippingTail tail) { final int leftmostSafeVariantPosition = Math.max(read.getSoftStart(), refCoord); return getReadCoordinateForReferenceCoordinate(read.getSoftStart(), read.getCigar(), leftmostSafeVariantPosition, tail, false); }
@Requires({"read != null", "read.getAlignmentStart() != -1", "read.getCigar() != null"}) public AlignmentStateMachine(final GATKSAMRecord read) { this.read = read; this.cigar = read.getCigar(); this.nCigarElements = cigar.numCigarElements(); initializeAsLeftEdge(); }
/** * Checks whether or not the read has any cigar element that is not H or S * * @param read the read * @return true if it has any M, I or D, false otherwise */ public static boolean readHasNonClippedBases(GATKSAMRecord read) { for (CigarElement cigarElement : read.getCigar().getCigarElements()) if (cigarElement.getOperator() != CigarOperator.SOFT_CLIP && cigarElement.getOperator() != CigarOperator.HARD_CLIP) return true; return false; }
@Override protected boolean isUsableRead(final GATKSAMRecord read, final int refLoc) { return super.isUsableRead(read, refLoc) && read.getSoftStart() + read.getCigar().getReadLength() > refLoc; } }
@Override protected Double getElementForPileupElement(final PileupElement p) { final int offset = AlignmentUtils.calcAlignmentByteArrayOffset(p.getRead().getCigar(), p, 0, 0); return (double)AnnotationUtils.getFinalVariantReadPosition(p.getRead(), offset); }
@Test(dataProvider = "ClippedReadLengthData", enabled = !DEBUG) public void testHardClipReadLengthIsRight(final int originalReadLength, final int nToClip) { GATKSAMRecord read = ReadClipperTestUtils.makeReadFromCigar(originalReadLength + "M", 0); read.getReadLength(); // provoke the caching of the read length final int expectedReadLength = originalReadLength - nToClip; GATKSAMRecord clipped = ReadClipper.hardClipByReadCoordinates(read, 0, nToClip - 1); Assert.assertEquals(clipped.getReadLength(), expectedReadLength, String.format("Clipped read length %d with cigar %s not equal to the expected read length %d after clipping %d bases from the left from a %d bp read with cigar %s", clipped.getReadLength(), clipped.getCigar(), expectedReadLength, nToClip, read.getReadLength(), read.getCigar())); }
public CigarCounter(GATKSAMRecord read) { CigarOperator[] operators = CigarOperator.values(); counter = new HashMap<CigarOperator, Integer>(operators.length); for (CigarOperator op : operators) counter.put(op, 0); for (CigarElement cigarElement : read.getCigar().getCigarElements()) counter.put(cigarElement.getOperator(), counter.get(cigarElement.getOperator()) + cigarElement.getLength()); }
public Integer map(ReferenceContext ref, GATKSAMRecord read, RefMetaDataTracker metaDataTracker) { // we can not deal with screwy records if ( read.getReadUnmappedFlag() || read.getCigar().numCigarElements() == 0 ) { emit(read); return 0; } // move existing indels (for 1 indel reads only) to leftmost position within identical sequence int numBlocks = AlignmentUtils.getNumAlignmentBlocks(read); if ( numBlocks == 2 ) { Cigar newCigar = AlignmentUtils.leftAlignIndel(IndelRealigner.unclipCigar(read.getCigar()), ref.getBases(), read.getReadBases(), 0, 0, true); newCigar = IndelRealigner.reclipCigar(newCigar, read); read.setCigar(newCigar); } emit(read); return 1; }
@Override public GATKSAMRecord apply(final GATKSAMRecord read) { if(read == null) throw new UserException.BadInput("try to transform a null GATKSAMRecord"); final Cigar originalCigar = read.getCigar(); if (originalCigar.isValid(read.getReadName(),-1) != null) throw new UserException.BadInput("try to transform a read with non-valid cigar string: readName: "+read.getReadName()+" Cigar String: "+originalCigar); read.setCigar(refactorNDNtoN(originalCigar)); return read; }
private Cigar makeExpectedCigar1(final Cigar originalCigar, final char indelOp, final int indelStart, final int indelSize, final int readLength) { if ( indelSize == 0 || indelStart < 17 || indelStart > (26 - (indelOp == 'D' ? indelSize : 0)) ) return originalCigar; final GATKSAMRecord read = ArtificialSAMUtils.createArtificialRead(header, "myRead", 0, 1, readLength); read.setCigarString(buildTestCigarString(indelOp, 0, 16, indelSize, readLength)); return read.getCigar(); }
public Pair<Long, Long> map(ReferenceContext ref, GATKSAMRecord read, RefMetaDataTracker tracker) { List<CigarElement> cigarElements = read.getCigar().getCigarElements(); CigarElement lastElement = null; for (CigarElement element : cigarElements) { if (element.getOperator() != CigarOperator.HARD_CLIP) lastElement = element; } if (lastElement == null) throw new UserException.MalformedBAM(read, "read does not have any bases, it's all hard clips"); long endsInIndel = lastElement.getOperator() == CigarOperator.INSERTION || lastElement.getOperator() == CigarOperator.DELETION? 1 : 0; long endsInSC = lastElement.getOperator() == CigarOperator.SOFT_CLIP ? 1 : 0; return new Pair<Long, Long>(endsInIndel, endsInSC); }
private Cigar makeExpectedCigar2(final Cigar originalCigar, final char indelOp, final int indelStart, final int indelSize, final int readLength) { if ( indelStart < 17 || indelStart > (26 - (indelOp == 'D' ? indelSize : 0)) ) return originalCigar; final GATKSAMRecord read = ArtificialSAMUtils.createArtificialRead(header, "myRead", 0, 1, readLength); if ( indelOp == 'I' && (indelSize == 1 || indelSize == 3) && indelStart % 2 == 1 ) read.setCigarString(buildTestCigarString(indelOp, 0, Math.max(indelStart - indelSize, 16), indelSize, readLength)); else if ( (indelSize == 2 || indelSize == 4) && (indelOp == 'D' || indelStart % 2 == 0) ) read.setCigarString(buildTestCigarString(indelOp, 0, 16, indelSize, readLength)); else return originalCigar; return read.getCigar(); }
private Cigar makeExpectedCigar3(final Cigar originalCigar, final char indelOp, final int indelStart, final int indelSize, final int readLength) { if ( indelStart < 17 || indelStart > (28 - (indelOp == 'D' ? indelSize : 0)) ) return originalCigar; final GATKSAMRecord read = ArtificialSAMUtils.createArtificialRead(header, "myRead", 0, 1, readLength); if ( indelSize == 3 && (indelOp == 'D' || indelStart % 3 == 1) ) read.setCigarString(buildTestCigarString(indelOp, 0, 16, indelSize, readLength)); else if ( (indelOp == 'I' && indelSize == 4 && indelStart % 3 == 2) || (indelOp == 'I' && indelSize == 2 && indelStart % 3 == 0) || (indelOp == 'I' && indelSize == 1 && indelStart < 28 && indelStart % 3 == 2) ) read.setCigarString(buildTestCigarString(indelOp, 0, Math.max(indelStart - indelSize, 16), indelSize, readLength)); else return originalCigar; return read.getCigar(); }
private List<PileupElement> makeReads( final int n, final int mapq, final String op ) { final int readLength = 3; final List<PileupElement> elts = new LinkedList<PileupElement>(); for ( int i = 0; i < n; i++ ) { GATKSAMRecord read = ArtificialSAMUtils.createArtificialRead(header, "read", 0, 1, readLength); read.setReadBases(Utils.dupBytes((byte) 'A', readLength)); read.setBaseQualities(Utils.dupBytes((byte) 30, readLength)); read.setCigarString("1M1" + op + "1M"); read.setMappingQuality(mapq); final int baseOffset = op.equals("M") ? 1 : 0; final CigarElement cigarElement = read.getCigar().getCigarElement(1); elts.add(new PileupElement(read, baseOffset, cigarElement, 1, 0)); } return elts; }