Skip to content

feat: Add bam::Record::set_cigar #477

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
wants to merge 3 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
61 changes: 61 additions & 0 deletions src/bam/mod.rs
Original file line number Diff line number Diff line change
Expand Up @@ -1933,6 +1933,67 @@ CCCCCCCCCCCCCCCCCCC"[..],
assert_eq!(rec.qname(), b"r0");
}

#[test]
fn test_set_cigar() {
let (names, _, seqs, quals, cigars) = gold();

assert!(names[0] != names[1]);

for i in 0..names.len() {
let mut rec = record::Record::new();
rec.set(names[i], Some(&cigars[i]), seqs[i], quals[i]);
rec.push_aux(b"NM", Aux::I32(15)).unwrap();

assert_eq!(rec.qname(), names[i]);
assert_eq!(*rec.cigar(), cigars[i]);
assert_eq!(rec.seq().as_bytes(), seqs[i]);
assert_eq!(rec.qual(), quals[i]);
assert_eq!(rec.aux(b"NM").unwrap(), Aux::I32(15));

// boring cigar
let new_cigar = CigarString(vec![Cigar::Match(rec.seq_len() as u32)]);
assert_ne!(*rec.cigar(), new_cigar);
rec.set_cigar(Some(&new_cigar));
assert_eq!(*rec.cigar(), new_cigar);

assert_eq!(rec.qname(), names[i]);
assert_eq!(rec.seq().as_bytes(), seqs[i]);
assert_eq!(rec.qual(), quals[i]);
assert_eq!(rec.aux(b"NM").unwrap(), Aux::I32(15));

// bizarre cigar
let new_cigar = (0..rec.seq_len())
.map(|i| {
if i % 2 == 0 {
Cigar::Match(1)
} else {
Cigar::Ins(1)
}
})
.collect::<Vec<_>>();
let new_cigar = CigarString(new_cigar);
assert_ne!(*rec.cigar(), new_cigar);
rec.set_cigar(Some(&new_cigar));
assert_eq!(*rec.cigar(), new_cigar);

assert_eq!(rec.qname(), names[i]);
assert_eq!(rec.seq().as_bytes(), seqs[i]);
assert_eq!(rec.qual(), quals[i]);
assert_eq!(rec.aux(b"NM").unwrap(), Aux::I32(15));

// empty cigar
let new_cigar = CigarString(Vec::new());
assert_ne!(*rec.cigar(), new_cigar);
rec.set_cigar(None);
assert_eq!(*rec.cigar(), new_cigar);

assert_eq!(rec.qname(), names[i]);
assert_eq!(rec.seq().as_bytes(), seqs[i]);
assert_eq!(rec.qual(), quals[i]);
assert_eq!(rec.aux(b"NM").unwrap(), Aux::I32(15));
}
}

#[test]
fn test_remove_aux() {
let mut bam = Reader::from_path(Path::new("test/test.bam")).expect("Error opening file.");
Expand Down
56 changes: 56 additions & 0 deletions src/bam/record.rs
Original file line number Diff line number Diff line change
Expand Up @@ -473,6 +473,62 @@ impl Record {
self.inner_mut().core.l_extranul = extranul as u8;
}

/// Replace current cigar with a new one.
pub fn set_cigar(&mut self, new_cigar: Option<&CigarString>) {
self.cigar = None;

let qname_data_len = self.qname_capacity();
let old_cigar_data_len = self.cigar_len() * 4;

// Length of data after cigar
let other_data_len = self.inner_mut().l_data - (qname_data_len + old_cigar_data_len) as i32;

let new_cigar_len = match new_cigar {
Some(x) => x.len(),
None => 0,
};
let new_cigar_data_len = new_cigar_len * 4;

if new_cigar_data_len < old_cigar_data_len {
self.inner_mut().l_data -= (old_cigar_data_len - new_cigar_data_len) as i32;
} else if new_cigar_data_len > old_cigar_data_len {
self.inner_mut().l_data += (new_cigar_data_len - old_cigar_data_len) as i32;

// Reallocate if necessary
if (self.inner().m_data as i32) < self.inner().l_data {
// Verbosity due to lexical borrowing
let l_data = self.inner().l_data;
self.realloc_var_data(l_data as usize);
}
}

if new_cigar_data_len != old_cigar_data_len {
// Move other data to new location
unsafe {
::libc::memmove(
self.inner.data.add(qname_data_len + new_cigar_data_len) as *mut ::libc::c_void,
self.inner.data.add(qname_data_len + old_cigar_data_len) as *mut ::libc::c_void,
other_data_len as usize,
);
}
}

// Copy cigar data
if let Some(cigar_string) = new_cigar {
let cigar_data = unsafe {
#[allow(clippy::cast_ptr_alignment)]
slice::from_raw_parts_mut(
self.inner.data.add(qname_data_len) as *mut u32,
cigar_string.len(),
)
};
for (i, c) in cigar_string.iter().enumerate() {
cigar_data[i] = c.encode();
}
}
self.inner_mut().core.n_cigar = new_cigar_len as u32;
}

fn realloc_var_data(&mut self, new_len: usize) {
// pad request
let new_len = new_len as u32;
Expand Down
2 changes: 1 addition & 1 deletion src/lib.rs
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@
//! ```
//!
//! We can reproduce that with Rust-Htslib. Reading BAM files and printing the header
//! to the the screen is as easy as
//! to the screen is as easy as
//!
//! ```
//! use rust_htslib::{bam, bam::Read};
Expand Down
Loading