Skip to content

Commit

Permalink
v1.1.0. Added --command option to CLI. Added sorting by POS. Fixed ba…
Browse files Browse the repository at this point in the history 
…ckspace-key segfault. Better drawing of mods when zoomed out. Other small bug fixes
  • Loading branch information
@kcleal
kcleal committed Sep 11, 2024
1 parent 348b96a commit 041efe4
Show file tree
Hide file tree
Showing 21 changed files with 653 additions and 193 deletions.
2 changes: 1 addition & 1 deletion .github/workflows/main.yml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@ on:
pull_request:

env:
version: 1.0.3
version: 1.1.0

jobs:
mingw:
Expand Down
1 change: 1 addition & 0 deletions .gitignore
View file
Original file line number Diff line number Diff line change
Expand Up @@ -41,3 +41,4 @@ libgw.so
/src/plot_commands.o
/src/ideogram.o
/.gw_session.ini
/include/glad.o
6 changes: 3 additions & 3 deletions .gw.ini
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@ saccer3=https://github.com/kcleal/ref_genomes/releases/download/v0.1.0/sacCer3.f
[general]
theme=dark
dimensions=1366x800
indel_length=10
indel_length=15
ylim=50
coverage=true
log2_cov=false
Expand All @@ -40,14 +40,14 @@ font=Menlo
font_size=14
sv_arcs=true
mods=false
mods_qual_threshold=50
mods_qual_threshold=150
session_file=

[view_thresholds]
soft_clip=20000
small_indel=100000
snp=500000
mod= 1000000
mod=1000000
edge_highlights=1000000
variant_distance=100000
low_memory=1500000
Expand Down
2 changes: 1 addition & 1 deletion deps/gw.desktop
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
[Desktop Entry]
Encoding=UTF-8
Version=1.0.3
Version=1.1.0
Type=Application
Terminal=true
Exec=bash -c "/usr/bin/gw"
Expand Down
129 changes: 92 additions & 37 deletions src/drawing.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -831,53 +831,107 @@ namespace Drawing {
void drawMods(SkCanvas *canvas, SkRect &rect, const Themes::BaseTheme &theme, const Utils::Region *region,
const Segs::Align &align,
float width, float xScaling, float xOffset, float mmPosOffset, float yScaledOffset,
float pH, int l_qseq, float monitorScale) { //SkPaint& fc5mc, SkPaint& fc5hmc, SkPaint& fcOther) {
float pH, int l_qseq, float monitorScale, bool as_dots) { //SkPaint& fc5mc, SkPaint& fc5hmc, SkPaint& fcOther) {
if (align.any_mods.empty()) {
return;
}

float precalculated_xOffset_mmPosOffset = xOffset + mmPosOffset + (0.5 * xScaling) - 2;
float precalculated_xOffset_mmPosOffset, h, top, middle, bottom, w;
if (as_dots) {
precalculated_xOffset_mmPosOffset = xOffset + mmPosOffset + (0.5 * xScaling);// - monitorScale;
top = yScaledOffset + (pH * 0.3333);
middle = yScaledOffset + pH - (pH * 0.5);
bottom = yScaledOffset + pH - (pH * 0.3333);
} else {
precalculated_xOffset_mmPosOffset = xOffset + mmPosOffset;
h = pH * 0.25;
top = yScaledOffset + h;
middle = top + h;
bottom = top + h;
w = std::fmax(monitorScale, xScaling);
}

auto mod_it = align.any_mods.begin();
auto mod_end = align.any_mods.end();

float top = yScaledOffset + (pH * 0.3333);
float middle = yScaledOffset + pH - (pH * 0.5);
float bottom = yScaledOffset + pH - (pH * 0.3333);
for (const auto& blk : align.blocks) {
if ((int)blk.end < region->start) {
continue;
} else if ((int)blk.start >= region->end) {
return;
}
int idx_start = blk.seq_index;
int idx_end = blk.seq_index + (blk.end - blk.start);
while (mod_it != mod_end && mod_it->index < idx_start) {
++mod_it;
}
while (mod_it != mod_end && mod_it->index < idx_end) {
float x = ((((int)blk.start + (int)mod_it->index - idx_start) - region->start) * xScaling) + precalculated_xOffset_mmPosOffset;
if (x < 0) {
auto mc_paint = &theme.ModPaints[0];
auto hmc_paint = &theme.ModPaints[1];
auto other_paint = &theme.ModPaints[2];

if (as_dots) {
for (const auto& blk : align.blocks) {
if ((int)blk.end < region->start) {
continue;
} else if ((int)blk.start >= region->end) {
return;
}
int idx_start = blk.seq_index;
int idx_end = blk.seq_index + (blk.end - blk.start);
while (mod_it != mod_end && mod_it->index < idx_start) {
++mod_it;
}
while (mod_it != mod_end && mod_it->index < idx_end) {
float x = ((((int)blk.start + (int)mod_it->index - idx_start) - region->start) * xScaling) + precalculated_xOffset_mmPosOffset;
if (x < 0) {
++mod_it;
continue;
}
int n_mods = mod_it->n_mods;
for (size_t j=0; j < (size_t)n_mods; ++j) {
switch (mod_it->mods[j]) {
case 'm': // 5mC
canvas->drawPoint(x, top, (*mc_paint)[ mod_it->quals[j] % 4 ]);
break;
case 'h': // 5hmC
canvas->drawPoint(x, bottom, (*hmc_paint)[ mod_it->quals[j] % 4 ]);
break;
default:
canvas->drawPoint(x, middle, (*other_paint)[ mod_it->quals[j] % 4 ]);
break;
}
}
++mod_it;
}
}
} else {
for (const auto& blk : align.blocks) {
if ((int)blk.end < region->start) {
continue;
} else if ((int)blk.start >= region->end) {
return;
}
int n_mods = mod_it->n_mods;
for (size_t j=0; j < (size_t)n_mods; ++j) {
switch (mod_it->mods[j]) {
case 'm': // 5mC
canvas->drawPoint(x, (n_mods == 1) ? middle : top, theme.ModPaints[0][ mod_it->quals[j] % 4 ]);
break;
case 'h': // 5hmC
canvas->drawPoint(x, (n_mods == 1) ? middle : bottom, theme.ModPaints[1][ mod_it->quals[j] % 4 ]);
break;
default:
canvas->drawPoint(x, middle, theme.ModPaints[2][ mod_it->quals[j] % 4 ]);
break;
int idx_start = blk.seq_index;
int idx_end = blk.seq_index + (blk.end - blk.start);
while (mod_it != mod_end && mod_it->index < idx_start) {
++mod_it;
}
while (mod_it != mod_end && mod_it->index < idx_end) {
float x = ((((int)blk.start + (int)mod_it->index - idx_start) - region->start) * xScaling) + precalculated_xOffset_mmPosOffset;
if (x < 0) {
++mod_it;
continue;
}
int n_mods = mod_it->n_mods;
for (size_t j=0; j < (size_t)n_mods; ++j) {
switch (mod_it->mods[j]) {
case 'm': // 5mC
rect.setXYWH(x, top, w, h);
canvas->drawRect(rect, (*mc_paint)[ mod_it->quals[j] % 4 ]);
break;
case 'h': // 5hmC
rect.setXYWH(x, bottom, w, h);
canvas->drawRect(rect, (*hmc_paint)[ mod_it->quals[j] % 4 ]);
break;
default:
rect.setXYWH(x, middle, w, h);
canvas->drawRect(rect, (*other_paint)[ mod_it->quals[j] % 4 ]);
break;
}
}
++mod_it;
}
++mod_it;
}
}

}

void drawCollection(const Themes::IniOptions &opts, Segs::ReadCollection &cl,
Expand Down Expand Up @@ -1077,10 +1131,6 @@ namespace Drawing {
drawMismatchesNoMD(canvas, rect, theme, cl.region, a, (float) width, xScaling, xOffset, mmPosOffset,
yScaledOffset, pH, l_qseq, mm_vector, cl.collection_processed);
}
if (opts.parse_mods && regionLen <= opts.mod_threshold) {
drawMods(canvas, rect, theme, cl.region, a, (float) width, xScaling, xOffset, mmPosOffset,
yScaledOffset, pH, l_qseq, monitorScale); //, theme.fc5mc, theme.fc5hmc, theme.fcOther);
}

// add insertions
if (!a.any_ins.empty()) {
Expand Down Expand Up @@ -1144,6 +1194,11 @@ namespace Drawing {
}
}
}
// Add modifications
if (opts.parse_mods && regionLen <= opts.mod_threshold) {
drawMods(canvas, rect, theme, cl.region, a, (float) width, xScaling, xOffset, mmPosOffset,
yScaledOffset, pH, l_qseq, monitorScale, regionLen <= 2000);
}
}

// draw text deletions + insertions
Expand Down
54 changes: 40 additions & 14 deletions src/hts_funcs.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -157,13 +157,22 @@ namespace HGW {
}
}

// char getRefBaseAtPos(const int sortReadsBy, Utils::Region *region) {
// if (sortReadsBy != 3) {
// return '\0';
// }
// region->setRefBaseAtPos(); // set lazily
// return region->refBaseAtPos;
// }

void collectReadsAndCoverage(Segs::ReadCollection &col, htsFile *b, sam_hdr_t *hdr_ptr,
hts_idx_t *index, int threads, Utils::Region *region,
bool coverage, std::vector<Parse::Parser> &filters, BS::thread_pool &pool,
const int parse_mods_threshold, int sortReadsBy) {
const int parse_mods_threshold) {

bam1_t *src;
hts_itr_t *iter_q;

int tid = sam_hdr_name2tid(hdr_ptr, region->chrom.c_str());
std::vector<Segs::Align>& readQueue = col.readQueue;
if (region->end - region->start < 1000000) {
Expand Down Expand Up @@ -193,7 +202,7 @@ namespace HGW {
readQueue.pop_back();
}

Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold, sortReadsBy == 2);
Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold);

if (!filters.empty()) {
applyFilters(filters, readQueue, hdr_ptr, col.bamIdx, col.regionIdx);
Expand Down Expand Up @@ -357,6 +366,12 @@ namespace HGW {
hts_itr_t *iter_q;
int tid = sam_hdr_name2tid(hdr_ptr, region->chrom.c_str());
std::vector<Segs::Align>& readQueue = col.readQueue;

if (col.levelsStart.empty()) {
col.levelsStart.resize(opts.ylim + col.vScroll, 1215752191);
col.levelsEnd.resize(opts.ylim + col.vScroll, 0);
}

if (!readQueue.empty()) {
for (auto &item: readQueue) {
bam_destroy1(item.delegate);
Expand Down Expand Up @@ -387,7 +402,8 @@ namespace HGW {
if (j < BATCH) {
continue;
}
Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold, false);
// No specialised sorting here
Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold);
if (filter) {
applyFilters_noDelete(filters, readQueue, hdr_ptr, col.bamIdx, col.regionIdx);
}
Expand All @@ -399,7 +415,10 @@ namespace HGW {
}
}
}
Segs::findY(col, readQueue, opts.link_op, opts, false, 0);
assert (opts.link_op == 0);

Segs::findYNoSortForward(readQueue, col.levelsStart, col.levelsEnd, col.vScroll);

Drawing::drawCollection(opts, col, canvas, trackY, yScaling, fonts, opts.link_op, refSpace, pointSlop, textDrop, pH, monitorScale);

for (int i=0; i < BATCH; ++ i) {
Expand All @@ -411,7 +430,7 @@ namespace HGW {
if (j < BATCH) {
readQueue.erase(readQueue.begin() + j, readQueue.end());
if (!readQueue.empty()) {
Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold, false);
Segs::init_parallel(readQueue, threads, pool, parse_mods_threshold);
if (!filters.empty()) {
applyFilters_noDelete(filters, readQueue, hdr_ptr, col.bamIdx, col.regionIdx);
}
Expand All @@ -423,7 +442,7 @@ namespace HGW {
}
}
}
Segs::findY(col, readQueue, opts.link_op, opts, false, 0);
Segs::findYNoSortForward(readQueue, col.levelsStart, col.levelsEnd, col.vScroll);
Drawing::drawCollection(opts, col, canvas, trackY, yScaling, fonts, opts.link_op, refSpace, pointSlop, textDrop, pH, monitorScale);
for (int i=0; i < BATCH; ++ i) {
Segs::align_clear(&readQueue[i]);
Expand All @@ -442,6 +461,10 @@ namespace HGW {
bam1_t *src;
hts_itr_t *iter_q;
int tid = sam_hdr_name2tid(hdr_ptr, region->chrom.c_str());
if (col.levelsStart.empty()) {
col.levelsStart.resize(opts.ylim + col.vScroll, 1215752191);
col.levelsEnd.resize(opts.ylim + col.vScroll, 0);
}
std::vector<Segs::Align>& readQueue = col.readQueue;
if (!readQueue.empty()) {
for (auto &item: readQueue) {
Expand All @@ -463,7 +486,7 @@ namespace HGW {
if (src->core.flag & 4 || src->core.n_cigar == 0) {
continue;
}
Segs::align_init(&readQueue.back(), parse_mods_threshold, false);
Segs::align_init(&readQueue.back(), parse_mods_threshold);
if (filter) {
applyFilters_noDelete(filters, readQueue, hdr_ptr, col.bamIdx, col.regionIdx);
if (readQueue.back().y == -2) {
Expand All @@ -475,7 +498,7 @@ namespace HGW {
int l_arr = (int)col.covArr.size() - 1;
Segs::addToCovArray(col.covArr, readQueue.back(), region->start, region->end, l_arr);
}
Segs::findY(col, readQueue, opts.link_op, opts, false, 0);
Segs::alignFindYForward(readQueue.back(), col.levelsStart, col.levelsEnd, col.vScroll);
Drawing::drawCollection(opts, col, canvas, trackY, yScaling, fonts, opts.link_op, refSpace, pointSlop, textDrop, pH, monitorScale);
Segs::align_clear(&readQueue.back());
}
Expand Down Expand Up @@ -539,15 +562,16 @@ namespace HGW {
*samMaxY = (maxY > *samMaxY || opts.tlen_yscale) ? maxY : *samMaxY;
}

void refreshLinked(std::vector<Segs::ReadCollection> &collections, Themes::IniOptions &opts, int *samMaxY, int sortReadsBy) {
void refreshLinked(std::vector<Segs::ReadCollection> &collections, std::vector<Utils::Region> &regions, Themes::IniOptions &opts, int *samMaxY) {
for (auto &cl : collections) {
refreshLinkedCollection(cl, opts, samMaxY, sortReadsBy);
assert (cl.regionIdx < regions.size());
refreshLinkedCollection(cl, opts, samMaxY, regions[cl.regionIdx].getSortOption());
}
}

void appendReadsAndCoverage(Segs::ReadCollection &col, htsFile *b, sam_hdr_t *hdr_ptr,
hts_idx_t *index, Themes::IniOptions &opts, bool coverage, bool left, int *samMaxY,
std::vector<Parse::Parser> &filters, BS::thread_pool &pool, int sortReadsBy) {
std::vector<Parse::Parser> &filters, BS::thread_pool &pool, Utils::Region &reg) {
bam1_t *src;
hts_itr_t *iter_q;
std::vector<Segs::Align>& readQueue = col.readQueue;
Expand All @@ -559,6 +583,7 @@ namespace HGW {
}
int lastPos;
const int parse_mods_threshold = (opts.parse_mods) ? 50 : 0;

if (!readQueue.empty()) {
if (left) {
lastPos = readQueue.front().pos; // + 1;
Expand Down Expand Up @@ -699,14 +724,15 @@ namespace HGW {
}

if (!newReads.empty()) {
Segs::init_parallel(newReads, opts.threads, pool, parse_mods_threshold, sortReadsBy == 2);
Segs::init_parallel(newReads, opts.threads, pool, parse_mods_threshold);
if (!filters.empty()) {
applyFilters(filters, newReads, hdr_ptr, col.bamIdx, col.regionIdx);
}

bool findYall = false;
int sort_state = Segs::getSortCodes(newReads, opts.threads, pool, region);
if (col.vScroll == 0 && opts.link_op == 0) { // only new reads need findY, otherwise, reset all below
int maxY = Segs::findY(col, newReads, opts.link_op, opts, left, sortReadsBy);
int maxY = Segs::findY(col, newReads, opts.link_op, opts, left, sort_state);
if (maxY > *samMaxY) {
*samMaxY = maxY;
}
Expand All @@ -721,7 +747,7 @@ namespace HGW {
col.readQueue = newReads;
}
if (findYall) {
refreshLinkedCollection(col, opts, samMaxY, sortReadsBy);
refreshLinkedCollection(col, opts, samMaxY, sort_state);
}
if (opts.link_op > 0) {
// move of data will invalidate some pointers, so reset
Expand Down
Loading

0 comments on commit 041efe4

Please sign in to comment.