Commit 300e75f6 authored by Mathieu Giraud's avatar Mathieu Giraud

core/windowExtractor.cpp, core/segment.cpp: better '-u' debug

Now 'out_unsegmented' causes are output in segment.cpp, displaying affectation for each germline
parent 03204136
......@@ -246,7 +246,7 @@ KmerSegmenter::~KmerSegmenter() {
// delete kaa;
}
KmerMultiSegmenter::KmerMultiSegmenter(Sequence seq, MultiGermline *multigermline)
KmerMultiSegmenter::KmerMultiSegmenter(Sequence seq, MultiGermline *multigermline, ostream *out_unsegmented)
{
int best_score = 0 ;
......@@ -257,6 +257,18 @@ KmerMultiSegmenter::KmerMultiSegmenter(Sequence seq, MultiGermline *multigermlin
KmerSegmenter kseg(seq, germline);
if (out_unsegmented)
{
// Debug, display k-mer affectation and segmentation result for this germline
*out_unsegmented << "#"
<< left << setw(4) << kseg.segmented_germline->code << " "
<< left << setw(20) << segmented_mesg[kseg.getSegmentationStatus()] << " ";
if (kseg.getSegmentationStatus() != UNSEG_TOO_SHORT)
*out_unsegmented << kseg.getKmerAffectAnalyser()->toString();
*out_unsegmented << endl ;
}
if (!best_score)
the_kseg = kseg;
......
......@@ -179,7 +179,7 @@ class KmerMultiSegmenter
* @param seq: An object read from a FASTA/FASTQ file
* @param multigermline: the multigermline
*/
KmerMultiSegmenter(Sequence seq, MultiGermline *multigermline);
KmerMultiSegmenter(Sequence seq, MultiGermline *multigermline, ostream *out_unsegmented);
~KmerMultiSegmenter();
KmerSegmenter the_kseg;
......
......@@ -20,8 +20,12 @@ WindowsStorage *WindowExtractor::extract(OnlineFasta *reads, MultiGermline *mult
while (reads->hasNext()) {
reads->next();
nb_reads++;
if (out_unsegmented) {
*out_unsegmented << reads->getSequence();
}
KmerMultiSegmenter kmseg(reads->getSequence(), multigermline);
KmerMultiSegmenter kmseg(reads->getSequence(), multigermline, out_unsegmented);
KmerSegmenter seg = kmseg.the_kseg ;
int read_length = seg.getSequence().sequence.length();
......@@ -43,18 +47,16 @@ WindowsStorage *WindowExtractor::extract(OnlineFasta *reads, MultiGermline *mult
if (out_segmented) {
*out_segmented << seg ; // KmerSegmenter output (V/N/J)
if (out_unsegmented)
if (out_unsegmented) {
*out_segmented << seg.getKmerAffectAnalyser()->toString() << endl;
*out_unsegmented << "#>" << reads->getSequence().label << " segmented on " << seg.segmented_germline->code << endl << endl;
}
}
nb_reads_germline[seg.system]++;
} else if (out_unsegmented) {
*out_unsegmented << reads->getSequence();
*out_unsegmented << "#" << segmented_mesg[seg.getSegmentationStatus()] << " " << seg.segmented_germline->code << endl;
if (seg.getSegmentationStatus() != UNSEG_TOO_SHORT) {
*out_unsegmented << seg.getKmerAffectAnalyser()->toString() << endl;
}
*out_unsegmented << "#>" << reads->getSequence().label << " not segmented" << endl << endl;
}
}
return windowsStorage;
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment