From 237bbdf363c9e42ee24e2fd63106dccf20d9bf2f Mon Sep 17 00:00:00 2001 From: Bo Li Date: Wed, 8 Feb 2012 18:07:26 -0600 Subject: [PATCH] The order of @SQ tags in SAM/BAM files can be arbitrary now --- BamConverter.h | 4 +- BamWriter.h | 27 ++++-------- Buffer.h | 6 +-- EM.cpp | 2 +- README.md | 6 ++- SamParser.h | 57 +++++++++++-------------- Transcripts.h | 48 ++++++++++++++++++++- WHAT_IS_NEW | 8 ++++ convert-sam-for-rsem | 99 ++++++-------------------------------------- makefile | 6 +-- parseIt.cpp | 14 +++---- utils.h | 22 ++++------ 12 files changed, 127 insertions(+), 172 deletions(-) diff --git a/BamConverter.h b/BamConverter.h index e7253ba..af984ac 100644 --- a/BamConverter.h +++ b/BamConverter.h @@ -50,6 +50,8 @@ BamConverter::BamConverter(const char* inpF, const char* outF, const char* chr_l in = samopen(inpF, "rb", NULL); assert(in != 0); + transcripts.buildMappings(in->header->n_targets, in->header->target_name); + bam_header_t *out_header = sam_header_read2(chr_list); refmap.clear(); for (int i = 0; i < out_header->n_targets; i++) { @@ -93,7 +95,7 @@ void BamConverter::process() { // at least one segment is not properly mapped if ((b->core.flag & 0x0004) || (isPaired && (b2->core.flag & 0x0004))) continue; - const Transcript& transcript = transcripts.getTranscriptAt(b->core.tid + 1); + const Transcript& transcript = transcripts.getTranscriptViaEid(b->core.tid + 1); convert(b, transcript); if (isPaired) { diff --git a/BamWriter.h b/BamWriter.h index ff11397..e014f8d 100644 --- a/BamWriter.h +++ b/BamWriter.h @@ -47,23 +47,13 @@ BamWriter::BamWriter(char inpType, const char* inpF, const char* fn_list, const } assert(in != 0); + //build mappings from external sid to internal sid + transcripts.buildMappings(in->header->n_targets, in->header->target_name); + //generate output's header bam_header_t *out_header = bam_header_dwt(in->header); - - if (out_header->n_targets != transcripts.getM()) { - fprintf(stderr, "Number of reference sequences recorded in the header is not correct! The header contains %d sequences while there should be %d sequences\n", out_header->n_targets, transcripts.getM()); - exit(-1); - } - for (int i = 0; i < out_header->n_targets; i++) { - const Transcript& transcript = transcripts.getTranscriptAt(i + 1); - if (out_header->target_name[i] != transcript.getTranscriptID()) { - fprintf(stderr, "Reference sequence %d's name recorded in the header is not correct! \n", i); - fprintf(stderr, "Name in the header: %s\n", out_header->target_name[i]); - fprintf(stderr, "Should be: %s\n", transcript.getTranscriptID().c_str()); - exit(-1); - } - out_header->target_len[i] = transcript.getLength(); // transcript length without poly(A) tail + out_header->target_len[i] = transcripts.getTranscriptViaEid(i + 1).getLength(); // transcript length without poly(A) tail } std::ostringstream strout; @@ -99,7 +89,7 @@ void BamWriter::work(HitWrapper wrapper) { hit = wrapper.getNextHit(); assert(hit != NULL); - assert(b->core.tid + 1 == hit->getSid()); + assert(transcripts.getInternalSid(b->core.tid + 1) == hit->getSid()); convert(b, hit->getConPrb()); if (b->core.qual > 0) samwrite(out, b); // output only when MAPQ > 0 } @@ -135,8 +125,8 @@ void BamWriter::work(HitWrapper wrapper) { hit = wrapper.getNextHit(); assert(hit != NULL); - assert(b->core.tid + 1 == hit->getSid()); - assert(b2->core.tid + 1 == hit->getSid()); + assert(transcripts.getInternalSid(b->core.tid + 1) == hit->getSid()); + assert(transcripts.getInternalSid(b2->core.tid + 1) == hit->getSid()); convert(b, hit->getConPrb()); convert(b2, hit->getConPrb()); @@ -159,8 +149,7 @@ void BamWriter::work(HitWrapper wrapper) { } void BamWriter::convert(bam1_t *b, double prb) { - int sid = b->core.tid + 1; - const Transcript& transcript = transcripts.getTranscriptAt(sid); + const Transcript& transcript = transcripts.getTranscriptViaEid(b->core.tid + 1); int pos = b->core.pos; int readlen = b->core.l_qseq; diff --git a/Buffer.h b/Buffer.h index 1ce28bc..5017796 100644 --- a/Buffer.h +++ b/Buffer.h @@ -62,15 +62,15 @@ private: std::streampos gap2 = std::streampos(nSamples - to) * FLOATSIZE; float *p = NULL; - ftmpOut.seekp(0, std::ios_base::beg); + ftmpOut.seekp(0, std::ios::beg); for (int i = 0; i < cvlen; i++) { p = buffer + i; - ftmpOut.seekp(gap1, std::ios_base::cur); + ftmpOut.seekp(gap1, std::ios::cur); for (int j = fr; j < to; j++) { ftmpOut.write((char*)p, FLOATSIZE); p += cvlen; } - ftmpOut.seekp(gap2, std::ios_base::cur); + ftmpOut.seekp(gap2, std::ios::cur); } cpos = 0; diff --git a/EM.cpp b/EM.cpp index b2493b2..dfdd0e4 100644 --- a/EM.cpp +++ b/EM.cpp @@ -748,7 +748,7 @@ int main(int argc, char* argv[]) { time_t b = time(NULL); - printTimeUsed(a, b); + printTimeUsed(a, b, "EM.cpp"); return 0; } diff --git a/README.md b/README.md index bb66d54..77c0693 100644 --- a/README.md +++ b/README.md @@ -109,8 +109,10 @@ and provide the SAM or BAM file as an argument. When using an alternative aligner, you may also want to provide the '--no-bowtie' option to 'rsem-prepare-reference' so that the Bowtie indices are not built. -Some aligners' (other than Bowtie) output might need to be converted -so that RSEM can use. For conversion, please run +RSEM requires all alignments of the same read group together. For +paired-end reads, RSEM also requires the two mates of any alignment be +adjacent. If the alternative aligner does not satisfy the first +requirement, you can use 'convert-sam-for-rsem' for conversion. Please run convert-sam-for-rsem --help diff --git a/SamParser.h b/SamParser.h index 425593a..c5a5035 100644 --- a/SamParser.h +++ b/SamParser.h @@ -12,9 +12,7 @@ #include "sam/sam.h" #include "utils.h" - -#include "RefSeq.h" -#include "Refs.h" +#include "my_assert.h" #include "SingleRead.h" #include "SingleReadQ.h" @@ -23,9 +21,11 @@ #include "SingleHit.h" #include "PairedEndHit.h" +#include "Transcripts.h" + class SamParser { public: - SamParser(char, const char*, Refs&, const char* = 0); + SamParser(char, const char*, Transcripts&, const char* = 0); ~SamParser(); /** @@ -51,6 +51,7 @@ private: bam_header_t *header; bam1_t *b, *b2; + Transcripts& transcripts; //tag used by aligner static char rtTag[STRLEN]; @@ -79,33 +80,23 @@ private: char SamParser::rtTag[STRLEN] = ""; // default : no tag, thus no Type 2 reads // aux, if not 0, points to the file name of fn_list -SamParser::SamParser(char inpType, const char* inpF, Refs& refs, const char* aux) { +SamParser::SamParser(char inpType, const char* inpF, Transcripts& transcripts, const char* aux) + : transcripts(transcripts) +{ switch(inpType) { case 'b': sam_in = samopen(inpF, "rb", aux); break; case 's': sam_in = samopen(inpF, "r", aux); break; default: assert(false); } - if (sam_in == 0) { fprintf(stderr, "Cannot open %s! It may not exist.\n", inpF); exit(-1); } - header = sam_in->header; - if (header == 0) { fprintf(stderr, "Fail to parse sam header!\n"); exit(-1); } - - // Check if the reference used for aligner is the transcript set RSEM generated - if (refs.getM() != header->n_targets) { - fprintf(stderr, "Number of transcripts does not match! Please align reads against the transcript set and use RSEM generated reference for your aligner!\n"); - exit(-1); - } - for (int i = 0; i < header->n_targets; i++) { - const RefSeq& refseq = refs.getRef(i + 1); - // If update int to long, chance the (int) conversion - if (refseq.getName().compare(header->target_name[i]) != 0 || refseq.getTotLen() != (int)header->target_len[i]) { - fprintf(stderr, "Transcript information does not match! Please align reads against the transcript set and use RSEM generated reference for your aligner!\n"); - exit(-1); - } - } - - b = bam_init1(); - b2 = bam_init1(); + general_assert(sam_in != 0, "Cannot open " + cstrtos(inpF) + "! It may not exist."); + header = sam_in->header; + general_assert(header != 0, "Fail to parse sam header!"); + + transcripts.buildMappings(header->n_targets, header->target_name); + + b = bam_init1(); + b2 = bam_init1(); } SamParser::~SamParser() { @@ -137,10 +128,10 @@ int SamParser::parseNext(SingleRead& read, SingleHit& hit) { if (!check(b)) { fprintf(stderr, "RSEM does not support gapped alignments, sorry!\n"); exit(-1); } if (getDir(b) > 0) { - hit = SingleHit(b->core.tid + 1, b->core.pos); + hit = SingleHit(transcripts.getInternalSid(b->core.tid + 1), b->core.pos); } else { - hit = SingleHit(-(b->core.tid + 1), header->target_len[b->core.tid] - b->core.pos - b->core.l_qseq); + hit = SingleHit(-transcripts.getInternalSid(b->core.tid + 1), header->target_len[b->core.tid] - b->core.pos - b->core.l_qseq); } } @@ -168,10 +159,10 @@ int SamParser::parseNext(SingleReadQ& read, SingleHit& hit) { if (!check(b)) { fprintf(stderr, "RSEM does not support gapped alignments, sorry!\n"); exit(-1); } if (getDir(b) > 0) { - hit = SingleHit(b->core.tid + 1, b->core.pos); + hit = SingleHit(transcripts.getInternalSid(b->core.tid + 1), b->core.pos); } else { - hit = SingleHit(-(b->core.tid + 1), header->target_len[b->core.tid] - b->core.pos - b->core.l_qseq); + hit = SingleHit(-transcripts.getInternalSid(b->core.tid + 1), header->target_len[b->core.tid] - b->core.pos - b->core.l_qseq); } } @@ -220,10 +211,10 @@ int SamParser::parseNext(PairedEndRead& read, PairedEndHit& hit) { } //assert(mp1->core.tid == mp2->core.tid); if (getDir(mp1) > 0) { - hit = PairedEndHit(mp1->core.tid + 1, mp1->core.pos, mp2->core.pos + mp2->core.l_qseq - mp1->core.pos); + hit = PairedEndHit(transcripts.getInternalSid(mp1->core.tid + 1), mp1->core.pos, mp2->core.pos + mp2->core.l_qseq - mp1->core.pos); } else { - hit = PairedEndHit(-(mp1->core.tid + 1), header->target_len[mp1->core.tid] - mp1->core.pos - mp1->core.l_qseq, mp1->core.pos + mp1->core.l_qseq - mp2->core.pos); + hit = PairedEndHit(-transcripts.getInternalSid(mp1->core.tid + 1), header->target_len[mp1->core.tid] - mp1->core.pos - mp1->core.l_qseq, mp1->core.pos + mp1->core.l_qseq - mp2->core.pos); } } @@ -271,10 +262,10 @@ int SamParser::parseNext(PairedEndReadQ& read, PairedEndHit& hit) { } //assert(mp1->core.tid == mp2->core.tid); if (getDir(mp1) > 0) { - hit = PairedEndHit(mp1->core.tid + 1, mp1->core.pos, mp2->core.pos + mp2->core.l_qseq - mp1->core.pos); + hit = PairedEndHit(transcripts.getInternalSid(mp1->core.tid + 1), mp1->core.pos, mp2->core.pos + mp2->core.l_qseq - mp1->core.pos); } else { - hit = PairedEndHit(-(mp1->core.tid + 1), header->target_len[mp1->core.tid] - mp1->core.pos - mp1->core.l_qseq, mp1->core.pos + mp1->core.l_qseq - mp2->core.pos); + hit = PairedEndHit(-transcripts.getInternalSid(mp1->core.tid + 1), header->target_len[mp1->core.tid] - mp1->core.pos - mp1->core.l_qseq, mp1->core.pos + mp1->core.l_qseq - mp2->core.pos); } } diff --git a/Transcripts.h b/Transcripts.h index 2750e0c..2637d35 100644 --- a/Transcripts.h +++ b/Transcripts.h @@ -10,7 +10,10 @@ #include #include #include +#include +#include +#include "my_assert.h" #include "Transcript.h" @@ -20,6 +23,8 @@ public: M = 0; this->type = type; transcripts.clear(); transcripts.push_back(Transcript()); + + e2i.clear(); i2e.clear(); } int getM() { return M; } @@ -42,9 +47,23 @@ public: void readFrom(const char*); void writeTo(const char*); + //Eid: external sid + int getInternalSid(int eid) { + assert(eid > 0 && eid <= M); + return e2i[eid]; + } + + const Transcript& getTranscriptViaEid(int eid) { + return transcripts[getInternalSid(eid)]; + } + + void buildMappings(int, char**); + private: int M, type; // type 0 from genome , 1 standalone transcriptome std::vector transcripts; + + std::vector e2i, i2e; // external sid to internal sid, internal sid to external sid }; void Transcripts::readFrom(const char* inpF) { @@ -55,8 +74,7 @@ void Transcripts::readFrom(const char* inpF) { fin>>M>>type; getline(fin, line); - transcripts.clear(); - transcripts.resize(M + 1); + transcripts.assign(M + 1, Transcript()); for (int i = 1; i <= M; i++) { transcripts[i].read(fin); } @@ -72,4 +90,30 @@ void Transcripts::writeTo(const char* outF) { fout.close(); } +void Transcripts::buildMappings(int n_targets, char** target_name) { + std::map dict; + std::map::iterator iter; + + general_assert(n_targets == M, "Number of transcripts does not match! Please check if the reads are aligned to a transcript set (instead of a genome)!"); + + dict.clear(); + for (int i = 1; i <= M; i++) { + const std::string& tid = transcripts[i].getTranscriptID(); + iter = dict.find(tid); + assert(iter == dict.end()); + dict[tid] = i; + } + + e2i.assign(M + 1, 0); + i2e.assign(M + 1, 0); + for (int i = 0; i < n_targets; i++) { + iter = dict.find(std::string(target_name[i])); + general_assert(iter != dict.end(), "RSEM can not recognize transcript " + cstrtos(target_name[i]) + "!"); + general_assert(iter->second > 0, "Reference sequence name " + cstrtos(target_name[i]) + " is duplicated!"); + e2i[i + 1] = iter->second; + i2e[iter->second] = i + 1; + iter->second = -1; + } +} + #endif /* TRANSCRIPTS_H_ */ diff --git a/WHAT_IS_NEW b/WHAT_IS_NEW index 08f9786..9ec79a7 100644 --- a/WHAT_IS_NEW +++ b/WHAT_IS_NEW @@ -1,3 +1,11 @@ +RSEM v1.1.17 + +- Fixed a bug related to parallezation of credibility intervals calculation +- Added --no-bam-output option to rsem-calculate-expression +- The order of @SQ tags in SAM/BAM files can be arbitrary now + +-------------------------------------------------------------------------------------------- + RSEM v1.1.16 - Added --time option to show time consumed by each phase diff --git a/convert-sam-for-rsem b/convert-sam-for-rsem index 68a42d7..348354b 100755 --- a/convert-sam-for-rsem +++ b/convert-sam-for-rsem @@ -4,9 +4,7 @@ use Getopt::Long; use Pod::Usage; use strict; -my $standard_output = "&STDOUT"; - -my $out_file = $standard_output; +my $out_file = ""; my @tmp_dirs = (); my $help = 0; @@ -16,86 +14,19 @@ GetOptions("o=s" => \$out_file, pod2usage(-verbose => 2) if ($help == 1); -pod2usage(-msg => "Invalid number of arguments!", -exitval => 2, -verbose => 2) if (scalar(@ARGV) != 2); +pod2usage(-msg => "Invalid number of arguments!", -exitval => 2, -verbose => 2) if (scalar(@ARGV) != 1); my $command; -my (@fields, @header) = (); -my $M; - -# Load fields -my @lines = (); -my $type; - -open(INPUT, "$ARGV[0].ti"); -@lines = ; -chomp(@lines); -close(INPUT); - -@fields = (); -($M, $type) = split(/ /, $lines[0]); -for (my $i = 0; $i < $M; $i++) { - push(@fields, "SN:$lines[$i * 6 + 1]"); -} - - -# Reorder header -my $line; - -open(INPUT, $ARGV[1]); -@header = (); -while (($line = ) && substr($line, 0, 1) eq '@') { - chomp($line); - push(@header, $line); -} -close(INPUT); - -my $n = scalar(@header); -if ($n > 0) { - my %hash = (); - my @ktable = (); - - my $tid = 0; - - for (my $i = 0; $i < $n; $i++) { - my @arr = split(/\t/, $header[$i]); - if ($arr[0] ne "\@SQ") { push(@ktable, ""); next; } - my $hasSN = 0; - foreach my $key (@arr) { - if (substr($key, 0, 3) eq "SN:") { - $hash{$key} = $i; - $hasSN = 1; - last; - } - } - if (!$hasSN) { print STDERR "\"$header[$i]\" does not have a SN tag!\n"; exit(-1); } - push(@ktable, $fields[$tid++]); - } - - if ($tid != $M) { print STDERR "Number of \@SQ lines is not correct!\n"; exit(-1); } - open(OUTPUT, ">$out_file"); - for (my $i = 0; $i < $n; $i++) { - if ($ktable[$i] eq "") { print OUTPUT $header[$i]; } - else { print OUTPUT $header[$hash{$ktable[$i]}]; } - print OUTPUT "\n"; - } - close(OUTPUT); -} - - -# extract alignment section -$command = "grep ^[^@] $ARGV[1] > $ARGV[1].__temp"; +# grep header section +$command = "grep ^@ $ARGV[0]"; +if ($out_file ne "") { $command .= " > $out_file"; } &runCommand($command); -# sort and output the alignment section -$command = "sort -k 1,1 -s"; +# sort alignment section +$command = "grep ^[^@] $ARGV[0] | sort -k 1,1 -s"; if (scalar(@tmp_dirs) > 0) { $" = " -T "; $command .= " -T @tmp_dirs"; } -$command .= " $ARGV[1].__temp"; -if ($out_file ne $standard_output) { $command .= " >> $out_file"; } -&runCommand($command); - -# delete temporary files -$command = "rm -f $ARGV[1].__temp"; +if ($out_file ne "") { $command .= " >> $out_file"; } &runCommand($command); # finish @@ -125,7 +56,7 @@ convert-sam-for-rsem =over - convert-sam-for-rsem [options] reference_name input_sam + convert-sam-for-rsem [options] input_sam =back @@ -133,10 +64,6 @@ convert-sam-for-rsem =over -=item B - -The name of the reference used. This should be the same name used by 'rsem-prepare-reference'. - =item B The SAM file (*.sam) generated by user's aligner. If the aligner produces a BAM file, please use samtools to convert it to a SAM file (with header information). @@ -165,20 +92,20 @@ Show help information. This program converts the SAM file generated by user's aligner into a SAM file which RSEM can process. However, users should make sure their aligners use 'reference_name.idx.fa' generated by 'rsem-prepare-reference' as their references. In addition, their aligners should output header information and make two mates of the same alignment adjacent to each other for paired-end data. This program will output the converted file into standard output by default for the purpose of piping. By setting '-o' option, users can make the converted file written into disk. -Note: You do not need to run this script if Bowtie (not Bowtie 2) is used, or the order of @SQ lines is the same as 'reference_name.idx.fa' and the alignment lines of a same read group together and the mates of the same alignment are adjacent each other for paired-end reads. +Note: You do not need to run this script if Bowtie (not Bowtie 2) is used, or the alignment lines of a same read group together and the mates of the same alignment are adjacent each other for paired-end reads. Note: This program can only recognize SAM files. See ARGUMENTS section. =head1 EXAMPLES -Suppose reference_name and input_sam are set to '/ref/mouse_125' and 'input.sam'. +Suppose input_sam is set to 'input.sam'. 1) Output to standard output and gzip the output to 'input_for_rsem.sam.gz': - convert-sam-for-rsem /ref/mouse_125 input.sam | gzip > input_for_rsem.sam.gz + convert-sam-for-rsem input.sam | gzip > input_for_rsem.sam.gz 2) Output to 'input_for_rsem.sam' directly: - convert-sam-for-rsem /ref/mouse_125 input.sam -o input_for_rsem.sam + convert-sam-for-rsem input.sam -o input_for_rsem.sam =cut diff --git a/makefile b/makefile index 45ab35d..cca54be 100644 --- a/makefile +++ b/makefile @@ -11,7 +11,7 @@ build-sam : Transcript.h : utils.h -Transcripts.h : Transcript.h +Transcripts.h : my_assert.h Transcript.h rsem-extract-reference-transcripts : utils.h GTFItem.h Transcript.h Transcripts.h extractRef.cpp $(CC) -Wall -O3 extractRef.cpp -o rsem-extract-reference-transcripts @@ -47,13 +47,13 @@ PairedEndHit.h : SingleHit.h HitContainer.h : GroupInfo.h -SamParser.h : sam/sam.h sam/bam.h utils.h SingleRead.h SingleReadQ.h PairedEndRead.h PairedEndReadQ.h SingleHit.h PairedEndHit.h RefSeq.h Refs.h +SamParser.h : sam/sam.h sam/bam.h utils.h my_assert.h SingleRead.h SingleReadQ.h PairedEndRead.h PairedEndReadQ.h SingleHit.h PairedEndHit.h Transcripts.h rsem-parse-alignments : parseIt.o sam/libbam.a $(CC) -o rsem-parse-alignments parseIt.o sam/libbam.a -lz -parseIt.o : utils.h GroupInfo.h Read.h SingleRead.h SingleReadQ.h PairedEndRead.h PairedEndReadQ.h SingleHit.h PairedEndHit.h HitContainer.h SamParser.h RefSeq.h Refs.h sam/sam.h sam/bam.h parseIt.cpp +parseIt.o : utils.h GroupInfo.h Read.h SingleRead.h SingleReadQ.h PairedEndRead.h PairedEndReadQ.h SingleHit.h PairedEndHit.h HitContainer.h SamParser.h Transcripts.h sam/sam.h sam/bam.h parseIt.cpp $(CC) $(COFLAGS) parseIt.cpp diff --git a/parseIt.cpp b/parseIt.cpp index e05585d..8620432 100644 --- a/parseIt.cpp +++ b/parseIt.cpp @@ -13,9 +13,7 @@ #include "utils.h" #include "GroupInfo.h" - -#include "RefSeq.h" -#include "Refs.h" +#include "Transcripts.h" #include "SingleRead.h" #include "SingleReadQ.h" @@ -34,12 +32,12 @@ int N[3]; // note, N = N0 + N1 + N2 , but may not be equal to the total number o int nHits; // # of hits int nUnique, nMulti, nIsoMulti; char fn_list[STRLEN]; -char refF[STRLEN], groupF[STRLEN]; +char groupF[STRLEN], tiF[STRLEN]; char imdName[STRLEN]; char datF[STRLEN], cntF[STRLEN]; -Refs refs; GroupInfo gi; +Transcripts transcripts; SamParser *parser; ofstream hit_out; @@ -55,7 +53,7 @@ void init(const char* imdName, char alignFType, const char* alignF) { char* aux = 0; if (strcmp(fn_list, "")) aux = fn_list; - parser = new SamParser(alignFType, alignF, refs, aux); + parser = new SamParser(alignFType, alignF, transcripts, aux); memset(cat, 0, sizeof(cat)); memset(readOutFs, 0, sizeof(readOutFs)); @@ -192,10 +190,10 @@ int main(int argc, char* argv[]) { verbose = !quiet; - sprintf(refF, "%s.seq", argv[1]); - refs.loadRefs(refF, 1); sprintf(groupF, "%s.grp", argv[1]); gi.load(groupF); + sprintf(tiF, "%s.ti", argv[1]); + transcripts.readFrom(tiF); sprintf(imdName, "%s.temp/%s", argv[2], argv[3]); sprintf(datF, "%s.dat", imdName); diff --git a/utils.h b/utils.h index 278e95e..e3f8bcb 100644 --- a/utils.h +++ b/utils.h @@ -120,20 +120,6 @@ inline std::string cleanStr(const std::string& str) { return (fr <= to ? str.substr(fr, to - fr + 1) : ""); } -void printTimeUsed(const time_t& a, const time_t& b, const char* filename = "") { - int hh = (b - a) / 3600; - int mm = (b - a) % 3600 / 60; - int ss = (b - a) % 60; - - printf("Time Used : %d h %02d m %02d s\n", hh, mm, ss); - - if (strcmp(filename, "")) { - FILE *fo = fopen(filename, "w"); - fprintf(fo, "Time Used : %d h %02d m %02d s\n", hh, mm, ss); - fclose(fo); - } -} - void genReadFileNames(const char* readFN, int tagType, int read_type, int& s, char readFs[][STRLEN]){ const char tags[3][STRLEN] = {"un", "alignable", "max"}; char suffix[STRLEN]; @@ -156,4 +142,12 @@ void genReadFileNames(const char* readFN, int tagType, int read_type, int& s, ch } } +void printTimeUsed(const time_t& a, const time_t& b, const char* program_name) { + int hh = (b - a) / 3600; + int mm = (b - a) % 3600 / 60; + int ss = (b - a) % 60; + + printf("Time Used for %s : %d h %02d m %02d s\n", program_name, hh, mm, ss); +} + #endif -- 2.39.2