Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Roundup16 #6

Open
wants to merge 3 commits into
base: dev
Choose a base branch
from
Open

Roundup16 #6

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
6fa0790
Add debug prints
joydddd Oct 25, 2022
24c0ad9
Debug print, sse
joydddd Oct 25, 2022
4d9fa93
Locate and <fix> extd2_sse initialization. Is this an error or intent…
joydddd Oct 25, 2022
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
43 changes: 26 additions & 17 deletions ksw2_extd2.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -264,11 +264,13 @@ int ksw_update_diag (
y2[t-1] = y2_new;
H[t-1] = H_new;
#ifdef DEBUG
// if (t >= t_st && t <= t_en+1)
// fprintf(align_debug_file,
// "t %d, s %d, u %d, v %d, x %d, y %d, x2 %d, y2
// %d H %d\n", t - 1, sc[t-1], u_new, v_new, x_new,
// y_new, x2_new, y2_new, H_new);
if (r == 129){
// if (t >= t_st && t <= t_en+1)
// fprintf(align_debug_file,
// "t %d, s %d, u %d, v %d, x %d, y %d, x2 %d, y2 %d H %d\n",
// t - 1, sc[t-1], u_new, v_new, x_new,
// y_new, x2_new, y2_new, H_new);
}
#endif
}

Expand Down Expand Up @@ -333,7 +335,7 @@ int ksw_update_diag (
rmax[t] = r;
}
#ifdef DEBUG
// fprintf(align_debug_file, "%d|%d ", prev_H, H_new);
// fprintf(align_debug_file, "p %x ", p[t]);
#endif
}
// update KZ matrix
Expand Down Expand Up @@ -577,22 +579,29 @@ void ksw_extd2_cpp(
// DEBUG: debug output
#ifdef DEBUG
fprintf(align_debug_file, "#%d (st=%d en=%d) ", r, st, en);
// for (int t = get_t(st0, r, n_col); t <= get_t(en0, r, n_col); ++t) {
for (int t = t_st; t <= t_en; ++t) {
// if (!align_debug_file) {
// align_debug_file = fopen("debug/test_sample_debug.output", "w+");
// }
fprintf(align_debug_file, "%d ", H[t]);
// fprintf(align_debug_file, "(%d,%d,%d|%d,%d,%d,%d,%d,%d,%d,%d,%d,0x%x)\n",
// r, t, get_i(t, r, n_col), ((int8_t *)u)[t], ((int8_t *)v)[t],
// fprintf(align_debug_file, "%d ", H[t]);
fprintf(align_debug_file, "[%d %d v-1 %d s %d]%d ", v[t], u[t], v[t-1], sc[t], H[t]);
// fprintf(align_debug_file,
// "(%d,%d,%d|%d,%d,%d,%d,%d,%d,%d,%d,%d,0x%x)\n",
// r, t, get_i(t, r, n_col), ((int8_t *)u)[t], ((int8_t
// *)v)[t],
// ((int8_t *)x)[t], ((int8_t *)y)[t], ((int8_t *)x2)[t],
// ((int8_t *)y2)[t], ((int32_t *)H)[t], ((int32_t *)Hmax)[t],
// ((int *)rmax)[t], p[r*min_n_col - t_st + t]); // for debugging
// ((int8_t *)y2)[t], ((int32_t *)H)[t], ((int32_t
// *)Hmax)[t],
// ((int *)rmax)[t], p[r*min_n_col - t_st + t]); // for
// debugging
// if (!align_score_file) {
// align_score_file = fopen("debug/test_sample_score.output", "w+");
// fprintf(align_score_file, "(r, t | u, v, x, y)\n");
// align_score_file = fopen("debug/test_sample_score.output",
// "w+"); fprintf(align_score_file, "(r, t | u, v, x, y)\n");
// }
// fprintf(align_score_file, "(%d,%d|%d,%d,%d,%d)",
// r, t-t_st+st, ((int8_t*)u)[t], ((int8_t*)v)[t], ((int8_t*)x)[t],
// fprintf(align_score_file, "(%d,%d|%d,%d,%d,%d)",
// r, t-t_st+st, ((int8_t*)u)[t], ((int8_t*)v)[t],
// ((int8_t*)x)[t],
// ((int8_t*)y)[t]); // for debugging
}
// fprintf(align_score_file, "\n");
Expand Down Expand Up @@ -622,8 +631,8 @@ void ksw_extd2_cpp(
// ez->score);
for (int i = 0; i < qlen + tlen - 1; i++) {
int len = i + 1 < qlen + tlen - i - 1 ? i + 1 : qlen + tlen - i - 1;
len = min_n_col < len ? min_n_col : len;
fprintf(align_score_file, "#%d ", i);
len = w < len ? w : len;
fprintf(align_score_file, "#%d (%d) ", i, len);
for (int j = 0; j < len; j++)
fprintf(align_score_file, "%x ", ((uint8_t *)p)[i * min_n_col + j]);
fprintf(align_score_file, "\n");
Expand Down
44 changes: 31 additions & 13 deletions ksw2_extd2_sse.c
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -158,13 +158,21 @@ void ksw_extd2_sse(void *km, int qlen, const uint8_t *query, int tlen, const uin
st0 = st, en0 = en;
st = st / 16 * 16, en = (en + 16) / 16 * 16 - 1;
// set boundary conditions
if (st > 0) {
if (st - 1 >= last_st && st - 1 <= last_en) {
x1 = x8[st - 1], x21 = x28[st - 1], v1 = v8[st - 1]; // (r-1,s-1) calculated in the last round
} else {
if (st0 > 0) {
if (st0 - 1 >= last_st && st0 - 1 <= last_en) {
if (r == 160)
printf("this is r == 160, last_st = %d, st0 = %d\n",
last_st, st0);
x1 = x8[st0 - 1], x21 = x28[st0 - 1],
v1 = v8[st0 - 1]; // (r-1,s-1) calculated in the last round
} else {
x1 = -q - e, x21 = -q2 - e2;
v1 = -q - e;
}
// v1 = -q - e;
Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For entries not in the band (but not outside of the tlen x qlen boundary), extd2_see initialize it to -q-e, instead of
v1 = r == 0? -q - e : r < long_thres? -e : r == long_thres? long_diff : -e2;
as described in the Minimap2 paper.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not sure if this is intentional, or a bug of the original minimap2 implementation.

v1 = r == 0 ? -q - e
: r < long_thres ? -e
: r == long_thres ? long_diff
: -e2;
}
} else {
x1 = -q - e, x21 = -q2 - e2;
v1 = r == 0? -q - e : r < long_thres? -e : r == long_thres? long_diff : -e2;
Expand Down Expand Up @@ -393,12 +401,18 @@ void ksw_extd2_sse(void *km, int qlen, const uint8_t *query, int tlen, const uin
if (r == qlen + tlen - 2 && en0 == tlen - 1)
ez->score = H0;
}
last_st = st, last_en = en;
last_st = st0, last_en = en0;
#ifdef DEBUG
fprintf(align_debug_file, "#%d (st=%d en=%d) ", r, st0, en0);
// for (t = st0; t <= en0; ++t) printf("(%d,%d)\t(%d,%d,%d,%d)\t%d\n", r, t, ((int8_t*)u)[t], ((int8_t*)v)[t], ((int8_t*)x)[t], ((int8_t*)y)[t], H[t]); // for debugging
// for (t = st0; t <= en0; ++t)
// printf("(%d,%d)\t(%d,%d,%d,%d)\t%d\n", r, t, ((int8_t*)u)[t],
// ((int8_t*)v)[t], ((int8_t*)x)[t], ((int8_t*)y)[t], H[t]); // for
// debugging
for (t = st0; t <= en0; ++t) {
fprintf(align_debug_file, "%d ", H[t]);
// fprintf(align_debug_file, "%d ", H[t]);
fprintf(align_debug_file, "[%d %d v1 %d s %d p %x]%d ", v8[t],
u8[t], v1, ((int8_t *)s)[t],
((uint8_t *)p)[r * n_col_ * 16 - st + t], H[t]);
}
fprintf(align_debug_file, "\n");
#endif // DEBUG
Expand All @@ -409,10 +423,14 @@ void ksw_extd2_sse(void *km, int qlen, const uint8_t *query, int tlen, const uin
n_col = (n_col < w + 1 ? n_col : w + 1);
for (int i = 0; i < qlen + tlen - 1; i++) {
int len = i + 1 < qlen + tlen - i - 1 ? i + 1 : qlen + tlen - i - 1;
len = n_col < len ? n_col : len;
fprintf(align_score_file, "#%d ", i);
for (int j = 0; j < len; j++) {
fprintf(align_score_file, "%x ", ((uint8_t*)p)[i * n_col_*16 + j]);
len = w < len ? w : len;
int st = 0;
if (st < i - qlen + 1) st = i - qlen + 1;
if (st < (i - wr + 1) >> 1) st = (i - wr + 1) >> 1; // take the ceil
fprintf(align_score_file, "#%d (%d) ", i, len, off[i] - st);
for (int j = st - off[i]; j < len + st - off[i]; j++) {
fprintf(align_score_file, "%x ",
((uint8_t *)p)[i * n_col_ * 16 + j]);
}
fprintf(align_score_file, "\n");
}
Expand Down