-
Notifications
You must be signed in to change notification settings - Fork 30
Description
I think this may be related to #54 but I'm getting this on a x86-64.
While working on getting profile-profile alignment to work in PyFAMSA I noticed that I was starting to see some assertion failures but only in debug mode. I traced them back to the indexing in CProfile::GetGapStats:
Lines 360 to 376 in 5a326d5
| void CProfile::GetGapStats(vector<size_t> &stats) | |
| { | |
| stats.clear(); | |
| stats.resize(width + 1, data.size()); | |
| vector<size_t> stats_term; | |
| for(auto &p : data) | |
| { | |
| size_t seq_pos = 0; | |
| for(size_t i = 0; i <= p->size; ++i) | |
| { | |
| --stats[seq_pos]; | |
| seq_pos += p->n_gaps[i] + 1; | |
| } | |
| } | |
| } |
I added explicit assertions here:
for(size_t i = 0; i <= p->size; ++
{
assert(seq_pos < stats.size());
--stats[seq_pos];
seq_pos += p->n_gaps[i] + 1;
}which now trigger even in the FAMSA binary with the provided input:
$ ./bin/famsa x.afa y.afa out.afa
FAMSA (Fast and Accurate Multiple Sequence Alignment)
version 2.4.1-45c9b2b-dirty (2025-05-09)
S. Deorowicz, A. Debudaj-Grabysz, A. Gudys
famsa: src/core/profile.cpp:372: void CProfile::GetGapStats(std::vector<long unsigned int>&): Assertion `seq_pos < stats.size()' failed.
fish: Job 1, './bin/famsa x.afa y.afa out.afa' terminated by signal SIGABRT (Abort)However, this only occurs in the refinement stage (which is why I think it may be linked to #54); the same command succeeds with -refine_mode off.
Proposed fix
In CFAMSA::RefineMostEmptyAndFullColumn, the CProfile::GetGapStats function is called once here:
Lines 23 to 29 in 5a326d5
| size_t size = profile_to_refine->data.front()->gapped_size; | |
| size_t card = profile_to_refine->data.size(); | |
| dest_prof_id.clear(); | |
| if (!valid_gap_stats) | |
| profile_to_refine->GetGapStats(gap_stats); |
However, it seems like the profile_to_refine->width may not have been properly set at that stage, I'm not sure why. By changing this part to:
size_t size = profile_to_refine->data.front()->gapped_size;
size_t card = profile_to_refine->data.size();
dest_prof_id.clear();
if (!valid_gap_stats)
{
profile_to_refine->width = size;
profile_to_refine->GetGapStats(gap_stats);
}The assertion doesn't trigger anymore and there is not segfault in my Python code. I'm not entirely sure however whether this is the correct place for the fix or if the fix should actually be applied upstream somewhere in the profile-profile code in profile_seq.cpp