diff --git a/src/bcf/header.rs b/src/bcf/header.rs index b4eb1f043..94b9f3f57 100644 --- a/src/bcf/header.rs +++ b/src/bcf/header.rs @@ -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="#; + /// 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 { let c_str = ffi::CString::new(name).unwrap(); unsafe { diff --git a/src/bcf/record.rs b/src/bcf/record.rs index 85e1adf6e..1625362e6 100644 --- a/src/bcf/record.rs +++ b/src/bcf/record.rs @@ -127,8 +127,22 @@ impl<'a, T: 'a + fmt::Debug + fmt::Display, B: Borrow + '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, @@ -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 = alleles .iter()