@@ -2048,8 +2048,7 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
20482048
20492049 // Fill this in with the indexes of pairs of alignments we will output
20502050 // each alignment is stored as <fragment index, alignment index> into alignments
2051- // fragment_index should be the same for both ends, unless one was rescued.
2052- // Note that for a failed rescue, we still include the "pair" here that has one read unmapped.
2051+ // fragment_index should be the same for both ends, unless one was rescued
20532052 vector<std::array<read_alignment_index_t , 2 >> paired_alignments;
20542053 paired_alignments.reserve (alignments.size ());
20552054
@@ -2062,11 +2061,8 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
20622061 vector<std::array<vector<vector<size_t >>, 2 >> alignment_groups (alignments.size ());
20632062
20642063 // Grab all the scores in order for MAPQ computation.
2065- // These correspond 1 to 1 with paired_alignments.
20662064 vector<double > paired_scores;
20672065 paired_scores.reserve (alignments.size ());
2068- // Record the fragment distances, which are 1 to 1 with paired_alignments
2069- // and feed into MAPQ capping.
20702066 vector<int64_t > fragment_distances;
20712067 fragment_distances.reserve (alignments.size ());
20722068
@@ -2206,6 +2202,7 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
22062202 unpaired_scores[r].reserve (unpaired_alignments.size ());
22072203 }
22082204
2205+
22092206 array<vector<read_alignment_index_t >, 2 > supplementaries;
22102207 if (!unpaired_alignments.empty ()) {
22112208 // If we found some clusters that had no pair in a fragment cluster
@@ -2356,15 +2353,14 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
23562353 attempt_rescue (mapped_aln, rescued_aln, minimizers_by_read[1 - index.read ], index.read == 0 );
23572354
23582355 bool properly_paired = false ;
2359- int64_t fragment_dist;
23602356 double score;
23612357
23622358 if (rescued_aln.path ().mapping_size () != 0 ) {
23632359 // If we actually found an alignment
23642360
23652361 // Compute the distance
2366- fragment_dist = index.read == 0 ? distance_between (mapped_aln, rescued_aln)
2367- : distance_between (rescued_aln, mapped_aln);
2362+ int64_t fragment_dist = index.read == 0 ? distance_between (mapped_aln, rescued_aln)
2363+ : distance_between (rescued_aln, mapped_aln);
23682364
23692365 // Use it to make a pair score
23702366 score = score_alignment_pair (mapped_aln, rescued_aln, fragment_dist);
@@ -2377,13 +2373,13 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
23772373 set_annotation (mapped_aln, " fragment_length" , distance_to_annotation (fragment_dist));
23782374 set_annotation (rescued_aln, " fragment_length" , distance_to_annotation (fragment_dist));
23792375
2376+ // Track it for the distribution
2377+ fragment_distances.emplace_back (fragment_dist);
23802378 } else {
23812379 // If there's no rescue result, the score of the pair is the score of the one actually-aligned read
23822380 score = mapped_aln.score ();
2383- // And the fragment distance is unreachable. But don't add an annotation for it.
2384- fragment_dist = std::numeric_limits<int64_t >::max ();
23852381 }
2386-
2382+
23872383 // Add the annotations that don't always need a fragment distance.
23882384 set_annotation (mapped_aln, " rescuer" , true );
23892385 set_annotation (rescued_aln, " rescued" , true );
@@ -2411,7 +2407,7 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
24112407 paired_alignments.back ().at (r).check_for_read_in (r, alignments);
24122408 }
24132409#endif
2414- fragment_distances. emplace_back (fragment_dist);
2410+
24152411 paired_scores.emplace_back (score);
24162412 pair_types.push_back (index.read == 0 ? rescued_from_first : rescued_from_second);
24172413 better_cluster_count_by_pairs.emplace_back (better_cluster_count[mapped_index.fragment ]);
@@ -2451,8 +2447,7 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
24512447 }
24522448
24532449 if (find_supplementaries) {
2454- supplementaries = std::move (identify_supplementary_alignments (alignments, paired_alignments, paired_scores,
2455- fragment_distances, pair_types, better_cluster_count_by_pairs,
2450+ supplementaries = std::move (identify_supplementary_alignments (alignments, paired_alignments, paired_scores, pair_types,
24562451 unpaired_alignments, attempted_rescue_from, funnels));
24572452 }
24582453 }
@@ -2463,16 +2458,6 @@ pair<vector<Alignment>, vector<Alignment>> MinimizerMapper::map_paired(Alignment
24632458 funnels[r].stage (" winner" );
24642459 }
24652460 }
2466-
2467- // Make sure we haven't dropped any parts of our multi-vector per-pair
2468- // alignment data records.
2469- //
2470- // TODO: Remove the possibility to get this wrong by making these all
2471- // members of a struct.
2472- crash_unless (paired_scores.size () == paired_alignments.size ());
2473- crash_unless (fragment_distances.size () == paired_alignments.size ());
2474- crash_unless (pair_types.size () == paired_alignments.size ());
2475- crash_unless (better_cluster_count_by_pairs.size () == paired_alignments.size ());
24762461
24772462 // Fill this in with the alignments we will output.
24782463 std::array<vector<Alignment>, 2 > mappings;
@@ -3565,9 +3550,7 @@ array<vector<MinimizerMapper::read_alignment_index_t>, 2>
35653550MinimizerMapper::identify_supplementary_alignments (vector<std::array<vector<Alignment>, 2 >>& alignments,
35663551 vector<std::array<read_alignment_index_t , 2 >>& paired_alignments,
35673552 vector<double >& paired_scores,
3568- vector<int64_t >& fragment_distances,
35693553 vector<MinimizerMapper::PairType>& pair_types,
3570- vector<size_t >& better_cluster_count_by_pairs,
35713554 const vector<alignment_index_t >& unpaired_alignments,
35723555 const vector<bool >& attempted_rescue_from,
35733556 array<Funnel, 2 >& funnels) const {
@@ -3718,9 +3701,7 @@ MinimizerMapper::identify_supplementary_alignments(vector<std::array<vector<Alig
37183701 // Shift into the front of the vector
37193702 paired_alignments[i - s] = paired_alignments[i];
37203703 paired_scores[i - s] = paired_scores[i];
3721- fragment_distances[i - s] = fragment_distances[i];
37223704 pair_types[i - s] = pair_types[i];
3723- better_cluster_count_by_pairs[i - s] = better_cluster_count_by_pairs[s];
37243705
37253706 }
37263707 }
@@ -3746,9 +3727,7 @@ MinimizerMapper::identify_supplementary_alignments(vector<std::array<vector<Alig
37463727
37473728 paired_alignments.resize (paired_alignments.size () - s);
37483729 paired_scores.resize (paired_alignments.size ());
3749- fragment_distances.resize (paired_alignments.size ());
37503730 pair_types.resize (paired_alignments.size ());
3751- better_cluster_count_by_pairs.resize (paired_alignments.size ());
37523731 }
37533732 else {
37543733 if (show_work) {
0 commit comments