Skip to content

Commit

Permalink
add examples to docs for name2rid and set_alleles (rust-bio#300)
Browse files Browse the repository at this point in the history
  • Loading branch information
mbhall88 authored Apr 29, 2021
1 parent 4ee5dd1 commit bf9e2de
Show file tree
Hide file tree
Showing 2 changed files with 55 additions and 3 deletions.
19 changes: 19 additions & 0 deletions src/bcf/header.rs
Original file line number Diff line number Diff line change
Expand Up @@ -313,6 +313,25 @@ impl HeaderView {
}
}

/// Retrieve the (internal) chromosome identifier
/// # Examples
/// ```rust
/// use rust_htslib::bcf::header::Header;
/// use rust_htslib::bcf::{Format, Writer};
///
/// let mut header = Header::new();
/// let contig_field = br#"##contig=<ID=foo,length=10>"#;
/// header.push_record(contig_field);
/// let mut vcf = Writer::from_stdout(&header, true, Format::VCF).unwrap();
/// let header_view = vcf.header();
/// let rid = header_view.name2rid(b"foo").unwrap();
/// assert_eq!(rid, 0);
/// // try and retrieve a contig not in the header
/// let result = header_view.name2rid(b"bar");
/// assert!(result.is_err())
/// ```
/// # Errors
/// If `name` does not match a chromosome currently in the VCF header, returns [`Error::BcfUnknownContig`]
pub fn name2rid(&self, name: &[u8]) -> Result<u32> {
let c_str = ffi::CString::new(name).unwrap();
unsafe {
Expand Down
39 changes: 36 additions & 3 deletions src/bcf/record.rs
Original file line number Diff line number Diff line change
Expand Up @@ -127,8 +127,22 @@ impl<'a, T: 'a + fmt::Debug + fmt::Display, B: Borrow<Buffer> + 'a> fmt::Display
}
}

/// A BCF record.
/// New records can be created by the `empty_record` methods of `bcf::Reader` and `bcf::Writer`.
/// A VCF/BCF record.
/// New records can be created by the `empty_record` methods of [`bcf::Reader`](crate::bcf::Reader)
/// and [`bcf::Writer`](crate::bcf::Writer).
/// # Example
/// ```rust
/// use rust_htslib::bcf::{Format, Writer};
/// use rust_htslib::bcf::header::Header;
///
/// // Create minimal VCF header with a single sample
/// let mut header = Header::new();
/// header.push_sample("sample".as_bytes());
///
/// // Write uncompressed VCF to stdout with above header and get an empty record
/// let mut vcf = Writer::from_stdout(&header, true, Format::VCF).unwrap();
/// let mut record = vcf.empty_record();
/// ```
#[derive(Debug)]
pub struct Record {
pub inner: *mut htslib::bcf1_t,
Expand Down Expand Up @@ -378,7 +392,26 @@ impl Record {
.collect()
}

/// Set alleles.
/// Set alleles. The first allele is the reference allele.
///
/// # Example
/// ```rust
/// # use rust_htslib::bcf::{Format, Writer};
/// # use rust_htslib::bcf::header::Header;
/// #
/// # // Create minimal VCF header with a single sample
/// # let mut header = Header::new();
/// # header.push_sample("sample".as_bytes());
/// #
/// # // Write uncompressed VCF to stdout with above header and get an empty record
/// # let mut vcf = Writer::from_stdout(&header, true, Format::VCF).unwrap();
/// # let mut record = vcf.empty_record();
/// assert_eq!(record.allele_count(), 0);
///
/// let alleles: &[&[u8]] = &[b"A", b"TG"];
/// record.set_alleles(alleles).unwrap();
/// assert_eq!(record.allele_count(), 2)
/// ```
pub fn set_alleles(&mut self, alleles: &[&[u8]]) -> Result<()> {
let cstrings: Vec<ffi::CString> = alleles
.iter()
Expand Down

0 comments on commit bf9e2de

Please sign in to comment.