diff --git a/java/src/org/broadinstitute/sting/gatk/contexts/StratifiedAlignmentContext.java b/java/src/org/broadinstitute/sting/gatk/contexts/StratifiedAlignmentContext.java index e14fec436..c4199eda6 100755 --- a/java/src/org/broadinstitute/sting/gatk/contexts/StratifiedAlignmentContext.java +++ b/java/src/org/broadinstitute/sting/gatk/contexts/StratifiedAlignmentContext.java @@ -38,7 +38,7 @@ import java.util.*; * User: ebanks * Modified: chartl (split by read group) */ -public class StratifiedAlignmentContext { +public class StratifiedAlignmentContext { // Definitions: // COMPLETE = full alignment context @@ -48,98 +48,35 @@ public class StratifiedAlignmentContext { public enum StratifiedContextType { COMPLETE, FORWARD, REVERSE } private GenomeLoc loc; - private AlignmentContext[] contexts = new AlignmentContext[StratifiedContextType.values().length]; - private boolean isExtended = false; // tells whether this alignment context is an extended event context + private RBP basePileup = null; - // todo -- why are you storing reads separately each time? There's a ReadBackedPileup object that's supposed to handle this -// private ArrayList[] reads = new ArrayList[StratifiedContextType.values().length]; -// private ArrayList[] offsets = new ArrayList[StratifiedContextType.values().length]; - - private ArrayList[] pileupElems = new ArrayList[StratifiedContextType.values().length]; // // accessors // public GenomeLoc getLocation() { return loc; } -// public ArrayList getReads(StratifiedContextType type) { return reads[type.ordinal()]; } -// public ArrayList getOffsets(StratifiedContextType type) { return offsets[type.ordinal()]; } - - public ArrayList getPileupElements(StratifiedContextType type) { - return pileupElems[type.ordinal()]; - } - -// public ArrayList getExtendedPileupElements(StratifiedContextType type) { -// if ( ! isExtended ) throw new StingException("Extended read backed pileups requested from StratifiedAlignmentContext that holds simple pileups"); -// -// return (ArrayList)(pileupElems[type.ordinal()]); -// } public StratifiedAlignmentContext(GenomeLoc loc) { - this(loc,false); + this(loc,null); } - public StratifiedAlignmentContext(GenomeLoc loc, boolean isExtended) { + public StratifiedAlignmentContext(GenomeLoc loc, RBP pileup) { this.loc = loc; - this.isExtended = isExtended; - for ( int i = 0; i < StratifiedContextType.values().length; i++) { - if ( isExtended ) pileupElems[i] = new ArrayList(); - else pileupElems[i] = new ArrayList(); - } + this.basePileup = pileup; } public AlignmentContext getContext(StratifiedContextType type) { - int index = type.ordinal(); - if ( contexts[index] == null ) { - if ( isExtended ) { - contexts[index] = new AlignmentContext(loc, new ReadBackedExtendedEventPileupImpl(loc,(List)getPileupElements(type))); - } else { - contexts[index] = new AlignmentContext(loc, new ReadBackedPileupImpl(loc,(List)getPileupElements(type))); - } + switch(type) { + case COMPLETE: + return new AlignmentContext(loc,basePileup); + case FORWARD: + return new AlignmentContext(loc,basePileup.getPositiveStrandPileup()); + case REVERSE: + return new AlignmentContext(loc,basePileup.getNegativeStrandPileup()); + default: + throw new StingException("Unable to get alignment context for type = " + type); } - return contexts[index]; } - public void add(SAMRecord read, int offset) { - if ( isExtended ) throw new StingException("Can not add read/offset without event type specified to the context holding extended events"); - if ( read.getReadNegativeStrandFlag() ) { - pileupElems[StratifiedContextType.REVERSE.ordinal()].add((PE)new PileupElement(read,offset)); - } else { - pileupElems[StratifiedContextType.FORWARD.ordinal()].add((PE)new PileupElement(read,offset)); - } - pileupElems[StratifiedContextType.COMPLETE.ordinal()].add((PE)new PileupElement(read,offset)); - } - - public void add(PE p) { -// if ( isExtended ) throw new StingException("Can not add simple pileup element to the context holding extended events"); - SAMRecord read = p.getRead(); - if ( read.getReadNegativeStrandFlag() ) { - pileupElems[StratifiedContextType.REVERSE.ordinal()].add(p); - } else { - pileupElems[StratifiedContextType.FORWARD.ordinal()].add(p); - } - pileupElems[StratifiedContextType.COMPLETE.ordinal()].add(p); - } - - public void add(SAMRecord read, int offset, int length, byte [] bases) { - if ( ! isExtended ) throw new StingException("Can not add read/offset with event type specified to the context holding simple events"); - if ( read.getReadNegativeStrandFlag() ) { - pileupElems[StratifiedContextType.REVERSE.ordinal()].add((PE)new ExtendedEventPileupElement(read,offset,length,bases)); - } else { - pileupElems[StratifiedContextType.FORWARD.ordinal()].add((PE)new ExtendedEventPileupElement(read,offset,length,bases)); - } - pileupElems[StratifiedContextType.COMPLETE.ordinal()].add((PE)new ExtendedEventPileupElement(read,offset,length,bases)); - } - -// public void add(ExtendedEventPileupElement p) { -// if ( ! isExtended ) throw new StingException("Can not add extended pileup element to the context holding simple events"); -// SAMRecord read = p.getRead(); -// if ( read.getReadNegativeStrandFlag() ) { -// pileupElems[StratifiedContextType.REVERSE.ordinal()].add(p); -// } else { -// pileupElems[StratifiedContextType.FORWARD.ordinal()].add(p); -// } -// pileupElems[StratifiedContextType.COMPLETE.ordinal()].add(p); -// } - /** * Splits the given AlignmentContext into a StratifiedAlignmentContext per sample. * @@ -148,8 +85,8 @@ public class StratifiedAlignmentContext { * @return a Map of sample name to StratifiedAlignmentContext * **/ - public static Map splitContextBySample(ReadBackedPileup pileup) { - return splitContextBySample(pileup, null, null); + public static Map splitContextBySample(RBP pileup) { + return splitContextBySample(pileup, null); } /** @@ -157,89 +94,28 @@ public class StratifiedAlignmentContext { * * @param pileup the original pileup * @param assumedSingleSample if not null, any read without a readgroup will be given this sample name - * @param collapseToThisSample if not null, all reads will be assigned this read group regardless of their actual read group * * @return a Map of sample name to StratifiedAlignmentContext * **/ - public static Map splitContextBySample(ReadBackedPileup pileup, String assumedSingleSample, String collapseToThisSample) { + public static Map splitContextBySample(RBP pileup, String assumedSingleSample) { - HashMap contexts = new HashMap(); GenomeLoc loc = pileup.getLocation(); - - for (PileupElement p : pileup ) - addToContext(contexts, p, loc, assumedSingleSample, collapseToThisSample); - - return contexts; - } - - /** - * Splits the given AlignmentContext into a StratifiedAlignmentContext per sample. - * - * @param pileup the original pileup - * - * @return a Map of sample name to StratifiedAlignmentContext - * - **/ - public static Map splitContextBySample(ReadBackedExtendedEventPileup pileup) { - return splitContextBySample(pileup, null, null); - } - - /** - * Splits the given AlignmentContext into a StratifiedAlignmentContext per sample. - * - * @param pileup the original pileup - * @param assumedSingleSample if not null, any read without a readgroup will be given this sample name - * @param collapseToThisSample if not null, all reads will be assigned this read group regardless of their actual read group - * - * @return a Map of sample name to StratifiedAlignmentContext - * - **/ - public static Map splitContextBySample(ReadBackedExtendedEventPileup pileup, String assumedSingleSample, String collapseToThisSample) { - HashMap contexts = new HashMap(); - GenomeLoc loc = pileup.getLocation(); - for (PileupElement p : pileup ) - addToContext(contexts, p, loc, assumedSingleSample, collapseToThisSample,true); + for(String sampleName: pileup.getSamples()) { + RBP pileupBySample = (RBP)pileup.getPileupForSample(sampleName); - return contexts; - } - - private static void addToContext(HashMap contexts, PileupElement p, GenomeLoc loc, String assumedSingleSample, String collapseToThisSample) { - addToContext(contexts, p, loc, assumedSingleSample, collapseToThisSample, false); - } - - private static void addToContext(HashMap contexts, PileupElement p, GenomeLoc loc, String assumedSingleSample, String collapseToThisSample, boolean isExtended) { - - // get the read - SAMRecord read = p.getRead(); - - // find the sample - String sample; - if ( collapseToThisSample != null ) { - sample = collapseToThisSample; - } else { - SAMReadGroupRecord readGroup = read.getReadGroup(); - if ( readGroup == null ) { - if ( assumedSingleSample == null ) - throw new StingException("Missing read group for read " + read.getReadName()); - sample = assumedSingleSample; - } else { - sample = readGroup.getSample(); + if(sampleName != null) + contexts.put(sampleName,new StratifiedAlignmentContext(loc,pileupBySample)); + else { + if(assumedSingleSample == null && pileupBySample.size() > 0) + throw new StingException("Missing read group for read " + pileupBySample.iterator().next().getRead()); + contexts.put(assumedSingleSample,new StratifiedAlignmentContext(loc,pileupBySample)); } } - // create a new context object if this is the first time we're seeing a read for this sample - StratifiedAlignmentContext myContext = contexts.get(sample); - if ( myContext == null ) { - myContext = new StratifiedAlignmentContext(loc,isExtended); - contexts.put(sample, myContext); - } - - // add the read to this sample's context - // note that bad bases are added to the context (for DoC calculations later) - myContext.add(p); + return contexts; } /** @@ -250,57 +126,43 @@ public class StratifiedAlignmentContext { * TODO - support for collapsing or assuming read groups if they are missing * **/ - public static Map splitContextByReadGroup(ReadBackedPileup pileup) { - HashMap contexts = new HashMap(); - - for ( PileupElement p : pileup ) { - SAMRecord read = p.getRead(); - - SAMReadGroupRecord readGroup = read.getReadGroup(); - if ( readGroup == null ) { - throw new StingException("Missing read group for read " + read.getReadName()); - } - - String group = readGroup.getReadGroupId(); - - StratifiedAlignmentContext myContext = contexts.get(group); - - if ( myContext == null ) { - myContext = new StratifiedAlignmentContext(pileup.getLocation()); - contexts.put(group,myContext); - } - - myContext.add(p); - } - - return contexts; + public static Map> splitContextByReadGroup(RBP pileup) { + HashMap> contexts = new HashMap>(); + for(String readGroupId: pileup.getReadGroups()) + contexts.put(readGroupId,new StratifiedAlignmentContext(pileup.getLocation(),(RBP)pileup.getPileupForReadGroup(readGroupId))); + return contexts; } - public static AlignmentContext joinContexts(Collection contexts, StratifiedContextType type) { - ArrayList pe = new ArrayList(); - - if ( contexts.size() == 0 ) - throw new StingException("BUG: joinContexts requires at least one context to join"); - - - Iterator it = contexts.iterator(); - StratifiedAlignmentContext context = it.next(); - boolean isExtended = context.isExtended; - GenomeLoc loc = context.getLocation(); - pe.addAll(context.getPileupElements(type)); - - while ( it.hasNext()) { - context = it.next(); - if ( ! loc.equals( context.getLocation() ) ) - throw new StingException("Illegal attempt to join contexts from different genomic locations"); - if ( context.isExtended != isExtended ) + public static AlignmentContext joinContexts(Collection contexts) { + + // validation + GenomeLoc loc = contexts.iterator().next().getLocation(); + boolean isExtended = contexts.iterator().next().basePileup instanceof ReadBackedExtendedEventPileup; + for(StratifiedAlignmentContext context: contexts) { + if(!loc.equals(context.getLocation())) + throw new StingException("Illegal attempt to join contexts from different genomic locations"); + if(isExtended != (context.basePileup instanceof ReadBackedExtendedEventPileup)) throw new StingException("Illegal attempt to join simple and extended contexts"); - pe.addAll(context.getPileupElements(type)); } - // dirty trick below. generics do not allow to cast pe (ArrayList) directly to ArrayList, - // so we first cast to "? extends" wildcard, then to what we actually need. - if ( isExtended ) return new AlignmentContext(loc, new ReadBackedExtendedEventPileupImpl(loc,(List)pe) ); - else return new AlignmentContext(loc, new ReadBackedPileupImpl(loc,(List)pe)); + AlignmentContext jointContext; + if(isExtended) { + List pe = new ArrayList(); + for(StratifiedAlignmentContext context: contexts) { + for(PileupElement pileupElement: context.basePileup) + pe.add((ExtendedEventPileupElement)pileupElement); + } + jointContext = new AlignmentContext(loc, new ReadBackedExtendedEventPileupImpl(loc,pe)); + } + else { + List pe = new ArrayList(); + for(StratifiedAlignmentContext context: contexts) { + for(PileupElement pileupElement: context.basePileup) + pe.add(pileupElement); + } + jointContext = new AlignmentContext(loc, new ReadBackedPileupImpl(loc,pe)); + } + + return jointContext; } } \ No newline at end of file diff --git a/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java b/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java index 67deb7b1e..3b452d703 100755 --- a/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java +++ b/java/src/org/broadinstitute/sting/gatk/iterators/LocusIteratorByState.java @@ -359,8 +359,8 @@ public class LocusIteratorByState extends LocusIterator { // In this case, the subsequent call to next() will emit the normal pileup at the current base // and shift the position. if (readInfo.generateExtendedEvents() && hasExtendedEvents) { - Map> fullExtendedEventPileup = - new HashMap>(); + Map fullExtendedEventPileup = + new HashMap(); SAMRecordState our1stState = readStates.getFirst(); // get current location on the reference and decrement it by 1: the indels we just stepped over @@ -422,7 +422,7 @@ public class LocusIteratorByState extends LocusIterator { nextAlignmentContext = new AlignmentContext(loc, new ReadBackedExtendedEventPileupImpl(loc, fullExtendedEventPileup)); } else { GenomeLoc location = getLocation(); - Map> fullPileup = new HashMap>(); + Map fullPileup = new HashMap(); // todo -- performance problem -- should be lazy, really for(String sampleName: sampleNames) { diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/HaplotypeScore.java b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/HaplotypeScore.java index 159da7025..e0e636084 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/HaplotypeScore.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/HaplotypeScore.java @@ -52,7 +52,7 @@ public class HaplotypeScore implements InfoFieldAnnotation, StandardAnnotation { if ( !vc.isBiallelic() || !vc.isSNP() || stratifiedContexts.size() == 0 ) // size 0 means that call was made by someone else and we have no data here return null; - AlignmentContext context = StratifiedAlignmentContext.joinContexts(stratifiedContexts.values(), StratifiedAlignmentContext.StratifiedContextType.COMPLETE); + AlignmentContext context = StratifiedAlignmentContext.joinContexts(stratifiedContexts.values()); int contextWingSize = Math.min(((int)ref.getWindow().size() - 1)/2, MIN_CONTEXT_WING_SIZE); int contextSize = contextWingSize * 2 + 1; diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/VariantAnnotator.java b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/VariantAnnotator.java index 1bc88e948..01f6ffe08 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/annotator/VariantAnnotator.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/annotator/VariantAnnotator.java @@ -200,9 +200,9 @@ public class VariantAnnotator extends RodWalker { Map stratifiedContexts; if ( BaseUtils.simpleBaseToBaseIndex(ref.getBase()) != -1 ) { if ( ! context.hasExtendedEventPileup() ) { - stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(context.getBasePileup(), ASSUME_SINGLE_SAMPLE, null); + stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(context.getBasePileup(), ASSUME_SINGLE_SAMPLE); } else { - stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(context.getExtendedEventPileup(), ASSUME_SINGLE_SAMPLE, null); + stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(context.getExtendedEventPileup(), ASSUME_SINGLE_SAMPLE); } if ( stratifiedContexts != null ) { annotatedVCs = engine.annotateContext(tracker, ref, stratifiedContexts, vc); diff --git a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java index 2a8ee1679..3e79b7094 100755 --- a/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java +++ b/java/src/org/broadinstitute/sting/gatk/walkers/genotyper/UnifiedGenotyperEngine.java @@ -179,7 +179,7 @@ public class UnifiedGenotyperEngine { return null; // stratify the AlignmentContext and cut by sample - Map stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(pileup, UAC.ASSUME_SINGLE_SAMPLE, null); + Map stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(pileup, UAC.ASSUME_SINGLE_SAMPLE); if ( stratifiedContexts == null ) return null; @@ -206,7 +206,7 @@ public class UnifiedGenotyperEngine { return null; // stratify the AlignmentContext and cut by sample - Map stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(pileup, UAC.ASSUME_SINGLE_SAMPLE, null); + Map stratifiedContexts = StratifiedAlignmentContext.splitContextBySample(pileup, UAC.ASSUME_SINGLE_SAMPLE); if ( stratifiedContexts == null ) return null; diff --git a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/AlleleBalanceHistogramWalker.java b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/AlleleBalanceHistogramWalker.java index 778f4c7c5..2846be7f9 100644 --- a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/AlleleBalanceHistogramWalker.java +++ b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/AlleleBalanceHistogramWalker.java @@ -68,7 +68,7 @@ public class AlleleBalanceHistogramWalker extends LocusWalker } private HashMap getAlleleBalanceBySample(VCFRecord vcf, ReferenceContext ref, AlignmentContext context) { - Map sampleContext = StratifiedAlignmentContext.splitContextBySample(context.getBasePileup(),null,null); + Map sampleContext = StratifiedAlignmentContext.splitContextBySample(context.getBasePileup(),null); HashMap balances = new HashMap(); System.out.println("----- "+ref.getLocus()+" -----"); int returnedBalances = 0; diff --git a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/MendelianViolationClassifier.java b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/MendelianViolationClassifier.java index 93c518f48..e5b18455b 100644 --- a/java/src/org/broadinstitute/sting/oneoffprojects/walkers/MendelianViolationClassifier.java +++ b/java/src/org/broadinstitute/sting/oneoffprojects/walkers/MendelianViolationClassifier.java @@ -533,11 +533,11 @@ public class MendelianViolationClassifier extends LocusWalker context) { + private Double getAlleleProportion(Allele a, StratifiedAlignmentContext context) { int numParental = 0; int total = 0; if ( context != null ) { - for ( PileupElement e : context.getPileupElements(StratifiedAlignmentContext.StratifiedContextType.COMPLETE)) { + for ( PileupElement e : context.getContext(StratifiedAlignmentContext.StratifiedContextType.COMPLETE).getBasePileup()) { if ( e.getQual() >= 10 && e.getMappingQual() >= 10 ) { total++; if ( e.getBase() == a.getBases()[0]) { diff --git a/java/src/org/broadinstitute/sting/utils/pileup/AbstractReadBackedPileup.java b/java/src/org/broadinstitute/sting/utils/pileup/AbstractReadBackedPileup.java index 37127b0c8..3e1088a3f 100644 --- a/java/src/org/broadinstitute/sting/utils/pileup/AbstractReadBackedPileup.java +++ b/java/src/org/broadinstitute/sting/utils/pileup/AbstractReadBackedPileup.java @@ -39,7 +39,7 @@ import net.sf.samtools.SAMRecord; * @author mhanna * @version 0.1 */ -public abstract class AbstractReadBackedPileup implements ReadBackedPileup { +public abstract class AbstractReadBackedPileup,PE extends PileupElement> implements ReadBackedPileup { protected final GenomeLoc loc; protected final PileupElementTracker pileupElementTracker; @@ -82,7 +82,7 @@ public abstract class AbstractReadBackedPileup pileup ) { + public AbstractReadBackedPileup(GenomeLoc loc, List pileup) { if ( loc == null ) throw new StingException("Illegal null genomeloc in ReadBackedPileup"); if ( pileup == null ) throw new StingException("Illegal null pileup in ReadBackedPileup"); @@ -115,10 +115,10 @@ public abstract class AbstractReadBackedPileup> pileupsBySample) { + protected AbstractReadBackedPileup(GenomeLoc loc, Map> pileupsBySample) { this.loc = loc; PerSamplePileupElementTracker tracker = new PerSamplePileupElementTracker(); - for(Map.Entry> pileupEntry: pileupsBySample.entrySet()) { + for(Map.Entry> pileupEntry: pileupsBySample.entrySet()) { tracker.addElements(pileupEntry.getKey(),pileupEntry.getValue().pileupElementTracker); addPileupToCumulativeStats(pileupEntry.getValue()); } @@ -192,7 +192,7 @@ public abstract class AbstractReadBackedPileup pileupElementTracker); + protected abstract AbstractReadBackedPileup createNewPileup(GenomeLoc loc, PileupElementTracker pileupElementTracker); protected abstract PE createNewPileupElement(SAMRecord read, int offset); // -------------------------------------------------------- @@ -217,10 +217,10 @@ public abstract class AbstractReadBackedPileup perSampleElements = tracker.getElements(sampleName); - AbstractReadBackedPileup pileup = (AbstractReadBackedPileup)createNewPileup(loc,perSampleElements).getPileupWithoutDeletions(); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getPileupWithoutDeletions(); filteredTracker.addElements(sampleName,pileup.pileupElementTracker); } - return createNewPileup(loc,filteredTracker); + return (RBP)createNewPileup(loc,filteredTracker); } else { @@ -232,7 +232,7 @@ public abstract class AbstractReadBackedPileup perSampleElements = tracker.getElements(sampleName); - AbstractReadBackedPileup pileup = (AbstractReadBackedPileup)createNewPileup(loc,perSampleElements).getOverlappingFragmentFilteredPileup(); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getOverlappingFragmentFilteredPileup(); filteredTracker.addElements(sampleName,pileup.pileupElementTracker); } - return createNewPileup(loc,filteredTracker); + return (RBP)createNewPileup(loc,filteredTracker); } else { Map filteredPileup = new HashMap(); @@ -288,7 +288,7 @@ public abstract class AbstractReadBackedPileup perSampleElements = tracker.getElements(sampleName); - AbstractReadBackedPileup pileup = (AbstractReadBackedPileup)createNewPileup(loc,perSampleElements).getPileupWithoutMappingQualityZeroReads(); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getPileupWithoutMappingQualityZeroReads(); filteredTracker.addElements(sampleName,pileup.pileupElementTracker); } - return createNewPileup(loc,filteredTracker); + return (RBP)createNewPileup(loc,filteredTracker); } else { @@ -324,13 +324,67 @@ public abstract class AbstractReadBackedPileup tracker = (PerSamplePileupElementTracker)pileupElementTracker; + PerSamplePileupElementTracker filteredTracker = new PerSamplePileupElementTracker(); + + for(String sampleName: tracker.getSamples()) { + PileupElementTracker perSampleElements = tracker.getElements(sampleName); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getPositiveStrandPileup(); + filteredTracker.addElements(sampleName,pileup.pileupElementTracker); + } + return (RBP)createNewPileup(loc,filteredTracker); + } + else { + UnifiedPileupElementTracker tracker = (UnifiedPileupElementTracker)pileupElementTracker; + UnifiedPileupElementTracker filteredTracker = new UnifiedPileupElementTracker(); + + for ( PE p : tracker ) { + if ( !p.getRead().getReadNegativeStrandFlag() ) { + filteredTracker.add(p); + } + } + return (RBP)createNewPileup(loc, filteredTracker); + } + } + + /** + * Gets the pileup consisting of only reads on the negative strand. + * @return A read-backed pileup consisting only of reads on the negative strand. + */ + public RBP getNegativeStrandPileup() { + if(pileupElementTracker instanceof PerSamplePileupElementTracker) { + PerSamplePileupElementTracker tracker = (PerSamplePileupElementTracker)pileupElementTracker; + PerSamplePileupElementTracker filteredTracker = new PerSamplePileupElementTracker(); + + for(String sampleName: tracker.getSamples()) { + PileupElementTracker perSampleElements = tracker.getElements(sampleName); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getNegativeStrandPileup(); + filteredTracker.addElements(sampleName,pileup.pileupElementTracker); + } + return (RBP)createNewPileup(loc,filteredTracker); + } + else { + UnifiedPileupElementTracker tracker = (UnifiedPileupElementTracker)pileupElementTracker; + UnifiedPileupElementTracker filteredTracker = new UnifiedPileupElementTracker(); + + for ( PE p : tracker ) { + if ( p.getRead().getReadNegativeStrandFlag() ) { + filteredTracker.add(p); + } + } + return (RBP)createNewPileup(loc, filteredTracker); + } + } + /** Returns subset of this pileup that contains only bases with quality >= minBaseQ, coming from * reads with mapping qualities >= minMapQ. This method allocates and returns a new instance of ReadBackedPileup. * @param minBaseQ @@ -345,11 +399,11 @@ public abstract class AbstractReadBackedPileup perSampleElements = tracker.getElements(sampleName); - AbstractReadBackedPileup pileup = (AbstractReadBackedPileup)createNewPileup(loc,perSampleElements).getBaseAndMappingFilteredPileup(minBaseQ,minMapQ); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements).getBaseAndMappingFilteredPileup(minBaseQ,minMapQ); filteredTracker.addElements(sampleName,pileup.pileupElementTracker); } - return createNewPileup(loc,filteredTracker); + return (RBP)createNewPileup(loc,filteredTracker); } else { UnifiedPileupElementTracker filteredTracker = new UnifiedPileupElementTracker(); @@ -360,7 +414,7 @@ public abstract class AbstractReadBackedPileup getReadGroups() { + Set readGroups = new HashSet(); + for(PileupElement pileupElement: this) + readGroups.add(pileupElement.getRead().getReadGroup().getReadGroupId()); + return readGroups; + } + + /** + * Gets the pileup for a given read group. Horrendously inefficient at this point. + * @param readGroupId Identifier for the read group. + * @return A read-backed pileup containing only the reads in the given read group. + */ + @Override + public RBP getPileupForReadGroup(String readGroupId) { + UnifiedPileupElementTracker filteredTracker = new UnifiedPileupElementTracker(); + for(PileupElement pileupElement: this) + filteredTracker.add((PE)pileupElement); + return (RBP)createNewPileup(loc,filteredTracker); + } + public Collection getSamples() { if(pileupElementTracker instanceof PerSamplePileupElementTracker) { PerSamplePileupElementTracker tracker = (PerSamplePileupElementTracker)pileupElementTracker; @@ -435,14 +514,14 @@ public abstract class AbstractReadBackedPileup pileup = (AbstractReadBackedPileup)createNewPileup(loc,perSampleElements); + AbstractReadBackedPileup pileup = createNewPileup(loc,perSampleElements); filteredTracker.addElements(sampleName,pileup.pileupElementTracker); } current++; } - return createNewPileup(loc,filteredTracker); + return (RBP)createNewPileup(loc,filteredTracker); } else { UnifiedPileupElementTracker tracker = (UnifiedPileupElementTracker)pileupElementTracker; @@ -455,7 +534,7 @@ public abstract class AbstractReadBackedPileup tracker = (PerSamplePileupElementTracker)pileupElementTracker; - return createNewPileup(loc,tracker.getElements(sampleName)); + return (RBP)createNewPileup(loc,tracker.getElements(sampleName)); } else { UnifiedPileupElementTracker filteredTracker = new UnifiedPileupElementTracker(); @@ -478,7 +557,7 @@ public abstract class AbstractReadBackedPileup0 ? createNewPileup(loc,filteredTracker) : null; + return filteredTracker.size()>0 ? (RBP)createNewPileup(loc,filteredTracker) : null; } } diff --git a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileup.java b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileup.java index 5bbfaa836..e85896755 100644 --- a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileup.java +++ b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileup.java @@ -67,6 +67,19 @@ public interface ReadBackedExtendedEventPileup extends ReadBackedPileup { * @return */ public ReadBackedExtendedEventPileup getPileupWithoutMappingQualityZeroReads(); + + + /** + * Gets the pileup consisting of only reads on the positive strand. + * @return A read-backed pileup consisting only of reads on the positive strand. + */ + public ReadBackedExtendedEventPileup getPositiveStrandPileup(); + + /** + * Gets the pileup consisting of only reads on the negative strand. + * @return A read-backed pileup consisting only of reads on the negative strand. + */ + public ReadBackedExtendedEventPileup getNegativeStrandPileup(); /** Returns subset of this pileup that contains only bases with quality >= minBaseQ, coming from * reads with mapping qualities >= minMapQ. This method allocates and returns a new instance of ReadBackedPileup. diff --git a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileupImpl.java b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileupImpl.java index 53457c473..0fc8e21df 100644 --- a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileupImpl.java +++ b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedExtendedEventPileupImpl.java @@ -31,7 +31,7 @@ import java.util.*; import net.sf.samtools.SAMRecord; -public class ReadBackedExtendedEventPileupImpl extends AbstractReadBackedPileup implements ReadBackedExtendedEventPileup { +public class ReadBackedExtendedEventPileupImpl extends AbstractReadBackedPileup implements ReadBackedExtendedEventPileup { private int nInsertions; private int maxDeletionLength; // cached value of the length of the longest deletion observed at the site /** @@ -59,7 +59,7 @@ public class ReadBackedExtendedEventPileupImpl extends AbstractReadBackedPileup< this.nInsertions = nInsertions; } - public ReadBackedExtendedEventPileupImpl(GenomeLoc loc, Map> pileupElementsBySample) { + public ReadBackedExtendedEventPileupImpl(GenomeLoc loc, Map pileupElementsBySample) { super(loc,pileupElementsBySample); } @@ -86,7 +86,7 @@ public class ReadBackedExtendedEventPileupImpl extends AbstractReadBackedPileup< } @Override - protected void addPileupToCumulativeStats(AbstractReadBackedPileup pileup) { + protected void addPileupToCumulativeStats(AbstractReadBackedPileup pileup) { super.addPileupToCumulativeStats(pileup); ReadBackedExtendedEventPileup extendedEventPileup = ((ReadBackedExtendedEventPileup)pileup); this.nInsertions += extendedEventPileup.getNumberOfInsertions(); @@ -94,7 +94,7 @@ public class ReadBackedExtendedEventPileupImpl extends AbstractReadBackedPileup< } @Override - protected ReadBackedExtendedEventPileup createNewPileup(GenomeLoc loc, PileupElementTracker tracker) { + protected ReadBackedExtendedEventPileupImpl createNewPileup(GenomeLoc loc, PileupElementTracker tracker) { return new ReadBackedExtendedEventPileupImpl(loc,tracker); } diff --git a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileup.java b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileup.java index 74fb7b6d9..4fce8f074 100644 --- a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileup.java +++ b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileup.java @@ -66,6 +66,18 @@ public interface ReadBackedPileup extends Iterable { */ public ReadBackedPileup getPileupWithoutMappingQualityZeroReads(); + /** + * Gets the pileup consisting of only reads on the positive strand. + * @return A read-backed pileup consisting only of reads on the positive strand. + */ + public ReadBackedPileup getPositiveStrandPileup(); + + /** + * Gets the pileup consisting of only reads on the negative strand. + * @return A read-backed pileup consisting only of reads on the negative strand. + */ + public ReadBackedPileup getNegativeStrandPileup(); + /** Returns subset of this pileup that contains only bases with quality >= minBaseQ, coming from * reads with mapping qualities >= minMapQ. This method allocates and returns a new instance of ReadBackedPileup. * @param minBaseQ @@ -102,6 +114,19 @@ public interface ReadBackedPileup extends Iterable { */ public boolean hasPileupBeenDownsampled(); + /** + * Gets a collection of all the read groups represented in this pileup. + * @return A collection of all the read group ids represented in this pileup. + */ + public Collection getReadGroups(); + + /** + * Gets all the reads associated with a given read group. + * @param readGroupId Identifier for the read group. + * @return A pileup containing only the reads in the given read group. + */ + public ReadBackedPileup getPileupForReadGroup(String readGroupId); + /** * Gets a collection of all the samples stored in this pileup. * @return Collection of samples in this pileup. diff --git a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileupImpl.java b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileupImpl.java index b1d871bee..03db78320 100644 --- a/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileupImpl.java +++ b/java/src/org/broadinstitute/sting/utils/pileup/ReadBackedPileupImpl.java @@ -29,7 +29,7 @@ import net.sf.samtools.SAMRecord; import java.util.List; import java.util.Map; -public class ReadBackedPileupImpl extends AbstractReadBackedPileup implements ReadBackedPileup { +public class ReadBackedPileupImpl extends AbstractReadBackedPileup implements ReadBackedPileup { public ReadBackedPileupImpl(GenomeLoc loc) { super(loc); @@ -47,7 +47,7 @@ public class ReadBackedPileupImpl extends AbstractReadBackedPileup> pileupElementsBySample) { + public ReadBackedPileupImpl(GenomeLoc loc, Map pileupElementsBySample) { super(loc,pileupElementsBySample); } @@ -65,7 +65,7 @@ public class ReadBackedPileupImpl extends AbstractReadBackedPileup tracker) { + protected ReadBackedPileupImpl createNewPileup(GenomeLoc loc, PileupElementTracker tracker) { return new ReadBackedPileupImpl(loc,tracker); }