Skip to content

Commit ad2bd9f

Browse files
committed
Check maple alignments
Add check that mutations in maple alignments do not have the same nucleotide as the reference sequence (this can cause segmentation faults later on).
1 parent fa99239 commit ad2bd9f

File tree

1 file changed

+12
-4
lines changed

1 file changed

+12
-4
lines changed

alignment/alignment.cpp

+12-4
Original file line numberDiff line numberDiff line change
@@ -910,16 +910,24 @@ void cmaple::Alignment::readMaple(std::istream& aln_stream) {
910910
"sequence length (" +
911911
convertPosTypeToString(static_cast<PositionType>(ref_seq.size())) + ")!");
912912
}
913-
} else if (cmaple::verbose_mode >= cmaple::VB_MED) {
914-
outWarning("Ignoring <Length> of " + tmp +
915-
". <Length> is only appliable for 'N' or '-'.");
913+
} else if (cmaple::verbose_mode >= cmaple::VB_MED) {
914+
outWarning("Ignoring <Length> of " + tmp +
915+
". <Length> is only appliable for 'N' or '-'.");
916+
}
916917
}
917-
}
918918

919919
// add a new mutation into mutations
920920
if (state == TYPE_N || state == TYPE_DEL) {
921921
mutations.emplace_back(state, pos - 1, length);
922922
} else {
923+
StateType refState = ref_seq[pos - 1];
924+
if(refState == state)
925+
{
926+
throw std::logic_error(
927+
"Mutation at position " + convertPosTypeToString(pos) +
928+
" in sequence " + seq_name +
929+
" is equal to reference nucleotide. Check reference and alignment are correct.");
930+
}
923931
mutations.emplace_back(state, pos - 1);
924932
}
925933
}

0 commit comments

Comments
 (0)