Skip to content

Commit

Permalink
[INFRA] Create associated cpp files for header files and move non-tem…
Browse files Browse the repository at this point in the history 
…plate functions into them.

Signed-off-by: Lydia Buntrock <[email protected]>
  • Loading branch information
@Irallia
Irallia committed Feb 22, 2021
1 parent 320d233 commit ccb013d
Show file tree
Hide file tree
Showing 9 changed files with 299 additions and 242 deletions.
134 changes: 10 additions & 124 deletions include/structures/aligned_segment.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,8 @@

#include <seqan3/alphabet/cigar/cigar.hpp>

#include "structures/breakend.hpp" // for strand

/*! \brief Read segment aligned to the reference genome (part of chimeric/split-aligned read). Contains information
* parsed from the SA tag of an alignment in the SAM/BAM file.
*
Expand All @@ -19,128 +21,19 @@ struct AlignedSegment
int32_t mapq;
std::vector<seqan3::cigar> cig;

int32_t get_reference_start() const
{
return pos;
}
int32_t get_reference_start() const;

int32_t get_reference_end() const
{
int32_t current_pos = pos;
for (auto [element_length, element_operation] : cig)
{
switch(element_operation.to_char())
{
case 'M':
current_pos += element_length;
break;
case 'D':
current_pos += element_length;
break;
case 'N':
current_pos += element_length;
break;
case 'X':
current_pos += element_length;
break;
case '=':
current_pos += element_length;
break;
}
}
return current_pos;
}
int32_t get_reference_end() const;

int32_t get_left_soft_clip() const
{
int32_t left_soft_clip = 0;
for (auto [element_length, element_operation] : cig)
{
if(element_operation.to_char() == 'S')
{
left_soft_clip += element_length;
}
else if(element_operation.to_char() == 'M' ||
element_operation.to_char() == '=' ||
element_operation.to_char() == 'X' ||
element_operation.to_char() == 'I')
{
break;
}
}
return left_soft_clip;
}
int32_t get_left_soft_clip() const;

int32_t get_right_soft_clip() const
{
int32_t right_soft_clip = 0;
for (auto [element_length, element_operation] : std::views::reverse(cig))
{
if(element_operation.to_char() == 'S')
{
right_soft_clip += element_length;
}
else if(element_operation.to_char() == 'M' ||
element_operation.to_char() == '=' ||
element_operation.to_char() == 'X' ||
element_operation.to_char() == 'I')
{
break;
}
}
return right_soft_clip;
}
int32_t get_right_soft_clip() const;

int32_t get_query_start() const
{
if (orientation == strand::forward)
{
return get_left_soft_clip();
}
else
{
return get_right_soft_clip();
}
}
int32_t get_query_start() const;

int32_t get_query_length() const
{
int32_t current_length = 0;
for (auto [element_length, element_operation] : cig)
{
switch(element_operation.to_char())
{
case 'M':
current_length += element_length;
break;
case 'S':
current_length += element_length;
break;
case 'I':
current_length += element_length;
break;
case 'X':
current_length += element_length;
break;
case '=':
current_length += element_length;
break;
}
}
return current_length;
}
int32_t get_query_length() const;

int32_t get_query_end() const
{
if (orientation == strand::forward)
{
return get_query_length() - get_right_soft_clip();
}
else
{
return get_query_length() - get_left_soft_clip();
}
}
int32_t get_query_end() const;
};

template <typename stream_t>
Expand All @@ -154,11 +47,4 @@ inline stream_t operator<<(stream_t && stream, AlignedSegment const & a)
return stream;
}

inline bool operator<(const AlignedSegment & lhs, const AlignedSegment & rhs)
{
return lhs.get_query_start() != rhs.get_query_start()
? lhs.get_query_start() < rhs.get_query_start()
: lhs.get_query_end() != rhs.get_query_end()
? lhs.get_query_end() < rhs.get_query_end()
: lhs.mapq < rhs.mapq;
}
bool operator<(const AlignedSegment & lhs, const AlignedSegment & rhs);
21 changes: 4 additions & 17 deletions include/structures/breakend.hpp
View file
Original file line number Diff line number Diff line change
@@ -1,5 +1,7 @@
#pragma once

#include <string>

enum struct sequence_type : uint8_t
{
reference,
Expand Down Expand Up @@ -42,21 +44,6 @@ inline stream_t operator<<(stream_t && stream, Breakend const & b)
return stream;
}

inline bool operator<(const Breakend & lhs, const Breakend & rhs)
{
return lhs.seq_type != rhs.seq_type
? lhs.seq_type < rhs.seq_type
: lhs.seq_name != rhs.seq_name
? lhs.seq_name < rhs.seq_name
: lhs.orientation != rhs.orientation
? lhs.orientation < rhs.orientation
: lhs.position < rhs.position;
}
bool operator<(const Breakend & lhs, const Breakend & rhs);

inline bool operator==(const Breakend & lhs, const Breakend & rhs)
{
return (lhs.seq_name == rhs.seq_name) &&
(lhs.position == rhs.position) &&
(lhs.orientation == rhs.orientation) &&
(lhs.seq_type == rhs.seq_type);
}
bool operator==(const Breakend & lhs, const Breakend & rhs);
79 changes: 6 additions & 73 deletions include/structures/cluster.hpp
View file
Original file line number Diff line number Diff line change
@@ -1,8 +1,8 @@
#pragma once

#include <cmath>
#include <stdexcept>
#include <utility>
#include <vector>

#include "structures/junction.hpp" // for class Junction

class Cluster
{
Expand All @@ -25,78 +25,11 @@ class Cluster
{

}
size_t get_cluster_size() const
{
return members.size();
}
size_t get_cluster_size() const;

Breakend get_average_mate1() const
{
sequence_type seq_type{};
std::string seq_name{};
uint64_t sum_positions = 0;
strand orientation{};
// Iterate through members of the cluster
for (uint32_t i = 0; i < members.size(); ++i)
{
Breakend mate1 = members[i].get_mate1();
if (i == 0)
{
seq_type = mate1.seq_type;
seq_name = mate1.seq_name;
orientation = mate1.orientation;
}
else
{
// Make sure that all members of the cluster have matching sequence types, names and orientations
if (mate1.seq_type != seq_type ||
mate1.seq_name != seq_name ||
mate1.orientation != orientation)
{
throw std::runtime_error("Junctions with incompatible breakends were clustered together (different seq_type, seq_name or orientation).");
}
}
// Add up breakend positions aross all members
sum_positions += mate1.position;
}
int32_t average_position = std::round(static_cast<double>(sum_positions) / members.size());
Breakend average_breakend{seq_name, average_position, orientation, seq_type};
return average_breakend;
}
Breakend get_average_mate1() const;

Breakend get_average_mate2() const
{
sequence_type seq_type{};
std::string seq_name{};
uint64_t sum_positions = 0;
strand orientation{};
// Iterate through members of the cluster
for (uint32_t i = 0; i < members.size(); ++i)
{
Breakend mate2 = members[i].get_mate2();
if (i == 0)
{
seq_type = mate2.seq_type;
seq_name = mate2.seq_name;
orientation = mate2.orientation;
}
else
{
// Make sure that all members of the cluster have matching sequence types, names and orientations
if (mate2.seq_type != seq_type ||
mate2.seq_name != seq_name ||
mate2.orientation != orientation)
{
throw std::runtime_error("Junctions with incompatible breakends were clustered together (different seq_type, seq_name or orientation).");
}
}
// Add up breakend positions aross all members
sum_positions += mate2.position;
}
int32_t average_position = std::round(static_cast<double>(sum_positions) / members.size());
Breakend average_breakend{seq_name, average_position, orientation, seq_type};
return average_breakend;
}
Breakend get_average_mate2() const;
};

template <typename stream_t>
Expand Down
35 changes: 8 additions & 27 deletions include/structures/junction.hpp
View file
Original file line number Diff line number Diff line change
@@ -1,8 +1,6 @@
#pragma once

#include <utility>

#include "breakend.hpp"
#include "structures/breakend.hpp"

class Junction
{
Expand Down Expand Up @@ -36,19 +34,12 @@ class Junction
this->mate2.flip_orientation();
}
}
Breakend get_mate1() const
{
return mate1;
}
Breakend get_mate2() const
{
return mate2;
}

std::string get_read_name() const
{
return read_name;
}
Breakend get_mate1() const;

Breakend get_mate2() const;

std::string get_read_name() const;
};

template <typename stream_t>
Expand All @@ -60,22 +51,12 @@ inline stream_t operator<<(stream_t && stream, Junction const & junc)
return stream;
}

inline bool operator<(const Junction & lhs, const Junction & rhs)
{
return lhs.get_mate1() < rhs.get_mate1()
? true
: rhs.get_mate1() < lhs.get_mate1()
? false
: lhs.get_mate2() < rhs.get_mate2();
}
bool operator<(const Junction & lhs, const Junction & rhs);

/*! \brief A junction is equal to another, if their mates are equal to each other. The read_name and supporting_reads
* are allowed to be unequal, because more than one read could support the same junction.
*
* \param lhs left side junction
* \param rhs right side junction
*/
inline bool operator==(const Junction & lhs, const Junction & rhs)
{
return (lhs.get_mate1() == rhs.get_mate1()) && (lhs.get_mate2() == rhs.get_mate2());
}
bool operator==(const Junction & lhs, const Junction & rhs);
6 changes: 5 additions & 1 deletion src/CMakeLists.txt
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,11 @@ cmake_minimum_required (VERSION 3.8)

# An object library (without main) to be used in multiple targets.
add_library ("${PROJECT_NAME}_lib" STATIC detect_breakends/junction_detection.cpp
find_deletions/deletion_finding_and_printing.cpp)
find_deletions/deletion_finding_and_printing.cpp
structures/aligned_segment.cpp
structures/breakend.cpp
structures/cluster.cpp
structures/junction.cpp)

target_link_libraries ("${PROJECT_NAME}_lib" PUBLIC seqan3::seqan3)
target_include_directories ("${PROJECT_NAME}_lib" PUBLIC ../include)
Expand Down
Loading

0 comments on commit ccb013d

Please sign in to comment.