X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=wiggle.cpp;h=c0cbdccc09ef8880b37ff5a01582d5e2869186f7;hb=543f20967e693ef9150eddc94e7869ed04bad962;hp=19f52b480a00dcd991b30b583fce2ce0177a628b;hpb=946f9a6adb2a82048c8453d44693cd3838d32939;p=rsem.git diff --git a/wiggle.cpp b/wiggle.cpp index 19f52b4..c0cbdcc 100644 --- a/wiggle.cpp +++ b/wiggle.cpp @@ -1,11 +1,13 @@ #include #include #include +#include #include #include "sam/bam.h" #include "sam/sam.h" +#include "utils.h" #include "wiggle.h" void add_bam_record_to_wiggle(const bam1_t *b, Wiggle& wiggle) { @@ -36,29 +38,42 @@ void build_wiggles(const std::string& bam_filename, WiggleProcessor& processor) { samfile_t *bam_in = samopen(bam_filename.c_str(), "rb", NULL); if (bam_in == 0) { fprintf(stderr, "Cannot open %s!\n", bam_filename.c_str()); exit(-1); } - //assert(bam_in != 0); + + bam_header_t *header = bam_in->header; + bool *used = new bool[header->n_targets]; + memset(used, 0, sizeof(bool) * header->n_targets); int cur_tid = -1; //current tid; - int cnt = 0; + HIT_INT_TYPE cnt; bam1_t *b = bam_init1(); Wiggle wiggle; while (samread(bam_in, b) >= 0) { if (b->core.flag & 0x0004) continue; if (b->core.tid != cur_tid) { - if (cur_tid >= 0) processor.process(wiggle); + if (cur_tid >= 0) { used[cur_tid] = true; processor.process(wiggle); } cur_tid = b->core.tid; - wiggle.name = bam_in->header->target_name[cur_tid]; - wiggle.read_depth.assign(bam_in->header->target_len[cur_tid], 0.0); + wiggle.name = header->target_name[cur_tid]; + wiggle.length = header->target_len[cur_tid]; + wiggle.read_depth.assign(wiggle.length, 0.0); } add_bam_record_to_wiggle(b, wiggle); ++cnt; - if (cnt % 1000000 == 0) fprintf(stderr, "%d FIN\n", cnt); + if (cnt % 1000000 == 0) std::cout<< cnt<< std::endl; } - if (cur_tid >= 0) processor.process(wiggle); + if (cur_tid >= 0) { used[cur_tid] = true; processor.process(wiggle); } + + for (int32_t i = 0; i < header->n_targets; i++) + if (!used[i]) { + wiggle.name = header->target_name[i]; + wiggle.length = header->target_len[i]; + wiggle.read_depth.clear(); + processor.process(wiggle); + } samclose(bam_in); bam_destroy1(b); + delete[] used; } UCSCWiggleTrackWriter::UCSCWiggleTrackWriter(const std::string& output_filename, @@ -75,9 +90,11 @@ UCSCWiggleTrackWriter::~UCSCWiggleTrackWriter() { void UCSCWiggleTrackWriter::process(const Wiggle& wiggle) { int sp, ep; + + if (wiggle.read_depth.empty()) return; sp = ep = -1; - for (size_t i = 0; i < wiggle.read_depth.size(); i++) { + for (size_t i = 0; i < wiggle.length; i++) { if (wiggle.read_depth[i] > 0) { ep = i; } @@ -102,9 +119,13 @@ ReadDepthWriter::ReadDepthWriter(std::ostream& stream) } void ReadDepthWriter::process(const Wiggle& wiggle) { + stream_ << wiggle.name << '\t' - << wiggle.read_depth.size() << '\t'; - for (size_t i = 0; i < wiggle.read_depth.size(); ++i) { + << wiggle.length << '\t'; + + if (wiggle.read_depth.empty()) { stream_ << "NA\n"; return; } + + for (size_t i = 0; i < wiggle.length; ++i) { if (i > 0) stream_ << ' '; stream_ << wiggle.read_depth[i]; }