diff --git a/lib_tsalign/src/a_star_aligner/alignment_geometry.rs b/lib_tsalign/src/a_star_aligner/alignment_geometry.rs index a37a82cb..6680bb94 100644 --- a/lib_tsalign/src/a_star_aligner/alignment_geometry.rs +++ b/lib_tsalign/src/a_star_aligner/alignment_geometry.rs @@ -52,41 +52,55 @@ impl AlignmentRange { self.offset.query..self.limit.query } - pub fn move_offsets_left(&self) -> Self { - Self::new_offset_limit( + #[must_use] + pub fn move_offsets_left(&self) -> Option { + Some(Self::new_offset_limit( AlignmentCoordinates::new( - self.offset.reference.checked_sub(1).unwrap(), - self.offset.query.checked_sub(1).unwrap(), + self.offset.reference.checked_sub(1)?, + self.offset.query.checked_sub(1)?, ), self.limit, - ) + )) } - pub fn move_offsets_right(&self) -> Self { - assert!(self.offset.reference < self.limit.reference); - assert!(self.offset.query < self.limit.query); + #[must_use] + pub fn move_offsets_right(&self) -> Option { + let new_ref_offset = self.offset.reference.checked_add(1)?; + let new_qry_offset = self.offset.query.checked_add(1)?; + + if new_ref_offset > self.limit.reference || new_qry_offset > self.limit.query { + return None; + } - Self::new_offset_limit( - AlignmentCoordinates::new(self.offset.reference + 1, self.offset.query + 1), + Some(Self::new_offset_limit( + AlignmentCoordinates::new(new_ref_offset, new_qry_offset), self.limit, - ) + )) } - pub fn move_limits_left(&self) -> Self { - assert!(self.offset.reference < self.limit.reference); - assert!(self.offset.query < self.limit.query); + #[must_use] + pub fn move_limits_left(&self) -> Option { + let new_ref_limit = self.limit.reference.checked_sub(1)?; + let new_qry_limit = self.limit.query.checked_sub(1)?; - Self::new_offset_limit( + if new_ref_limit < self.offset.reference || new_qry_limit > self.offset.query { + return None; + } + Some(Self::new_offset_limit( self.offset, - AlignmentCoordinates::new(self.limit.reference - 1, self.limit.query - 1), - ) + AlignmentCoordinates::new(new_ref_limit, new_qry_limit), + )) } - pub fn move_limits_right(&self) -> Self { - Self::new_offset_limit( + #[must_use] + pub fn move_limits_right(&self) -> Option { + Some(Self::new_offset_limit( self.offset, - AlignmentCoordinates::new(self.limit.reference + 1, self.limit.query + 1), - ) + AlignmentCoordinates::new( + self.limit.reference.checked_add(1)?, + self.limit.query.checked_add(1)?, + ), + )) } } diff --git a/lib_tsalign/src/a_star_aligner/alignment_result.rs b/lib_tsalign/src/a_star_aligner/alignment_result.rs index 57818fda..f1bfc48a 100644 --- a/lib_tsalign/src/a_star_aligner/alignment_result.rs +++ b/lib_tsalign/src/a_star_aligner/alignment_result.rs @@ -237,6 +237,7 @@ impl AlignmentResult> AlignmentResult { + #[allow(clippy::too_many_lines)] pub fn extend_beyond_range_with_equal_cost< AlphabetType: Alphabet, SubsequenceType: GenomeSequence + ?Sized, @@ -248,6 +249,7 @@ impl> config: &TemplateSwitchConfig, ) { let Some(range) = range else { + trace!("No range given to extend beyond."); return; }; let Self::WithTarget { @@ -255,6 +257,7 @@ impl> statistics, } = self else { + trace!("There is no alignment, therefore we cannot postprocess it."); return; }; if config.left_flank_length > 0 || config.right_flank_length > 0 { @@ -263,7 +266,7 @@ impl> } // Compute cost before extending. - let initial_cost = alignment.compute_cost( + let mut current_cost = alignment.compute_cost( reference, query, range.reference_offset(), @@ -271,18 +274,18 @@ impl> config, ); let alignment_cost = (statistics.cost.round().raw() as u64).into(); - assert_eq!( - initial_cost, + debug_assert_eq!( + current_cost, alignment_cost, - "computed cost {initial_cost} != alignment cost {alignment_cost}; {}", + "computed cost {current_cost} != alignment cost {alignment_cost}; {}", alignment.cigar() ); // Extend left with equal cost. - while range.reference_offset() > 0 && range.query_offset() > 0 { + while let Some(new_range) = range.move_offsets_left() { // Determine extension alignment type. - let reference_character = reference[range.reference_offset() - 1].clone(); - let query_character = query[range.query_offset() - 1].clone(); + let reference_character = reference[new_range.reference_offset()].clone(); + let query_character = query[new_range.query_offset()].clone(); let alignment_type = if reference_character == query_character { AlignmentType::PrimaryMatch } else { @@ -290,15 +293,16 @@ impl> }; // Insert alignment type at beginning of alignment. - if alignment.inner_mut()[0].1 == alignment_type { + if alignment + .inner_mut() + .first() + .is_some_and(|(_, ty)| *ty == alignment_type) + { alignment.inner_mut()[0].0 += 1; } else { alignment.inner_mut().insert(0, (1, alignment_type)); } - // Update ranges. - let new_range = range.move_offsets_left(); - // Compute cost. let new_cost = alignment.compute_cost( reference, @@ -308,25 +312,29 @@ impl> config, ); - if new_cost > initial_cost { + if new_cost > current_cost { // If cost is not equal, revert alignment and break. alignment.inner_mut()[0].0 -= 1; if alignment.inner_mut()[0].0 == 0 { alignment.inner_mut().remove(0); } break; - } else { - // If cost is equal, commit range and continue. - *range = new_range; - assert_eq!(new_cost, initial_cost); } + current_cost = new_cost; + *range = new_range; } // Extend right with equal cost. while range.reference_limit() < reference.len() && range.query_limit() < query.len() { + // Update ranges. + let Some(new_range) = range.move_limits_right() else { + // Never occurs due to loop invariant + break; + }; + // Determine extension alignment type. - let reference_character = reference[range.reference_limit()].clone(); - let query_character = query[range.query_limit()].clone(); + let reference_character = reference[new_range.reference_limit() - 1].clone(); + let query_character = query[new_range.query_limit() - 1].clone(); let alignment_type = if reference_character == query_character { AlignmentType::PrimaryMatch } else { @@ -334,15 +342,19 @@ impl> }; // Insert alignment type at end of alignment. - if alignment.inner_mut().last().unwrap().1 == alignment_type { - alignment.inner_mut().last_mut().unwrap().0 += 1; + if alignment + .inner_mut() + .last() + .is_some_and(|(_, ty)| *ty == alignment_type) + { + let Some((n, _)) = alignment.inner_mut().last_mut() else { + unreachable!("Due to if condition") + }; + *n += 1; } else { alignment.inner_mut().push((1, alignment_type)); } - // Update ranges. - let new_range = range.move_limits_right(); - // Compute cost. let new_cost = alignment.compute_cost( reference, @@ -352,24 +364,26 @@ impl> config, ); - if new_cost > initial_cost { + if new_cost > current_cost { // If cost is not equal, revert alignment and break. - alignment.inner_mut().last_mut().unwrap().0 -= 1; - if alignment.inner_mut().last().unwrap().0 == 0 { + let Some((n, _)) = alignment.inner_mut().last_mut() else { + unreachable!("As we have just added at least one element") + }; + *n -= 1; + if *n == 0 { alignment.inner_mut().pop(); } - break; - } else { - // If cost is equal, commit range and continue. - *range = new_range; - assert_eq!(new_cost, initial_cost); } + current_cost = new_cost; + *range = new_range; } + // Update statistics from updated range statistics.reference_offset = range.reference_offset(); statistics.query_offset = range.query_offset(); } + #[allow(clippy::too_many_lines)] pub fn compute_ts_equal_cost_ranges< AlphabetType: Alphabet, SubsequenceType: GenomeSequence + ?Sized, @@ -385,6 +399,7 @@ impl> statistics, } = self else { + trace!("There is no alignment, therefore we cannot postprocess it."); return; }; if config.left_flank_length > 0 || config.right_flank_length > 0 { @@ -392,163 +407,161 @@ impl> return; } - let reference_offset = range - .as_ref() - .map(|range| range.reference_offset()) - .unwrap_or(0); + let reference_offset = range.as_ref().map_or( + 0, + super::alignment_geometry::AlignmentRange::reference_offset, + ); let query_offset = range .as_ref() - .map(|range| range.query_offset()) - .unwrap_or(0); + .map_or(0, super::alignment_geometry::AlignmentRange::query_offset); for i in 0..alignment.inner_mut().len() { let alignment_type = alignment.inner_mut()[i].1; - match alignment_type { - super::template_switch_distance::AlignmentType::TemplateSwitchEntrance { - mut equal_cost_range, - .. - } => { - equal_cost_range.min_start = 0; - equal_cost_range.max_start = 0; - equal_cost_range.min_end = 0; - equal_cost_range.max_end = 0; - - let initial_cost = alignment.compute_cost( + if let super::template_switch_distance::AlignmentType::TemplateSwitchEntrance { + mut equal_cost_range, + .. + } = alignment_type + { + equal_cost_range.min_start = 0; + equal_cost_range.max_start = 0; + equal_cost_range.min_end = 0; + equal_cost_range.max_end = 0; + + // Keep track of the current cost, which might get less while extending but never increase + let mut current_cost = alignment.compute_cost( + reference, + query, + reference_offset, + query_offset, + config, + ); + debug_assert_eq!(current_cost, (statistics.cost.round().raw() as u64).into()); + + // Move start to the left + { + let mut min_start_alignment = alignment.clone(); + let mut i = i; + + while min_start_alignment.move_template_switch_start_backwards( reference, query, reference_offset, query_offset, - config, - ); - assert_eq!(initial_cost, (statistics.cost.round().raw() as u64).into()); - - { - let mut min_start_alignment = alignment.clone(); - let mut i = i; - - while min_start_alignment.move_template_switch_start_backwards( + &mut i, + ) { + let new_cost = min_start_alignment.compute_cost( reference, query, reference_offset, query_offset, - &mut i, - ) { - let new_cost = min_start_alignment.compute_cost( - reference, - query, - reference_offset, - query_offset, - config, + config, + ); + if new_cost > current_cost { + trace!( + "Stopping moving TS start backwards because cost {new_cost} > {current_cost} with alignment {min_start_alignment}" ); - if new_cost > initial_cost { - trace!( - "Stopping moving TS start backwards because cost {new_cost} > {initial_cost} with alignment {min_start_alignment}" - ); - break; - } else { - assert_eq!(new_cost, initial_cost); - equal_cost_range.min_start -= 1; - } + break; } + current_cost = new_cost; + equal_cost_range.min_start -= 1; } + } - { - let mut max_start_alignment = alignment.clone(); - let mut i = i; + // Move start to the right + { + let mut max_start_alignment = alignment.clone(); + let mut i = i; - while max_start_alignment.move_template_switch_start_forwards( + while max_start_alignment.move_template_switch_start_forwards( + reference, + query, + reference_offset, + query_offset, + &mut i, + ) { + let new_cost = max_start_alignment.compute_cost( reference, query, reference_offset, query_offset, - &mut i, - ) { - let new_cost = max_start_alignment.compute_cost( - reference, - query, - reference_offset, - query_offset, - config, + config, + ); + if new_cost > current_cost { + trace!( + "Stopping moving TS start forwards because cost {new_cost} > {current_cost} with alignment {max_start_alignment}" ); - if new_cost > initial_cost { - trace!( - "Stopping moving TS start forwards because cost {new_cost} > {initial_cost} with alignment {max_start_alignment}" - ); - break; - } else { - assert_eq!(new_cost, initial_cost); - equal_cost_range.max_start += 1; - } + break; } + current_cost = new_cost; + equal_cost_range.max_start += 1; } + } - { - let mut min_end_alignment = alignment.clone(); - while min_end_alignment.move_template_switch_end_backwards( + // Move end to the left + { + let mut min_end_alignment = alignment.clone(); + while min_end_alignment.move_template_switch_end_backwards( + reference, + query, + reference_offset, + query_offset, + i, + ) { + let new_cost = min_end_alignment.compute_cost( reference, query, reference_offset, query_offset, - i, - ) { - let new_cost = min_end_alignment.compute_cost( - reference, - query, - reference_offset, - query_offset, - config, + config, + ); + if new_cost > current_cost { + trace!( + "Stopping moving TS end backwards because cost {new_cost} > {current_cost} with alignment {min_end_alignment}" ); - if new_cost > initial_cost { - trace!( - "Stopping moving TS end backwards because cost {new_cost} > {initial_cost} with alignment {min_end_alignment}" - ); - break; - } else { - assert_eq!(new_cost, initial_cost); - equal_cost_range.min_end -= 1; - } + break; } + current_cost = new_cost; + equal_cost_range.min_end -= 1; } + } - { - let mut max_end_alignment = alignment.clone(); - while max_end_alignment.move_template_switch_end_forwards( + // Move end to the right + { + let mut max_end_alignment = alignment.clone(); + while max_end_alignment.move_template_switch_end_forwards( + reference, + query, + reference_offset, + query_offset, + i, + ) { + let new_cost = max_end_alignment.compute_cost( reference, query, reference_offset, query_offset, - i, - ) { - let new_cost = max_end_alignment.compute_cost( - reference, - query, - reference_offset, - query_offset, - config, + config, + ); + if new_cost > current_cost { + trace!( + "Stopping moving TS end forwards because cost {new_cost} > {current_cost} with alignment {max_end_alignment}" ); - if new_cost > initial_cost { - trace!( - "Stopping moving TS end forwards because cost {new_cost} > {initial_cost} with alignment {max_end_alignment}" - ); - break; - } else { - assert_eq!(new_cost, initial_cost); - equal_cost_range.max_end += 1; - } + break; } + current_cost = new_cost; + equal_cost_range.max_end += 1; } - - let super::template_switch_distance::AlignmentType::TemplateSwitchEntrance { - equal_cost_range: alignment_equal_cost_range, - .. - } = &mut alignment.inner_mut()[i].1 - else { - unreachable!() - }; - *alignment_equal_cost_range = equal_cost_range; } - _ => { /* Do nothing. */ } + + let super::template_switch_distance::AlignmentType::TemplateSwitchEntrance { + equal_cost_range: alignment_equal_cost_range, + .. + } = &mut alignment.inner_mut()[i].1 + else { + unreachable!() + }; + *alignment_equal_cost_range = equal_cost_range; } } } diff --git a/lib_tsalign/src/a_star_aligner/alignment_result/alignment/template_switch_specifics.rs b/lib_tsalign/src/a_star_aligner/alignment_result/alignment/template_switch_specifics.rs index 363eb8bd..4cf16da1 100644 --- a/lib_tsalign/src/a_star_aligner/alignment_result/alignment/template_switch_specifics.rs +++ b/lib_tsalign/src/a_star_aligner/alignment_result/alignment/template_switch_specifics.rs @@ -3,6 +3,7 @@ use compact_genome::interface::{ sequence::GenomeSequence, }; use generic_a_star::cost::AStarCost; +use log::error; use crate::{ a_star_aligner::{ @@ -20,10 +21,12 @@ impl Alignment { /// Moves the start of a template switch one character pair to the left if possible, moving the character pair into the template switch. /// /// Only works if the alignment preceding the template switch is a match or substitution, and not a flank. - /// In this case it returns true, and otherwise false. + /// + /// In case of success, it returns true, and otherwise false. /// /// Compact index identifies the start of the template switch in terms of the compact representation of the alignment. /// See [`Self::iter_compact`]. + #[allow(clippy::too_many_lines)] pub fn move_template_switch_start_backwards< AlphabetType: Alphabet, SubsequenceType: GenomeSequence + ?Sized, @@ -35,24 +38,32 @@ impl Alignment { query_offset: usize, compact_index: &mut usize, ) -> bool { - let AlignmentType::TemplateSwitchEntrance { - first_offset, - primary, - secondary, - direction, - .. - } = self.alignment[*compact_index].1 + let Some(&( + _, + AlignmentType::TemplateSwitchEntrance { + first_offset, + primary, + secondary, + direction, + .. + }, + )) = self.alignment.get(*compact_index) else { - panic!() + error!( + "(Programming) Error: compact index {compact_index} does not point to a template switch entrance" + ); + return false; }; - if matches!( - self.alignment.get(compact_index.saturating_sub(1)), - Some(( - _, - AlignmentType::PrimaryMatch | AlignmentType::PrimarySubstitution - )) - ) { + if *compact_index > 0 + && matches!( + self.alignment.get(*compact_index - 1), + Some(( + _, + AlignmentType::PrimaryMatch | AlignmentType::PrimarySubstitution + )) + ) + { // Compute TS inner first indices. let mut stream = AlignmentStream::new_with_offset(reference_offset, query_offset); stream.push_all(self.iter_compact_cloned().take(*compact_index)); @@ -60,34 +71,39 @@ impl Alignment { TemplateSwitchPrimary::Reference => stream.head_coordinates().reference(), TemplateSwitchPrimary::Query => stream.head_coordinates().query(), }; - let ts_inner_secondary_index = usize::try_from( - isize::try_from(match secondary { - TemplateSwitchSecondary::Reference => stream.head_coordinates().reference(), - TemplateSwitchSecondary::Query => stream.head_coordinates().query(), - }) - .unwrap() - .checked_add(first_offset) - .unwrap(), - ) - .unwrap(); + if ts_inner_primary_index == 0 { + // We cannot extend more backwards. + return false; + } + let Some(ts_inner_secondary_index) = isize::try_from(match secondary { + TemplateSwitchSecondary::Reference => stream.head_coordinates().reference(), + TemplateSwitchSecondary::Query => stream.head_coordinates().query(), + }) + .ok() + .and_then(|i| usize::try_from(i.checked_add(first_offset)?).ok()) else { + error!( + "(Programming) Error: finding inner secondary index -- integer math does not check out" + ); + return false; + }; // Check if indices can be moved while staying in bounds. match direction { - TemplateSwitchDirection::Forward => { - if ts_inner_secondary_index == 0 { - return false; - } - } - TemplateSwitchDirection::Reverse => { - if ts_inner_secondary_index >= secondary.get(reference, query).len() { - return false; - } + TemplateSwitchDirection::Forward if ts_inner_secondary_index == 0 => return false, + TemplateSwitchDirection::Reverse + if ts_inner_secondary_index >= secondary.get(reference, query).len() => + { + return false; } + _ => {} } // Remove one match or substitution from before the TS. let multiplicity = &mut self.alignment[*compact_index - 1].0; - assert!(*multiplicity > 0); + if *multiplicity == 0 { + error!("Invalid input alignment! Cannot have multiplicity of 0."); + return false; + } *multiplicity -= 1; // Remove the alignment entry if it has zero multiplicity. @@ -97,9 +113,11 @@ impl Alignment { } // Check if the new inner pair is a match or a substitution. + // ts_inner_primary_index > 0, otherwise we would've returned false earlier let primary_char = primary.get(reference, query)[ts_inner_primary_index - 1].clone(); let secondary_char = match direction { TemplateSwitchDirection::Forward => { + // ts_inner_secondary_index > 0, otherwise we would've returned false earlier secondary.get(reference, query)[ts_inner_secondary_index - 1].clone() } TemplateSwitchDirection::Reverse => { @@ -113,7 +131,11 @@ impl Alignment { }; // Insert new inner alignment. - if self.alignment[*compact_index + 1].1 == inner_alignment_type { + if self + .alignment + .get(*compact_index + 1) + .is_some_and(|(_, ty)| *ty == inner_alignment_type) + { self.alignment[*compact_index + 1].0 += 1; } else { self.alignment @@ -126,7 +148,7 @@ impl Alignment { let AlignmentType::TemplateSwitchEntrance { first_offset, .. } = &mut self.alignment[*compact_index].1 else { - unreachable!(); + unreachable!("reborrow of already borrowed element, definitely exists"); }; *first_offset += 2; } @@ -137,7 +159,7 @@ impl Alignment { .iter_mut() .find(|(_, alignment_type)| alignment_type.is_template_switch_exit()) else { - unreachable!(); + unreachable!("There should be a TS exit.."); }; *anti_primary_gap += 1; @@ -170,22 +192,30 @@ impl Alignment { let AlignmentType::TemplateSwitchEntrance { direction, .. } = self.alignment[*compact_index].1 else { - panic!() + error!( + "Error: compact index {compact_index} does not point to a template switch entrance" + ); + return false; }; // Assert that no flanks are involved. - assert!( - self.alignment + if *compact_index != 0 + && !self + .alignment .get(*compact_index - 1) - .map(|(_, alignment_type)| !matches!( - alignment_type, - AlignmentType::PrimaryFlankDeletion - | AlignmentType::PrimaryFlankInsertion - | AlignmentType::PrimaryFlankSubstitution - | AlignmentType::PrimaryFlankMatch - )) - .unwrap_or(true) - ); + .is_none_or(|(_, alignment_type)| { + !matches!( + alignment_type, + AlignmentType::PrimaryFlankDeletion + | AlignmentType::PrimaryFlankInsertion + | AlignmentType::PrimaryFlankSubstitution + | AlignmentType::PrimaryFlankMatch + ) + }) + { + error!("We did not want flanks to get involved..."); + return false; + } if let Some((_, AlignmentType::SecondaryMatch | AlignmentType::SecondarySubstitution)) = self.alignment.get(*compact_index + 1) @@ -202,8 +232,12 @@ impl Alignment { } // Remove one match or substitution from inside the TS. + // compact_index + 1 is valid since the index points to the TS entrance, and that is certainly not the last entry let multiplicity = &mut self.alignment[*compact_index + 1].0; - assert!(*multiplicity > 0); + if *multiplicity == 0 { + error!("Invalid input alignment! Cannot have multiplicity of 0."); + return false; + } *multiplicity -= 1; // Remove the alignment entry if it has zero multiplicity. @@ -221,7 +255,9 @@ impl Alignment { }; // Insert new outer alignment. - if self.alignment.get(*compact_index - 1).map(|t| t.1) == Some(outer_alignment_type) { + if *compact_index != 0 + && self.alignment.get(*compact_index - 1).map(|t| t.1) == Some(outer_alignment_type) + { self.alignment[*compact_index - 1].0 += 1; } else { self.alignment @@ -235,7 +271,7 @@ impl Alignment { let AlignmentType::TemplateSwitchEntrance { first_offset, .. } = &mut self.alignment[*compact_index].1 else { - unreachable!(); + unreachable!("Merely a reborrow"); }; *first_offset -= 2; } @@ -263,6 +299,7 @@ impl Alignment { /// /// Compact index identifies the start of the template switch in terms of the compact representation of the alignment. /// See [`Self::iter_compact`]. + #[allow(clippy::too_many_lines)] pub fn move_template_switch_end_forwards< AlphabetType: Alphabet, SubsequenceType: GenomeSequence + ?Sized, @@ -282,19 +319,22 @@ impl Alignment { .. } = self.alignment[compact_index].1 else { - panic!() + error!( + "Error: compact index {compact_index} does not point to a template switch entrance" + ); + return false; }; - let mut exit_index = compact_index - + self - .alignment - .iter() - .skip(compact_index) - .enumerate() - .find(|(_, (_, alignment_type))| { + let Some((exit_index_offset, _)) = + self.alignment.iter().skip(compact_index).enumerate().find( + |(_, (_, alignment_type))| { matches!(alignment_type, AlignmentType::TemplateSwitchExit { .. }) - }) - .unwrap() - .0; + }, + ) + else { + error!("There should be a TS exit after the TS entrance"); + return false; + }; + let mut exit_index = compact_index + exit_index_offset; let inner_secondary_length = self .alignment .iter() @@ -327,21 +367,26 @@ impl Alignment { TemplateSwitchPrimary::Reference => stream.head_coordinates().reference(), TemplateSwitchPrimary::Query => stream.head_coordinates().query(), }; - let ts_inner_secondary_index = usize::try_from( - isize::try_from(match secondary { - TemplateSwitchSecondary::Reference => stream.tail_coordinates().reference(), - TemplateSwitchSecondary::Query => stream.tail_coordinates().query(), - }) - .unwrap() - .checked_add(first_offset) - .unwrap(), - ) - .unwrap(); + + let Some(ts_inner_secondary_index) = isize::try_from(match secondary { + TemplateSwitchSecondary::Reference => stream.tail_coordinates().reference(), + TemplateSwitchSecondary::Query => stream.tail_coordinates().query(), + }) + .ok() + .and_then(|i| usize::try_from(i.checked_add(first_offset)?).ok()) else { + error!("Error finding inner secondary index -- integer math does not check out"); + return false; + }; let ts_inner_secondary_index = match direction { TemplateSwitchDirection::Forward => { - let ts_inner_secondary_index = ts_inner_secondary_index - .checked_add(inner_secondary_length) - .unwrap(); + let Some(ts_inner_secondary_index) = + ts_inner_secondary_index.checked_add(inner_secondary_length) + else { + error!( + "ts_inner_secondary_index {ts_inner_secondary_index} should be at least {inner_secondary_length} away from any boundary" + ); + return false; + }; // Check if indices can be moved while staying in bounds. if ts_inner_secondary_index >= secondary.get(reference, query).len() { return false; @@ -349,9 +394,14 @@ impl Alignment { ts_inner_secondary_index } TemplateSwitchDirection::Reverse => { - let ts_inner_secondary_index = ts_inner_secondary_index - .checked_sub(inner_secondary_length) - .unwrap(); + let Some(ts_inner_secondary_index) = + ts_inner_secondary_index.checked_sub(inner_secondary_length) + else { + error!( + "ts_inner_secondary_index {ts_inner_secondary_index} should be at least {inner_secondary_length} away from any boundary" + ); + return false; + }; // Check if indices can be moved while staying in bounds. if ts_inner_secondary_index == 0 { return false; @@ -362,7 +412,10 @@ impl Alignment { // Remove one match or substitution from after the TS. let multiplicity = &mut self.alignment[exit_index + 1].0; - assert!(*multiplicity > 0); + if *multiplicity == 0 { + error!("Invalid input alignment! Cannot have multiplicity of 0."); + return false; + } *multiplicity -= 1; // Remove the alignment entry if it has zero multiplicity. @@ -377,6 +430,7 @@ impl Alignment { secondary.get(reference, query)[ts_inner_secondary_index].clone() } TemplateSwitchDirection::Reverse => { + // ts_inner_secondary_index > 0 since we checked that earlier and would have returned otherwise. secondary.get(reference, query)[ts_inner_secondary_index - 1].complement() } }; @@ -398,7 +452,7 @@ impl Alignment { let AlignmentType::TemplateSwitchExit { anti_primary_gap } = &mut self.alignment[exit_index].1 else { - unreachable!(); + unreachable!("Merely a reborrow"); }; *anti_primary_gap += 1; @@ -428,35 +482,44 @@ impl Alignment { query_offset: usize, compact_index: usize, ) -> bool { - assert!(matches!( + if !matches!( self.alignment[compact_index].1, AlignmentType::TemplateSwitchEntrance { .. } - )); - let mut exit_index = compact_index - + self - .alignment - .iter() - .skip(compact_index) - .enumerate() - .find(|(_, (_, alignment_type))| { + ) { + error!( + "Error: compact index {compact_index} does not point to a template switch entrance" + ); + return false; + } + let Some((exit_index_offset, _)) = + self.alignment.iter().skip(compact_index).enumerate().find( + |(_, (_, alignment_type))| { matches!(alignment_type, AlignmentType::TemplateSwitchExit { .. }) - }) - .unwrap() - .0; + }, + ) + else { + error!("There must be a TS exit after the TS entrance!"); + return false; + }; + let mut exit_index = compact_index + exit_index_offset; // Assert that no flanks are involved. - assert!( - self.alignment - .get(exit_index + 1) - .map(|(_, alignment_type)| !matches!( + if !self + .alignment + .get(exit_index + 1) + .is_none_or(|(_, alignment_type)| { + !matches!( alignment_type, AlignmentType::PrimaryFlankDeletion | AlignmentType::PrimaryFlankInsertion | AlignmentType::PrimaryFlankSubstitution | AlignmentType::PrimaryFlankMatch - )) - .unwrap_or(true) - ); + ) + }) + { + error!("We have some unexpected flanks!"); + return false; + } if let Some((_, AlignmentType::SecondaryMatch | AlignmentType::SecondarySubstitution)) = self.alignment.get(exit_index - 1) @@ -474,7 +537,10 @@ impl Alignment { // Remove one match or substitution from inside the TS. let multiplicity = &mut self.alignment[exit_index - 1].0; - assert!(*multiplicity > 0); + if *multiplicity == 0 { + error!("Invalid input alignment! Cannot have multiplicity of 0."); + return false; + } *multiplicity -= 1; // Remove the alignment entry if it has zero multiplicity. @@ -506,7 +572,7 @@ impl Alignment { .iter_mut() .find(|(_, alignment_type)| alignment_type.is_template_switch_exit()) else { - unreachable!(); + unreachable!("Just a reborrow"); }; *anti_primary_gap -= 1;