|
|
@@ -13,6 +13,17 @@ fn main() {
|
|
|
.replace("\n", "")
|
|
|
.replace(" ", "");
|
|
|
|
|
|
+ let sequence = "
|
|
|
+ CATTATATTCTTTTTTACAGTGCCACGGGTATATTTAGTAAACATTGAAATGTGTGGTTG
|
|
|
+ GAATTGTTCATCAGTGATATTTGTTTATGTTCTGAAGCTATTGCTGTA
|
|
|
+ GACCAACATGGAGTAAACAGAAAATAATTGGGACACAGGGGCAATAAAACTCATATATTA
|
|
|
+ TTTTCCCATGACACTCAGAGGACTCTTCCCAATGCTGAGACATACATAGCGGGTCACGGC
|
|
|
+ TGGGTGTTTTTGGACAAAGGCTTAATGCACTCCAACCTAGGTCCAAGAACTACAACGGCT
|
|
|
+ GCTGCCGGCCGGGCCCGCTCCCTGGCCTACCTTGCACAAGAAGTCGACATCGTAGAAGGC
|
|
|
+ GGACAGAAGTGTGGTCGACATGTTTCTGGAA"
|
|
|
+ .replace("\n", "")
|
|
|
+ .replace(" ", "");
|
|
|
+
|
|
|
println!("sequence len: {}", sequence.len());
|
|
|
println!(">sequence\n{}", sequence);
|
|
|
|
|
|
@@ -37,8 +48,10 @@ fn main() {
|
|
|
println!("{:?}", e.cigar.expand());
|
|
|
}
|
|
|
|
|
|
- println!("-------------------------------")
|
|
|
+ println!("-------------------------------");
|
|
|
+
|
|
|
});
|
|
|
+ make_aln(&sequence, &aln);
|
|
|
}
|
|
|
|
|
|
pub fn bwa_aln(fasta_ref_path: &str, sequence: &str) -> Vec<Sam> {
|
|
|
@@ -174,16 +187,61 @@ fn switch_base(c:char) -> char {
|
|
|
}
|
|
|
|
|
|
fn make_aln(seq_query: &str, alignements: &Vec<Sam>) {
|
|
|
- let consumes_query = vec!["M", "I", "S", "=", "X", "H"];
|
|
|
- let consumes_ref = vec!["M", "D", "N", "=", "X"];
|
|
|
+ let n_line = 80;
|
|
|
+ let consumes_query = vec![b'M', b'I', b'S', b'=', b'X', b'H'];
|
|
|
+ // let consumes_ref = vec![b'M', b'D', b'N', b'=', b'X'];
|
|
|
+ let consumes_ref = vec![b'M', b'D', b'I', b'N', b'=', b'X'];
|
|
|
+ let no_bar_with = vec![b'I', b'H', b'S'];
|
|
|
|
|
|
let mut stack: Vec<Vec<u8>> = Vec::new();
|
|
|
stack.push(seq_query.as_bytes().to_vec());
|
|
|
|
|
|
for alignement in alignements {
|
|
|
- // for pos in 0..seq_query.len() {
|
|
|
-
|
|
|
- // }
|
|
|
+ let mut ref_seq = alignement.seq.clone().as_bytes().to_vec();
|
|
|
+ let mut cigar = alignement.cigar.expand().as_bytes().to_vec();
|
|
|
+ if alignement.is_reverse() {
|
|
|
+ ref_seq = revcomp(&alignement.seq).as_bytes().to_vec();
|
|
|
+ cigar = cigar.into_iter().rev().collect();
|
|
|
+ }
|
|
|
+ let mut matches: Vec<u8> = Vec::new();
|
|
|
+ let mut ref_matches: Vec<u8> = Vec::new();
|
|
|
+ let mut ref_seq = ref_seq.iter();
|
|
|
+ for pos in 0..cigar.len() {
|
|
|
+ let op = *cigar.get(pos).unwrap();
|
|
|
+ if op == b'M' {
|
|
|
+ matches.push(b'|');
|
|
|
+ } else {
|
|
|
+ if no_bar_with.contains(&op) { matches.push(b' '); }
|
|
|
+ }
|
|
|
+ if !consumes_ref.contains(&op) {
|
|
|
+ ref_matches.push(b' ');
|
|
|
+ } else {
|
|
|
+ ref_matches.push(*ref_seq.next().unwrap());
|
|
|
+ }
|
|
|
+ }
|
|
|
+ stack.push(matches);
|
|
|
+ stack.push(ref_matches);
|
|
|
+ }
|
|
|
|
|
|
+ let max_len = stack.iter().map(|e| e.len()).max().unwrap();
|
|
|
+
|
|
|
+ let mut lines: Vec<Vec<String>> = Vec::new();
|
|
|
+ let mut acc_lines: Vec<String> = vec![String::new(); stack.len()];
|
|
|
+ for p in 1..(max_len + 1) {
|
|
|
+ let pos = p - 1;
|
|
|
+ for (i, l) in stack.iter().enumerate() {
|
|
|
+ if let Some(al) = l.get(pos) {
|
|
|
+ acc_lines[i].push(*al as char);
|
|
|
+ } else {
|
|
|
+ acc_lines[i].push(' ');
|
|
|
+ }
|
|
|
+ }
|
|
|
+
|
|
|
+ if p % n_line == 0 {
|
|
|
+ lines.push(acc_lines.clone());
|
|
|
+ acc_lines.iter().for_each(|e| println!("{}", e));
|
|
|
+ acc_lines = vec![String::new(); stack.len()];
|
|
|
+ }
|
|
|
}
|
|
|
+
|
|
|
}
|