zzh
/
gatk-3.8
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

synchronizing latest changes 

git-svn-id: file:///humgen/gsa-scr1/gsa-engineering/svn_contents/trunk@212 348d0f76-0448-11de-a6fe-93d51630548a
This commit is contained in:
asivache 2009-03-27 14:15:44 +00:00
parent 230c1ad161
commit c6d9848d08
3 changed files with 33 additions and 35 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

38
java/src/org/broadinstitute/sting/playground/indels/IndelInspector.java
View File

@ -26,7 +26,7 @@ public class IndelInspector extends CommandLineProgram {
@Option(shortName="V",doc="Verbosity level: SILENT, PILESUMMARY, ALIGNMENTS", optional=true) public String VERBOSITY_LEVEL; @Option(shortName="V",doc="Verbosity level: SILENT, PILESUMMARY, ALIGNMENTS", optional=true) public String VERBOSITY_LEVEL;
@Option(doc="Output file (sam or bam) for non-indel related reads and indel reads that were not improved") public String OUT1; @Option(doc="Output file (sam or bam) for non-indel related reads and indel reads that were not improved") public String OUT1;
@Option(doc="Output file (sam or bam) for improved (realigned) indel related reads") public String OUT2; @Option(doc="Output file (sam or bam) for improved (realigned) indel related reads") public String OUT2;
@Option(doc="[paranoid] Output \"control\" file (sam or bam): all reads picked and processed by this tool will be also saved, unmodified, into this file", optional=true) public String OUTC; @Option(doc="[paranoid] If true, all reads that would be otherwise picked and processed by this tool will be saved, unmodified, into OUT1", optional=true) public boolean CONTROL_RUN;
@Option(doc="Error counting mode: MM - count mismatches only, ERR - count errors (arachne style), MG - count mismatches and gaps as one error each") public String ERR_MODE; @Option(doc="Error counting mode: MM - count mismatches only, ERR - count errors (arachne style), MG - count mismatches and gaps as one error each") public String ERR_MODE;
@Option(doc="Maximum number of errors allowed (see ERR_MODE)") public Integer MAX_ERRS; @Option(doc="Maximum number of errors allowed (see ERR_MODE)") public Integer MAX_ERRS;
// @Option(shortName="R", doc="Reference fasta or fasta.gz file") public File REF_FILE; // @Option(shortName="R", doc="Reference fasta or fasta.gz file") public File REF_FILE;
@ -63,23 +63,26 @@ public class IndelInspector extends CommandLineProgram {
IndelRecordPileCollector col = null; IndelRecordPileCollector col = null;
PassThroughWriter ptWriter = new PassThroughWriter(OUT1,samReader.getFileHeader()); PassThroughWriter ptWriter = new PassThroughWriter(OUT1,samReader.getFileHeader());
PileBuilder pileBuilder = new PileBuilder(OUT2,samReader.getFileHeader(),ptWriter); PileBuilder pileBuilder = null;
if ( ! CONTROL_RUN ) pileBuilder = new PileBuilder(OUT2,samReader.getFileHeader(),ptWriter);
SAMFileWriter controlWriter = null;
if ( OUTC != null ) controlWriter = new SAMFileWriterFactory().makeSAMOrBAMWriter(samReader.getFileHeader(),false,new File(OUTC));
try { try {
col = new IndelRecordPileCollector(ptWriter, pileBuilder, controlWriter ); if ( CONTROL_RUN ) col = new IndelRecordPileCollector(ptWriter, new DiscardingPileReceiver() );
else col = new IndelRecordPileCollector(ptWriter, pileBuilder );
} catch(Exception e) { System.err.println(e.getMessage()); } } catch(Exception e) { System.err.println(e.getMessage()); }
if ( col == null ) return 1; if ( col == null ) return 1;
if ( VERBOSITY_LEVEL == null ) VERBOSITY_LEVEL = new String("SILENT"); col.setControlRun(CONTROL_RUN);
if ( VERBOSITY_LEVEL.toUpperCase().equals("SILENT")) pileBuilder.setVerbosity(pileBuilder.SILENT);
else if ( VERBOSITY_LEVEL.toUpperCase().equals("PILESUMMARY") ) pileBuilder.setVerbosity(pileBuilder.PILESUMMARY); if ( ! CONTROL_RUN ) {
else if ( VERBOSITY_LEVEL.toUpperCase().equals("ALIGNMENTS") ) pileBuilder.setVerbosity(pileBuilder.ALIGNMENTS); if ( VERBOSITY_LEVEL == null ) VERBOSITY_LEVEL = new String("SILENT");
else { if ( VERBOSITY_LEVEL.toUpperCase().equals("SILENT")) pileBuilder.setVerbosity(pileBuilder.SILENT);
System.out.println("Unrecognized VERBOSITY_LEVEL setting."); else if ( VERBOSITY_LEVEL.toUpperCase().equals("PILESUMMARY") ) pileBuilder.setVerbosity(pileBuilder.PILESUMMARY);
return 1; else if ( VERBOSITY_LEVEL.toUpperCase().equals("ALIGNMENTS") ) pileBuilder.setVerbosity(pileBuilder.ALIGNMENTS);
else {
System.out.println("Unrecognized VERBOSITY_LEVEL setting.");
return 1;
}
} }
String cur_contig = null; String cur_contig = null;
@ -97,7 +100,7 @@ public class IndelInspector extends CommandLineProgram {
contig_seq = reference.get(r.getReferenceIndex()); contig_seq = reference.get(r.getReferenceIndex());
String refstr = new String(contig_seq.getBases()); String refstr = new String(contig_seq.getBases());
col.setReferenceSequence(refstr); col.setReferenceSequence(refstr);
pileBuilder.setReferenceSequence(refstr); if (!CONTROL_RUN) pileBuilder.setReferenceSequence(refstr);
System.out.println("loaded contig "+cur_contig+" (index="+r.getReferenceIndex()+"); length="+contig_seq.getBases().length+" tst="+contig_seq.toString()); System.out.println("loaded contig "+cur_contig+" (index="+r.getReferenceIndex()+"); length="+contig_seq.getBases().length+" tst="+contig_seq.toString());
} }
} }
@ -137,13 +140,14 @@ public class IndelInspector extends CommandLineProgram {
} }
pileBuilder.printStats(); if ( ! CONTROL_RUN ) {
pileBuilder.printStats();
pileBuilder.close();
}
System.out.println("done."); System.out.println("done.");
col.printLengthHistograms(); col.printLengthHistograms();
samReader.close(); samReader.close();
pileBuilder.close();
ptWriter.close(); ptWriter.close();
if ( controlWriter != null ) controlWriter.close();
return 0; return 0;
} }

23
java/src/org/broadinstitute/sting/playground/indels/IndelRecordPileCollector.java
View File

@ -89,7 +89,7 @@ public class IndelRecordPileCollector implements RecordReceiver {
private RecordReceiver defaultReceiver; // we will send there records that do not overlap with regions of interest private RecordReceiver defaultReceiver; // we will send there records that do not overlap with regions of interest
private RecordPileReceiver indelPileReceiver; // piles over indel regions will be sent there private RecordPileReceiver indelPileReceiver; // piles over indel regions will be sent there
private SAMFileWriter controlWriter; private boolean controlRun = false;
private String referenceSequence; private String referenceSequence;
@ -99,12 +99,8 @@ public class IndelRecordPileCollector implements RecordReceiver {
//+" Bndries="+mIndelRegionStart +":"+ mIndelRegionStop; //+" Bndries="+mIndelRegionStart +":"+ mIndelRegionStop;
} }
public IndelRecordPileCollector(RecordReceiver rr, RecordPileReceiver rp) throws java.io.IOException {
this(rr,rp,null);
}
public IndelRecordPileCollector(RecordReceiver rr, RecordPileReceiver rp) throws java.io.IOException {
public IndelRecordPileCollector(RecordReceiver rr, RecordPileReceiver rp, SAMFileWriter cw) throws java.io.IOException {
mRecordPile = new LinkedList<SAMRecord>(); mRecordPile = new LinkedList<SAMRecord>();
mAllIndels = new TreeSet<CountedObject<Indel> >( mAllIndels = new TreeSet<CountedObject<Indel> >(
new CountedObjectComparatorAdapter<Indel>(new IntervalComparator())); new CountedObjectComparatorAdapter<Indel>(new IntervalComparator()));
@ -120,7 +116,6 @@ public class IndelRecordPileCollector implements RecordReceiver {
defaultReceiver = rr; defaultReceiver = rr;
indelPileReceiver = rp; indelPileReceiver = rp;
referenceSequence = null; referenceSequence = null;
controlWriter = cw;
setWaitState(); setWaitState();
} }
@ -136,6 +131,8 @@ public class IndelRecordPileCollector implements RecordReceiver {
mState = WAIT_STATE; // got to do this if we were in avoid_region state mState = WAIT_STATE; // got to do this if we were in avoid_region state
} }
public void setControlRun(boolean c) { controlRun = c; }
public void setReferenceSequence(String contig) { public void setReferenceSequence(String contig) {
referenceSequence = contig; referenceSequence = contig;
} }
@ -154,7 +151,6 @@ public class IndelRecordPileCollector implements RecordReceiver {
SAMRecord r = i.next(); SAMRecord r = i.next();
if ( r.getAlignmentEnd() <= pos ) { if ( r.getAlignmentEnd() <= pos ) {
defaultReceiver.receive(r); defaultReceiver.receive(r);
if ( controlWriter != null ) controlWriter.addAlignment(r);
i.remove(); i.remove();
} else break; } else break;
} }
@ -170,7 +166,6 @@ public class IndelRecordPileCollector implements RecordReceiver {
SAMRecord r = i.next(); SAMRecord r = i.next();
if ( r.getAlignmentStart() >= pos ) { if ( r.getAlignmentStart() >= pos ) {
defaultReceiver.receive(r); defaultReceiver.receive(r);
if ( controlWriter != null ) controlWriter.addAlignment(r);
i.remove(); i.remove();
} else break; } else break;
} }
@ -204,7 +199,12 @@ public class IndelRecordPileCollector implements RecordReceiver {
public void receive(final SAMRecord r) throws RuntimeException { public void receive(final SAMRecord r) throws RuntimeException {
if ( r.getReadUnmappedFlag() ) return; // read did not align, nothing to do if ( r.getReadUnmappedFlag() ) return; // read did not align, nothing to do
if ( controlRun ) {
defaultReceiver.receive(r);
return;
}
int currContig = r.getReferenceIndex(); int currContig = r.getReferenceIndex();
int currPos = r.getAlignmentStart(); int currPos = r.getAlignmentStart();
@ -262,7 +262,6 @@ public class IndelRecordPileCollector implements RecordReceiver {
// no indels or avoiding indels in bad region: send all records to defaultReceiver and clear the pile // no indels or avoiding indels in bad region: send all records to defaultReceiver and clear the pile
for ( SAMRecord r : mRecordPile ) { for ( SAMRecord r : mRecordPile ) {
defaultReceiver.receive(r); defaultReceiver.receive(r);
if ( controlWriter != null ) controlWriter.addAlignment(r);
} }
setWaitState(); setWaitState();
return; return;
@ -338,11 +337,9 @@ public class IndelRecordPileCollector implements RecordReceiver {
System.out.print(finalPile.size() + " reads in the pile;") ; System.out.print(finalPile.size() + " reads in the pile;") ;
System.out.println(formatRange(finalTrain)); System.out.println(formatRange(finalTrain));
indelPileReceiver.receive(finalPile); indelPileReceiver.receive(finalPile);
if ( controlWriter != null ) for ( SAMRecord r : finalPile ) controlWriter.addAlignment(r);
} else { } else {
for ( SAMRecord r : finalPile ) { for ( SAMRecord r : finalPile ) {
defaultReceiver.receive(r); defaultReceiver.receive(r);
controlWriter.addAlignment(r);
} }
} }
finalPile.clear(); finalPile.clear();

7
java/src/org/broadinstitute/sting/playground/indels/PileBuilder.java
View File

@ -336,8 +336,8 @@ public class PileBuilder implements RecordPileReceiver {
failedPileReceiver.receive(r); // nothing to do, send failed piles directly for writing failedPileReceiver.receive(r); // nothing to do, send failed piles directly for writing
continue; continue;
} }
// we improved stuff!! let's reset the alignment parameters!
// we improved stuff!! let's reset the alignment parameters!
int cons_offset = ma.getOffsetWrtConsensus(id); // offset of the read 'id' wrt multiple alignment's full consensus seq int cons_offset = ma.getOffsetWrtConsensus(id); // offset of the read 'id' wrt multiple alignment's full consensus seq
// offset of the realigned read r on the reference // offset of the realigned read r on the reference
@ -351,13 +351,10 @@ public class PileBuilder implements RecordPileReceiver {
r.setAttribute("NM",new Integer(AlignmentUtils.numMismatches(r,referenceSequence))); r.setAttribute("NM",new Integer(AlignmentUtils.numMismatches(r,referenceSequence)));
if ( r.getAlignmentStart() == 713655 ) {
System.out.println("!!!----> "+r.format());
System.out.println("!!!----> "+AlignmentUtils.toString(cig) +" --- " +AlignmentUtils.toString(r.getCigar()));
}
// System.out.println("writing " + id); // System.out.println("writing " + id);
samWriter.addAlignment(r); samWriter.addAlignment(r);
} }
} }